2. 测序数据格式与质量控制

好了,假设你已经拿到了测序公司的原始数据。打开一看,全是 .fastq 后缀的文件。这东西长什么样?我当年第一次看到 FastQ 文件时,愣了半天——全是字母和符号,怎么看都不像生物数据。

说白了,FastQ 就是测序仪吐出来的「毛坯房」数据。咱们得先验验货,看看质量怎么样,才能决定下一步怎么处理。这一章我就带你搞清楚三件事:FastQ 格式到底存了什么、怎么用 FastQC 做质量检查,以及如何用 MultiQC 把几十个样本的报告汇总到一块。

2.1 FastQ 格式:每四行一个「生命单位」

每个测序 read 在 FastQ 文件里占刚好四行。我记得刚入行时有个同事非说那四行是随机文本,结果跑完比对才发现全是低质量数据,白算了一整天。

@EAS139:136:FC706VJ:2:2104:15343:197393 1:Y:18:ATCACG
CTGTTCCTGCTGAACTTCAGGGTGAGTGCGCCAGCCCTCCACACACAGCACCC
+
IIIIIIIIIIIIIIIIIIIIIIIIIJIIIIIIIIIJIJJJJJJJJJJJJJJJ

这四行的含义如下:

行号标识内容说明
第1行@序列标识符(仪器ID、flowcell编号、坐标、barcode等信息)
第2行序列碱基序列(A/T/C/G/N)
第3行+分隔符,有时后面重复序列标识符
第4行质量值与第2行等长的ASCII字符串,每个字符代表一个碱基的测序质量

关键在第4行。质量值的编码方式叫 Phred+33,每个 ASCII 字符的数值减去 33 就是 Phred 质量分数 Q。Q 值计算公式很简单:

Q = -10 * log10(P)
其中P是碱基测序错误的概率
经验换算: Q30 相当于错误率 1/1000(99.9% 准确率),Q20 是 1/100(99% 准确率)。我一般要求 RNA-seq 数据至少 Q30 比例在 80% 以上,低于这个的样本我会直接找测序公司沟通。

ASCII 字符对照表里,从 '!'(ASCII 33,Q=0)到 '~'(ASCII 126,Q=93),覆盖了测序仪的整个质量范围。注意,Illumina 1.8+ 版本统一用 Phred+33 编码。如果你遇到老版本数据(比如 1.3 或 1.5),编码方式不同,需要用工具转换。我之前接过一个合作项目,对方给了 1.5 格式的 FastQ,我直接用 FastQC 跑完才发现碱基质量全标错了——嗯,从此我养成了先看版本号的习惯。

2.2 为什么要做质量控制?

测序过程不是完美的。你想想看,测序仪在那么小的纳米孔里读取荧光信号,信号衰减、光学衍射、试剂沉淀、cluster 重叠……哪个环节出点小问题,最终都会反映到 FastQ 数据里。

常见问题:

  • 接头污染:片段比测序读长短,测到了接头序列。RNA-seq 中 poly(A) 富集后的数据经常在 3' 端出现接头残留。
  • 低质量碱基:读段末端质量衰减——Illumina 平台的常见现象。
  • GC 含量异常:可能有污染 RNA(比如 rRNA 没有去除干净)或 PCR 扩增偏好。
  • 重复序列过高:可能来自过度 PCR 扩增,或者 mRNA 本身高表达导致的生物学重复。
  • N 碱基过多:说明测序质量差或 DNA 模板有问题。

如果不做质量控制就强行往下走,低质量 read 会污染比对结果,假阳性差异表达基因的数量会噌噌往上冒。我曾经帮人审过一个项目,对方说自己数据不好用,我一看 FastQC 报告,接头污染率超过 60%,直接把比对率拉低到 30% 以下。

2.3 FastQC:一个文件一个报告的「显微镜」

FastQC 是 Babraham Bioinformatics 开发的工具,专门用来评估高通量测序数据质量。它会对每个 FastQ 文件跑出 11 项核心指标,生成一个直观的 HTML 报告。

用法极其简单:

fastqc sample1.fastq.gz -o ./fastqc_reports/ -t 8
  • sample1.fastq.gz:输入文件,支持 gzip 压缩
  • -o:输出目录
  • -t:线程数,多核环境加速

输出两个文件:sample1_fastqc.html(可视化报告)和 sample1_fastqc.zip(原始数据)。

报告的核心模块包括:

模块名称评估内容常见红色警报原因
Per base sequence quality每个碱基位置的箱线图质量分布5'端质量低(酶偏好)、3'端严重衰减
Per tile sequence qualityflowcell 每个 tile 的质量热图某个 tile 整体偏红→仪器问题
Per sequence quality scores每条 read 的平均质量分布大量 read 平均 Q<20
Per base sequence content每个位置 A/T/C/G 比例前几个碱基 AT/GC 分离→随机引物偏好
Per sequence GC content整体 GC 含量分布双峰→可能有两种物种的 DNA 混合
Per base N content每个位置 N 的比例末端 N 升高→质量衰减
Sequence Length Distributionread 长度分布长短不一→片段化或接头残留
Sequence Duplication Levels序列复制程度重复比例过高→PCR bias 或真实高表达
Overrepresented sequences过度呈现的序列接头序列、rRNA、线粒体 DNA
Adapter Content接头含量累积曲线3'端接头指数上升→接头需要去除
Kmer Content富集的短 k-mer接头或引物二聚体
别被「红色叉号」吓到。 FastQC 的通过/警告/失败阈值设定得很保守。比如 Per base sequence quality 如果只在最后几个碱基出现黄色警告,通常不影响整体分析。我个人经验:只要前 80% 的碱基质量中位数不低于 Q30,末端的 10-20 bp 可以果断剪掉。

2.4 MultiQC:多份报告的「集控中心」

当你手上有 50 个样本时,一个个打开 FastQC 的 HTML 页面就很痛苦了。MultiQC 就是来解决这个问题的——它能把 FastQC 和其他多种生物信息学工具(如 STAR、Kallisto、Salmon、featureCounts 等)的输出汇总到一个交互式 HTML 报告里。

安装很简单:

pip install multiqc

用法也相当直接:

multiqc ./fastqc_reports/ -o ./multiqc_report/ -n my_rnaseq_qc

Multiqc 会自动扫描目录下的 *_fastqc.zip 文件,提取数据,然后生成一个统一的报告。我强烈建议在每个 RNA-seq 项目的 QC 阶段都生成一份 MultiQC 报告,原因有三:

  • 横向对比:把所有样本的 GC 含量曲线堆叠在一起,一眼就能看出哪个样本异常。
  • 分组统计:按实验组/对照组着色,快速定位 batch effect 是否严重。
  • 汇总表格:MultiQC 会自动生成样本的质量参数表格,方便你决定保留哪些样本。

我的实战经验: 在正式跑差异分析之前,我会先看 MultiQC 的「Per sequence GC content」和「Adapter Content」两个面板。如果某个样本的 GC 分布出现双峰,或者接头含量曲线在 3' 端快速上扬,我马上就知道——这个样本跑完差异分析后大概率是假阳性大户。

2.5 质量控制的核心流程

下面这张图是我每次做 RNA-seq QC 时遵循的基本流程:

RNA-seq 质量控制核心流程 原始 FastQ 数据 FastQC 逐个样本评估 判断是否通过 QC 阈值 碱基质量 | GC 含量 | 接头含量 | 重复率 跳过/重新测序 数据清洗(Trimming) 直接进入下游 MultiQC 汇总报告 高质量数据 🔍 关键检查点 N 比例 < 5% | GC 无双峰

流程说白了就是三步走:

  1. 跑 FastQC:对所有样本逐一生成质量报告
  2. 看 MultiQC:汇总后横向对比,定位问题样本
  3. 做决策:确定哪些样本需要清洗、哪些直接过、哪些干脆放弃
一个实用技巧: 我习惯在 MultiQC 报告里把「Per base sequence quality」和「Adapter Content」两个面板并排截图,直接贴到项目报告里。这样合作者一眼就能看到数据质量概貌,比写一堆文字更直观。

2.6 实战中的避坑指南

做了这么多年 RNA-seq QC,有几个坑我掉进去过不止一次,分享给你:

  • 不要只看单个 FastQC 报告。我曾经在某个项目里只看了一个样本的 FastQC,觉得质量很好,结果跑完差异分析发现几个样本对比时假阳性特别高。后来用 MultiQC 一看,有个样本的重复序列率高达 90%,明显是 PCR 过度扩增。单样本报告根本发现不了这个问题。
  • 看清编码格式。前面提过,老版本的 Illumina(如 GAII)用 Phred+64 编码。如果你直接用 Phred+33 的工具去读,所有碱基质量都会被误判为超级高。我遇到过有人用错了工具版本,QC 报告全是绿色通过,但实际上数据质量一塌糊涂。
  • GC 双峰不一定是污染。RNA-seq 数据中,如果某个基因异乎寻常地高表达(比如 Actin 或 GAPDH),GC 含量可能会拉偏分布。再加上 poly(A) 富集后的偏差,GC 曲线可能会偏。所以看到双峰先别急着定罪,可以局部看看是哪个 reads 引起的。
  • 接头检测阈值要调。FastQC 默认识别接头用的是一个固定库。如果你用的是特殊接头(比如带 UMI 或 barcode 的定制接头),FastQC 可能识别不出来。这种情况下可以用 --adapters 参数手动指定接头序列文件。