ComBat 的数学原理：位置尺度调整模型与经验贝叶斯估计
嗯，这一章是很多人觉得“劝退”的地方——数学原理。但我想说，你真不用怕它。我在早年做多中心转录组整合时，一开始也是硬着头皮调包，结果换了个数据集就崩了。后来才明白，不懂原理，你调参就是瞎蒙。
位置尺度调整模型：说白了就是“对齐均值、拉齐方差”
ComBat 的核心思想，其实很直观。每个批次都有自己的“风格”——有的批次整体表达量偏高，有的批次数据波动特别大。我们要做的就是：把每个批次的数据，拉到同一个“基准”上。
数学上，对于基因 \( g \) 在样本 \( j \)（属于批次 \( i \)）中的表达值 \( Y_{ijg} \)，ComBat 假设：
Y_{ijg} = α_g + Xβ_g + γ_{ig} + δ_{ig} ε_{ijg}这里：
- α_g：基因 \( g \) 的全局基线表达。可以理解为所有批次统一后的“理想水平”。
- Xβ_g：生物学协变量的影响（比如疾病组 vs 对照组）。这部分要保留，不能校掉。
- γ_{ig}：批次 \( i \) 的加法偏移——这是“位置参数”。说白了，就是批次导致整体上移或下移了多少。
- δ_{ig}：批次 \( i \) 的缩放因子——这是“尺度参数”。数据太散或太紧，就靠它拉齐。
- ε_{ijg}：随机误差，假设服从 \( N(0, σ_g^2) \)。
那么校正后的表达值就是：
Y*_{ijg} = (Y_{ijg} - α_g - Xβ_g - γ_{ig}) / δ_{ig} + α_g + Xβ_g你看，减掉批次效应，再除以批次标准差，然后加回全局基线和生物学效应。简单说：去批次、保生物学。
实战要点：我在项目中遇到过一个问题——如果某个批次样本量特别小（比如只有 2 个样本），位置参数 γ 很容易被极端值带偏。这时候可以考虑用后面讲的经验贝叶斯来“收缩”它，不要只用简单平均。
经验贝叶斯估计：用数据本身来“借力”
位置尺度模型看起来简单，但问题在于：每个基因在每批次里都有独立的 γ 和 δ，参数太多了！如果批次少还好，批次一多（比如 10 个批次，2 万个基因），那就是 40 万个参数要估。直接估计会严重过拟合。
这时经验贝叶斯就派上用场了。它假设每个基因的批次参数不是孤立的，而是来自一个共同的“先验分布”。我用个不严谨但好懂的说法：
- 先看所有基因的 γ 和 δ 在批次内的分布，得到一个“整体印象”。
- 对于单个基因，如果它的 γ 离全体均值很远（可能是噪声），就把它往均值方向“拉”一点——这叫收缩。
- 如果它的 γ 很稳定，就保留原始估计。
具体步骤分为两步：
- 先验估计：对每个基因 \( g \)，用批次内数据粗糙地估计 \( \hat{γ}_{ig} \) 和 \( \hat{δ}_{ig}^2 \)。然后假设这些估计值服从某种分布（例如 γ 为正态分布，δ² 为逆伽马分布），通过矩估计或极大似然得到先验的超参数。
- 后验估计：对每个基因 \( g \)，结合先验和当前基因的数据，计算后验均值和后验方差，作为最终的 γ 和 δ 估计。这里可以用迭代法（EM 算法）或直接公式。
个人习惯：我通常用 ComBat 时，默认用非参数先验（non-parametric），它对分布假设更宽松，鲁棒性更好。只有在样本量非常充足（每批次 > 30）时才考虑参数先验。
整体流程与关键参数
我把整个流程浓缩成一张图，你一看就明白ComBat是怎么把“数据输入”变成“校正输出”的：
注意：我曾经在 TCGA 数据上跑 ComBat，忘记把临床协变量放入 X 矩阵，结果把肿瘤 vs 正常的差异全部校掉了。后来检查发现，校正后的数据两组之间完全没差异…… 所以请务必确认你的生物学变量在模型中。
核心参数速查
| 参数 | 含义 | 典型取值 |
|---|---|---|
| prior.plots | 是否绘制先验密度图（诊断用） | FALSE（不画，或 TRUE 调试） |
| par.prior | 是否使用参数先验 | FALSE（推荐非参数，更稳健） |
| mean.only | 是否只校均值（不校方差） | FALSE（一般要同时校） |
调参小经验：如果你数据里某些批次方差差异特别大，可以试试把 mean.only=FALSE 但把非参数先验的迭代次数调高到 5 次（默认 3）。我试过，效果改善挺明显。
好了，数学原理就到这。你只需要记住两件事：位置尺度模型是“对齐+缩放”，经验贝叶斯是“借力防过拟合”。下次你在代码里看到 ComBat(dat, batch)，心里就知道它在背后干了什么。