单细胞转录组基础分析五:细胞再聚类
单细胞数据分析中,一般需要对可以细分的细胞再聚类,比如本次分析中的T细胞群体可以细分为Naive T cells、CD8+ T cells、Treg cells、Tmemory cells等。细胞子集的提取使用subset函数,主要依据meta.data的分类项目选择,也可以自定义细胞barcode提取。
提取细胞子集
library(Seurat)
library(monocle)
library(tidyverse)
library(patchwork)
rm(list=ls())
dir.create("subcluster")
scRNA <- readRDS("scRNA.rds")
##提取细胞子集
Cells.sub <- subset(scRNA@meta.data, celltype=="T_cells")
scRNAsub <- subset(scRNA, cells=row.names(Cells.sub))
重新降维聚类
因为再聚类的细胞之间差异比较小,所以聚类函数FindClusters()控制分辨率的参数建议调高到resolution = 0.9。
##PCA降维
scRNAsub <- FindVariableFeatures(scRNAsub, selection.method = "vst", nfeatures = 2000)
scale.genes <- rownames(scRNAsub)
scRNAsub <- ScaleData(scRNAsub, features = scale.genes)
scRNAsub <- RunPCA(scRNAsub, features = VariableFeatures(scRNAsub))
ElbowPlot(scRNAsub, ndims=20, reduction="pca")
pc.num=1:10
##细胞聚类
scRNAsub <- FindNeighbors(scRNAsub, dims = pc.num)
scRNAsub <- FindClusters(scRNAsub, resolution = 0.9)
table(scRNAsub@meta.data$seurat_clusters)
metadata <- scRNAsub@meta.data
cell_cluster <- data.frame(cell_ID=rownames(metadata), cluster_ID=metadata$seurat_clusters)
write.csv(cell_cluster,'subcluster/cell_cluster.csv',row.names = F)
##非线性降维
#tSNE
scRNAsub = RunTSNE(scRNAsub, dims = pc.num)
embed_tsne <- Embeddings(scRNAsub, 'tsne')
write.csv(embed_tsne,'subcluster/embed_tsne.csv')
plot1 = DimPlot(scRNAsub, reduction = "tsne")
ggsave("subcluster/tSNE.pdf", plot = plot1, width = 8, height = 7)
ggsave("subcluster/tSNE.png", plot = plot1, width = 8, height = 7)
#UMAP
scRNAsub <- RunUMAP(scRNAsub, dims = pc.num)
embed_umap <- Embeddings(scRNAsub, 'umap')
write.csv(embed_umap,'subcluster/embed_umap.csv')
plot2 = DimPlot(scRNAsub, reduction = "umap")
ggsave("subcluster/UMAP.pdf", plot = plot2, width = 8, height = 7)
ggsave("subcluster/UMAP.png", plot = plot2, width = 8, height = 7)
#合并tSNE与UMAP
plotc <- plot1+plot2+ plot_layout(guides = 'collect')
ggsave("subcluster/tSNE_UMAP.pdf", plot = plotc, width = 10, height = 5)
ggsave("subcluster/tSNE_UMAP.png", plot = plotc, width = 10, height = 5)
Cluster差异分析
diff.wilcox = FindAllMarkers(scRNAsub)
all.markers = diff.wilcox %>% select(gene, everything()) %>% subset(p_val<0.05)
top10 = all.markers %>% group_by(cluster) %>% top_n(n = 10, wt = avg_logFC)
write.csv(all.markers, "subcluster/diff_genes_wilcox.csv", row.names = F)
write.csv(top10, "subcluster/top10_diff_genes_wilcox.csv", row.names = F)
SingleR细胞鉴定
Subcluster的细胞同样可以使用SingleR鉴定细胞类型。使用的时候注意调整参考数据库和分类标签,以便鉴定结果更有针对性。上节使用SingleR时使用的参考数据库是人类主要细胞图谱(HumanPrimaryCellAtlasData),采用分类标签是主分类标签(label.main);这次建议使用人类免疫细胞数据(MonacoImmuneData),分类标签采用精细分类标签(label.fine)。希望详细了解SingleR的朋友可以到github看看:
https://github.com/dviraran/singler
##细胞类型鉴定
library(SingleR)
refdata <- MonacoImmuneData()
testdata <- GetAssayData(scRNAsub, slot="data")
clusters <- scRNAsub@meta.data$seurat_clusters
cellpred <- SingleR(test = testdata, ref = refdata, labels = refdata$label.fine,
method = "cluster", clusters = clusters,
assay.type.test = "logcounts", assay.type.ref = "logcounts")
celltype = data.frame(ClusterID=rownames(cellpred), celltype=cellpred$labels, stringsAsFactors = F)
write.csv(celltype,"subcluster/celltype_singleR.csv",row.names = F)
scRNAsub@meta.data$celltype = "NA"
for(i in 1:nrow(celltype)){
scRNAsub@meta.data[which(scRNAsub@meta.data$seurat_clusters == celltype$ClusterID[i]),'celltype'] <- celltype$celltype[i]}
p1 = DimPlot(scRNAsub, group.by="celltype", label=T, label.size=5, reduction='tsne')
p2 = DimPlot(scRNAsub, group.by="celltype", label=T, label.size=5, reduction='umap')
p3 = plotc <- p1+p2+ plot_layout(guides = 'collect')
ggsave("subcluster/tSNE_celltype.pdf", p1, width=7 ,height=6)
ggsave("subcluster/UMAP_celltype.pdf", p2, width=7 ,height=6)
ggsave("subcluster/celltype.pdf", p3, width=10 ,height=5)
ggsave("subcluster/celltype.png", p3, width=10 ,height=5)
获取帮助
本教程的目的在于把常用的单细胞分析流程串起来,适合有一定R语言基础的朋友参考。分析原理和代码我没有详细解释,官网有详细的教程和权威的解释,翻译好的中文教程也有多个版本,有兴趣的可以搜索一下。如果您学习本教程有一定困难,可以分享此篇文章到朋友圈,截图微信发给Kinesin(文末二维码),我会抽时间组群直播讲解单细胞数据分析的全过程。本专题所用的软件、数据、代码脚本和中间结果,我打包放在了百度云上,需要的朋友添加Kinesin微信索取。