3. 原始数据预处理

测序仪下机的那一瞬间，你拿到的 FASTQ 文件其实挺「脏」的。里面不光有我们想要的 cDNA 片段，还混着接头（Adapter）、测序引物，以及一堆越往后读质量越差的碱基。我早年刚跑 RNA-seq 时，觉得这步无所谓，直接拿 raw data 去比对——结果比对率连 60% 都不到，差异分析更是乱七八糟。后来才明白，预处理这步没做干净，后面所有结论都得打问号。

这一章，我就把 Trimmomatic 和 Trim Galore 这两个主流工具怎么选、怎么配、怎么避坑，一次性给你讲透。

核心目标：从 raw reads 中切除接头序列、修剪低质量碱基、过滤过短 reads，最终得到可用于下游比对和定量分析的 clean reads。
原始数据预处理核心流程 原始数据预处理流程 —— 从 raw reads 到 clean reads Raw FASTQ Trimmomatic（PE / SE） Trim Galore（自动检测接头） 切除接头 · 修剪低质量 · 过滤短 reads Clean FASTQ（QC 达标） ⚠️ 关键质控指标 Q30 占比 > 85% 接头残留 < 0.1%

接头污染是怎么来的？为什么非切不可？

RNA-seq 建库时，我们需要把 cDNA 片段两头接上测序接头，才能上机簇生成和测序。但片段长度分布是随机的——有些片段比读长短，测序仪一读就「越界」了，把接头序列也读进去。这就是所谓的 接头污染。

你想想看，如果比对时 reads 尾巴上挂着一截接头，比对软件会怎么处理？要么比对不上，要么比对到基因组上错误的位置。我见过一个极端案例：某批数据 FastQC 报告里接头污染率高达 12%，没处理就直接跑比对，结果基因定量结果整整偏了一倍。所以，接头切除是预处理的第一个硬性门槛。

⚠️ 千万不要凭肉眼判断「我的数据好像没接头」。 很多接头序列在测序仪上已经被部分读取，肉眼看不出来，但 FastQC 的 Overrepresented sequences 模块可以检测。我每次拿到新数据，第一件事就是跑 FastQC，看接头污染的比例。

低质量序列为什么必须修剪？

测序质量从 5' 到 3' 是逐渐下降的——这是 Illumina 平台的物理特性决定的。3' 端的碱基越往后，错误率越高。如果你不修剪这些低质量尾巴，它们会直接导致：

  • 假阳性 SNP / 突变：低质量碱基被错误识别为变异位点。
  • 比对率下降：低质量区域比对软件无法匹配。
  • 定量偏差：低质量 reads 被错误分配到多个基因上。

我习惯用 滑动窗口修剪策略——从 5' 端开始，每次取一个窗口（比如 4 bp），计算平均质量，如果低于阈值就从窗口位置直接切掉。这样既保住了高质量的前段，又能干净地切掉后段「垃圾」。

Trimmomatic：工业级的老牌工具

Trimmomatic 是我用得最多的工具。它支持双端（PE）和单端（SE），参数设计很工程化，适合批量处理。来看看它的核心参数：

java -jar trimmomatic-0.39.jar PE \
  -threads 8 \
  -phred33 \
  input_R1.fastq.gz input_R2.fastq.gz \
  output_R1_paired.fastq.gz output_R1_unpaired.fastq.gz \
  output_R2_paired.fastq.gz output_R2_unpaired.fastq.gz \
  ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 \
  LEADING:3 \
  TRAILING:3 \
  SLIDINGWINDOW:4:15 \
  MINLEN:36

我把这个命令拆解开，一个个参数说清楚：

参数含义我常用的值
ILLUMINACLIP切除接头：fa文件:允许错配数:比对分值阈值:palindrome分值2:30:10
LEADING切除开头的低质量碱基（低于该值直接切）3
TRAILING切除结尾的低质量碱基3
SLIDINGWINDOW滑动窗口修剪：窗口大小:平均质量阈值4:15
MINLEN过滤掉修剪后长度小于该值的 reads36
💡 palindrome 模式是什么？ 双端测序时，如果两端的 reads 在接头处发生 overlap，Trimmomatic 会用 palindrome 模式检测并对齐，直接切除整个接头 dimer。我遇到过一批建库失败的样本，接头 dimer 占比超过 30%，全靠这个模式才救回来。

Trim Galore：新手友好，自动检测接头

Trim Galore 本质上是 Cutadapt + FastQC 的封装。它最大的优点是能自动检测接头类型——你甚至不用提供接头序列文件。它会在运行前先扫描一小部分 reads，推测接头序列，然后再全量切除。

trim_galore --paired \
  --quality 20 \
  --phred33 \
  --length 36 \
  --stringency 3 \
  --cores 8 \
  --output_dir ./trimmed \
  sample_R1.fastq.gz sample_R2.fastq.gz

参数解读：

  • --quality 20：切除质量低于 Q20 的尾部碱基（相当于 Trimmomatic 的 TRAILING）。
  • --length 36：过滤掉短于 36 bp 的 reads。
  • --stringency 3：接头匹配的严格度——值越高越保守。我一般用 3，遇到接头残留多的样本会降到 1。
  • --cores：多核加速，我通常开 8 核。
⚠️ Trim Galore 的自动检测不是万能的。 如果建库用了比较冷门的接头，或者接头被严重降解，自动检测可能会漏掉。我建议每次运行后都检查一下 Trim Galore 的 report 文件，看看检测到了什么接头序列，确认是否和建库试剂盒一致。

Trimmomatic vs Trim Galore：怎么选？

对比维度TrimmomaticTrim Galore
接头检测方式需手动提供接头 FASTA自动检测 + 手动指定均可
滑动窗口修剪原生支持，参数精确不支持（仅尾部修剪）
双端 pairing 保留强（输出 paired + unpaired）强（自动同步两端）
运行速度略快（Java 多线程）中等（单线程 + 自动检测）
适合场景批量生产、参数精细控制快速上手、自动化流程

我个人习惯是：建库试剂盒固定的项目用 Trimmomatic，因为参数可精确调优；公共数据或来源不明的数据用 Trim Galore，让它自动识别接头，省去手动查接头序列的麻烦。

避坑指南：我踩过的 3 个坑

  1. 接头序列文件版本不对。 Trimmomatic 自带的 TruSeq3-PE.fa 对应的是老版 TruSeq 试剂盒。如果建库用了 TruSeq RNA CD 或其它改版，接头序列可能不同。我有一批数据用了 TruSeq3-PE.fa，结果残留率一直降不下来——后来换成 TruSeq2-PE.fa 才解决。所以，每次都要确认试剂盒版本。
  2. 双端 reads 不同步。 处理双端数据时，如果一端的 read 被过滤掉了，另一端也要同步丢弃。Trimmomatic 的 PE 模式和 Trim Galore 的 --paired 都会自动处理，但如果你自己写脚本切接头，很容易搞乱 pairing 顺序。我建议直接用成熟工具，别自己造轮子。
  3. 过度修剪导致 reads 太短。 有一回我把 SLIDINGWINDOW:4:20 设得太严，结果 40% 的 reads 被切到 36 bp 以下丢掉了。后来我改成 4:15，保留了更多有效数据。阈值不是越高越好，要根据数据质量灵活调整。

预处理后的质检：你怎么知道做干净了？

跑完预处理，一定要用 FastQC 再检查一遍。我重点关注三个指标：

  • Per base sequence quality：修剪后大部分区域应该保持在绿色区（Q30 以上）。
  • Adapter Content：接头含量应接近 0%。如果还有 0.5% 以上的残留，说明切除参数太松。
  • Sequence Length Distribution：长度分布应该是一个光滑的峰形，没有异常短 read 的拖尾。
一个经验值：双端 150 bp 的 reads，经过合理预处理后，paired reads 保留率应该在 85% - 95%。如果低于 80%，说明要么数据质量太差，要么参数设得太严。我有一批来自降解样本的数据，保留率只有 62%——但那是样本问题，不是工具问题。

预处理这步，说白了就是「脏活累活」，但也是整个 RNA-seq 流程里最值得花时间打磨的环节。你把这一步做扎实了，后面的比对、定量、差异分析都会顺很多。工具本身不复杂，复杂的是根据你的数据特点去调参数——而这些经验，只能靠一次次跑完质检报告、再回头调整参数来积累。

← 上一章 📖 返回目录 下一章 →