▎第 2 章：BLAST 结果格式解析

BLAST 跑完了，拿到结果的那一刻，很多人第一反应是「好大一坨文字」。其实它背后藏着三种截然不同的格式——默认的 pairwise、简洁的表格格式、以及结构化的 XML。我从第一次接触 BLAST 就被默认输出搞得有点懵，后来才搞清楚每种格式到底长什么样、里面每个字段在说什么。

这一章我带你把三种格式从头到尾拆一遍。我会告诉你每个字段的含义、什么时候该用哪种、以及我在实际项目中遇到的小坑。

3.1 默认 pairwise 格式（带 -outfmt 0）

这是 BLAST 最原始的输出方式。你什么参数都不加，结果就是 pairwise 格式。它用一个「人类可读」的方式展示每条比对，从 query 序列开始，再到 subject 序列，一行一行对齐给你看。

长什么样？

BLASTP 2.15.0+
Query= seq_001
Length=345
Subject= sp|P12345|MYG_HUMAN
Length=154
 Score =  245 bits (625),  Expect = 2e-67
 Identities = 125/154 (81%), Positives = 134/154 (87%), Gaps = 0/154 (0%)
 Query  1  MVHLTPEEKSAVTALWGKVNLDEVGGETLGRLLVVYPWTQRFFESFGDLS  50
            MVHLTPEEKSAVTALWGKV--DEVGGETLGRLLVVYPWTQRFFESFGDLS
 Sbjct  1  MVHLTPEEKSAVTALWGKVNVDEVGGETLGRLLVVYPWTQRFFESFGDLS  50
 ...

我习惯用 pairwise 来做「人工审查」。当你怀疑某个比对结果是不是有异常时，把原始排列拿出来一行一行看，比任何统计数字都直观。但项目量大时，用它做批量分析就太慢了。

3.2 表格格式 outfmt 6 / 7 / 10

这是做大规模分析的首选格式。你指定 -outfmt 6，BLAST 会输出制表符分隔的行，每条比对一行。它特别适合直接喂给 awk、perl、R 和 Python。

-outfmt 7 和 -outfmt 10 本质相同——7 加了个注释行（以 # 开头），10 用逗号分隔而不是制表符。

格式示例（outfmt 6）：

seq_001	sp|P12345|MYG_HUMAN	81.17	154	29	0	1	154	1	154	2e-67	245.0
seq_002	sp|Q9Y6K1|CYC_HUMAN	72.50	120	33	0	2	121	1	120	1e-42	162.0

12 个字段按顺序分别为：

列号	字段名	含义
1	qseqid	Query 序列 ID
2	sseqid	Subject 序列 ID
3	pident	百分比序列相似度
4	length	比对长度（碱基/氨基酸数）
5	mismatch	错配数
6	gapopen	gap 开放数
7	qstart	Query 起始位置（1-based）
8	qend	Query 结束位置
9	sstart	Subject 起始位置
10	send	Subject 结束位置
11	evalue	E-value（科学计数法）
12	bitscore	比特分

我一般这么用：

# 只用 awk 就能过滤出 E-value < 1e-20 且相似度 > 90% 的比对
awk '$11 < 1e-20 && $3 > 90' blast_output.tab

3.3 XML 格式（-outfmt 5）

XML 格式是所有格式里信息最完整的。它把每条比对、每个高分段、每个 alignment 都变成结构化的 XML 节点。你可以用任何支持 XML 的语言或库来解析，比如 Python 的 ElementTree 或者 Java 的 DOM/SAX。

结构概览（简化示例）：

<BlastOutput>
  <BlastOutput_program>blastp</BlastOutput_program>
  <BlastOutput_version>BLASTP 2.15.0+</BlastOutput_version>
  <BlastOutput_db>swissprot</BlastOutput_db>
  <BlastOutput_iterations>
    <Iteration>
      <Iteration_iter-num>1</Iteration_iter-num>
      <Iteration_hits>
        <Hit>
          <Hit_num>1</Hit_num>
          <Hit_id>sp|P12345|MYG_HUMAN</Hit_id>
          <Hit_len>154</Hit_len>
          <Hit_hsps>
            <Hsp>
              <Hsp_score>245.0</Hsp_score>
              <Hsp_evalue>2e-67</Hsp_evalue>
              <Hsp_identity>125</Hsp_identity>
              <Hsp_positive>134</Hsp_positive>
              <Hsp_gaps>0</Hsp_gaps>
              <Hsp_qseq>MVHLTPEEKSAVTALWGKVNLDEVGGETLGRLLVVYPWTQRFFESFGDLS</Hsp_qseq>
              <Hsp_hseq>MVHLTPEEKSAVTALWGKVNVDEVGGETLGRLLVVYPWTQRFFESFGDLS</Hsp_hseq>
              <Hsp_midline>                                </Hsp_midline>
            </Hsp>
          </Hit_hsps>
        </Hit>
      </Iteration_hits>
    </Iteration>
  </BlastOutput_iterations>
</BlastOutput>

XML 格式的优点是：你可以写脚本一次性提取你想要的所有字段，而不用去解析人类可读的 pairwise 格式。我自己的习惯是：当需要生成报告或建立数据库时，就用 XML；如果只是日常检查、快速统计，表格格式就够了。

3.4 三种格式对比一览

我整理了一个表，方便你一眼看出区别：

特性	Pairwise（outfmt 0）	表格（outfmt 6/7/10）	XML（outfmt 5）
可读性	最好	一般（需要配合工具）	差（需要解析器）
字段覆盖面	全（包含比对序列）	只含统计指标	全（包含比对序列）
解析难度	低（人工）	中（可用 awk 等）	高（需要 XML 库）
批量处理效率	低	高	中
适合场景	人工审查单条比对	批量过滤、排序、统计	程序化处理与数据入库
是否包含序列排列	是	否	是

3.5 知识体系图：三种格式的适用流程

下面这张图帮你理解：拿到结果后，如何根据你的目标选择格式。图中用流程串联了 pairwise、表格格式和 XML 的典型使用场景。

3.6 格式转换示例：从 pairwise 到表格，或 XML 到表格

有时候你已经有了 pairwise 输出的结果，但需要转成表格格式做批量分析。别慌，我自己也遇到过。最常用的做法是：重新跑一次 BLAST，指定你想要的格式。但如果你不想重新跑，也可以用正则抓取字段。

用 awk 从 pairwise 中提取 E-value 和 Score：

awk '/Score = / { split($0, a, "Expect = "); split(a[2], b, ","); print $4, $5, b[1] }' blast_output.txt

但说实话，这种做法容易出错，特别是 pairwise 格式在不同 BLAST 版本间有细微差异。我建议你直接：

• 保留原 query 和 subject 序列文件
• 用 -outfmt 6 重跑一次
• 把结果和 pairwise 做个交叉验证

这样既快又准确，还能顺便检验一下原始结果的一致性。