单细胞测序—2次分群
Seurat里的FindClusters函数设置的resolution数值越大,分群的数量就越多,但是当单细胞数量太多的时候,会遇到resolution再变大,分群的数量也不再增加的情况。一次分群分不开时就会需要二次分群。
这里的示例数据seu.obj.Rdata是GSE218208降维聚类分群的结果,参照单细胞测序—GSE218208(流程简化)
代码语言:r复制rm(list = ls())
library(Seurat)
library(dplyr)
load("../2.GSE218208/seu.obj.Rdata")
p1 = DimPlot(seu.obj, reduction = "umap",label=T) NoLegend()
p1
1 二次分群
这里以树突细胞(DC)为例进行二次分群,想要切换别的细胞类型直接修改下面的my_sub即可。
核心就是提取感兴趣的亚群的细胞,后面就是标准流程和可视化了,没有区别
代码语言:r复制my_sub = "DC"
sub.cells <- subset(seu.obj, idents = my_sub)
f = "obj.Rdata"
if(!file.exists(f)){
sub.cells = sub.cells %>%
NormalizeData() %>%
FindVariableFeatures() %>%
ScaleData(features = rownames(.)) %>%
RunPCA(features = VariableFeatures(.)) %>%
FindNeighbors(dims = 1:15) %>%
FindClusters(resolution = 0.5) %>%
RunUMAP(dims = 1:15)
save(sub.cells,file = f)
}
load(f)
DimPlot(sub.cells, reduction = 'umap',label = T) NoLegend()
注:
subset(seu.obj, idents = my_sub)
从 seu.obj 对象中提取所有标记为 "DC" 的细胞,创建一个新的 Seurat对象 sub.cells。这个对象包含了所有的树突细胞。
2 Marker基因及其可视化
代码语言:r复制sub.cells.markers <- FindAllMarkers(sub.cells, only.pos = TRUE,
min.pct = 0.25, logfc.threshold = 0.25)
top10 <- sub.cells.markers %>%
group_by(cluster) %>%
top_n(n = 10, wt = avg_log2FC) %>%
pull(gene);top10
## [1] "JCHAIN" "IGKC" "MZB1" "PACSIN1" "WNT10A" "MAP1A" "VASH2"
## [8] "NIBAN3" "SMPD3" "TNFRSF4" "LYZ" "TIMP1" "GPAT3" "ITGAX"
## [15] "SAMSN1" "OLR1" "FPR3" "EREG" "FCN1" "AKAP12"
VlnPlot(sub.cells, features = top10)
RidgePlot(sub.cells, features = top10)
FeaturePlot(sub.cells, features = top10)
DotPlot(sub.cells,features = top10) RotatedAxis()
DoHeatmap(sub.cells, features = top10) NoLegend()
3 放回原有的seurat对象中
上面的umap图是感兴趣的单独的展示,也可以把它放回原来的seurat对象里。
代码语言:r复制sub.cells@meta.data$celltype = paste0("M",sub.cells$seurat_clusters)
#整合二次分群结果到原始对象
seu.obj$celltype = as.character(Idents(seu.obj))
seu.obj$celltype = ifelse(seu.obj$celltype==my_sub,
sub.cells$celltype[match(colnames(seu.obj),colnames(sub.cells))],
seu.obj$celltype)
Idents(seu.obj) = seu.obj$celltype
p2 = DimPlot(seu.obj,label = T) NoLegend()
p1 p2
对比二次分群前的结果,可以看到DC被进一步划分为M1,M0两群。
注:
sub.cells@meta.data$celltype
sub.cells@meta.data$celltype:在 sub.cells 的元数据中创建一个新列 celltype,用于存储每个细胞的细胞类型信息。paste0("M", sub.cells$seurat_clusters):为每个细胞分配一个标识符,前缀 "M" 表示这些细胞来自于二次分群的树突细胞群体,后面跟随的是其聚类编号。
seu.obj$celltype = as.character(Idents(seu.obj))
将原始Seurat对象 seu.obj 中每个细胞的群体标识复制到新的 celltype 列中。这一步确保了所有细胞都有默认的群体标识。
ifelse
seu.obj$celltype == my_sub:检查 seu.obj 中的细胞是否属于原始树突细胞群体(即 my_sub)。sub.cells$celltypematch(colnames(seu.obj), colnames(sub.cells)):如果细胞属于树突细胞群体,则使用 sub.cells 中的 celltype 信息替换其 seu.obj 中的 celltype**。**match 函数用于对齐 seu.obj 和 sub.cells中的细胞。否则,保持 seu.obj$celltype 的原始值不变。
Idents(seu.obj) = seu.obj$celltype
Idents(seu.obj) = seu.obj$celltype:更新 seu.obj中的身份信息,使其与 celltype列保持一致。这一步使得 seu.obj中的每个细胞都有更新的标识,包含了二次分群的结果。