Cutadapt与Trimmomatic：去接头与数据清洗

fastp 适合标准 RNA-seq 质控场景，但当接头序列复杂（自定义 adapter、多重 adapter）或需要精细控制裁剪策略时，Cutadapt 和 Trimmomatic 是更专业的选择。本文从算法原理、命令行实操和输出解读三个维度对比这两个工具，帮助选型。

实测环境：Debian 13，Cutadapt 4.9，Trimmomatic 0.39。

1. 宿主环境与安装#

1
$ cat /etc/os-release | head -2
2
PRETTY_NAME="Debian GNU/Linux 13 (trixie)"
3
NAME="Debian GNU/Linux"
4

5
$ uname -m && free -h | head -2
6
x86_64
7
Mem:  7.7Gi
8

9
# 安装
10
mamba install -c bioconda cutadapt trimmomatic -y
11

12
$ cutadapt --version
13
4.9
14

15
$ trimmomatic -version
16
0.39

2. 原理对比——它们怎么找到接头#

2.1 Cutadapt：半全局比对 + 容错匹配#

Cutadapt 的核心算法是半全局比对（semi-global alignment）。给定一条 read 序列 $R$ 和接头序列 $A$ ，它在 $R$ 中寻找 $A$ 的最佳匹配位置：

$\text{score} = \max_{i,j} \sum_{k=1}^{|A|} s(R[i+k], A[j+k])$

其中 $s(a,b)$ 是碱基匹配得分（match=+1, mismatch=-1, gap=-2）。

关键：这个比对允许错配——默认允许接头长度 10% 的错误率。接头在测序过程中可能出现测序错误，Cutadapt 通过设定最大错误率来处理这个问题：

$\text{error\_rate} = \frac{\text{mismatches}}{\text{adapter\_length}} \leq 0.1$

2.2 Trimmomatic：滑动窗口 + k-mer 哈希#

Trimmomatic 用更简单直接的方法——种子匹配 + 滑动窗口延伸：

取接头前 16bp 作为种子（seed）
在 read 中找种子精确匹配的位置
从匹配位置向两端延伸，允许少量错配
用滑动窗口质量过滤（4bp 窗口平均 Q<20 就切）

Trimmomatic 的滑动窗口公式：

$Q_{window} = \frac{1}{k} \sum_{i=1}^{k} Q_i$

其中 $k$ 默认 = 4， $Q_i$ 是每个碱基的 Phred 质量分数。如果 $Q_{window} < 20$ ，从该位置切开，保留前面部分。

2.3 核心区别表#

特性	Cutadapt	Trimmomatic
接头匹配	半全局比对，允许错配	种子匹配，精确
质量过滤	需要额外参数（不内置窗口过滤）	内置滑动窗口 SLIDINGWINDOW
双端处理	自动同步配对	需要手动配对模式
速度	快（单遍扫描）	较慢（多步串联）
灵活性	极高（几乎可切任何位置）	中等（预设步骤有限）
输出	单文件 + 统计 JSON	多文件（配对/未配对）

3. Cutadapt 实操——精细到每个碱基#

3.1 基础去接头#

1
# 单端数据：去掉3'接头并过滤短reads
2
cutadapt -a AGATCGGAAGAGCACACGTCTGAACTCCAGTCA \
3
    -m 36 \
4
    -o sample_clean.fastq.gz \
5
    sample_raw.fastq.gz

参数解释：

-a：3’接头序列。3’接头连接在 read 的 3’端
-m 36：丢弃处理后长度 < 36bp 的 reads
-o：输出文件

3.2 双端数据——务必用 -A 指定 R2 接头#

1
# 双端：-a 是 R1 接头，-A 是 R2 接头
2
cutadapt \
3
    -a AGATCGGAAGAGCACACGTCTGAACTCCAGTCA \
4
    -A AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGT \
5
    -m 36 \
6
    -o sample_clean_R1.fastq.gz \
7
    -p sample_clean_R2.fastq.gz \
8
    sample_raw_R1.fastq.gz \
9
    sample_raw_R2.fastq.gz

很多人漏掉 -A，只给 -a 去 R1 的接头。结果 R2 上带着接头继续往下跑比对，比对率狂掉。

3.3 质量修剪——-q 和 —quality-cutoff#

1
# 从3'端修剪低质量碱基
2
cutadapt \
3
    -a AGATCGGAAGAGCACACGTCTGAACTCCAGTCA \
4
    -q 20 \            # 从3'端修剪 <Q20 的碱基
5
    --quality-base=33 \ # Phred+33（Sanger/Illumina 1.8+）
6
    -m 36 \
7
    -o sample_clean.fastq.gz \
8
    sample_raw.fastq.gz

-q 20 的行为是从 3’端往前扫描，遇到第一个 Q<20 的碱基就从那里切断。这比 Trimmomatic 的滑动窗口更激进——一个低质量碱基就触发修剪。

3.4 切除固定长度的碱基#

1
# 切除 read 开头 10bp（去除随机引物/UMI）
2
cutadapt -u 10 -o output.fastq.gz input.fastq.gz
3

4
# 切除末尾 5bp
5
cutadapt -u -5 -o output.fastq.gz input.fastq.gz
6

7
# 保留 read 的前 75bp（WGS 低深度测序常用）
8
cutadapt --length 75 -o output.fastq.gz input.fastq.gz

3.5 批量处理脚本#

1
#!/bin/bash
2
set -euo pipefail
3

4
ADAPTER_R1="AGATCGGAAGAGCACACGTCTGAACTCCAGTCA"
5
ADAPTER_R2="AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGT"
6
THREADS=8
7

8
mkdir -p clean_data cutadapt_reports
9

10
for r1 in raw_data/*_R1.fastq.gz; do
11
    sample=$(basename "$r1" _R1.fastq.gz)
12
    r2="${r1/_R1/_R2}"
13

14
    echo "[$(date +%H:%M:%S)] Processing: $sample"
15

16
    cutadapt \
17
        -a "$ADAPTER_R1" -A "$ADAPTER_R2" \
18
        -q 20 --quality-base=33 \
19
        -m 36 \
20
        --max-n 0.1 \
21
        --trim-n \
22
        -j "$THREADS" \
23
        -o "clean_data/${sample}_R1.fastq.gz" \
24
        -p "clean_data/${sample}_R2.fastq.gz" \
25
        "$r1" "$r2" \
26
        > "cutadapt_reports/${sample}_report.txt"
27

28
    echo "  Done: $sample"
29
done

额外参数：

--max-n 0.1：含 N 超过 10% 的 read 丢弃
--trim-n：把两端 N 切掉
-j 8：8 线程并行

3.6 读懂 Cutadapt 输出统计#

1
=== Summary ===
2
Total read pairs processed:        10,000,000
3
  Read 1 with adapter:               850,000 (8.5%)
4
  Read 2 with adapter:               920,000 (9.2%)
5
Pairs that were too short:           120,000 (1.2%)
6
Pairs written (passing filters):   9,880,000 (98.8%)

重点关注：

“Read with adapter” 比例：>10% 说明接头残留严重，检查文库插入片段是否太短
“too short”：太高（>5%）说明你 -m 设大了，或者数据本身片段短（如 miRNA-seq）
保留率 = “Pairs written” / “Total read pairs”，通常期望 >95%

4. Trimmomatic 实操——四步流水线#

Trimmomatic 的设计理念是固定步骤的流水线。它一次完成：去接头 → 滑动窗口质量修剪 → 最小长度过滤。步骤之间用 : 分隔，参数用空格分隔（非常不直观）。

4.1 基础用法（双端）#

1
trimmomatic PE \
2
    -phred33 \
3
    -threads 8 \
4
    sample_raw_R1.fastq.gz sample_raw_R2.fastq.gz \
5
    sample_paired_R1.fastq.gz sample_unpaired_R1.fastq.gz \
6
    sample_paired_R2.fastq.gz sample_unpaired_R2.fastq.gz \
7
    ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 \
8
    SLIDINGWINDOW:4:20 \
9
    MINLEN:36

输出四个文件：

_paired_R1/R2：配对保留的（这是你后续分析要用的）
_unpaired_R1/R2：配对中只剩一条的（一般丢弃或特殊处理）

4.2 步骤详解#

1
trimmomatic PE \
2
    sample_R1.fq.gz sample_R2.fq.gz \
3
    out_paired_R1.fq.gz out_unpaired_R1.fq.gz \
4
    out_paired_R2.fq.gz out_unpaired_R2.fq.gz \
5
    ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 \
6
    SLIDINGWINDOW:4:20 \
7
    LEADING:3 \
8
    TRAILING:3 \
9
    CROP:150 \
10
    MINLEN:36

步骤	格式	含义
ILLUMINACLIP	`接头文件:种子错配数:回文剪辑阈值:简单剪辑阈值`	2:30:10 是推荐值
SLIDINGWINDOW	`窗口大小:平均质量阈值`	4:20 即 4bp 窗口平均 Q≥20
LEADING	`Q值`	切掉开头低于此 Q 的碱基
TRAILING	`Q值`	切掉末尾低于此 Q 的碱基
CROP	`保留长度`	把所有 read 切成统一长度
MINLEN	`最小长度`	丢弃低于此长度的 read

4.3 ILLUMINACLIP 参数详解#

1
ILLUMINACLIP:TruSeq3-PE.fa:2:30:10
2
              ^           ^ ^  ^ ^
3
              |           | |  | +-- simpleClipThreshold（简单模式阈值）
4
              |           | |  +---- palindromeClipThreshold（回文模式阈值）
5
              |           | +------- seedMismatches（种子允许错配数）
6
              |           +--------- 接头FASTA文件路径

seedMismatches=2：种子匹配允许 2 个错配。设 1 更严格但可能漏掉测序错误的接头
palindromeClipThreshold=30：回文模式——当 R1 和 R2 的接头互相比对上时（即插入片段短于 read 长度）
simpleClipThreshold=10：简单模式——接头直接在 read 上比对上的最低分数

Trimmomatic 自带接头文件在 $CONDA_PREFIX/share/trimmomatic/adapters/ 下。常见的有 TruSeq3-PE.fa（双端）、TruSeq3-SE.fa（单端）、NexteraPE-PE.fa（Nextera 转座酶文库）。

4.4 单端数据#

1
trimmomatic SE \
2
    -phred33 \
3
    -threads 8 \
4
    sample_raw.fastq.gz \
5
    sample_clean.fastq.gz \
6
    ILLUMINACLIP:TruSeq3-SE.fa:2:30:10 \
7
    SLIDINGWINDOW:4:20 \
8
    MINLEN:36

5. 速度与质量实测对比#

用同一份模拟数据（1000 万对 150bp PE reads，含 8.5% 接头污染）对比：

1
# Cutadapt（双端）
2
time cutadapt \
3
    -a AGATCGGAAGAGCACACGTCTGAACTCCAGTCA \
4
    -A AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGT \
5
    -q 20 -m 36 -j 8 \
6
    -o cutadapt_R1.fq.gz -p cutadapt_R2.fq.gz \
7
    raw_R1.fq.gz raw_R2.fq.gz
8
# real 2m18s，保留率 96.3%
9

10
# Trimmomatic（双端）
11
time trimmomatic PE -phred33 -threads 8 \
12
    raw_R1.fq.gz raw_R2.fq.gz \
13
    trim_R1_paired.fq.gz trim_R1_unpaired.fq.gz \
14
    trim_R2_paired.fq.gz trim_R2_unpaired.fq.gz \
15
    ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 \
16
    SLIDINGWINDOW:4:20 MINLEN:36
17
# real 5m42s，保留率 95.8%

指标	Cutadapt	Trimmomatic
运行时间	2m18s	5m42s
保留率	96.3%	95.8%
输出文件数	2	4
JSON报告	✓	✗

结论： Cutadapt 更快、保留率略高、输出更简洁。但 Trimmomatic 提供滑动窗口质量修剪，对低质量数据更友好。

6. 选择指南#

1
# 用 Cutadapt 的场景：
2
# 1. 只需要去接头，不需要质量修剪（用 fastp 做完质量再跑）
3
# 2. 需要切除固定位置碱基（UMI、barcode）
4
# 3. 需要精确的接头统计报告
5
# 4. 处理大量样本需要速度
6
cutadapt -a ADAPTER_SEQ -m 36 -o clean.fq.gz raw.fq.gz
7

8
# 用 Trimmomatic 的场景：
9
# 1. 需要去接头 + 滑动窗口 + 质量修剪一条龙
10
# 2. 数据质量参差，需要激进的质量过滤
11
# 3. 实验室有固定的 Trimmomatic 参数传统
12
trimmomatic PE ... ILLUMINACLIP:adapter.fa:2:30:10 SLIDINGWINDOW:4:20 MINLEN:36
13

14
# 最推荐的组合：fastp（质量+接头）→ 后续流程
15
# fastp 内部用了类似 cutadapt 的算法，同时有滑动窗口
16
fastp -i raw_R1.fq.gz -I raw_R2.fq.gz \
17
    -o clean_R1.fq.gz -O clean_R2.fq.gz \
18
    --detect_adapter_for_pe \
19
    --qualified_quality_phred 20 \
20
    --length_required 36 \
21
    --cut_window_size 4 --cut_mean_quality 20

7. 踩坑记录#

坑1：Cutadapt 的 -a 和 -g 用反了#

1
# 错误：3'接头用了 -g（5'接头参数）
2
cutadapt -g AGATCGGAAGAGC... -o clean.fq.gz raw.fq.gz
3
# 结果：接头根本切不掉，还浪费了计算时间
4

5
# 正确：3'接头用 -a，5'接头用 -g
6
cutadapt -a AGATCGGAAGAGC... -o clean.fq.gz raw.fq.gz

记忆口诀： a = adapter at 3’ end, g = germ at 5’ end（g 是 5’端的意思）。或者——R1 的测序引物在 5’ 端用 -g，常见的接头污染在 3’ 端用 -a。

坑2：Trimmomatic 输出文件名写反了#

1
# 正确顺序：
2
trimmomatic PE \
3
    input_R1.fq.gz input_R2.fq.gz \
4
    output_paired_R1.fq.gz output_unpaired_R1.fq.gz \
5
    output_paired_R2.fq.gz output_unpaired_R2.fq.gz \
6
    ...
7
#  R1 paired → R1 unpaired → R2 paired → R2 unpaired

这是我踩过最多次的坑——把 paired_R2 的内容写到了 paired_R1 的文件位置。然后下游 STAR 比对报 “inconsistent pairs”，查了三个小时回到 Trimmomatic 这一步。

检查方法：

1
# 对比行数（FASTQ 每 4 行一条 read）
2
zcat output_paired_R1.fq.gz | wc -l
3
zcat output_paired_R2.fq.gz | wc -l
4
# 两个数字必须相等！

坑3：Cutadapt 忘了加 -m 导致零长度 reads#

症状：接头比 read 本身还长（比如 miRNA-seq 的 22bp reads 配 33bp 的接头），Cutadapt 把整个 read 当接头切了，输出 0bp 的空序列。

解决： 永远加 -m 参数。

1
cutadapt -a ADAPTER -m 18 -o clean.fq.gz raw.fq.gz
2
# miRNA 至少保留 18bp（成熟 miRNA 约 22bp）

坑4：Trimmomatic 的 SLIDINGWINDOW 窗口大小设成 1#

1
# 错误：每个碱基独立判断
2
SLIDINGWINDOW:1:20
3

4
# 正确：4bp 窗口平均
5
SLIDINGWINDOW:4:20

窗口设成 1，等于看到任何一个 Q<20 的碱基就切。Illumina 测序中偶尔出现单一低质量碱基是正常的，用 4bp 窗口能平滑这个波动。窗口设太大（比如 10bp）又会太保守，低质量尾巴切不干净。

坑5：忘记检查 Trimmomatic 的接头文件路径#

1
# 错误：相对路径在脚本换了目录后找不到
2
ILLUMINACLIP:TruSeq3-PE.fa:2:30:10
3

4
# 解决：用绝对路径
5
ADAPTER_FILE="${CONDA_PREFIX}/share/trimmomatic/adapters/TruSeq3-PE.fa"
6
trimmomatic PE ... ILLUMINACLIP:${ADAPTER_FILE}:2:30:10 ...

最好在脚本开头检查文件存在：

1
ADAPTER_FILE="${CONDA_PREFIX}/share/trimmomatic/adapters/TruSeq3-PE.fa"
2
if [[ ! -f "$ADAPTER_FILE" ]]; then
3
    echo "ERROR: Adapter file not found: $ADAPTER_FILE"
4
    exit 1
5
fi

本文于 2025-04-17 在 Debian 13 + Cutadapt 4.9 + Trimmomatic 0.39 上实测。