单细胞分析入门(3)—— 细胞注释:给每个细胞群贴上“身份标签”
乔粒:基于上篇内容:我们进行了pbmc数据的读入、质控、降维,教程:
单细胞分析入门(1)pbmc数据读入、Seurat象创建与基础质控,
单细胞分析入门(2)-PCA降维+聚类:锁定细胞类群的关键两步
本篇我们来进行细胞注释及绘图。
三、细胞注释
9.寻找差异表达特征(簇生物标志物)
寻找每个细胞群(Cluster)的差异表达基因(Marker Genes)。通过找出这些Marker基因,推断出每个细胞群究竟是哪种免疫细胞(比如T细胞、B细胞或单核细胞)。
①将第2号细胞群单独拎出来,与数据集里的所有其他细胞群(背景)进行基因表达对比。找出2群相对于全部的群,哪些基因显著高表达或低表达。
②将5号群单独拎出来,与第0群和第3群进行对比。当我们发现第5群与第0、3群在UMAP图上离得很近、容易混淆时,用这种精准对比能找出它们之间微小的差异基因,从而区分开这两个相似的亚型。
③自动循环遍历所有的细胞群,分别找出每个群相对于其余群的高表达基因(正向标记)。avg_log2FC代表对数2的倍数变化,数值大于1意味着该基因在这个群里的平均表达量是其他群的2倍以上,过滤掉那些差异微弱的基因,只留下高可信度的强特异性标记物。
④单独针对第0号群,设置 logfc.threshold = 0.25,放宽了基因表达的差异起步门槛。 test.use = “roc”。放弃了默认的Wilcoxon秩和检验,改用ROC曲线分析,输出的核心指标不再是“倍数变化”,而是“分类能力”或“预测能力”。数值范围在0-1之间,越接近1,用该基因的表达量来区分“第0群”和“非第0群”的准确率越高。当只想知道“哪个基因最能代表这个群的独特性”时,用ROC检验。
#===================== 9. 寻找差异表达特征(簇生物标志物)====================================# 单独拎出簇2的所有标记,与数据集里的所有其他细胞群(背景)进行基因表达对比cluster2.markers<-FindMarkers(pbmc,ident.1=2)head(cluster2.markers,n=5)#找到区分第5个集群与第0和第3个集群的所有标记比对cluster5.markers<-FindMarkers(pbmc,ident.1=5,ident.2=c(0,3))head(cluster5.markers,n=5)# 找到每个聚类相对于所有剩余细胞的标记,仅报告正向结果# 1个pbmc.markers<-FindAllMarkers(pbmc,only.pos=TRUE)pbmc.markers%>%group_by(cluster)%>%dplyr::filter(avg_log2FC>1)# 范围为0 - 随机,1 - 完美cluster0.markers<-FindMarkers(pbmc,ident.1=0,logfc.threshold=0.25,test.use="roc",only.pos=TRUE)可视化:
#可视化: 我们包含多种用于可视化标记表达的工具。 (显示表达式概率分布 跨群集)、# 以及(可视化特征) tSNE或PCA图上的表达)是我们最常用的可视化方法。# VlnPlot()FeaturePlot()RidgePlot()CellScatter()DotPlot()VlnPlot(pbmc,features=c("MS4A1","CD79A"))# 绘制原始数据计数VlnPlot(pbmc,features=c("NKG7","PF4"),slot="counts",log=TRUE)FeaturePlot(pbmc,features=c("MS4A1","GNLY","CD3E","CD14","FCER1A","FCGR3A","LYZ","PPBP","CD8A"))# DoHeatmap()生成给定细胞的表达热图 以及功能。在这里,我们绘制的是前10个标记(或全部) 每个簇的标记值均为小于10)。pbmc.markers%>%group_by(cluster)%>%dplyr::filter(avg_log2FC>1)%>%slice_head(n=10)%>%ungroup()->top10 DoHeatmap(pbmc,features=top10$gene)+NoLegend()①小提琴图:针对所有细胞群,画出这两个基因的表达量密度分布。这两个基因是B细胞的特异性标志物。如果显示,在某个特定的群(比如第3群)里,这两个基因的表达两特别高,而在其他群里几乎为零,那就说明这个群大概率就是B细胞,聚类成功。②小提琴图:直接画原始的UMI计数(原始转录本数量),Y轴取对数。
目的:看清在极端高表达的情况下,原始数据的真实分布情况,不受标准化算法的平滑影响。NKG7代表NK细胞或杀伤性T细胞,PF4则是血小板的标志物。
③表达量的空间定位(特征图 FeaturePlot),
把基因表达量直接投射到之前跑出来的UMAP二维平面图上。每个点代表一个细胞,点的颜色深浅代表该基因的表达高低。
这些基因的高亮区域在UMAP图上互不重叠,各自占据特定的区域。这就能从空间分布上说明,聚类算法确实把不同谱系的细胞在物理层面分开了。
MS4A1 和 CD79A 是 B 细胞。
CD3E 和 CD8A 是 T 细胞(特别是CD8阳性T细胞)。
GNLY 是 NK 细胞或杀伤性T细胞。
CD14 和 FCGR3A 是单核细胞(其中CD14是经典单核,FCGR3A是非经典单核或CD16阳性单核)。
FCER1A 通常是浆细胞样树突状细胞。
LYZ 在髓系细胞(单核/巨噬)中高表达。PPBP 是血小板的强烈标志物。
④把每个细胞群最特异的前10个高表达基因(pbmc.markers),全部画在一张热图上。
10.注释每簇
关于每个聚类簇细胞的注释,我们可以参考已经发表的文献,也可以参考一些专门用于细胞注释的网站:https://mp.weixin.qq.com/s/6JgBFoQygCfzUEdIOBngeg,其中CellMarker3.0专注于人类、小鼠两大模式生物,是目前最全面的细胞标记基因权威资源库之一。PanglaoDB,Single Cell Portal、Single Cell Expression Atlas - EMBL-EBI、scRNASeqDB、Human cell atlas、HCA、DISCO、CancerSEA都是不错的细胞注释的网站,可以参考。
不过幸运的是,本次我们使用的数据有已知的细胞类型相匹配,我们就直接定义了。
# =========================10.匹配簇名称===========================================new.cluster.ids<-c("Naive CD4 T","CD14+ Mono","Memory CD4 T","B","CD8 T","FCGR3A+ Mono","NK","DC","Platelet")names(new.cluster.ids)<-levels(pbmc)pbmc<-RenameIdents(pbmc,new.cluster.ids)DimPlot(pbmc,reduction="umap",label=TRUE,pt.size=0.5)+NoLegend()library(ggplot2)plot<-DimPlot(pbmc,reduction="umap",label=TRUE,label.size=4.5)+xlab("UMAP 1")+ylab("UMAP 2")+theme(axis.title=element_text(size=18),legend.text=element_text(size=18))+guides(colour=guide_legend(override.aes=list(size=10)))plot ggsave(filename='./pbmc3k_umap.jpg',height=7,width=12,plot=plot,quality=50)#保存注释图片saveRDS(pbmc,file="./pbmc3k_final.rds")#保存rds文件留待后用在图中直观得显示每簇的名称。
参考:https://satijalab.org/seurat/articles/pbmc3k_tutorial.html#set-up-the-seurat-object
好了,到这里使用pbmc数据进行单细胞分析的流程到此基本结束啦!我们下期见!
乔粒科研工坊:以解牛之法析生信,观微雀之形览科研。欢迎评论区交流!
