Enter your search text in the box above

Select a result to preview

bcftools实战教程

#garden/🌿Budding

bcftools 实战教程

官方文档：https://samtools.github.io/bcftools/bcftools.html
当前版本：1.22（2025-06）

一、基本概念

二、索引

# 建索引（默认 CSI，推荐；TBI 不支持超长染色体）
bcftools index file.vcf.gz
bcftools index -c file.vcf.gz   # 强制 CSI
bcftools index -t file.vcf.gz   # 强制 TBI（tabix 兼容）

# 同时建索引（写文件时）
bcftools view -Oz -W file.vcf > out.vcf.gz   # -W 自动建 CSI

三、位点数量统计 ⚡ 奇技淫巧

# 【最快】直接读索引，不扫描文件，瞬间返回
bcftools index -n file.vcf.gz

# 按染色体分别统计（输出：chrom / n_samples / n_records）
bcftools index -s file.vcf.gz

# 常规方法（需扫描全文件）
bcftools view -H file.vcf.gz | wc -l

# 更快的扫描方式（跳过基因型字段解析）
bcftools query -f '\n' file.vcf.gz | wc -l

⚠️ index -n 统计的是行数（VCF记录数），多等位位点一行算一个。通常这就是你想要的。

四、格式转换

# VCF → BGZ
bgzip -c file.vcf > file.vcf.gz && bcftools index file.vcf.gz

# BGZ VCF → BCF
bcftools view -Ob file.vcf.gz -o file.bcf && bcftools index file.bcf

# BCF → VCF
bcftools view -Ov file.bcf -o file.vcf

# 去除基因型，只保留位点信息（SITES-only）
bcftools view -G file.vcf.gz -Oz -o sites.vcf.gz

五、子集提取（view）

# 按样本提取
bcftools view -s sample1,sample2 file.vcf.gz -Oz -o subset.vcf.gz

# 排除样本（^前缀）
bcftools view -s ^sample1,sample2 file.vcf.gz

# 按区域提取
bcftools view -r chr1,chr2 file.vcf.gz
bcftools view -r chr1:1000000-2000000 file.vcf.gz
bcftools view -R regions.bed file.vcf.gz   # 从文件读取区域

# 按类型过滤
bcftools view --type snps file.vcf.gz      # 只保留 SNP
bcftools view --type indels file.vcf.gz    # 只保留 indel

六、过滤（filter / view -i/-e）

# 保留 PASS 位点
bcftools view -f PASS file.vcf.gz

# 按表达式过滤（-i 包含，-e 排除）
bcftools view -i 'QUAL>=30 && DP>=10' file.vcf.gz
bcftools view -e 'F_MISSING > 0.1' file.vcf.gz    # 排除缺失率>10%位点
bcftools view -i 'MAF>=0.01' file.vcf.gz           # MAF 过滤
bcftools view -i 'INFO/AC>=2' file.vcf.gz

# 软过滤（写入 FILTER 列，不删除位点）
bcftools filter -s LowQual -e 'QUAL<30' file.vcf.gz

# 常用组合：MAF + 缺失率过滤
bcftools view -i 'MAF>=0.01 && F_MISSING<=0.1' file.vcf.gz -Oz -o filtered.vcf.gz

⚠️ INFO 字段如 AC/AN 默认不随样本子集更新，需先 bcftools view -s ... | bcftools view（view 是例外，会更新部分 tag）

七、信息查询（query）

# 提取 CHROM POS REF ALT
bcftools query -f '%CHROM\t%POS\t%REF\t%ALT\n' file.vcf.gz

# 提取带样本基因型
bcftools query -f '%CHROM\t%POS[\t%GT]\n' file.vcf.gz

# 提取特定 INFO 字段
bcftools query -f '%CHROM\t%POS\t%INFO/AF\n' file.vcf.gz

# 只输出特定区域 + 特定样本的 GT
bcftools query -r chr1:1-1000000 -s sample1 -f '[%GT\n]' file.vcf.gz

# 查看 header 中的样本列表
bcftools query -l file.vcf.gz

# 统计样本数量
bcftools query -l file.vcf.gz | wc -l

八、填充/计算 INFO 标签（+fill-tags 插件）

# 填充 AC, AN, AF, MAF, NS 等常用标签
bcftools +fill-tags file.vcf.gz -- -t AC,AN,AF,MAF,NS

# 计算所有可用标签
bcftools +fill-tags file.vcf.gz -- -t all

# 管道中使用（VCF 缺少 AC/AN 时必须）
bcftools +fill-tags file.vcf.gz -- -t AC,AN | bcftools query -f '%INFO/AC\t%INFO/AN\n'

⚠️ 原始 VCF 缺少 AC/AN/AF 时，直接 query 这些字段会报错或返回空，必须先 +fill-tags

九、合并与拆分

# 合并同一样本集的多个 VCF（按染色体分割后合并）
bcftools concat -f vcf_list.txt -Oz -o merged.vcf.gz --naive   # --naive 最快，要求格式完全一致

# 合并不同样本集（merge，位点取并集）
bcftools merge file1.vcf.gz file2.vcf.gz -Oz -o merged.vcf.gz --missing-to-ref

# 多等位位点拆分为双等位（biallelic）
bcftools norm -m -any file.vcf.gz -Oz -o biallelic.vcf.gz

# 合并多等位位点
bcftools norm -m +any file.vcf.gz -Oz -o multiallelic.vcf.gz

# Indel 标准化
bcftools norm -f ref.fa file.vcf.gz -Oz -o normalized.vcf.gz

十、集合操作（isec）

# 两个 VCF 的交集
bcftools isec file1.vcf.gz file2.vcf.gz -p output_dir

# 只保留两个文件共有的位点
bcftools isec -n=2 file1.vcf.gz file2.vcf.gz -p out/

# 只保留 file1 独有的位点（差集）
bcftools isec -C file1.vcf.gz file2.vcf.gz   # 输出到 stdout

十一、统计（stats）

# 基本统计
bcftools stats file.vcf.gz > stats.txt
grep "^SN" stats.txt   # 摘要统计

# 分样本统计
bcftools stats -s sample1,sample2 file.vcf.gz

# 可视化（需 matplotlib）
plot-vcfstats stats.txt -p output_prefix/

十二、注释（annotate）

# 用 dbSNP rsID 注释 ID 列
bcftools annotate -a dbsnp.vcf.gz -c ID file.vcf.gz

# 添加 INFO 字段
bcftools annotate -a annot.tsv.gz -c CHROM,POS,INFO/MyTag -h header.txt file.vcf.gz

# 删除所有 INFO 和 FORMAT（只保留 GT）
bcftools annotate -x INFO,^FORMAT/GT file.vcf.gz

# 重命名染色体（如 1 → chr1）
bcftools annotate --rename-chrs chr_name_map.txt file.vcf.gz

十三、ROH 检测

bcftools roh --AF-tag AF file.vcf.gz -o roh.txt
# 需要 AF 标签，建议先 +fill-tags

十四、管道组合技巧

# 过滤 + 提取 + 压缩 一步完成
bcftools view -i 'MAF>=0.01' -Ou file.vcf.gz | \
  bcftools norm -m -any -Ou | \
  bcftools view -Oz -o out.vcf.gz -W

# 多染色体并行处理（SLURM array）
for i in $(seq 1 29) X; do
  sbatch --wrap="bcftools view -r chr${i} input.vcf.gz -Oz -o chr${i}.vcf.gz -W"
done

# 快速检查 VCF 是否有问题
bcftools view file.vcf.gz | head -100 | bcftools stats | grep "^SN"

十五、常见陷阱

问题 原因 解决
INFO/AC 查询结果错误 样本子集后 AC 未更新 用 +fill-tags 重新计算
-r 报错 文件未建索引 bcftools index
染色体名不匹配 1 vs chr1 --rename-chrs
concat 后位点重复 区域有重叠 加 --remove-duplicates
内存爆炸 merge 时样本太多 分批 merge
F_MISSING 报错 老版本不支持 升级到 ≥1.14

参考