2021-05-14 scRNA 基础分析：细胞类型鉴定

细胞完成聚类和降维可视化之后，细胞群的异质性已经体现出来了，那么各个cluster的细胞如何鉴定细胞类型呢？通用的方法有两种：一种是通过细胞类型特异性表达的marker基因识别，就好像流式细胞仪用特定的抗体筛选细胞一样；第二种是建立已知细胞类型的转录谱数据库，将未知细胞类型的转录谱数据与之比较相似性，就知道它最有可能是哪种细胞了。第一种方法需要人工收集marker基因比对各个cluster的显著高表达基因综合分析，第二种方法可以使用SingleR包自动识别细胞类型。建议两种方法结合起来进行细胞鉴定。

#单细胞转录组基础分析四：细胞类型鉴定
library(Seurat)
library(tidyverse)
library(patchwork)
rm(list=ls())
scRNA1 <- readRDS("scRNA1.rds")

marker基因鉴定

通过细胞类型特异性表达的marker基因识别（常用），就好像流式细胞仪用特定的抗体筛选细胞一样。
最关键的是要找到各个cluster的显著高表达的基因（cluster间基因表达差异分析），主要是利用Seurat包的FindAllMarkers函数，通过设置其test.use参数，从而选择不同的method

####细胞类型的注释一般有三种方法.1、利用marker基因查找网站进行注释  2、使用singler进行注释 3、根据已有的生物学知识或者文献，按照dotplot来注释。
##现在使用方法一寻找marker基因使用网站注释 找marker基因有以下方法三选一，建议第一种
#默认wilcox方法
colnames(all.markers)[3] <- 'avg_logFC'
diff.wilcox = FindAllMarkers(scRNA1)
#大概1-2min
all.markers = diff.wilcox %>% select(gene, everything()) %>% subset(p_val<0.05)
#先把基因这一列放在第一列，然后选取p值小于0.05的行(结果行数不变，说明都挺好的)
top10 = all.markers %>% group_by(cluster) %>% top_n(n = 10, wt = avg_logFC)
#0-9共10个cluster，每个cluster选取top10高变基因
###将marker基因保存一下
write.csv(all.markers, "diff_genes_wilcox.csv", row.names = F)
###把top10marker基因保存一下
write.csv(top10, "top10_diff_genes_wilcox.csv", row.names = F)

##top10基因绘制热图
top10_genes = CaseMatch(search = as.vector(top10$gene), match = rownames(scRNA1)) 
plot1 = DoHeatmap(scRNA1, features = top10_genes, group.by = "seurat_clusters", group.bar = T, size = 4)
plot1

heatmap.png

选择基因作图展示

#挑选部分基因
select_genes <- c('LYZ','CD79A','CD8A','CD8B','GZMB','FCGR3A')
#vlnplot展示
p1 <- VlnPlot(scRNA1, features = select_genes, pt.size=0, ncol=2)
ggsave("cell_identify/selectgenes_VlnPlot.png", p1, width=6 ,height=8)
#featureplot展示
p2 <- FeaturePlot(scRNA1, features = select_genes, reduction = "tsne", label=T, ncol=2)
ggsave("cell_identify/selectgenes_FeaturePlot.png", p2, width=8 ,height=12)
p3=p1|p2

SingleR鉴定细胞类型

BiocManager::install("SingleR")
library(SingleR)
load("H:/SingleR_ref/ref_Human_all.RData")
refdata <- ref_Human_all
testdata <- GetAssayData(scRNA1, slot="data")
###把scRNA数据中的seurat_clusters提取出来，注意这里是因子类型的
clusters <- scRNA1@meta.data$seurat_clusters
###开始用singler分析
cellpred <- SingleR(test = testdata, ref = refdata, labels = refdata$label.main, 
                    method = "cluster", clusters = clusters, 
                    assay.type.test = "logcounts", assay.type.ref = "logcounts")
###制作细胞类型的注释文件
celltype = data.frame(ClusterID=rownames(cellpred), celltype=cellpred$labels, stringsAsFactors = FALSE)
###保存一下
write.csv(celltype,"celltype_singleR.csv",row.names = FALSE)


##把singler的注释写到metadata中 有两种方法
###方法一
scRNA1@meta.data$celltype = "NA"
for(i in 1:nrow(celltype)){
  scRNA1@meta.data[which(scRNA1@meta.data$seurat_clusters == celltype$ClusterID[i]),'celltype'] <- celltype$celltype[i]}
table(scRNA1@meta.data$celltype)
#B_cell Monocyte  NK_cell  T_cells 
#175      319       53      466 
###方法二：
celltype = data.frame(ClusterID=rownames(cellpred), celltype=cellpred$labels, stringsAsFactors = F) 
scRNA1@meta.data$singleR=celltype[match(clusters,celltype$ClusterID),'celltype']

鉴定结果展示

p1 = DimPlot(scRNA1, group.by="celltype", label=T, label.size=5, reduction='tsne')
p2 = DimPlot(scRNA1, group.by="celltype", label=T, label.size=5, reduction='umap')
p3 = plotc <- p1+p2+ plot_layout(guides = 'collect')

有了celltype后 重新选择基因作图展示

select_genes <- c('LYZ','CD79A','GZMB','FCGR3A')
#vlnplot展示
p1 <- VlnPlot(scRNA1, features = select_genes, pt.size=0, group.by="celltype", ncol=2)
#featureplot展示
p2 <- FeaturePlot(scRNA1, features = select_genes, reduction = "tsne", label=T, ncol=2)
p3=p1|p2

参考 小贝学生信 公众号：生信会客厅