单细胞聚类注释Seurat¶

一句话概述：Seurat v5是单细胞RNA-seq分析的标准工具，从质控、归一化、聚类到细胞类型注释提供完整流程，新版支持百万级细胞和多模态数据整合。

核心知识点速览¶

概念	白话解释
Seurat对象	装单细胞数据的"大箱子"，包含表达矩阵+元数据+降维结果
质控（QC）	过滤低质量细胞（死细胞、双细胞、空液滴）
nFeature	一个细胞检测到多少个基因（太少=空液滴，太多=双细胞）
nCount	一个细胞的总UMI数（读段总数）
percent.mt	线粒体基因占比（太高=死细胞/破损细胞）
HVG	高变基因，细胞间差异最大的基因（用于降维）
PCA	主成分分析，第一步降维
UMAP/tSNE	可视化降维，把高维数据画到2D
Leiden/Louvain	图聚类算法，给细胞自动分组
Resolution	聚类分辨率，越大cluster越多
FindMarkers	找每个cluster的标记基因
细胞注释	给cluster贴标签（如"T细胞"、"B细胞"）

一、标准分析流程¶

1.1 安装和数据加载¶

# 安装Seurat v5
install.packages("Seurat")  # 安装
library(Seurat)               # 加载

# 读取10x Genomics数据（cellranger输出）
data_dir <- "filtered_feature_bc_matrix/"  # 10x输出目录
counts <- Read10X(data.dir = data_dir)     # 读取三个文件

# 创建Seurat对象
seurat_obj <- CreateSeuratObject(
  counts = counts,            # 计数矩阵
  project = "PBMC3k",        # 项目名
  min.cells = 3,             # 基因至少在3个细胞中表达
  min.features = 200         # 细胞至少检测到200个基因
)
seurat_obj  # 查看对象摘要

1.2 质控过滤¶

# 计算线粒体基因比例
seurat_obj[["percent.mt"]] <- PercentageFeatureSet(
  seurat_obj,
  pattern = "^MT-"  # 人类线粒体基因以MT-开头（小鼠用mt-）
)

# 查看QC指标分布
VlnPlot(seurat_obj,
        features = c("nFeature_RNA", "nCount_RNA", "percent.mt"),
        ncol = 3, pt.size = 0)  # 小提琴图

# 散点图检查指标间关系
FeatureScatter(seurat_obj, feature1 = "nCount_RNA", feature2 = "nFeature_RNA")

# 过滤低质量细胞
seurat_obj <- subset(seurat_obj,
  nFeature_RNA > 200 &        # 至少200个基因
  nFeature_RNA < 5000 &        # 最多5000个基因（排除双细胞）
  percent.mt < 20              # 线粒体<20%（排除死细胞）
)
cat("过滤后细胞数:", ncol(seurat_obj), "\n")

1.3 归一化¶

# 方法1：LogNormalize（默认，快速）
seurat_obj <- NormalizeData(
  seurat_obj,
  normalization.method = "LogNormalize",  # log(x/total * 10000 + 1)
  scale.factor = 10000                     # 缩放因子
)

# 方法2：SCTransform（推荐，更准确）
# seurat_obj <- SCTransform(seurat_obj, verbose = FALSE)

1.4 高变基因选择¶

# 选择高变基因（HVG）
seurat_obj <- FindVariableFeatures(
  seurat_obj,
  selection.method = "vst",  # VST方法（推荐）
  nfeatures = 2000           # 选2000个高变基因
)

# 查看高变基因
top10 <- head(VariableFeatures(seurat_obj), 10)  # 前10个
VariableFeaturePlot(seurat_obj)  # 高变基因图

# 标注前10个基因
LabelPoints(VariableFeaturePlot(seurat_obj),
            points = top10, repel = TRUE)

1.5 缩放和PCA¶

# 缩放数据（零均值单位方差）
seurat_obj <- ScaleData(
  seurat_obj,
  features = rownames(seurat_obj)  # 对所有基因缩放
  # vars.to.regress = "percent.mt"  # 可选：回归掉线粒体效应
)

# PCA降维
seurat_obj <- RunPCA(
  seurat_obj,
  features = VariableFeatures(seurat_obj),  # 用高变基因
  npcs = 50                                  # 计算50个主成分
)

# 碎石图（确定用多少个PC）
ElbowPlot(seurat_obj, ndims = 50)
# 看"拐点"在哪里，通常取10-30个PC

1.6 聚类¶

# 构建KNN图
seurat_obj <- FindNeighbors(
  seurat_obj,
  dims = 1:20  # 用前20个PC
)

# Leiden聚类（Seurat v5默认）
seurat_obj <- FindClusters(
  seurat_obj,
  resolution = 0.5,       # 分辨率（0.1-2.0，越大cluster越多）
  algorithm = 4            # 4=Leiden（推荐），1=Louvain
)

# 查看聚类结果
table(Idents(seurat_obj))  # 每个cluster的细胞数

1.7 UMAP可视化¶

# UMAP降维
seurat_obj <- RunUMAP(
  seurat_obj,
  dims = 1:20  # 用前20个PC
)

# 画UMAP
DimPlot(seurat_obj,
        reduction = "umap",
        label = TRUE,            # 显示cluster编号
        label.size = 5,          # 标签大小
        pt.size = 0.5) +         # 点大小
  ggtitle("UMAP Clustering")

二、细胞类型注释¶

2.1 找标记基因（手动注释）¶

# 找每个cluster的标记基因
all_markers <- FindAllMarkers(
  seurat_obj,
  only.pos = TRUE,        # 只看上调的标记
  min.pct = 0.25,         # 至少在25%的细胞中表达
  logfc.threshold = 0.5   # log2FC > 0.5
)

# 每个cluster的Top5标记基因
top5 <- all_markers %>%
  group_by(cluster) %>%
  slice_max(avg_log2FC, n = 5)

# 画热图
DoHeatmap(seurat_obj, features = top5$gene) +
  scale_fill_viridis_c()

# 画气泡图
DotPlot(seurat_obj,
        features = unique(top5$gene)) +
  RotatedAxis()  # 旋转X轴标签

2.2 已知标记基因注释¶

# 用已知的细胞类型标记基因确认注释
# PBMC常见标记：
marker_genes <- list(
  "T cells" = c("CD3D", "CD3E", "IL7R"),
  "CD8 T" = c("CD8A", "CD8B"),
  "NK cells" = c("NKG7", "GNLY", "KLRD1"),
  "B cells" = c("MS4A1", "CD79A", "CD79B"),
  "Monocytes" = c("CD14", "LYZ", "CST3"),
  "DC" = c("FCER1A", "CD1C"),
  "Platelets" = c("PPBP", "PF4")
)

# 逐个检查
FeaturePlot(seurat_obj,
            features = c("CD3D", "CD14", "MS4A1", "NKG7"),
            ncol = 2)

# 气泡图检查所有标记
DotPlot(seurat_obj,
        features = unlist(marker_genes),
        group.by = "seurat_clusters") +
  RotatedAxis()

# 手动注释（根据标记基因表达模式）
new_ids <- c("CD4 T", "CD14 Mono", "B cells", "CD8 T",
             "NK cells", "CD16 Mono", "DC", "Platelets")
names(new_ids) <- levels(seurat_obj)
seurat_obj <- RenameIdents(seurat_obj, new_ids)

# 重新画UMAP
DimPlot(seurat_obj, label = TRUE, pt.size = 0.5) +
  ggtitle("Cell Type Annotation")

2.3 自动注释（SingleR）¶

# 安装SingleR
BiocManager::install("SingleR")     # 自动注释工具
BiocManager::install("celldex")     # 参考数据集

library(SingleR)
library(celldex)

# 获取参考数据集
ref <- celldex::HumanPrimaryCellAtlasData()  # 人类细胞图谱参考

# 运行自动注释
predictions <- SingleR(
  test = GetAssayData(seurat_obj, layer = "data"),  # 表达矩阵
  ref = ref,                                          # 参考数据
  labels = ref$label.main                             # 参考标签
)

# 添加注释到Seurat对象
seurat_obj$SingleR_labels <- predictions$labels

# 可视化
DimPlot(seurat_obj, group.by = "SingleR_labels",
        label = TRUE, repel = TRUE, pt.size = 0.5)

三、Seurat v5新特性¶

3.1 Layer-based数据组织¶

# Seurat v5用"层"管理多样本数据
# 同一个对象可以包含多个样本的数据

# 按样本分层
seurat_obj[["RNA"]] <- split(seurat_obj[["RNA"]],
                              f = seurat_obj$sample)  # 按样本拆分

# 查看层
Layers(seurat_obj)  # 列出所有层

# 每个层独立归一化
seurat_obj <- NormalizeData(seurat_obj)
seurat_obj <- FindVariableFeatures(seurat_obj)
seurat_obj <- ScaleData(seurat_obj)

3.2 整合分析¶

# Seurat v5整合（比v3/v4更简洁）
seurat_obj <- IntegrateLayers(
  seurat_obj,
  method = CCAIntegration,     # 整合方法（CCA/RPCA/Harmony）
  orig.reduction = "pca",
  new.reduction = "integrated.cca"
)

# 也可以用Harmony整合
seurat_obj <- IntegrateLayers(
  seurat_obj,
  method = HarmonyIntegration,
  orig.reduction = "pca",
  new.reduction = "harmony"
)

# 合并层（整合后需要）
seurat_obj <- JoinLayers(seurat_obj)

常见报错与解决¶

报错信息	原因	解决方案
`Not enough memory`	数据太大	用Seurat v5的sketch-based方法
`No variable features found`	没运行FindVariableFeatures	按顺序运行流程
`Error in FindNeighbors: dims`	PC数超过计算的数	检查RunPCA的npcs参数
`resolution: cluster太多/太少`	分辨率参数不合适	调整resolution（0.1-2.0）
`SingleR: no reference`	参考数据未下载	检查网络，确认celldex安装
`Read10X: file not found`	路径错误或文件不全	需要barcodes/features/matrix三个文件

速查表¶

# Seurat v5标准流程（8步）
Read10X() → CreateSeuratObject() → QC过滤
→ NormalizeData() → FindVariableFeatures()
→ ScaleData() → RunPCA() → ElbowPlot()
→ FindNeighbors() → FindClusters()
→ RunUMAP() → DimPlot()
→ FindAllMarkers() → 细胞注释

# 关键参数
min.cells: 3（基因最少出现在3个细胞）
min.features: 200（细胞最少200个基因）
nFeature上限: 2500-5000（排除双细胞）
percent.mt上限: 10-20%（排除死细胞）
nfeatures: 2000（高变基因数）
dims: 10-30（PCA维度数，看碎石图）
resolution: 0.3-1.5（聚类分辨率）

# 可视化函数
DimPlot()       — UMAP/tSNE散点图
FeaturePlot()   — 基因表达UMAP
VlnPlot()       — 小提琴图
DotPlot()       — 气泡图
DoHeatmap()     — 热图
FeatureScatter() — 散点图

# 注释方法
手动: FindAllMarkers() + 已知标记基因
自动: SingleR + celldex参考数据集
辅助: scType, CellTypist, Azimuth

面试高频问题¶

Q1：Seurat v5相比v4有哪些重要改进？¶

答：①Layer-based数据结构——多样本数据存在同一对象的不同"层"中，不再需要merge+split；②Sketch-based分析——支持百万级细胞，先用代表性子集快速分析，再映射回全量数据；③BPCells后端——支持on-disk存储，突破内存限制；④简化的整合接口——IntegrateLayers()统一了CCA/RPCA/Harmony等整合方法；⑤presto加速——差异分析速度大幅提升。

Q2：如何确定聚类的分辨率（resolution）？¶

答：没有固定标准，需要结合生物学知识。一般策略：①从多个resolution（0.1, 0.3, 0.5, 0.8, 1.0, 1.5）尝试，用clustree包可视化不同分辨率之间cluster的"继承关系"；②检查是否每个cluster都有独特的标记基因——如果两个cluster的标记基因几乎一样，说明resolution太高被过度分割了；③参考已知的细胞类型数量——如果已知PBMC应有7-8种细胞类型，resolution不应产生20+个cluster。

Q3：如何判断一个细胞是否是双细胞（doublet）？¶

答：双细胞是两个细胞被包在同一个液滴中，混合了两种细胞的RNA。判断方法：①nFeature异常高——基因数远超同类细胞的中位数；②nCount异常高——UMI总数远超平均；③表达混合标记——同时高表达两种不同细胞类型的标记基因；④使用专门工具——DoubletFinder、Scrublet等算法可以自动检测双细胞。通常在QC步骤用nFeature上限（如5000）过滤，或运行DoubletFinder后移除。

Q4：为什么要选择高变基因？选多少合适？¶

答：高变基因（HVG）是细胞间表达差异最大的基因，代表了细胞间的生物学差异（而非技术噪音）。选HVG的原因：①降低噪声——管家基因（所有细胞都一样表达）不包含细胞类型信息；②减少计算量——只用2000-3000个基因做PCA比用2万个快很多。通常选2000-3000个HVG。太少可能遗漏重要信号，太多会引入噪声。如果关注的稀有细胞类型标记不在HVG中，可以手动添加。

Q5：手动注释和自动注释各有什么优缺点？¶

答：手动注释（FindMarkers+已知标记）：优点是灵活、可解释性强、可以发现新细胞类型；缺点是耗时、需要专业知识、主观性强。自动注释（SingleR/Azimuth）：优点是快速、客观、可重复；缺点是依赖参考数据库质量、可能遗漏新类型、对稀有细胞注释不准。最佳实践是两者结合：先自动注释获得初步结果，再用标记基因手动验证和调整。