4.4 归一化与差异表达分析

4.3 节的实验流程以图像分析结束：对每个点读取前景与背景荧光，产出一张「基因 × 样本」的强度矩阵（4.1 节的表达矩阵）。但从原始强度到可信的生物学结论，中间还隔着三道统计关卡——背景校正（扣除点周围环形区域估计的非特异荧光）、归一化(normalization)（消除通道间与芯片间的系统偏差）与差异检验（在噪声中辨别真实的表达变化）。本节依次讨论这三道关卡：4.4.1 剖析偏差的来源并引入对数变换；4.4.2 梳理归一化方法的谱系；4.4.3 建立差异检验的线性模型框架；4.4.4 处理大规模多重检验；4.4.5 与 4.4.6 讨论结果的规范呈现与可重复性。这些方法是为微阵列发展的，但其中的统计思想——方差稳定、信息借用、假发现率控制——同样适用于 4.7 节的 RNA 测序数据。

4.4.1 系统偏差的来源与对数变换

微阵列测量的是荧光强度，而非 mRNA 分子数；从前者到后者，隔着多层与生物学无关的系统性因素。双色点样阵列（4.2 节）上最著名的干扰是染料偏倚(dye bias)：Cy3 与 Cy5 两种花菁染料的标记效率、量子产额与光电倍增管响应并不相同，经验上 Cy5 通道往往整体偏高，使未经校正的数据中所有基因的比值系统性地偏离 1。第二类是打印针效应(print-tip effect)：点样阵列由排成栅格的多个打印头逐点完成，不同针尖的点样量与点形态存在差异，同一打印头所印的点会呈现整体偏移。第三类是批次效应(batch effect)：单色平台的比较在「芯片之间」进行，杂交日期、试剂批次、洗涤条件与扫描仪增益的漂移都会造成整片偏移；在多芯片实验中，批次效应的解释方差有时甚至超过生物学因素本身——这正是 4.3 节强调随机化与区组化的原因。第四类偏差最隐蔽：偏倚的大小往往不是常数，而随信号强度非线性变化，低强度端的比值偏倚最为严重（背景扣除后信噪比恶化、两种染料的检测下限不同）。

这些偏差叠加在一个不利的统计性质之上：荧光强度的均值—方差关系。原始强度的标准差近似随均值增大（约呈 SD ∝ μ 的常数变异系数关系），高强度基因的绝对噪声远大于低强度基因；若直接在原始刻度上检验，高强度基因的微小比例变化会被误判为显著，而低强度基因的真实变化会被淹没。对数变换同时缓解两个问题：其一，由 delta 方法可知 SD ∝ μ 的数据在 log 刻度上标准差近似为常数（方差稳定）；其二，对比值取对数使上调与下调对称——2 倍上调是 +1，2 倍下调是 −1，而原始比值分别是 2 与 0.5，不对称且难以直接平均。对双色阵列，习惯上把每个点的红绿通道强度 R 与 G 转换为一对坐标：

4.4.2 归一化方法的谱系

归一化的共同前提是一个近乎哲学的假设：绝大多数基因不发生差异表达。在此假设下，点云整体偏离零线的部分应归咎于技术偏差而非生物学，可以安全扣除。最简单的实现是全局中位数归一化：把一张芯片（或一个通道）内全部 M 值的中位数平移到 0。它隐含着「偏倚与强度无关」的假定——一个全局常数加性偏移。上一小节已经说明这一假定常常不成立：染料偏倚在低强度端最强，随 A 增大而衰减，是强度依赖的非线性偏倚。此时全局平移会校正过度与校正不足并存：低强度端仍残留正偏，高强度端反被压成负偏。

校正强度依赖偏倚的标准做法是局部加权回归(locally weighted scatterplot smoothing, LOWESS/loess)：在 M-A 散点图上，以每个点为中心取一个滑动窗口，对窗口内的点做加权多项式回归（离中心越近权重越大，离群点降权），得到拟合曲线 M̂(A)；逐点计算残差 M′ = M − M̂(A)，即完成校正（Yang et al., 2002）。LOWESS 不预设曲线形状，因而能追踪任意形态的强度依赖偏倚。对点样阵列可按打印头分组分别拟合（print-tip LOWESS），把打印针效应一并吸收；片间还可再做一次尺度缩放，使各张芯片的 M 分布具有一致的离散程度。图 4.4-1 给出校正前后的对照。

与 LOWESS「校正均值偏倚」的思路互补的是方差稳定变换(variance stabilizing transformation)一族的 VSN 方法：它以参数化模型同时描述均值—方差关系并施加变换（如双曲反正弦 arcsinh），在低强度端行为类似 log、在高强度端趋近线性，把校正与方差稳定一步完成。单色寡核苷酸阵列（Affymetrix 平台）则发展出另一套流水线，其比较发生在芯片之间，归一化的对象是各芯片的整个强度分布：

4.4.3 差异表达的统计框架

归一化之后的表达矩阵（行 = 基因，列 = 样本）进入检验环节。最朴素的办法是对每个基因分别做双样本 t 检验：以两组均值之差除以合并标准误。它的脆弱之处在于方差估计不稳：每组只有 3 个生物学重复时，样本方差只有 4 个自由度，抽样波动很大；偶尔出现的极小样本方差会制造出巨大的 t 值，把噪声基因推上显著列表。换言之，小样本下不仅功效低，假阳性也集中出现在「碰巧方差小」的基因上——这些基因的生物学特征（如低表达、探针质量差）与显著性混杂，形成系统性偏差。

现代做法是把一次实验的全部因素写进一个线性模型统一估计：对基因 g，样本 j 的表达值 y_gj = x_j^Tβ_g + ε_gj，其中 x_j 是设计矩阵的第 j 行——分组、染料、批次、个体等因素都以 0/1 协变量的形式编码在内，关心的对比（如处理 − 对照）表现为参数向量 β_g 的一个线性组合（contrast）。与逐对 t 检验相比，线性模型能在一个模型内消化复杂设计（配对、区组、时间进程、2×2 因子），并直接检验交互项。经典的方差分析（ANOVA）分解是同一思想的加性形式：把表达值拆解为基因效应、样本（处理）效应、染料效应、阵列效应与基因×处理交互项之和，其中基因×处理交互项正是「这个基因是否响应处理」的度量；批次与个体作为额外因子进入分解后，其解释方差被显式剥离，不再污染差异估计（Gibson & Muse, 2009）。

线性模型解决了均值部分的参数化，方差部分的不稳则由经验贝叶斯收缩(empirical Bayes shrinkage)解决——这是 limma 框架的核心（Smyth, 2004）。其直觉是：一万个基因各有一份带 4 个自由度的方差估计，但它们共享同一个「方差总体」；可以把全基因组估计出的方差分布作为先验，将每个基因自己的方差向这个先验「收缩」，得到更稳定的后验方差，再构造检验统计量：

公式 (4.4-2) 的代价与收益都值得强调。收益是每个基因的方差都不再孤军奋战——哪怕只有两三个重复，只要全基因组提供几千份方差信息，检验就相当稳健，这使 limma 成为微阵列时代事实上的标准。代价是引入了「基因间方差同属一个分布族」的假设：若真实存在一小簇方差天然极小的基因（如高度稳定的看家基因），收缩会略微高估它们的方差、低估其显著性——一个可以接受、但应当知道的偏倚方向。同一时期发展的 SAM 方法以「接近零的阈值 + 重复重排」控制假发现率，是解决多重检验问题的另一条著名路径（Tusher et al., 2001）；调节 t 与 SAM 的比较研究表明两者在多数场景下结论相近，前者在建模灵活性上更强。

4.4.4 多重检验与假发现率

做完一万个基因的检验，接下来面对的问题是：如何读出这张 p 值清单。单基因时代 p < 0.05 的惯例在基因组尺度上彻底失效。做一个简单的算术：若对 10⁴ 个基因各做一次 α = 0.05 的检验，且所有基因实际上都不差异，则每个基因以 5% 的概率越线，期望出现 10⁴ × 0.05 = 500 个「显著」基因——全部为假阳性。这 500 个噪声基因与真实的差异基因混在同一张列表里，裸 p 值无法区分它们。控制这一问题的第一代标准是族错误率(family-wise error rate, FWER)——「至少犯一次假阳性错误」的概率；邦费罗尼校正把每个检验的阈值收紧为 α/m，便可将 FWER 压在 α 以下。但 m = 10⁴ 意味着阈值是 5×10⁻⁶：只有效应极强、方差极小的基因才可能通过，大量真基因被一同牺牲。对以「筛出候选清单、后续逐一验证」为目的的探索性筛选，更合适的误差度量是假发现率(false discovery rate, FDR)（Benjamini & Hochberg, 1995）：

Benjamini–Hochberg（BH）步进程序给出了一条在任意 p 值组合下把 FDR 控制在 q 以内的实用规则：

BH 程序在所有原假设为真时把 FDR 严控在 q 以下，但当真信号占比较高时它会明显过保守——分母 m 中包含了大量不可能为真的检验。Storey（2002）的q 值(q-value)从数据中估计真原假设的比例 π₀ 再以 m·π₀ 替换 m，在保持 FDR 控制的同时提高功效；q 值本身定义为「把该基因及更显著基因全部收入清单时所能达到的最小 FDR」，可视为多重检验校正后的 p 值，实践中直接报告每个基因的 q 值并按 q ≤ 0.05 取清单即可。无论采用哪种变体，共同的报告纪律是：给 q 值而非裸 p 值，给效应量而非仅给显著性。

4.4.5 结果呈现：火山图、效应量与富集预告

检验完成后，把每个基因的效应大小（log₂ 折叠变化）作为横轴、统计显著性（−log₁₀ p 或 −log₁₀ q）作为纵轴作图，即得表达分析最常用的图形——火山图(volcano plot)。多数基因聚在底部山脚，真实差异基因冲向两上角。图上的两条阈值线（竖线 |log₂FC| = 1、横线 FDR 5%）把平面切成几个语义不同的区域，解读时应同时使用两个维度：

火山图的价值正在于它把两个维度同时摆在眼前：p 值回答「变化可靠吗」，折叠变化回答「变化重要吗」。统计上极显著而倍数只有 1.05 的基因，与倍数达 8 倍但 p = 0.2 的基因，都不能直接进入结论——前者的效应量可能低于检测平台的技术重复性，后者多半是组内方差过大或个别离群样本所致。由此可以整理出差异表达结果的报告规范：对每个报告的基因给出点估计（log₂FC）与置信区间（或后验标准误），显著性一律用校正后的 q 值；对全表给出检验总数、采用的校正方法与阈值；缩略语与预处理参数（如 LOWESS 的窗口宽度、RMA 的版本）随文给出，使他人可以复现整条链路。置信区间的宽度本身就是重要信息：一个 FC = 2.0、95% 置信区间 [1.2, 3.3] 的基因，其可信程度远低于区间 [1.8, 2.2] 的基因，尽管点估计相同。

最后需要对下游分析做出预告。单个基因的清单只是第一步：差异基因往往数十上百个，逐条罗列难以形成机制图景。标准做法是把清单交给功能富集分析——检验差异基因是否在某个功能类别（如基因本体 GO 的分子功能、生物过程、细胞组分）中比例过高。其统计骨架正是 2.5 节注释部分与概率论中的超几何检验(hypergeometric test)：从 N 个基因的「总体」中抽出 n 个差异基因，某类别（共 K 个成员）中被抽中 k 个以上的概率。富集分析将作为功能注释的统计延伸在 2.5 节基因组注释的基础上于后续章节（如 6.4 节的多组学整合）系统展开；此处只需记住：富集分析同样面对多重检验（数千个类别），FDR 框架在这里再度登场。

4.4.6 可重复性实践

微阵列分析链条长、可调参数多，同一份原始数据经不同预处理可能得出不同清单——这一观察在 2000 年代引发了领域内对可重复性的系统反思，也留下了至今适用的实践规范。第一层是分析与设计的预先声明：在看到数据之前确定分组定义、比较组合、将使用的归一化方法与检验框架，即「预注册」或分析计划存档；事后挑选「效果最好的组合」会以不可量化的方式抬高假阳性。第二层是代码与数据的存档：原始图像、强度矩阵与处理脚本一并存入公共数据库（GEO、ArrayExpress），并以 MIAME（关于微阵列实验的最小信息量）标准记录样本、处理与阵列的元数据；只有原始数据可及，他人才可能用不同假设重走整条链路。第三层是面向读者的诊断图自查，其中最有用的是 p 值分布：在「多数基因无差异」的前提下，未校正 p 值的直方图应大体均匀分布在 [0, 1] 上、并在接近 0 的一端因真信号而隆起；若分布整体右偏（大量小 p 值堆积在 0 附近之外）或呈异常驼峰，提示归一化不足、批次效应或离群样本。相应的分位数—分位数图（Q-Q 图）把观测 p 值的分位数与均匀分布的期望分位数对照，理想情形是沿对角线的直线在右端抬升——偏离形态可据此定位问题环节。

第四层是独立验证。统计框架再精巧，也只是「与模型一致」，不等于「与生物学一致」。惯例的做法分两级：对筛选出的少量关键基因，用定量 PCR（4.1 节的 ΔΔCt 方法）在同批样本上验证方向与幅度；对整个结论（如一个表达签名能否区分表型），在独立队列上检验——不同批次、不同平台乃至不同实验室的样本。验证通过的关键基因进入后续机制研究，验证不通过则回溯诊断图与原始数据。需要指出，「用同一份微阵列数据再算一遍」不构成独立验证：共同的系统偏差会被一同复制。真正的独立要求测量过程与样本来源至少其一重新开始。本节的统计链路至此闭环：设计（4.3 节）限定噪声结构，归一化剥离技术偏差，线性模型与收缩稳定推断，FDR 校正清单污染，而独立验证为整个链条提供外部锚点。