遇到两个项目的比对质量差一截，很多人先查测序质量或比对参数。有时答案藏在参考序列里。我们给 BWA-MEM2 两条各长两千二百碱基的人工参考序列，其中五百碱基完全相同。从共享区取一条一百碱基的读段，面对完整参考时，它能全长匹配，定位质量却为零。另建一份索引，只保留第一条参考，同一条读段以同样的分数匹配到同一段序列，定位质量升到六十。下面用另外三条已知来源的读段做对照，判断这个变化究竟能证明什么。
比对得分、错配数与定位质量分别量什么
读取一行 SAM，先把三个问题分开。第一，这段序列在报告的位置上有多吻合。第二，参考中是否还有接近的候选位置。第三，写在输出里的这条记录是主比对、次要比对，还是未比对。第一项可以结合比对得分、编辑距离和比对操作检查。第二项需要结合定位质量、候选得分和重复区背景判断。第三项由标志位说明。碱基质量则来自测序环节，描述具体碱基的判读。这几项都用数字表示，含义和适用的位置却各自不同。SAM 格式规范与字段定义
本例控制原始碱基和质量字符串始终相同。所有查询都是单端，固定一条线程，同时要求报告找到的其他单端比对。这样，出现第二行时能直接查看它的坐标和得分。默认输出只保存部分候选时，表里少一行并不能证明参考只存在一个吻合位置。即便打开这个选项，启发式搜索也没有给出所有可能位置的数学枚举。定位质量是在当前参考和运行条件下对定位歧义的估计，参考中未收录的序列自然不会出现在竞争者清单里。BWA-MEM2 作者的安装与运行说明
实跑经验｜共享区内部与跨边界读段给出不同答案
安装 BWA-MEM2 2.2.1 后，我们独立构建双序列参考与单序列参考的索引。第一条读段完全落在共享区，第二、第三条分别落在两段独有区域，第四条跨越第一条参考的共享区边界。为了能用纸笔推算结果，参考由固定随机种子构造，刻意保留一个较长的相同片段。实际比对都重新从同一份 FASTQ 读取四条查询。下表位置按 SAM 的一基坐标解释，主比对表示程序选作主要报告的记录。
| 读段 | 参考同时含 A、B | 参考只含 A |
|---|---|---|
| 完全在共享区 | B 的 871 位为主，A 的 871 位为次，两者 MAPQ 0 | A 的 871 位为主，MAPQ 60 |
| A 独有区 | A 的 141 位，MAPQ 60 | 同一位置，MAPQ 60 |
| B 独有区 | B 的 141 位，MAPQ 60 | 未比对，FLAG 4 |
| A 与共享区交界 | A 上全长，另有 B 上 55 碱基的次比对 | A 上全长，MAPQ 60 |
在双序列参考里，共享区读段两条记录的比对得分都为一百，编辑距离都为零，主记录还给出了同为一百的次优得分。原始切片客观上取自第一条参考，本次主记录偏偏写在第二条。由于两处序列逐个碱基相同，单看这一百个碱基，读者无从恢复抽样时的标签。主记录在这里方便保存一个代表位置，把它当作来源鉴定会产生错误。去掉第二条参考后，主记录的比对得分仍为一百，次优得分字段变成零。这一次增加的定位质量来自可供比较的位置变少。
跨边界那条读段提供了更接近实际诊断的例子。完整参考中，第一条参考上得到一百碱基全长吻合，定位质量六十。第二条上只对上共享部分的五十五个碱基，次记录的操作串为 45H55M，得分五十五。这里的硬剪切意味着次记录的序列字段不包含被剪去的四十五个碱基，原始 FASTQ 中它们仍然存在。普通的匹配操作符允许匹配或错配，要确认这个例子真的逐碱基一致，还应查看编辑距离或直接对照序列。跨出共享区的那些碱基提供了足以区分两种位置的信息。BWA 作者的选项与输出说明
bwa-mem2 index both/ref.fa
bwa-mem2 index only_A/ref.fa
bwa-mem2 mem -t 1 -a both/ref.fa query.fq > both.sam
bwa-mem2 mem -t 1 -a only_A/ref.fa query.fq > only_A.sam
awk -F '\t' '$1=="shared" || $1=="junctionA"' both.sam
比较两批数据前，怎样把参考的差异查出来
先从比对结果的头部导出参考序列名称与长度，检查是否包含同一组常染色体、替代序列、未定位片段及实验所需的其他序列。名称和长度一致也可能隐藏内部碱基差异，因此应保留参考文件的校验值，并核对索引实际来自哪份文件。两个目录都叫人类参考，或者流程里都写了同一个版本简称，不足以说明内容一致。对高度同源基因、假基因和外源构建序列，哪一方缺少竞争序列会直接改变判读。本例的独有区读段也说明，缩小参考有时还会让原本匹配很好的证据彻底消失。
抽取若干可疑读段，带上读段名、两端信息、主次标志、定位质量、得分、编辑距离和操作串，用同一套参数对两个索引交叉比对。将变化按原因展开，比只比较两条比对率曲线更容易行动。若是大量重复区读段获得高分，应优先盘点消失的候选参考。若参考不变而读长缩短，检查原本能跨越重复边界的碱基是否在质控时被剪去。对于双端数据，还应看另一端落点、方向与片段长度能否提供额外限制，不能直接把这四条单端读段的输出推广到全样本。
下游设置定位质量门槛之前，建议把被保留和被剔除的目标位点分别抽样，检查高同源区域是否总被同一规则排除。研究目标若要求区分两个高度相似的拷贝，统计独有区覆盖和跨界读段数通常更有解释力。若问题允许合并报告某一基因家族，可另行制定家族汇总方法并交代归属不确定性。这样得到的计数才与研究问题相符，而某一次出现六十这个分数，并不自动使两个拷贝变得可区分。
参考与参数核对后，若要把同一份索引用于整批样本重比对，可在云生信的生物信息分析服务器中集中保存参考、命令和结果。
