单细胞分析入门（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文件留待后用