FastQC质检与fastp清洗：测序数据预处理

2.1 基础用法#

1
# 创建输出目录
2
mkdir -p qc_reports
3

4
# 对单个文件跑
5
fastqc sample_R1.fastq.gz -o qc_reports/
6

7
# 批量跑（-t 并行线程）
8
fastqc -t 8 -o qc_reports/ *.fastq.gz

输出两个文件：

• sample_R1_fastqc.html：可视化报告（浏览器打开）
• sample_R1_fastqc.zip：原始数据（给 MultiQC 用）

2.2 报告怎么看——关注这五个图#

打开 HTML 报告，左边导航有十几个模块。新手不用全看，先抓这五个最关键的：

1
Read 1 adapter: AGATCGGAAGAGCACACGTCTGAACTCCAGTCA
2
Read 2 adapter: AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGT

2.3 批量 FastQC 脚本模板#

1
#!/bin/bash
2
# fastqc_batch.sh - 批量质检脚本
3
# Debian 13 实测 2025-02-03
4

5
set -euo pipefail
6

7
DATA_DIR="./raw_data"
8
OUT_DIR="./qc_reports"
9
THREADS=8
10

11
mkdir -p "$OUT_DIR"
12

13
echo "=== FastQC Started: $(date) ==="
14
echo "Processing files in: $DATA_DIR"
15

16
# 一次性跑所有fastq.gz
17
fastqc -t "$THREADS" -o "$OUT_DIR" "${DATA_DIR}"/*.fastq.gz
18

19
echo "=== FastQC Finished: $(date) ==="
20
echo "Reports in: $OUT_DIR"
21
echo "Open *_fastqc.html in browser to view"

3.1 完整命令及参数解读#

1
fastp \
2
    -i sample_R1.fastq.gz \        # 输入read1
3
    -I sample_R2.fastq.gz \        # 输入read2
4
    -o clean_R1.fastq.gz \         # 输出read1
5
    -O clean_R2.fastq.gz \         # 输出read2
6
    -h fastp_report.html \         # HTML质控报告
7
    -j fastp_report.json \         # JSON数据（给程序读）
8
    --detect_adapter_for_pe \      # 自动检测双端接头
9
    --qualified_quality_phred 20 \ # Q≥20才算合格碱基
10
    --length_required 36 \         # 过滤后短于36bp就扔掉
11
    --cut_front \                  # 从5'端裁剪低质量
12
    --cut_tail \                   # 从3'端裁剪低质量
13
    --cut_window_size 4 \          # 滑动窗口4bp
14
    --cut_mean_quality 20 \        # 窗口平均Q<20就裁掉
15
    --correction \                 # PE overlap区域纠错
16
    --overlap_len_require 30 \     # overlap至少30bp
17
    -w 8                           # 8线程并行

3.2 滑动窗口裁剪原理#

fastp 的滑动窗口过滤是这个工具最精妙的设计之一。算法：

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

其中 $w=4$（窗口大小），$Q_i$ 是窗口中第 i 个碱基的质量值。如果 $Q_{window} < 20$，就从该窗口开始裁掉后面所有的。

为什么不用全局平均？ 因为局部质量暴跌才是问题——read 开头 Q40 结尾 Q5，全局平均可能还有 Q30 但你用不了那条 read。

3.3 参数调优经验#

1
# 场景1：数据质量较好（Q30>90%）
2
fastp --qualified_quality_phred 15 --length_required 50
3

4
# 场景2：数据质量较差（Q30<80%）
5
fastp --qualified_quality_phred 20 --length_required 36 \
6
      --cut_front --cut_tail --cut_window_size 4 --cut_mean_quality 20
7

8
# 场景3：怀疑有接头污染但接头序列不确定
9
fastp --detect_adapter_for_pe --adapter_sequence auto
10

11
# 场景4：只想看报告，不实际过滤（干跑）
12
fastp --disable_quality_filtering --disable_length_filtering \
13
      --disable_adapter_trimming

3.4 解读 fastp 输出#

运行完后终端会打印过滤统计：

1
Read1 before filtering:
2
total reads: 30000
3
total bases: 4500000
4
Q20 bases: 4500000(100%)
5
Q30 bases: 4500000(100%)
6

7
Read1 after filtering:
8
total reads: 29850
9
total bases: 4477500
10

11
Filtering result:
12
reads passed filter: 29850
13
reads failed due to low quality: 120
14
reads failed due to too short: 30
15
reads with adapter trimmed: 850
16
bases trimmed due to adapters: 2800
17

18
Duplication rate: 0%

重点关注这几个指标：

3.5 批量 fastp 脚本#

1
#!/bin/bash
2
# Debian 13 实测 2025-02-03
3

4
set -euo pipefail
5

6
INPUT_DIR="./raw_data"
7
OUTPUT_DIR="./clean_data"
8
REPORT_DIR="./qc_reports"
9
THREADS=8
10

11
mkdir -p "$OUTPUT_DIR" "$REPORT_DIR"
12

13
# 从R1文件名推导R2文件名
14
for r1 in "$INPUT_DIR"/*_R1.fastq.gz; do
15
    r2="${r1/_R1.fastq.gz/_R2.fastq.gz}"
16
    sample=$(basename "$r1" _R1.fastq.gz)
17

18
    echo "Processing: $sample"
19

20
    fastp \
21
        -i "$r1" -I "$r2" \
22
        -o "$OUTPUT_DIR/${sample}_clean_R1.fastq.gz" \
23
        -O "$OUTPUT_DIR/${sample}_clean_R2.fastq.gz" \
24
        -h "$REPORT_DIR/${sample}_fastp.html" \
25
        -j "$REPORT_DIR/${sample}_fastp.json" \
26
        --detect_adapter_for_pe \
27
        --qualified_quality_phred 20 \
28
        --length_required 36 \
29
        --cut_front --cut_tail \
30
        --cut_window_size 4 --cut_mean_quality 20 \
31
        -w "$THREADS"
32

33
    echo "  Done: $sample"
34
done
35

36
echo "=== All samples processed ==="

关键点： for r1 in ..._R1.fastq.gz 用 _R1 自动找对应的 _R2，这种命名规则是 Illumina 下机数据的默认格式。如果你的文件命名不标准，改掉 _R1 和 _R2 的替换规则。

坑1：FastQC 报 Can't allocate memory#

30GB 的 FASTQ 直接喂给 FastQC，8GB 内存的服务器直接 OOM。

解决： FastQC 默认读整个文件到内存。用 --nogroup 禁用一个内存密集型的画图功能：

1
fastqc --nogroup sample.fastq.gz

或者先用 seqtk sample 随机抽 100 万条跑个初检：

1
seqtk sample -s42 sample.fastq.gz 1000000 | gzip > sample_sub.fastq.gz
2
fastqc sample_sub.fastq.gz

坑2：fastp 输出 reads 数骤减#

症状：输入 5000 万 reads，输出只剩 3000 万。

排查顺序：

1. 先看 --length_required 是不是设太高（比如设 100，很多短 read 被误杀）
2. 再看 --qualified_quality_phred 是不是太严（设 30 而不是 20）
3. 最后看数据是不是真的差——跑个 FastQC 确认

经验数值： --length_required 36（最常用），--qualified_quality_phred 20。如果这个参数下还是 >20% 被过滤，数据确实有问题。

坑3：PE 和 SE 数据混用导致 fastp 报错#

1
# 错误做法：单端数据用了双端参数
2
fastp -i single_end.fastq.gz -I NOT_EXIST.fastq.gz  # 报错！
3

4
# 正确做法：单端只用 -i
5
fastp -i single_end.fastq.gz -o clean.fastq.gz

加个判断脚本：

1
if [[ -f "${r1/_R1/_R2}" ]]; then
2
    # 双端模式
3
    fastp -i "$r1" -I "$r2" ...
4
else
5
    # 单端模式
6
    fastp -i "$r1" ...
7
fi

坑4：MultiQC 报 No analysis results found#

症状：MultiQC 扫描目录输出了空报告。

原因：你给 MultiQC 的目录里没有它能识别的文件。MultiQC 认识的 FastQC 输出是 *_fastqc.zip，不是 *_fastqc.html。

1
# 确认 FastQC 输出文件存在
2
ls *_fastqc.zip
3

4
# 再跑 MultiQC
5
multiqc .

坑5：接头序列检测不到#

fastp 的 --detect_adapter_for_pe 不是万能的。新测序平台或定制文库的接头，可能不在 fastp 内置的接头库中。

解决： 手动指定接头序列：

1
fastp --adapter_sequence AGATCGGAAGAGCACACGTCTGAACTCCAGTCA \
2
      --adapter_sequence_r2 AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGT

或者用 cutadapt 专门去接头：

1
cutadapt -a AGATCGGAAGAGC... -A AGATCGGAAGAGCG... \
2
    -o clean_R1.fastq.gz -p clean_R2.fastq.gz \
3
    sample_R1.fastq.gz sample_R2.fastq.gz