1. 批次效应与批量校正简介

做了这么多年组学数据分析,我越来越觉得批次效应是最容易被忽视、但后果最严重的坑。你花三个月建库、两个月测序、又花一个月跑差异分析,结果发现显著基因全是批次间的差异——那种欲哭无泪的感觉,我经历过。

说白了,批次效应就是非生物学变量带来的系统偏差。样品采集时间不同、试剂盒批号不同、操作人员不同、测序仪器不同,都会引入这个“幽灵”。它不除掉,后面的分析全是扯淡。

什么是批次效应

想象一下你比较两组患者的 RNA-seq 数据:一组是去年做的,一组是今年做的。实验室换了新试剂盒,测序仪也升级了。结果差异基因列表里,一半是核糖体蛋白、一半是线粒体基因——典型的批次信号。

从统计角度看,批次效应是一种混杂因子。它的特点是:

  • 系统性:影响同一批次内所有样本
  • 非随机:与生物学分组可能相关或正交
  • 不可忽略:如果不处理,假阳性率可以飙到 80% 以上
维度批次效应特征生物学效应特征
来源实验操作、试剂、仪器基因表达、突变、甲基化
可重复性批次内一致,批次间漂移在不同条件下可验证
统计影响增加方差、混淆组间差异代表真实差异
校正方法SVA、ComBat、limma 等不需要校正
⚠️ 曾经踩过的坑:有一次我做多中心整合分析,三个中心的 RNA-seq 数据分别来自不同平台。我偷懒只做了 TMM 标准化就合并了,结果聚类图是按中心分三团,而不是按疾病分型。那一周都在加班重做——之后我再也不敢跳过批次校正。

为什么要校正批次效应

你可能会想:“我设计实验时把批次和分组平衡了行不行?”理论上可以,但现实中很困难。样本不是同时收集的,随诊数据逐年增加,公共数据库的合并更是家常便饭。

不校正的后果很直接:

  1. 假阳性差异基因——那些高表达的“基因”可能只是某个批次噪音大
  2. 假阴性遗漏——真实信号被批次方差淹没
  3. 聚类失真——PCA 第一主成分一定是批次,而不是生物学分组
  4. 后续模型失效——机器学习模型学到的是批次规律,不是生物学规律
💡 我的经验:拿到新数据后,第一件事不是差异分析,而是用 PCA 或 UMAP 按批次着色看看。如果样本按批分成几团,说明批次效应显著,必须校正。这个习惯帮我避免了至少三次“假阳性狂欢”。

批量校正的常见策略

校正方法大致分两类:基于模型的和基于对齐的。我个人更习惯先用无监督方法看看严重程度,再决定用哪种。

  • 均值中心化:简单但粗暴,假设批次只影响均值,不考虑方差变化。适用于跨平台微阵列数据。
  • ComBat:经验贝叶斯方法,同时调整均值和方差。我一般在 RNA-seq 和甲基化数据上用,效果稳定。
  • SVA(Surrogate Variable Analysis):用奇异值分解估计潜在变量,再把它们作为协变量纳入模型。适合你不知道具体批次来源时使用。
  • limma removeBatchEffect:线性模型移除批次效应,保留生物学对比。适合表达谱芯片的快速处理。
  • Harmony / scVI:单细胞数据专用,基于深度学习和迭代聚类对齐。

看到这么多方法,别慌。选哪个取决于你的数据特性和分析目标。我通常这样判断:

数据类型推荐方法原因
微阵列表达谱ComBat方差稳健,已知批次信息
RNA-seq 批量合并SVA + 差异模型未知隐变量也能捕捉
单细胞多批次Harmony / scVI高稀疏、异质性大
甲基化芯片ComBat(beta 值)Probe 水平批次效应明显
关键原则:校正不能引入人为的生物学信号。校正后的数据要能保留真实组间差异,同时消除批次间的分布差异。验证方法是用校正后的数据做 PCA,按批次着色——如果各批次混在一起,说明校正有效;如果还是分开,说明方法不对或批次效应太强。

SVA 与 ComBat 简介

这门课重点讲这两个方法,因为它们互补性强。SVA 擅长处理未知的“隐藏”批次,ComBat 擅长处理已知批次的显式校正。

SVA(Surrogate Variable Analysis)

我最初接触 SVA 是在关联分析中。你有一个表型 Y,一个大表达式矩阵 X,但总有些没测到的混杂因素(比如样本保存时间、RNA 质量)。SVA 会从残差中提取“替代变量”(surrogate variables),然后作为协变量塞进回归模型里。说白了,就是让数据自己告诉你还有什么没控制住。

library(sva)
# 假设 expression_matrix 是基因×样本,pheno 是样本表型
mod <- model.matrix(~pheno)
mod0 <- model.matrix(~1, data=pheno)
svobj <- sva(expression_matrix, mod, mod0)
# 将 svobj$sv 加入后续差异分析模型

ComBat(Combating Batch Effects)

ComBat 是 Johnson 等人 2007 年提出的。它假设批次效应表现为均值和方差偏移,用经验贝叶斯收缩估计每个批次的调整参数。这个方法在微阵列时代就是黄金标准,现在在 RNA-seq 和甲基化里也是主流。

library(sva)
# batch 向量指明每个样本所属批次
adjusted_data <- ComBat(dat = expression_matrix, batch = batch_vector)
💡 我自己的经验:ComBat 默认参数在大多数场景下工作很好。但如果批次间样本量差异极大(比如 A 批 100 个,B 批 3 个),建议设置 mean.only=TRUE,只调均值不调方差,防止小批次被过度校正。

好了,理论知识先到这里。下面的 SVG 图总结了整个批次校正流程,从发现问题到选择方法再到验证,一目了然。

批次效应校正全流程 原始多批次数据 PCA / UMAP 按批次着色 有批次效应? 无 → 直接下游分析 有,已知批次 → ComBat 有,未知隐藏 → SVA 再次PCA验证校正效果 虚线箭头表示校正后直接进入验证,实线箭头表示常规流程

图里我把两种场景都画了——知道批次就用 ComBat,不知道批次来源就用 SVA。校正完一定要回到 PCA 验证,这是闭环。很多初学者做完校正就跑了,结果发现批次效应还在——嗯,我自己也犯过一次这样的错。

好了,这一章先把批次效应的概念和校正概览讲清楚。接下来的章节,我会手把手带你过 ComBat 和 SVA 的具体代码、参数调优以及实战案例。记住一条原则:没有批次校正的多组学整合,都是在做玩具分析。