1. 为什么基因组比对值得单独写一篇实战指南

做生物信息分析的人，绕不开比对这一步。测序仪吐出来的 FASTQ 文件，本质上就是几千万到几亿条几十到一百多碱基的短字符串，而参考基因组是几个 G 的长字符串。把短串精准地"贴"回长串上，这个动作叫比对（alignment），它是变异检测、表达定量、宏基因组分析、进化研究几乎所有下游流程的第一道关口。BWA（Burrows-Wheeler Aligner）就是这道关口上被用得最多的工具之一，尤其是它的bwa mem算法，在二代短读长数据上几乎是事实标准。

这篇内容面向三类人：刚进实验室、第一次跑比对流程的新手；已经会用但总在索引、参数、报错上卡壳的中间层；以及需要把流程固化下来、给团队做标准化的人。我会从安装、建索引、比对、格式转换、质控到常见报错排查，把整条链路拆开讲清楚，重点放在"为什么这么选"和"踩过的坑"上，而不是把官方文档翻译一遍。你照着做，能跑通；你读完，能明白每一步在干什么。

需要先说明一点：BWA 本身只负责比对，输出 SAM 格式，后续的排序、去重、索引通常交给 SAMtools 或 Picard 这类工具。所以这篇指南里会自然带出这些配套工具，但核心始终围绕 BWA 展开。另外，本文所有示例数据、样本名、项目背景均为虚构代称，仅用于演示流程，不指向任何真实项目。

2. 安装与运行环境：别在第一步就埋雷

2.1 三种安装方式的取舍逻辑

BWA 的安装方式主要有三种：源码编译、包管理器安装、容器镜像。很多人上来就apt install bwa，图省事，但这里有个容易被忽略的问题——版本。不同版本的 BWA 在mem算法的默认行为和输出细节上可能有差异，尤其是 0.7.17 之后的一些修复，对结果有实际影响。如果你的分析需要可复现，版本必须锁死。

我个人的选择顺序是：生产流程用容器镜像，本地调试用源码编译，快速验证用包管理器。容器的好处是把 BWA、SAMtools、参考基因组版本全部固化在一个镜像里，换台机器结果一致；源码编译的好处是你能控制编译选项，比如针对特定 CPU 架构开启优化，跑起来更快。

源码编译的典型流程如下，注意make之前确认 gcc 和 zlib 开发库已就位：

# 下载源码包后解压进入目录 make # 编译完成后当前目录会生成可执行文件 bwa ./bwa # 看到用法输出即编译成功

编译完把bwa复制到PATH覆盖的目录，或者直接把所在目录加进环境变量。这里有个小细节：BWA 编译产物是单个可执行文件，没有复杂的依赖树，所以拷贝到别的同架构机器上通常直接能用，这点比很多 Python 工具省心。

2.2 参考基因组的准备与版本管理

比对之前必须有参考基因组 FASTA。这里第一个大坑是版本混乱。同一个物种往往有多个组装版本，染色体命名规则也不统一（有的用chr1，有的用1，有的还带NC_编号）。如果你用的注释文件（GTF/GFF）和 FASTA 版本对不上，后面要么报错，要么更糟——不报错但结果全错。

我的做法是给每个参考基因组建一个独立目录，目录名里带版本号和来源标识，并在目录里放一个README记录下载地址和校验值。FASTA 文件本身要做两件事：确认没有重复的序列名，确认序列名里没有空格和特殊字符。有些组装文件里序列名带空格，BWA 建索引时会把空格后的内容截断，导致后续按名字取序列时对不上。

提示：建索引前用grep ">" 参考.fa | sort | uniq -d检查有没有重复序列名，有的话必须先处理，否则索引行为不可预期。

2.3 硬件与系统层面的现实考量

BWA 的索引加载进内存后，人类基因组大约占 3 到 4 GB 内存，加上比对过程中的开销，建议机器至少 16 GB 内存起步，32 GB 更从容。bwa mem是多线程的，-t参数控制线程数，但线程不是越多越好——超过物理核心数后收益递减，还可能因为内存带宽争抢反而变慢。我实测下来，线程数设成物理核心数的 80% 左右比较稳，留一点余量给系统和其他进程。

磁盘方面，SAM 文件是纯文本，体积可能是 FASTQ 的好几倍，中间文件建议放高速盘，最终结果再归档。如果磁盘 IO 是瓶颈，比对速度会被拖得很明显，这时候考虑把参考索引放在本地 SSD 而不是网络存储上。

3. 建索引：一次做好，长期受益

3.1 索引到底建了什么

BWA 的索引不是单一文件，而是一组文件，常见后缀包括.amb、.ann、.bwt、.pac、.sa。它们分别承担不同职责：.bwt是 Burrows-Wheeler 变换后的核心结构，.sa是后缀数组，.pac是打包后的参考序列，.ann和.amb记录注释和歧义碱基信息。理解这一点很重要，因为很多人看到一堆文件不知道哪个能删——答案是都不能单独删，它们是一套。

建索引的命令很简单：

bwa index -a bwtsw 参考基因组.fa

-a参数指定建索引算法，bwtsw适合较大的基因组，is适合较小的基因组。选错的后果不是报错，而是索引构建极慢或者内存爆掉。经验规则：基因组大于 500 Mb 用bwtsw，小于这个量级用is更快。人类、小鼠这类大基因组毫无疑问用bwtsw。

3.2 建索引的时间与资源预估

建索引是个重活，人类基因组在普通服务器上大概要一到几个小时不等，取决于磁盘和 CPU。这个过程内存占用也不低，峰值可能到十几 GB。所以别在笔记本上随手跑，容易卡死。我的习惯是建索引单独安排，建完做一次校验，确认文件齐全、大小合理，然后归档保存，后续所有比对直接复用，不再重复建。

这里有个实操心得：建完索引后立刻做一次小规模比对测试。拿几条 reads 跑一下bwa mem，确认能正常输出 SAM，再投入大批量任务。我见过太多次索引文件损坏或者不完整，结果跑了半天才发现比对全失败，白白浪费机时。

3.3 索引的复用与迁移

索引文件可以跨机器复用，前提是 BWA 版本一致、参考序列完全一致。迁移时把整套索引文件一起拷走，不要只拷.bwt。如果换了 BWA 大版本，建议重建索引，因为内部格式可能有变。另外，索引文件对参考序列的路径不敏感，但.ann文件里记录了序列名和长度，所以参考序列改名后索引仍然可用，只是后续按名字取序列时要保证名字一致。

4. 比对实战：参数背后的逻辑

4.1 bwa mem 的基本调用与输出解读

最常用的比对命令长这样：

bwa mem -t 8 -R '@RG\tID:sample1\tSM:sample1\tPL:ILLUMINA' 参考索引前缀 read1.fq read2.fq > output.sam

-t是线程数，-R是读组信息，这个参数经常被忽略但非常关键。读组（read group）里的ID、SM（样本名）、PL（测序平台）会写进 SAM 头部，下游的变异检测、去重、多样本合并都依赖这些信息。如果SM写错或者多个样本用了同一个SM，后面合并时会把不同样本当成同一个，结果直接废掉。这是我在实际项目里见过的高频错误，而且往往到下游才暴露，排查成本很高。

输出的 SAM 是文本格式，每行一条比对记录，字段用制表符分隔。前 11 个字段是强制的，包括 read 名、flag、参考序列名、位置、比对质量（MAPQ）、CIGAR 字符串等。CIGAR 描述比对的具体形态，比如100M表示 100 个碱基全部匹配，50M10I40M表示中间有 10 个碱基插入。看懂 CIGAR 是排查比对问题的基本功。

4.2 关键参数的取舍与计算

bwa mem的参数不少，但真正需要根据数据调整的没几个。我挑几个最常动的说。

-k是最小种子长度，默认 19。短读长数据（比如 50 bp 以下）可以适当调小，长读长或者高错误率数据可以调大。这个参数影响敏感度和速度的平衡，调小更敏感但更慢。

-A是匹配得分，-B是错配罚分。默认-A 1 -B 4，这个比例决定了比对算法对错配的容忍度。对于高保真数据，默认值通常够用；对于错误率偏高的数据，可以适当调整，但不要随意改，改了要重新评估结果。

-M参数用于标记较短的 split hits 为 secondary，这个参数在需要兼容某些下游工具时有用。很多老流程会加-M，新流程里是否加要看下游工具的要求。

-Y使用软剪切（soft clipping）处理 supplementary alignment，对结构变异检测有帮助。默认行为在不同版本间有差异，需要明确。

参数调整的核心原则是：先跑默认，看结果，再针对性调。不要一上来就堆一堆参数，那样出了问题你根本不知道是哪个参数导致的。

4.3 双端数据的处理细节

双端测序的两个文件要按顺序传给bwa mem，read1 在前，read2 在后。BWA 会自动识别配对关系，在 SAM 的 flag 字段里标记。这里有个坑：如果两个文件的 read 顺序不一致（比如一个文件被单独过滤过），配对会错乱，表现为大量 read 变成单端比对或者 flag 异常。所以比对前要确认两个 FASTQ 的 read 数量和顺序严格对应。

另一个细节是插入片段大小（insert size）。bwa mem会在比对过程中估计插入片段分布，写进 SAM 头部。如果数据质量差或者插入片段分布异常，这个估计可能不准，影响后续的配对判断。可以在比对后检查 SAM 头部的IS行，看看估计的均值和标准差是否合理。

5. 比对后的处理链路：SAM 到可用结果

5.1 排序与格式转换

BWA 输出的 SAM 是按参考序列顺序排列的，但下游很多工具要求按坐标排序。排序用 SAMtools：

samtools sort -@ 8 -o sorted.bam output.sam samtools index sorted.bam

sort会生成 BAM（二进制格式，体积小、读取快），index生成.bai索引，供快速随机访问。这两步几乎是标配，别省。排序时-@指定线程，-m可以控制每个线程的内存，大数据量时调大-m能减少临时文件数量，加快速度。

5.2 去重与质量过滤

PCR duplicates 是测序数据里的常见问题，尤其是文库构建时扩增次数多的情况。去重可以用 SAMtools 的markdup或者 Picard 的MarkDuplicates。去重前必须先排序，因为去重算法依赖坐标顺序。

质量过滤主要看 MAPQ。MAPQ 为 0 表示该 read 在参考上有多处同样好的匹配，位置不可信；MAPQ 较高表示匹配唯一。很多下游分析会过滤掉 MAPQ 低于某个阈值的 read，阈值取多少取决于应用场景——变异检测通常要求较高，表达定量可以放宽。

5.3 比对质量评估

比对完不能直接往下走，要看一眼质量。samtools flagstat给出总体统计，包括总 read 数、比对率、配对率、properly paired 比例等。比对率低通常意味着参考基因组不对、污染、或者 reads 质量太差。properly paired 比例低可能是插入片段估计不准或者文库问题。

更细的评估可以用samtools stats，输出大量指标。我习惯关注几个关键数：比对率、唯一比对率、重复率、插入片段分布。这几个数正常，后面的分析才有意义。

6. 常见报错与排查实录

6.1 索引相关报错

最常见的报错是找不到索引文件，提示类似fail to locate the index files。原因通常是索引前缀写错，或者索引文件和参考序列不在同一目录。BWA 找索引是按前缀匹配的，你给的前缀必须能对应到那套.amb/.ann/.bwt/.pac/.sa文件。

另一个是索引不完整，比如只拷了部分文件。解决办法是重新拷全套，或者重建。还有一种情况是参考序列里有非法字符，建索引时报错，需要先清理序列。

6.2 内存与线程问题

大数据量比对时可能遇到内存不足。bwa mem的内存占用和参考基因组大小、线程数、数据量都有关。如果内存吃紧，先降线程数试试，因为每个线程都有额外开销。实在不够就分批处理，把 FASTQ 拆成几份分别比对再合并。

线程相关的报错比较少见，但有时会出现线程数设得过高导致系统卡死。这时候看系统负载，如果 load average 远超核心数，就是线程开太多了。

6.3 结果异常排查

比对率异常低是最让人头疼的。排查顺序建议是：先确认参考基因组版本对不对，再确认 reads 有没有被污染，然后看 FASTQ 质量（用 FastQC 之类工具），最后检查比对参数有没有被误改。我遇到过一次比对率只有百分之十几，查了半天发现是参考基因组用错了物种，这种低级错误反而最难想到。

SAM 文件异常大或者异常小也要警惕。异常大可能是重复序列太多或者参数导致大量多重比对；异常小可能是比对提前中断。看 SAM 头部的命令行记录（@PG行）能确认实际用的参数，这个信息在排查时非常有用。

6.4 常见问题速查表

7. 流程固化与效率提升的个人经验

把上面这些串起来，一条完整的比对链路是：准备参考与索引 → 检查 FASTQ 配对与质量 →bwa mem比对 → SAMtools 排序索引 → 去重与过滤 → 质量评估。这条链路我建议直接写成脚本或者流程文件，参数集中管理，避免每次手敲出错。

效率上，几个实测有效的点：索引放本地高速盘；线程数按物理核心的八成设；SAM 转 BAM 尽早做，别拿纯文本 SAM 到处传；中间文件及时清理，磁盘满了比什么都麻烦。还有一点，比对日志一定要留，bwa mem的标准错误输出里有运行时间、线程数、处理 read 数等信息，出问题时这些日志是重要线索。

最后分享一个我踩过的坑：早期做多样本项目时，图省事在循环里用同一个-R的SM值，结果下游合并样本时全混在一起，返工重跑花了两天。从那以后，读组信息我都是按样本单独生成，并且在流程里加了一步校验，确认每个样本的SM唯一。这种看起来琐碎的检查，恰恰是保证结果可信的关键。