4. 碱基质量分数校正失败：BaseRecalibrator报错'No known sites file found'或'Reference mismatch'

4.1 BQSR 是个好东西，但脾气不小

BaseRecalibrator 是 GATK 碱基质量分数校正的核心工具。说白了，它就是利用已知的 SNP 位点，去纠正测序仪给的质量分数，让后续的变异检测更准。

但很多新手一上来就被两个报错给卡住了："No known sites file found" 和 "Reference mismatch"。我刚开始带项目时也栽过跟头，后来才发现，这些问题大多出在数据准备阶段。你想想看，校正用的参考数据都没准备好，工具能不罢工吗？

4.2 "No known sites file found"——known sites 没找到

这个报错字面意思很清楚：已知位点文件没找到。但背后的原因，我总结下来主要有三种。

4.2.1 路径写错了

最常见的，就是 --known-sites 参数给的路径不对。绝对路径少了个斜杠，相对路径当前目录没对上，这些我都会在排查时先检查。

# 错误示例：路径写错了
gatk BaseRecalibrator \
   -R /data/reference/hg38.fasta \
   -I sample.bam \
   --known-sites known_sites.vcf.gz \        # 没写完整路径，当前目录又没有这个文件
   -O recal_data.table

# 正确示例：使用绝对路径
gatk BaseRecalibrator \
   -R /data/reference/hg38.fasta \
   -I sample.bam \
   --known-sites /data/known_sites/dbsnp_138.hg38.vcf.gz \
   -O recal_data.table

4.2.2 文件格式不对

BaseRecalibrator 要求 known sites 文件必须是 VCF 或 BCF 格式，而且最好带索引。你扔一个文本文件或者 gzip 压缩但没索引的 VCF 给它，它照样报错。

# 检查文件格式
file /data/known_sites/dbsnp_138.hg38.vcf.gz
# 输出应该是：gzip compressed data, ... VCF variant call format

# 检查索引文件
ls -la /data/known_sites/dbsnp_138.hg38.vcf.gz.tbi
# 如果没找到索引文件，先创建
gatk IndexFeatureFile -I /data/known_sites/dbsnp_138.hg38.vcf.gz

4.2.3 权限问题

这个比较隐蔽。尤其在集群上跑任务时，known sites 文件在某个用户目录下，其他用户读不了。我记得有一次帮同事排查，折腾了半天，最后发现是文件权限设成了 700。

# 检查文件权限
ls -l /data/known_sites/dbsnp_138.hg38.vcf.gz
# 确保其他用户有读权限（至少 644）
chmod 644 /data/known_sites/dbsnp_138.hg38.vcf.gz
chmod 644 /data/known_sites/dbsnp_138.hg38.vcf.gz.tbi

4.3 "Reference mismatch"——参考基因组对不上

这个报错更让人头疼，因为它涉及多个文件之间的版本一致性。BaseRecalibrator 会校验 BAM 头部记录的参考基因组与当前提供的参考是否一致。

4.3.1 BAM 头部与参考版本不一致

BAM 文件在比对时用了 hg19，而你给 BaseRecalibrator 传了 hg38，那肯定报错。这个问题我在项目中遇到过好几次，尤其是从公共数据库下载的数据，参考版本五花八门。

# 查看 BAM 头部的参考基因组信息
samtools view -H sample.bam | grep @SQ

# 输出示例（hg19）
@SQ SN:chr1 LN:249250621
@SQ SN:chr2 LN:243199373
...

# 参考基因组 dict 文件的 SN 行必须与 BAM 头部一致
cat /data/reference/hg19.dict | grep @SQ
@SQ SN:chr1 LN:249250621 UR:/data/reference/hg19.fasta M5:...

4.3.2 参考基因组 fasta 与 dict 不匹配

有时候参考基因组 fasta 文件和它的 .dict 字典文件对不上，比如 fasta 更新了但 dict 没重新生成。BaseRecalibrator 加载参考时会校验 dict，不一致就报错。

# 重新生成 dict 文件
samtools dict /data/reference/hg38.fasta > /data/reference/hg38.dict

# 或者用 GATK 的 CreateSequenceDictionary
gatk CreateSequenceDictionary \
   -R /data/reference/hg38.fasta \
   -O /data/reference/hg38.dict

4.4 排查思路速查表

我把两种报错的关键区别和排查方向整理了一下，方便你快速定位问题。

4.5 避坑指南与最佳实践

这些年被这两个报错折磨过不少次，我总结了几条经验，可以帮大家少走弯路。

• 统一参考基因组版本：整个项目从比对到变异检测，全程只用同一个参考基因组版本（比如 hg38 或 GRCh38）。不要比对用 hg19，校正用 hg38。
• known sites 文件集中管理：我建议你建一个固定的目录，比如 /data/known_sites/，把所有的已知位点文件按版本分类放好，写路径时直接复制粘贴，不容易出错。
• 写一个简单的校验脚本：在运行 BaseRecalibrator 之前，先跑一个检查脚本，确认文件存在、格式正确、参考基因组一致。花 5 分钟做这件事，能省下半天排查时间。

# 一个简单的校验脚本思路（bash 伪代码）
#!/bin/bash
KNOWN_SITES="/data/known_sites/dbsnp_138.hg38.vcf.gz"
BAM="sample.bam"
REF="/data/reference/hg38.fasta"

# 检查 known sites
if [ ! -f "$KNOWN_SITES" ]; then
    echo "ERROR: Known sites file not found!"
    exit 1
fi

# 检查 BAM 头部与参考 dict 的 SN 行是否一致
samtools view -H "$BAM" | grep @SQ | cut -f2 > bam_sn.txt
grep @SQ "$REF.dict" | cut -f2 > ref_sn.txt
if ! diff bam_sn.txt ref_sn.txt > /dev/null; then
    echo "ERROR: Reference mismatch between BAM and FASTA dict!"
    exit 1
fi

echo "All checks passed. Ready to run BaseRecalibrator."

其实这两个报错都不难解决，只要你在准备数据时多留个心眼，把文件路径、格式、参考版本这三样东西对齐，BaseRecalibrator 就能顺利跑起来。说到底，GATK 的报错信息已经写得很清楚了，关键是你愿不愿意静下心来看一眼日志，找到真正的原因。