差异表达分析概述
差异表达分析，说白了就是找出两组样本之间哪些基因的表达量不一样。你在RNA-seq项目里花大价钱做测序，核心目标往往就是这个。我个人习惯把RNA-seq分析分成两大块：数据预处理和统计推断。今天这章先跟大家聊聊整个流程，看看差异表达分析站在哪个位置，以及它到底能干啥。
RNA-seq数据分析全流程
先拆解一下完整的RNA-seq分析管道。我刚开始做这行的时候，总喜欢一股脑把所有数据丢进软件跑，结果经常踩坑。后来才明白，每一步都有自己的规矩，跳不得。
下面这张图是我认为最精简也最实用的流程框架。你们先有个整体印象，后面每一章我会把细节掰开揉碎了讲。
这张图里每个模块我都用不同颜色区分了。你会发现差异表达分析不是终点——它产出的基因列表通常要送到下游做可视化或者功能注释，不然很难讲出生物学故事。反过来，如果上游实验设计有缺陷，后面再怎么调参数也救不回来。我在项目里吃过这个亏：曾经有一批样本，批次效应没处理好，差异基因一筛几百个，结果全是假阳性。嗯，这个教训让我养成了先看实验设计再跑分析的习惯。
差异表达分析的目的
差异表达分析的目标很明确——在控制假阳性的前提下，识别出不同实验条件下表达量有统计显著差异的基因。你想想看，一个基因在对照组表达量是100，处理组飙升到500，这很可能意味着它在响应某种刺激。但问题来了：生物学重复之间的波动怎么处理？测序技术本身也有噪音。所以我们需要一个统计模型来帮忙判断，到底是真差异还是随机波动。
核心输出就是一张差异基因列表，包含每个基因的log₂ Fold Change和校正后的p值。我通常用 padj < 0.05 和 |log₂FC| > 1 作为筛选阈值，但这不是铁律。曾经我在一个低表达量的项目中，把阈值放宽到 padj < 0.1 才筛出足够多的候选基因做验证，结果实验验证了70%的阳性率——所以阈值这东西，要结合你的实验背景来定。
关键指标：差异表达分析本质上是假设检验的集合。对每个基因做一次检验，然后用BH方法（Benjamini-Hochberg）控制FDR。
为什么要控制多重假设检验？假如你一次性检验2万个基因，就算所有基因都没差异，按0.05的显著性水平也能找出约1000个“显著”的假阳性。这就是为什么我们不用原始p值，而用校正后的padj。
应用场景
差异表达分析的应用场景比你想象的要广。我把它归纳成几类，你们看看自己手头的项目属于哪一类。
| 场景类型 | 对比方式 | 典型例子 |
|---|---|---|
| 疾病 vs 正常 | 两组比较 | 肿瘤组织 vs 癌旁正常组织 |
| 处理 vs 对照 | 两组比较 | 药物处理细胞 vs 溶剂对照 |
| 时间序列 | 多组比较 | 发育过程中多个时间点 |
| 基因敲除 vs 野生型 | 两组比较 | CRISPR敲除某基因后的转录变化 |
| 条件多因素 | 复杂设计 | 处理×时间交互作用 |
我当时做的第一个RNA-seq项目是乳腺癌亚型比较。三阴性vs Luminal B，每组只有3个生物学重复，样本量很小。我一开始准备用edgeR，但后来发现数据里有一个样本的文库大小异常，质控阶段差点被直接剔除。我建议你们在分析前一定要做样本间的相关性热图，它能帮你发现潜在的离群样本。
我的小技巧：拿到原始表达矩阵后，先做一个简单的PCA或者MDS图。如果处理组和对照组的样本在图上能明显分开，那差异分析的结果通常会比较靠谱；如果混杂在一起，即使软件跑出几百个显著基因，也要打个问号。
差异分析工具的“三驾马车”
现在主流工具就三个：DESeq2、edgeR 和 limma-voom。其中DESeq2和edgeR都基于负二项分布假设，对文库大小和基因长度做了不同的处理方式。它们俩的输出格式和思路很接近，但细节上各有优劣。我从一个实际项目里对比过这两个工具，同样的数据DESeq2筛出了153个差异基因，edgeR筛出了178个，重叠部分有124个。那28个不重叠的基因到底信谁？后面我会专门开一章讲这个。
下面是两类工具在底层模型上的差异对照表，供你们快速参考：
| 特性 | DESeq2 | edgeR |
|---|---|---|
| 统计模型 | 负二项 + 拟合方差 | 负二项 + 经验贝叶斯 |
| 归一化方法 | 中位数比值（MED） | TMM（加权截尾均值） |
| 离散度估计 | 先验收缩 | 稳健离散度 |
| 处理小样本 | 较好（有收缩） | 需要谨慎 |
| 适用数据 | 有生物学重复 | 可处理无重复（需设置BCV） |
你会注意到DESeq2对离散度做了“先验收缩”，这在生物学重复少的时候特别有用。我在做小鼠组织样本分析时，每组只有2个重复，用edgeR跑出来结果方差很大，换成DESeq2之后的基因列表重复性就稳多了。但这不是说edgeR不好，它处理大样本量、多分组时速度更快，而且TMM归一化在某些场景下比DESeq2的中位数比值更鲁棒。具体怎么选，我会在后面用例数据带着你们跑一遍。
重要提醒：这两个工具要求输入原始counts（整数），千万不能把FPKM/RPKM或TPM直接丢进去！它们内部会自动做文库大小的归一化，用标准化后的表达值反而会破坏统计假设。我见过有人把FPKM当输入，跑完差异分析几乎找不到显著基因——因为归一化逻辑已经被破坏掉了。
在实际项目中，我往往会用两个工具交叉验证。如果只有个别基因在DESeq2和edgeR之间有冲突，我会去查原始counts和基因自身的表达量分布，再做决定。但如果两个工具的结果大相径庭，那我得先从数据质控找原因——批次效应？异常样本？基因过滤阈值不对？这些才是根子上的问题。
好了，这一章主要帮大家建立起完整的流程框架。记住一句话：差异分析不是点点鼠标就完事的，它需要你对上游数据和下游生物学背景都有手感。后面每一章我都会结合自己踩过的坑，手把手带你们用R跑一遍DESeq2和edgeR，然后把它们的结果画在一起对比。