各个单细胞亚群独立在两个分组做差异分析

前些天我们的单细胞学徒培养有小伙伴分享了文章；在 JCI Insight 2022  https://doi.org/10.1172/jci.insight.152616  ，里面对第一次降维聚类分群后的各个单细胞亚群独立在两个分组做差异分析 ，如下所示：

各个单细胞亚群独立在两个分组做差异分析

可以看到，每个单细胞亚群都有自己的差异分析火山图，会议上有人提问这个分析如何做。其实主要是大家可能是初次接触生物信息学就是单细胞数据处理，所以基础知识有点欠缺。它就是普通的表达量矩阵分析而已，我七八年前就写过系列笔记，公众号推文在：

我们这里以大家熟知的pbmc3k数据集为例。大家先安装这个数据集对应的包，并且对它进行降维聚类分群，参考前面的例子：人人都能学会的单细胞聚类分群注释

# 0.安装R包 ---- 
# InstallData("pbmc3k") 

library(SeuratData) #加载seurat数据集  
getOption('timeout')
options(timeout=10000)
#InstallData("pbmc3k")  
data("pbmc3k")  
sce <- pbmc3k.final   
library(Seurat)
table(Idents(sce))
DimPlot(sce,label = T)

这个时候，因为它pbmc3k数据集并没有分组，所以我们没办法做差异分析。不过我们可以简单的模拟一个分组。如下所示：

sce$celltype = Idents(sce)
sce$group = sample(1:2,ncol(sce),replace = T)
table(sce$celltype,sce$group )


# 如下所示，两个随机赋予的分组，每个分组里面的都是有这些不同单细胞亚群
  Naive CD4 T  349 348
  Memory CD4 T 231 252
  CD14+ Mono   224 256
  B            186 158
  CD8 T        135 136
  FCGR3A+ Mono  81  81
  NK            66  89
  DC            12  20
  Platelet       9   5

下面对它们进行批量差异分析：

Idents(sce) = paste0('c',sce$group )
table(Idents(sce))
degs = lapply(unique(sce$celltype), function(x){
  FindMarkers(sce[,sce$celltype==x],ident.1 = 'c1',
              ident.2 = 'c2')
})
x=degs[[1]]
do.call(rbind,lapply(degs, function(x){
  table(x$avg_log2FC > 0 )
}))

值得注意的是这个FindMarkers函数并不是最好的单细胞转录组表达量矩阵的差异分析方法，我这里仅仅是举例哦！

可以看到，如果是以为 avg_log2FC 标准，这个时候很容易得到假阳性的差异基因：

FindMarkers函数的结果

其实我们更应该关心的是 "pct.1"  和  "pct.2" 的差异，这个也是各个单细胞亚群特异性高表达量基因的金标准，不过跟我们这个时候的差异分析不太一样的需求哦，需要自己多思考。

另外，这个时候因为我们是随机赋值的两个分组，所以差异分析理论上就没有意义，所以"p_val_adj"  在这里基本上没有统计学显著性。

往期回顾

---

---

---

如果你对单细胞转录组研究感兴趣，但又不知道如何入门，也许你可以关注一下下面的课程

看完记得顺手点个“在看”哦！

生物 | 单细胞 | 转录组丨资料每天都精彩

长按扫码可关注