第二章 RNA-seq 差异表达分析简介
RNA-seq 的常规流程
说回 RNA-seq 分析，大家基本都走一条差不多的路。我个人习惯把整个流程拆成几个模块，这样心里有数。
从拿到测序仪的 raw data 开始，到最终拿到有生物学意义的基因列表，中间大概要过这么几关：
- 实验设计与建库：生物学重复的数目、测序深度这些，一开始就得想清楚。我见过不少项目，重复数不够，后期再怎么补统计方法也救不回来。
- 数据质控与预处理：主要是 FastQC 看质量，Trim Galore 去接头和低质量碱基。这步我建议别偷懒，脏数据进去，脏结果出来。
- 比对或定量：传统做法用 STAR 或 HISAT2 比对到参考基因组，现在流行直接做 quasi-mapping，比如 Salmon 或 kallisto。速度快很多，准确率其实也够用。
- 基因水平的表达定量：得到每个基因的 read count。这里有个坑——有些工具比如 featureCounts 要求你先有比对好的 BAM 文件，而 Salmon 这类工具可以直接输出 count。
- 差异表达分析：这就是咱们这节课的核心了。DESeq2 和 edgeR 都是在这个环节上场。
- 功能注释与富集分析：拿到差异基因列表，丢到 GO 或 KEGG 数据库里看看，哪些通路被激活了。
流程大致是这个样子。嗯，这里要注意——差异表达分析是整个流程的“脊梁骨”。前面所有步骤都是在为这个环节准备数据，后面所有生物学解释都依赖差异基因的列表质量。
为什么差异表达分析是核心？
说白了，一个 RNA-seq 项目花了大价钱，你总得找出几个「菜」来吧？
差异表达分析要回答的问题很简单：在两种或多种条件之间，哪些基因的表达水平发生了有统计学意义的变化？
但操作起来并不简单。你需要同时处理三个问题：
- 技术噪音：同一个样本跑两遍，你可能会得到略有不同的 count 数。这不是错误，而是测序本身固有的波动。
- 生物学变异：同一种处理下的不同个体，基因表达本身就有差异。这一点经常被忽视。
- 多重假设检验：一次分析通常检测 20000 个左右的基因。如果显著性水平设在 0.05，那么光凭随机因素就能找到 1000 个「显著」基因。必须校正。
我记得刚接触 RNA-seq 那会儿，有个项目对照组和实验组各只有 2 个生物学重复。用 limma 做出来的结果看起来很美，但后来补了 3 个重复之后，之前所谓的「显著基因」有一大半都消失了。从那以后，我对重复数和技术噪音的关系就特别敏感。
差异表达分析的核心目标：在控制假阳性的前提下，尽可能多地识别出真实差异表达的基因。这不是简单的「做个 t 检验就能交差」的事。
RNA-seq 数据的统计特性
为什么不能用普通的 t 检验？原因很直接——RNA-seq 的 read count 数据不符合正态分布。
你看一下原始 count 的分布就知道了：
- 大部分基因的 count 数很小，甚至为零
- 少数高表达基因的 count 数很高，形成长尾
- 方差随着均值增大而增大，且不是简单的正比关系
这种数据更适合用 负二项分布（Negative Binomial, NB）来建模。DESeq2 和 edgeR 都采用了这个框架。它们之间的差异主要体现在：
| 特性 | DESeq2 | edgeR |
|---|---|---|
| 方差估计策略 | 通过所有基因的信息来收缩单个基因的方差 | 使用经验贝叶斯方法估计分散度 |
| 标准化方式 | 利用几何均值的中位数比值 | TMM（加权截断均值）标准化 |
| 统计检验 | Wald 检验或 LRT | 精确检验（类似 Fisher 精确检验的 NB 版本） |
| 适用场景 | 小样本、有异常值时更稳健 | 大样本、速度快、早期经典方法 |
这两种方法我都用过。说实话，在大多数常规数据上，它们给出的结果高度一致。差异主要出现在样本量极小时，或者实验设计比较复杂（比如有多个批次效应需要校正）的时候。
差异表达分析的关键挑战
就算你选好了工具，还是有几个坑要注意。
多重假设检验校正：默认你用 BH 方法控制 FDR（False Discovery Rate）。FDR 设定在 0.05 意味着，你找到的差异基因里平均有 5% 是假的。如果你只看 p-value，不做校正，那基本就是在乱来。
我过去有个项目，客户自己用 Excel 按 p-value 筛了一遍基因，拿给我看的时候说「这些肯定都是真实差异的」。我一查原始数据，对照组 3 个样本里有一个样本的某个基因表达量是 0，其他 2 个样本都是正常水平，结果那个基因就被筛出来了。嗯，这其实是缺乏生物学重复和稳健性检验的典型案例。
实战建议：在运行 DESeq2 或 edgeR 前，先画一张 MA 图看看数据分布是否正常。如果散点在低表达区大量偏离基线，可能意味着数据质量或标准化有问题。我每次都会在正式跑差异分析前做这一步，能堵住不少坑。
还有一点——批次效应。RNA-seq 实验如果分多批做，不同批次之间可能存在系统性偏移。这种偏移如果不处理，会被差异表达分析误认为是处理效应。DESeq2 和 edgeR 都可以把批次信息作为协变量纳入模型，但前提是你得在实验设计阶段就保存好批次信息。
最后说下可视化。差异分析的结果光靠表格是不行的。火山图、热图、MA 图是标配。而这些图的呈现效果，恰恰取决于你选的差异表达方法输出的统计量——log2 Fold Change 和 adjusted p-value。DESeq2 的收缩效应（apeglm）处理后的 log2FC 会更美观，不会出现极端的飘移点。edgeR 默认不做收缩，做火山图时偶尔会看到一些高表达基因的 FC 值飞得很远。
本章知识结构
下面这张图帮你快速梳理本章的核心脉络。你仔细看，基本上就是 RNA-seq 从原始数据到差异基因列表的「流水线」。
上面这张图的核心逻辑就是：前面每一步的质量，最终都投射到差异表达分析的结果上。实验设计不好，数据预处理粗糙，再好的统计方法也救不了。