为什么要关注批次效应？我踩过的那些坑

做组学数据分析的，十有八九都会被批次效应坑过。我自己入行第三年就吃过一次大亏——两个批次的血样RNA-seq数据合并后，差异基因列表里居然有30%是假的，最后发现是试剂盒换批次了。从那时候起，我就把批次效应检测列为数据质控的第一步。

说白了，批次效应就是非生物学变量引入的系统性偏差。比如同一个人在两台机器上测出来的表达量能差两倍，你想想看，这种差异会直接淹没真实的生物学信号。

批次效应到底长什么样？

我习惯用PCA先看一眼。正常的数据，同一个生物学条件下的样本会聚成一团；有批次效应的话，同一批次的样本会先抱团，把生物学分组打散。下面这张图是我早年项目里的真实数据效果：

哪些场景最容易踩雷？

我梳理了一下，至少有五类常见场景，几乎每个做多批次数整合的人都会遇到：

场景	典型特征	风险等级
多批测序数据合并	不同时间/不同测序仪/不同试剂盒	极高
公共数据库整合	GEO / TCGA / SRA 来源混合	高
多中心临床样本	不同医院、不同操作人员、不同保存条件	高
长时间跨度实验	样本采集跨月/跨年	中高
技术重复与生物重复混淆	同一批制备的多次测量 vs 独立采样	中

为什么必须检测？不做会怎样

简单说，不检测批次效应，你的结论可能就是垃圾。我个人经历过：

• 假阳性率飙升——差异基因列表里充满批次相关的基因，比如线粒体基因、热休克蛋白，这些都是技术偏差的常见受害者。
• 聚类完全被批次的标签控制——你以为是发现了新的亚型，其实是第1批和第2批的试剂不一样。
• 生物标志物不可复现——换个批次验证，marker直接消失。

我早年发了一篇paper，审稿人要求用另一批独立数据验证——结果差异基因对不上，最后花了一个月重做实验才发现是批次效应。所以我现在的习惯是：任何合并数据集的第一步，先做batch detection。

检测批次效应的常用思路

其实检测方法不复杂，我按从粗到细列出来：

1. 可视化检查——PCA / t-SNE / UMAP，看样本点是否按批次聚集。
2. 定量指标——LME模型、PCA投影的批次占比、PCA的轮廓系数。
3. 统计检验——Kruskal-Wallis检验每个基因的表达是否在不同批次间有显著差异。
4. 机器学习方法——训练一个分类器，看能不能仅根据表达量准确预测批次标签。

这里我给一个最轻量的命令，用R语言的sva包就能快速跑起来：

library(sva)
# exprs: 基因表达矩阵 (行=基因, 列=样本)
# batch: 批次标签向量
# mod: 你关注的生物学协变量 (比如处理组/对照组)
mod <- model.matrix(~condition, data = pData)
batch <- pData$batch
adjusted_exprs <- ComBat(dat = exprs, batch = batch, mod = mod)
# 然后对比前后PCA图，看批次聚集是否被削弱

现实中如何区分批次效应与真实生物学差异？

这是个灵魂拷问。我的经验是：没有绝对完美的分离方法，只能结合实验设计来找线索。

• 看重复样本——同一个样本分到不同批次处理，如果它们表达谱的距离远大于同一批次的重复，那就是批次效应。
• 看已知的批次无关基因——比如管家基因、内参基因，如果在不同批次间出现系统性漂移，说明问题很大。
• 看批次变量的来源——如果是“第1批 vs 第2批”且刚好对应不同的RNA提取试剂盒，几乎板上钉钉。

我有时会用混杂PCA（CPCA）来做解耦，但更实用的还是把所有元数据（包括日期、操作员、试剂盒批号、温度湿度）全部记录下来——等你发现批次效应的时候，这些都是溯源的关键变量。

一个现实案例：多中心RNA-seq数据

曾经合作过一个项目，三个医院的结直肠癌样本要一起分析。我拿到数据后照例画了PCA：

# R代码 - 标准化并画PCA
exprs_scaled <- t(scale(t(exprs)))
pca <- prcomp(t(exprs_scaled))
plot(pca$x[,1], pca$x[,2], col = as.numeric(batch))  # 按医院批次着色

结果很明显——每个医院成一团，处理组和对照组混在一起。我用上面的ComBat调整后，PCA图才勉强按疾病分组。但这还不够，后续的差异分析里我特意加入了一个协变量“医院”，把残留的批次效应回归掉。最终只保留了那些在两个医院中方向一致、效应相近的基因——这才是真正靠谱的结果。