任意细胞亚群的差异分析
我们以 seurat 官方教程为例：
rm(list = ls())
library(Seurat)
# devtools::install_github('satijalab/seurat-data')
library(SeuratData)
library(ggplot2)
library(patchwork)
library(dplyr)
load(file = 'basic.sce.pbmc.Rdata')
DimPlot(pbmc, reduction = 'umap',
label = TRUE, pt.size = 0.5) + NoLegend()
sce=pbmc
如果你不知道basic.sce.pbmc.Rdata 这个文件如何得到的，麻烦自己去跑一下 可视化单细胞亚群的标记基因的5个方法，自己 save(pbmc,file = 'basic.sce.pbmc.Rdata') ，我们后面的教程都是依赖于这个 文件哦！
对各个细胞亚群找高表达量的标记基因
代码如下：
if (file.exists('sce.markers.all_10_celltype.Rdata')) {
load('sce.markers.all_10_celltype.Rdata')
}else {
sce.markers <- FindAllMarkers(object = sce, only.pos = TRUE,
min.pct = 0.25,
thresh.use = 0.25)
save(sce.markers,file = 'sce.markers.all_10_celltype.Rdata')
}
# 并且可视化它
head(sce.markers)
table(sce.markers$cluster)
# 首先挑选基因
kp=grepl('Mono',sce.markers$cluster)
table(kp)
cg_sce.markers = sce.markers [ kp ,]
# 然后挑选细胞
kp=grepl('Mono', Idents(sce ) )
table(kp)
sce=sce[,kp]
sce
table( Idents(sce ))
cg_sce.markers=cg_sce.markers[cg_sce.markers$avg_logFC>2,]
dim(cg_sce.markers)
DoHeatmap(subset(sce, downsample = 15),
unique(cg_sce.markers$gene),
slot = 'counts',
size=3 )
如下所示，
对指定的两个细胞亚群找差异
levels(Idents(sce))
markers_df <- FindMarkers(object = sce,
ident.1 = 'FCGR3A+ Mono',
ident.2 = 'CD14+ Mono',
#logfc.threshold = 0,
min.pct = 0.25)
head(markers_df)
cg_markers_df=markers_df[abs(markers_df$avg_logFC) >1,]
dim(cg_markers_df)
DoHeatmap(subset(sce, downsample = 15),
slot = 'counts',
unique(rownames(cg_markers_df)),size=3)
如下所示：
任意划分亚群再找差异
# drop-out
highCells= colnames(subset(x = sce, subset = FCGR3A > 1,
slot = 'counts'))
highORlow=ifelse(colnames(sce) %in% highCells,'high','low')
table(highORlow)
table(Idents(sce))
sce@meta.data$highORlow=highORlow
可以看到两个方法的划分亚群对比情况 ：
> table(highORlow)
highORlow
high low
160 482
> table(Idents(sce))
CD14+ Mono FCGR3A+ Mono
480 162
> table(Idents(sce),highORlow)
highORlow
high low
CD14+ Mono 15 465
FCGR3A+ Mono 145 17
然后再找差异：
markers <- FindMarkers(sce, ident.1 = "high",
group.by = 'highORlow' )
head(x = markers)
cg_markers=markers[abs(markers$avg_logFC) >1,]
dim(cg_markers)
DoHeatmap(subset(sce, downsample = 15),
rownames(cg_markers) ,
size=3 )
如下所示：