scCODE_scSTAR

专注于单细胞差异分析的算法

Hi，大家好，我是晨曦
今天这期推文的内容是晨曦在前段时间听邹老师讲座的时候发现邹老师团队开发出了很多scRNA-seq的处理算法，重点是代码编写的十分简洁，简直就是tidy党的福音，所以赶快mark了一下并且写下推文分享给大家

下面将以模块的方式介绍邹老师团队开发的一些算法

scCODE：获取差异基因的全方位策略

现在已经有几十种方法可以从单细胞数据中判定差异基因，但是不同的差异基因算法得到的差异基因是不同的，不同的单细胞数据需要使用不同的方法，scCODE可以针对不同的数据给出最优的个性化解决方案

简单来说就是scCODE可以告诉你哪些差异基因可能更有意义（排除不同算法的干扰）

那么我们直接开始进行代码实战

library(scCODE)

#下载该R包需要提前安装一些依赖包，建议上github上查看完整安装教程https://github.com/XZouProjects/scCODE

data<-data_sccode#demo数据#观测是Gene ID#变量是barcode ID#数据要求就是TPM、FPKM、counts，但是不能是log后的数据group<-group_sccode#分组信息（两两比较)#1为control；2为treatresults<-scCODE(data,group,light = TRUE,top_ranked=5)#light参数为只选择一些优质的差异算法，加速整个过程#top_ranked打分为只选择前5个打分最优质的算法?scCODEhead(results[["DE_results"]])#Gene_name就是Gene symbol#Detected_times：我们选择了top5的方法里有几个方法检测到了这个Gene#注意：Detected_times和logFC是没有关系的，logFC是多少就是多少不会改变#注意：Seurat对象中表达矩阵是稀疏矩阵，需要通过GetAssayData函数获取原始矩阵library(Matrix)data("pbmc_small")ct=GetAssayData(object = pbmc_small, assay = "RNA", slot = "counts") ct <- as.data.frame(ct)

晨曦解读

这个代码写的十分的简洁，可以让我们发现scRNA-seq中差异基因中更具有差异的那些基因（在排除算法的干扰下），更多内容可以参考这个地址：：https://github.com/XZouProjects/scCODE

scSTAR：首次提出单个细胞动态建模，深挖更多差异

目前来说基本上全部的单细胞数据挖掘方法基本上都是只能挖掘信号大于噪声的那部分数据，没有任何一种方法可以挖掘到信号小于噪声的数据

传统的单细胞静态建模，单细胞可以看成是一个反映系统动态过程的一个阶段，每一个细胞都是一个动态过程中静态的一个点，因此，传统单细胞数据处理思路都是尽可能真实地反应每一个细胞的静态状态，进而用一系列静态状态细胞重构动态特征

但是，由于噪声的存在，通过静态特征并不能100%还原信息，尤其是比较微弱的生物学信息

为了解决这个问题，提出了“单个细胞动态建模”的概念，把每一个细胞用动态特征描述而不是静态特征描述，这样就可以放大细胞间微小的差别，同时也能从根本上解决细胞间随机异质性的干扰

晨曦解读

总结来说就一句，这个算法可以挖掘出单细胞更多的差异

这个R包嵌套了一个这个团队自己开发的算法——OGFSC，目前，去除噪声主要依赖人工硬阈值，这样导致数据会受到“损伤”。这个算法利用统计优化方法，代替人工判断，构建基因滤除阈值，以达到尽可能滤除随机噪声污染基因，并保留数据信息的目的

那么我们直接开始进行代码实战

####0.准备工作####library(R.matlab) library(tsne) library(scSTAR) library(pls)library(RColorBrewer)####1.准备输入数据####extdir<-system.file("extdata",package="scSTAR")data=readMat(paste0(extdir,'/demo_data.mat'))data1 <- as.matrix((data$data1))#表达矩阵 data2 <- as.matrix((data$data2))#表达矩阵 geneList <- as.matrix(unlist(data$geneList))#Gene symbol矩阵####2.基因筛选#####这一步通过去除噪声污染的基因来减少随机噪声的干扰data <- cbind(data1, data2)#查看数据data[1:10,1:10]#可以很清楚的看到就是一个大型的表达矩阵idx_OGFSC = filter_gene(data, nBins = 20, plot_option = 1)$OGFSC_idx

#OGFSC是利用统计优化方法，代替人工判断，构建基因滤除阈值，以达到尽可能滤除随机噪声污染基因，并保留数据信息的目的

#nBins参数：广义线性模式构筑格子的数量，简单理解就是一个自动QC阈值的设定，越大剩下的基因越少

Ctr_filtered = data1[idx_OGFSC,]#对照组 Case_filtered = data2[idx_OGFSC,]#实验组 ####3.identify anchor cells using KNN####

anchorCells = findAnchors(Ctr_filtered, Case_filtered, 3)#因为这个算法主要是发现处理前后变化的细微marker，所有必须得有两种状态的表达矩阵

#这个3代表最近的细胞数量，默认即可####4.single-cell kinetics estimation####out = scSTAR(Ctr_filtered, Case_filtered, anchorCells, PLScomp = 4)

#评价单个细胞动力学的函数。当 PLScomp = 0时，可以自动估计 PLScomp 的个数。但是它也可以手动设置。建议使用者从3-5开始查找数值

Case_kinetics = out$Case_kinetics idx_DE_up = out$idx_DE_up#表达升高的gene ####5.data illustration (optional)#####可以把表达升高的Gene进行聚类K = 3 out2 = autoClustering(Case_kinetics, K, geneList[idx_OGFSC], idx_DE_up) clusterIdx = out2$clusterIdx#上调基因划分为三类markerGenes = out2$markerGenes[1]#第一类的上调基因

晨曦解读

总的来说，这套代码中需要修改的参数其实不多，而且大多数作者都建议默认或者小范围修改调试，所以我们可以更换自己的数据集然后走这个流程，这个包其实目的就是为了深度挖掘某个表达矩阵前后状态更深层次的marker基因，可以挖掘到上调/下调更具意义的基因，并且可以把上调/下调的基因再度进行聚类划分亚群，总的来说是一个挺方便的R包

小结

这两个R包都是针对scRNA-seq差异基因的再度挖掘，一个是为了解决哪个marker更有意义的问题，一个则是为了找到更细微的marker，所以如果我们有了scRNA-seq后也可以尝试这两种算法然后得到更加精确或者说是不容易注意到的new marker

那么本期推文到这里就结束了，这期又是“浪学”的一期，下一期我们再一次学点有意义的东西~

我是晨曦，我们下期再见

参考教程

1.https://github.com/XZouProjects/scCODE
2.https://github.com/XZouProjects/scSTAR