第四章 接头与低质量序列去除：Trimmomatic、Cutadapt 参数配置与实战

数据清洗这件事，看似简单，其实门道不少。我见过不少新手拿着 fastq 直接跑组装，结果出来的 contig 全是 N 和错误碱基。原因呢？多半是接头没切干净，或者低质量尾巴没去掉。

这一章，咱们就聊透两个主力工具——Trimmomatic 和 Cutadapt。我会从参数原理讲到实战坑点，最后给你一个我自己常用的配置模版。

4.1 为什么要切接头？为什么要去低质量？

宏基因组测序的读长通常只有 150bp，接头序列却占了几十个碱基。你想想看，如果接头不被切除，后续的比对和组装会怎样？接头序列会错误地比对到参考基因组上，导致假阳性的 SNPs。更麻烦的是，两个读长如果都有接头尾巴，拼接时会把两个本来不相关的片段连在一起——这就是嵌合体。

低质量碱基的影响更直接：质量值 Q20 以下的碱基错误率已经大于 1%，如果这些位置出现在 reads 的 3′ 端，会影响 k‑mer 的准确性。我个人习惯把质量阈值设在 Q20 以上，至少把尾巴剪掉。

一点经验：不要盲目相信默认参数。Trimmomatic 的 LEADING 默认值是 10，但很多人换成 15 效果反而更好。你得根据自己数据的质量分布来调。

4.2 Trimmomatic：参数配置与实战

Trimmomatic 是我最早接触的清洗工具，优点在多线程和滑动窗口去低质量。它的参数设计很直观：告诉它“从哪一端切、用多大窗口、平均质量阈值多少”。

4.2.1 常用参数速查表

参数 作用 典型值 我的建议
ILLUMINACLIP 切除接头及 adapter 序列 adapters.fa:2:30:10 保守用 2:30:10，严格用 2:40:15
LEADING 从 5′ 端切除低质量碱基 10 推荐 15，尤其对宏基因组数据
TRAILING 从 3′ 端切除低质量碱基 10 与 LEADING 保持一致，设 15
SLIDINGWINDOW 滑动窗口内平均质量低于阈值则切除 4:15 4:20 更严格，5:15 更温和
MINLEN 切除后 reads 最短长度，低于则丢弃 50 我一般用 60，太短没法用
TOPHRED33 强制转换成 Phred33 质量体系 无默认 建议加上，避免版本不兼容

4.2.2 单端与双端数据命令对比

双端数据要特别注意 PAIRED 输出模式：它会同时保留成对的 reads，丢掉没有配对的单端。我把两组典型命令写在下面，你直接复制改路径就行。

# 单端数据
java -jar trimmomatic-0.39.jar SE \
  -phred33 \
  sample_R1.fastq.gz sample_cleaned.fastq.gz \
  ILLUMINACLIP:adapters.fa:2:30:10 \
  LEADING:15 TRAILING:15 \
  SLIDINGWINDOW:4:20 MINLEN:60

# 双端数据
java -jar trimmomatic-0.39.jar PE \
  -phred33 \
  sample_R1.fastq.gz sample_R2.fastq.gz \
  sample_R1_paired.fastq.gz sample_R1_unpaired.fastq.gz \
  sample_R2_paired.fastq.gz sample_R2_unpaired.fastq.gz \
  ILLUMINACLIP:adapters.fa:2:30:10 \
  LEADING:15 TRAILING:15 \
  SLIDINGWINDOW:4:20 MINLEN:60
我曾经踩过的坑：ILLUMINACLIP 的第三个数字（palindromeClipThreshold）设得太低会误删大量 reads。有一次我把 30 改成了 20，结果 30% 的 reads 被判定为 adapter 二聚体切掉了……数据直接废了一半。所以保守起见，我不建议低于 30。

4.3 Cutadapt：灵活但需要手动处理双端

Cutadapt 的优势在于接头切除的灵活性——支持自定义 adapter 序列、支持多接头同时匹配、支持 5′ 和 3′ 端分别指定。但它没有内置的滑动窗口去低质量功能，需要配合其他工具（比如 seqtk）来实现。

4.3.1 核心参数说明

参数 作用 典型值 我的建议
-a 切除 3′ 端 adapter AGATCGGAAGAGC 直接拿文库试剂盒给出的序列
-g 切除 5′ 端 adapter CTGTCTCTTATA 一般用于双端，但在宏基因组中很少用
-e 允许的错配比例 0.1 0.1 已经很宽了，别超过 0.2
--minimum-length 切除后 reads 最短长度 50 与 Trimmomatic 保持一致，用 60
-q 从 3′ 端切除低质量碱基（类似 TRAILING） 15 我一般设 20，要求高一点
--nextseq-trim NextSeq 平台专用：切除 polyG 尾巴 20 如果有 NextSeq 数据必须加

4.3.2 实战命令示例

# 单端数据：切除 3′ adapter 并低质量尾部裁切
cutadapt -a AGATCGGAAGAGC \
  -q 20 \
  --minimum-length 60 \
  -o sample_cleaned.fastq.gz \
  sample.fastq.gz

# 双端数据：两端分别切除 adapter
cutadapt -a AGATCGGAAGAGC -A AGATCGGAAGAGC \
  -q 20 \
  --minimum-length 60 \
  -o sample_R1_cleaned.fastq.gz \
  -p sample_R2_cleaned.fastq.gz \
  sample_R1.fastq.gz sample_R2.fastq.gz

# 双端 + NextSeq polyG 去除（我常用的一套）
cutadapt -a AGATCGGAAGAGC -A AGATCGGAAGAGC \
  -q 20 \
  --nextseq-trim 20 \
  --minimum-length 60 \
  -o R1_clean.fastq.gz -p R2_clean.fastq.gz \
  R1_raw.fastq.gz R2_raw.fastq.gz
小窍门：Cutadapt 支持 --report 参数输出清洗报告，里面详细写了每种 adapter 匹配了多少条 reads，错配率是多少。我每次跑完都会看一眼，如果 adapter 匹配率超过 5%，说明文库制备可能有问题。

4.4 Trimmomatic vs Cutadapt：我该选哪个？

这是个经典问题。我自己的答案是：看数据特点。如果你手上是常规的 Illumina 双端 150bp，我推荐 Trimmomatic，因为它一次性搞定接头切除、滑动窗口去低质量和长度过滤，一个命令解决问题。如果你的数据来自不同平台（比如既有 HiSeq 又有 NextSeq），或者文库的 adapter 序列比较特殊，Cutadapt 的灵活性更占优势。

另外，在高通量场景下（比如几百个样本），Trimmomatic 的多线程支持更好；Cutadapt 虽然也支持多线程，但内存管理不如 Trimmomatic 稳定。我曾在 48 核的机器上用 Cutadapt 同时跑 20 个样本，有几次直接 OOM 了，换成 Trimmomatic 就没问题。

我的日常工作流：
  • 首选 Trimmomatic：双端 + 滑动窗口 + MINLEN=60
  • 次选 Cutadapt：当 adapter 序列不标准，或需要 NextSeq polyG 修剪时
  • 两者结合：先用 Cutadapt 切除特殊 adapter，再用 Trimmomatic 做质量过滤
接头与低质量去除流程决策图 原始 FASTQ 数据 adapter 标准? Trimmomatic ILLUMINACLIP + SW Cutadapt 自定义 adapter 质量过滤 & 长度筛选 LEADING/TRAILING + MINLEN 干净 FASTQ，可接组装 SW: 滑动窗口过滤 minlen: 推荐60bp

4.5 实战：完整参数配置模版

这里我把自己在项目里用了一年多的配置模版贴出来。它经过了 500+ 个宏基因组样本的检验，出问题的概率很低。

# Trimmomatic 双端清洗（推荐优先使用）
java -jar trimmomatic-0.39.jar PE \
  -phred33 -threads 16 \
  raw_R1.fastq.gz raw_R2.fastq.gz \
  clean_R1_paired.fastq.gz clean_R1_unpaired.fastq.gz \
  clean_R2_paired.fastq.gz clean_R2_unpaired.fastq.gz \
  ILLUMINACLIP:/path/to/TruSeq3-PE.fa:2:30:10 \
  LEADING:15 \
  TRAILING:15 \
  SLIDINGWINDOW:4:20 \
  MINLEN:60

# 如果你必须用 Cutadapt，我的版本如下
cutadapt -a AGATCGGAAGAGC -A AGATCGGAAGAGC \
  -q 20 --minimum-length 60 \
  --max-n 0.05 \
  --nextseq-trim 20 \
  -j 16 \
  -o clean_R1.fastq.gz -p clean_R2.fastq.gz \
  raw_R1.fastq.gz raw_R2.fastq.gz
避坑指南：我曾经犯过一个低级错误——在双端清洗时忘了指定 -p（Trimmomatic）或 -p（Cutadapt）参数，结果两条 reads 各跑各的，配对被破坏。那一次跑完后比对率直接掉了 20%。所以每次写完命令，我一定先检查参数是否成对出现。

4.6 清洗质量如何评估？

你切完了，怎么知道效果好不好？我习惯看两个指标：

  1. 保留率：清洗后剩下的 paired reads 数量 / 原始数量。一般 80% 以上就算正常。如果低于 70%，可能是参数太严格，或者数据本身质量太差。
  2. FastQC 前后对比：我每次清洗完都会跑一次 FastQC，重点看 per‑base sequence quality 和 adapter content 两个图。如果 adapter content 还有残留，说明 ILLUMINACLIP 或 Cutadapt 的参数需要调严格一些。
快速检查命令： 
# 统计清洗前后 reads 数
echo "原始 reads: $(zcat raw_R1.fastq.gz | wc -l)/4"
echo "清洗后 reads: $(zcat clean_R1_paired.fastq.gz | wc -l)/4"

# 保留率计算（Bash 一行）
echo "scale=2; $(zcat clean_R1_paired.fastq.gz | wc -l)/4 / $(zcat raw_R1.fastq.gz | wc -l)/4 * 100" | bc
← 上一章 📖 返回目录 下一章 →