什么是批次效应？
批次效应（Batch Effect）是指在不同批次实验中产生的非生物学差异。这些差异可能来源于：
🔬 技术因素
- • 不同测序平台（10x vs Smart-seq）
- • 不同试剂盒批次
- • 实验操作人员差异
📅 时间因素
- • 不同实验日期
- • 样本处理时间差异
- • 试剂保存时间差异
为什么必须去除批次效应？
批次效应会掩盖真实的生物学差异，导致错误的聚类结果和生物学结论。 例如，不同批次的相同细胞类型可能会被聚类成不同的亚群， 或者在差异分析中产生大量与批次相关的假阳性基因。
三种主流方法对比
| 方法 | 原理 | 优点 | 缺点 | 适用场景 |
|---|---|---|---|---|
| CCA | 典型相关分析，寻找锚点对齐 | 经典方法，结果可靠 | 内存消耗大，运行较慢 | 跨平台整合，追求准确性 |
| SCT | 基于SCTransform的标准化 | 整合效果最好 | 运行时间最长 | 数据质量参差不齐 |
| Harmony | 迭代校正低维嵌入 | 速度快，内存友好 | 对稀有细胞敏感度稍低 | 大规模数据整合 |
方法一：CCA（Canonical Correlation Analysis）
CCA是Seurat包最早推出的数据整合方法，通过寻找不同数据集之间的典型相关变量（锚点）来对齐数据。 这种方法适合于跨平台、跨物种的复杂数据整合。
核心概念
实战代码
# 准备工作
library(Seurat)
library(SeuratData)
library(patchwork)
# 加载示例数据
InstallData("ifnb")
ifnb <- LoadData("ifnb")
# 按刺激条件分割数据集
ifnb.list <- SplitObject(ifnb, split.by = "stim")
# 对每个数据集进行标准化和寻找高可变基因
ifnb.list <- lapply(X = ifnb.list, FUN = function(x) {
x <- NormalizeData(x)
x <- FindVariableFeatures(x, selection.method = "vst", nfeatures = 2000)
})
# 寻找锚点并进行整合
ifnb.anchors <- FindIntegrationAnchors(
object.list = ifnb.list,
dims = 1:30
)
ifnb.integrated <- IntegrateData(
anchorset = ifnb.anchors,
dims = 1:30
)性能提示
CCA方法会消耗大量内存，对于超过3个数据集的整合，建议使用reference参数指定参考数据集，
这可以大幅减少运行时间。例如：reference = c(1, 3)表示先整合第1和第3个数据集作为参考。
后续分析流程
# 设置使用整合后的数据
DefaultAssay(ifnb.integrated) <- "integrated"
# 标准流程（注意：不需要再运行ScaleData）
ifnb.integrated <- ScaleData(ifnb.integrated, verbose = F)
ifnb.integrated <- RunPCA(ifnb.integrated, npcs = 30, verbose = F)
ifnb.integrated <- RunUMAP(ifnb.integrated, reduction = "pca", dims = 1:30)
ifnb.integrated <- FindNeighbors(ifnb.integrated, reduction = "pca", dims = 1:30)
ifnb.integrated <- FindClusters(ifnb.integrated, resolution = 0.5)
# 可视化整合效果
p1 <- DimPlot(ifnb.integrated, reduction = "umap", group.by = "stim")
p2 <- DimPlot(ifnb.integrated, reduction = "umap",
group.by = "seurat_annotations", label = T) + NoLegend()
p1 + p2重要提示
聚类使用integrated，差异分析使用RNA：
整合后的数据仅用于聚类和可视化，后续的差异基因分析仍需使用原始RNA数据。
可以通过DefaultAssay(object) <- "RNA"切换数据源。
方法二：SCT（SCTransform-based Integration）
SCT方法基于SCTransform标准化技术，在数据整合的同时进行标准化。 相比CCA，SCT能够更彻底地去除批次效应，但运行时间也相应更长。
实战代码
# 准备工作
library(Seurat)
library(SeuratData)
# 加载并分割数据
ifnb <- LoadData("ifnb")
ifnb.list <- SplitObject(ifnb, split.by = "stim")
# 使用SCTransform进行标准化
ifnb.list <- lapply(X = ifnb.list, FUN = SCTransform)
# 选择整合特征
features <- SelectIntegrationFeatures(
object.list = ifnb.list,
nfeatures = 3000
)
# 准备整合
ifnb.list <- PrepSCTIntegration(
object.list = ifnb.list,
anchor.features = features
)
# 寻找锚点并整合
ifnb.anchors <- FindIntegrationAnchors(
object.list = ifnb.list,
normalization.method = "SCT",
anchor.features = features
)
ifnb.combined <- IntegrateData(
anchorset = ifnb.anchors,
normalization.method = "SCT"
)SCTransform的优势
SCTransform可以一次性完成标准化、高可变基因筛选和线性归一化三个步骤， 且能有效去除测序深度相关的影响。因此SCT整合方法在数据质量参差不齐时表现更好。
后续分析
# PCA和聚类
ifnb.combined <- RunPCA(ifnb.combined, verbose = F)
ifnb.combined <- RunUMAP(ifnb.combined, reduction = "pca", dims = 1:30)
# 可视化
p1 <- DimPlot(ifnb.combined, reduction = "umap", group.by = "stim")
p2 <- DimPlot(ifnb.combined, reduction = "umap",
group.by = "seurat_annotations", label = T)
p1 + p2注意事项
不要重复运行ScaleData： SCTransform已经包含了数据标准化步骤， 因此在整合后不需要再运行ScaleData，这与CCA方法不同。
方法三：Harmony
Harmony是一种基于迭代校正的批次效应去除方法，直接在PCA降维后的空间中进行校正。 相比CCA和SCT，Harmony具有明显的速度优势，特别适合大规模数据整合。
实战代码
# 安装Harmony（仅需一次）
install.packages("harmony")
# 加载包
library(harmony)
library(Seurat)
library(tidyverse)
# 准备数据（Harmony可以直接处理合并后的数据）
ifnb <- LoadData("ifnb")
ifnb.list <- SplitObject(ifnb, split.by = "stim")
# 合并数据
ifnb.merged <- merge(ifnb.list[[1]], ifnb.list[[2]])
# 标准流程（使用管道符简化代码）
ifnb.merged <- ifnb.merged %>%
NormalizeData() %>%
FindVariableFeatures() %>%
ScaleData() %>%
RunPCA(verbose = F)
# 运行Harmony校正
ifnb.harmony <- RunHarmony(
ifnb.merged,
group.by.vars = "stim",
max.iter.harmony = 20
)
# 使用Harmony校正后的数据进行聚类
ifnb.harmony <- RunUMAP(
ifnb.harmony,
reduction = "harmony",
dims = 1:30
)
ifnb.harmony <- FindNeighbors(
ifnb.harmony,
reduction = "harmony",
dims = 1:30
) %>%
FindClusters()
# 可视化
p1 <- DimPlot(ifnb.harmony, reduction = "umap", group.by = "stim")
p2 <- DimPlot(ifnb.harmony, reduction = "umap",
group.by = "seurat_annotations", label = T)
p1 + p2Harmony的参数说明
- •
group.by.vars: 指定批次变量，如"stim"、"tech"、"orig.ident" - •
max.iter.harmony: 最大迭代次数，默认为10 - •
lambda: 整合力度参数（0.5-2），值越小整合力度越大
高级话题
寻找保守Marker基因
在多数据集分析中，我们可能需要寻找在所有条件下都稳定表达的保守marker基因。 这些基因不受批次效应影响，是更可靠的细胞类型标记。
# 切换到原始RNA数据
DefaultAssay(ifnb.integrated) <- "RNA"
# 寻找保守marker基因
conserved.markers <- FindConservedMarkers(
ifnb.integrated,
ident.1 = 1, # 指定cluster
grouping.var = "stim", # 批次变量
verbose = F
)
# 查看结果
head(conserved.markers)结果解读
FindConservedMarkers会在每个分组中分别进行差异分析，
然后汇总结果。输出结果包含每个分组的p值、logFC等信息，
最后还有合并的统计量。选择在所有分组中都显著的基因作为保守marker。
指定参考数据集
当有多个数据集需要整合时，可以通过指定参考数据集来大幅减少运行时间。
# 选择部分数据集
pancreas.list <- pancreas.list[c("celseq", "celseq2", "smartseq2")]
# 指定第1和第3个数据集作为参考
pancreas.anchors <- FindIntegrationAnchors(
object.list = pancreas.list,
dims = 1:30,
reference = c(1, 3) # 关键参数
)
pancreas.integrated <- IntegrateData(
anchorset = pancreas.anchors,
dims = 1:30
)策略建议
对于5个以上的数据集，建议选择数据质量最好或样本量最大的2-3个数据集作为参考。 这样可以在保证整合质量的前提下，显著减少计算时间。
方法选择建议
追求准确性
选择CCA或SCT，适合跨平台整合或数据质量差异大的情况
追求速度
选择Harmony，适合大规模数据整合或快速探索性分析
不确定？
三种方法都试试，选择可视化效果最好的那个