测序质量分数基础

大家在做测序数据分析时，第一步肯定会看到FASTQ文件里那些奇怪的符号和数字。嗯，其实就是质量分数。很多人直接跳过这步，拿着数据就开始比对，结果下游分析跑出异常才回头查。我当年也犯过这种错，后来学乖了——质量分数不过关，后面所有分析都是白搭。

今天咱们就把这个最基础的东西掰开揉碎讲清楚。你想想看，测序仪是怎么知道它读出来的碱基是对是错的？说白了，它并不能百分之百确定，只能给一个“可信度”的量化指标。这个指标就是Phred质量分数。

什么是Phred质量分数

Phred质量分数这个概念，最早是1998年由华盛顿大学的Phil Green和他的团队提出的，用在Phred碱基识别软件里。我记得当时Sanger测序还是主流，那个软件能从荧光信号轨迹里直接算出碱基的出错的概率。后来二代测序兴起，这个体系就被完整继承了下来，一直用到现在。

Phred质量分数，简称Q值，是一个对数转换的数值。公式很简单：

Q = -10 × log10(P)

其中P是碱基识别错误的概率

举例：
如果P = 0.01（百分之一），则Q = -10 × log10(0.01) = 20
如果P = 0.001（千分之一），则Q = -10 × log10(0.001) = 30

为什么非要用对数？我个人理解是：错误概率的变化范围可以从10⁻¹到10⁻⁴甚至更低，线性数字实在太长了，用对数一次就能压缩到直观的两位数。你想想看，Q20和Q30只差10，但它们的错误概率差了整整一个数量级。

Q20与Q30在实战中的意义

我刚入行时觉得“99%准确率”已经挺不错了。后来做过一次人全外显子组测序，Q20以下的数据里找了二十几个假阳性SNP，还好没有直接发报告给医生，不然就出大问题了。从那以后我就把Q30当成了自己的“心理安全线”。

为什么会这样？你要知道测序仪每次跑一个lane就能产出几十亿个碱基，哪怕只有1%的错误率，按3亿碱基的基因组算，错误位点就有300万个。如果做低频突变检测，比如肿瘤组织里5%频率的突变，Q20的数据根本没法信任。

FASTQ文件里质量分数的编码方式

实际文件里看到的不是数字Q值，而是ASCII字符。这是因为当年设计FASTQ格式时，为了方便用一个字节存储，就把Q值加上33或64之后映射到可打印字符上。这里有个经典的坑——编码偏移。

最常见的两种编码：

Sanger格式: ASCII = Q + 33
Illumina 1.3+: ASCII = Q + 64

例如：
Q = 30
Sanger编码: 30 + 33 = 63 → ASCII字符 '?'
Illumina 1.3+: 30 + 64 = 94 → ASCII字符 '~'

注意，字符 '?' 和 '~' 完全不同

我在一次多批次数据分析项目中遇到过这个问题。之前同事把不同平台的FASTQ直接合并了，我一看质量分布图有点奇怪，排查了半小时才发现他们混用了33偏移和64偏移。后来我建议大家一起用Sanger / Phred+33格式，统一后再也不出这种幺蛾子了。

# 用 head 查看前4行，观察第4行的ASCII字符
head -n4 sample.fastq

# 用 FastQC 直接识别
fastqc sample.fastq

from Bio import SeqIO
records = SeqIO.parse("sample.fastq", "fastq")
for rec in records:
    print(rec.letter_annotations["phred_quality"][:10])
    break

质量分数与批次关联的直觉

为什么质量分数会和批次挂钩？你想想看，同一台测序仪，不同run之间试剂批次不同、芯片批次不同、操作人员不同，甚至实验室温度湿度都有差异。这些因素都会直接影响碱基识别的信噪比，最终体现为Q20 / Q30比例的波动。

我之前做一批200个肿瘤样本的WES，第一次跑完发现其中两板96孔板的Q30比例差了8%。一开始以为是操作失误，后来查了仪器日志才知道——那批芯片在运输过程中受到过冻融。从那以后，我每批数据进来都第一时间做“质量分数-批次”的关联分析，把异常批次提前标记出来，避免带到下游分析里混淆。