FASTQ 文件结构与质量控制：序列、质量值、Phred 编码详解

真正开始跑一个空间转录组 pipeline 之前，有一件事我劝你别偷懒——搞懂 FASTQ 文件。

这玩意儿看着就是个四行一组的文本块，很多人直接扔进上游工具就完事了。但我见过太多 case 了：因为质量值编码搞错，整个比对结果全是垃圾；因为 adapter 没切干净，导致聚类图上一堆莫名其妙的噪音。嗯，这块我们得说透。

FASTQ 的四行结构

每一份 FASTQ 记录都包含四行，这是硬性规定。缺一行或顺序乱了，下游工具直接报错。

@EAS139:136:FC706VJ:2:2104:15343:197393 1:Y:18:ATCACG
CTCAAGTTGTCAGCAATGGTGGTGCAATTGGCAGGATGAATATTA
+
B@@FFFFFHHHHHJJJJJJIIIIIIJJJJJJJJIJJJJJJJJJJJH
行号 标识符 内容
第1行 @ 序列标识 + 可选注释（如流道位置、barcode）
第2行 — 碱基序列（A/T/C/G/N）
第3行 + 分隔符，有时跟第1行相同的ID
第4行 — 质量值字符串，与第2行碱基一一对应
我的习惯：拿到新测序数据后先用 head -n 40 看一眼，确认第四行全是可打印字符，没有换行符乱入。有一次我从云盘下数据，第四行被中间截成两段——后来查出来是传输过程没加 -b 选项搞的。

质量值是什么？

每个碱基的质量值代表这个碱基被测序仪判定正确的概率。说白了就是测序仪自己给的“信心指数”。

质量值的数学定义很简单：
Q = -10 × log₁₀(P)
其中 P 是错误概率。比如 Q=30 表示错误概率 0.001，也就是一千个碱基里大概错一个。

Phred 编码——小心这个坑

质量值存到文件里不直接存数字，而是用 ASCII 字符表示。不同时代、不同测序平台用的编码偏移不一样。

Phred+33: ASCII码 = Q + 33  （Illumina 1.8+，当前主流）
Phred+64: ASCII码 = Q + 64  （早期 Illumina/Solexa，已淘汰）

举个例子，质量值 Q=30，在 Phred+33 下存为 ASCII 码 63（即字符 ?），而在 Phred+64 下存为 ASCII 码 94（即字符 ^）。

我最惨的一次踩坑：曾经拿到一批 2012 年的数据，FASTQ 头文件里没写明编码。我默认用了 Phred+33 去跑 cutadapt，结果 adapter 没切掉几个，后来一看质量值分布图——中位数才 Q10，明显是编码错位。改用 Phred+64 之后，一切正常。从那以后我养成了一个习惯：新数据先跑 fastqc 看 per base sequence quality 图，如果中位数低于 10 或者出现断崖式下降，十有八九是编码问题。

质量控制：从 FASTQ 到干净数据

我一般把 QC 分成三个关卡：

  1. 基础统计——总读段数、碱基数、GC含量、是否存在过表达序列（adapter）。
  2. 质量过滤——切除低质量尾部、去掉太短的读段、去除 N 碱基过多的读段。
  3. adapter 切除——尤其是空间转录组的 read1 包含 barcode/UMI，要特别注意不要误切。

常用的工具组合：

# 推荐流程：先用 FastQC 看报告
fastqc sample_R1.fastq.gz sample_R2.fastq.gz -o qc_reports/

# 然后用 fastp 一站式处理（质量过滤 + adapter 切除 + 统计报告）
fastp -i sample_R1.fastq.gz -I sample_R2.fastq.gz \
      -o clean_R1.fastq.gz -O clean_R2.fastq.gz \
      --detect_adapter_for_pe \
      --cut_front --cut_tail \
      --qualified_quality_phred 20 \
      --unqualified_percent_limit 40 \
      --length_required 50 \
      --json fastp_report.json

为什么我偏爱 fastp？

  • 一步到位，不需要先 adapter 识别再调用 cutadapt
  • 支持 per-read 质量裁切，比滑动窗口更敏感
  • 自动检测双端 adapter 方向
  • 输出 JSON/HTML 报告，方便集成到 pipeline

从质量控制到 Visium 聚类前的预处理

Visium 空间转录组是双端测序：R1 通常只读 28bp（包含 16bp barcode + 12bp UMI），R2 读 91bp（基因序列）。这意味着 adapter 位置是固定的，不像是 RNA-seq 那样 adapter 随机出现。所以你想想看，如果 R1 被错误切掉几个碱基，barcode 就全乱了，后续聚类就等于瞎猜。

特别注意：Visium 数据中 R1 的 barcode 和 UMI 紧贴着读段开头，质量裁切如果开启了 --cut_front，可能导致 barcode 被切掉前几个碱基。我个人建议在 fastp 命令中显式使用 --cut_front 0 来避免这个问题，或者单独处理 R1 和 R2。

下图展示了我每次做空间转录组 QC 时脑子里过的流程：

FASTQ 质量控制与清洗流程（空间转录组 Visium 视角） 原始 FASTQ（R1+R2） 第一关：基础统计 FastQC → 总读段、GC%、过表达序列 第二关：质量过滤 fastp → 裁切低质量尾、去短读、去N clean FASTQ R1 特殊检查 barcode/UMI 完整性 adapter 切除策略 固定位置切除，避免误切 barcode ⚠ 空间转录组特有注意事项 R1 仅 28bp，质量裁切 if --cut_front 必须关掉，否则 barcode 前 6bp 被切 → 聚类错位 建议先跑一个测试小样本，检查 barcode 分配率是否 > 80%

实战检查清单

我每次新到一个项目，都会拿下面这个清单过一遍，避免低级错误：

检查项 工具/方法 合格指标
编码类型确认 FastQC 的 per base sequence quality 图 中位数 Q 值落在 20-40 范围，无底层平台偏移
GC 含量异常 FastQC / fastp 报告 正态分布，无明显双峰或偏移
adapter 残留 FastQC 的 Overrepresented sequences adapter 占比 < 1%
Barcode 分配率 先跑 spaceranger mkfastq 看指标 分配率 > 75%（我的经验：低于 60% 要重做 QC）
双端读段配对完整性 wc -l 对比 R1/R2 行数 行数相等，R1 总读段 = R2 总读段
小技巧：快速判断编码方式还有一个土办法——看第四行字符范围。如果包含 ! 到 I 之间（ASCII 33-73），大概率是 Phred+33；如果包含很多 B 到 h 之间的字符（ASCII 66-104），可能是 Phred+64（旧数据）。我的做法是用 python 读一行看看最小 ASCII 码：
head -n 4 sample.fastq | tail -n 1 | od -An -tuC | awk '{min=$1; for(i=1;i<=NF;i++) if($i
如果最小值 33 就是普通 Phred+33；如果最小值在 59-64 之间，说明是旧版编码。

嗯，到这里 FASTQ 结构、质量值、Phred 编码这些基础你应该有底了。质量控制这块我强调再多也不为过——数据进 pipeline 之前多花 10 分钟检查，后面能省下好几天排查 bug 的时间。接下来我们继续往下走，看看这些干净的 FASTQ 数据如何和参考基因组对上号。

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