WhatsHap — 基于读段的单倍型分型工具¶

一句话说明¶

WhatsHap（v2.8）把变异检测出来的VCF中"杂合变异"归类到不同的染色体拷贝上——就像区分"这个变异来自父亲那条染色体还是母亲那条染色体"，对长读数据尤其强大，短读也支持。

安装与配置¶

# 方法1：conda安装（推荐，v2.8）
conda create -n whatshap python=3.10
conda activate whatshap
conda install -c bioconda whatshap   # 安装最新版（v2.8）

# 验证安装
whatshap --version       # 查看版本号（应显示2.8）
whatshap phase --help    # 分型命令帮助

# 方法2：pip安装
pip install whatshap     # Python生态安装

# 依赖：需要samtools（用于读取BAM）
conda install -c bioconda samtools

核心用法¶

基本单倍型分型（phase）¶

# 最基本用法：用比对reads对VCF中的变异进行单倍型分型
whatshap phase \
    --output phased.vcf.gz \           # 输出带单倍型信息的VCF
    variants.vcf.gz \                  # 输入VCF（已检测的变异）
    reads.bam                          # 输入比对BAM

# 多样本联合分型
whatshap phase \
    --output phased.vcf.gz \
    variants.vcf.gz \
    sample1.bam sample2.bam            # 可多个BAM

# 参考基因组（某些情况下需要，特别是CRAM输入）
whatshap phase \
    --output phased.vcf.gz \
    --reference reference.fasta \      # 参考基因组（CRAM输入必须）
    variants.vcf.gz \
    reads.cram                         # CRAM格式也支持

单倍型标签（haplotag）——给reads打标签¶

# 根据分型结果，给每条read打上HP标签（HP:i:1或HP:i:2）
whatshap haplotag \
    --output haplotagged.bam \         # 输出带HP标签的BAM
    --reference reference.fasta \      # 参考基因组
    phased.vcf.gz \                    # 已分型VCF
    reads.bam                          # 输入BAM

# 建立新BAM的索引
samtools index haplotagged.bam         # 索引标签后的BAM

单倍型统计分析¶

# 统计分型结果
whatshap stats \
    --tsv stats.tsv \                  # 输出TSV格式统计
    phased.vcf.gz                      # 已分型VCF

# 对比两种分型结果（评估质量）
whatshap compare \
    --sample SAMPLE_NAME \             # 样本名
    phased1.vcf.gz phased2.vcf.gz     # 对比两个分型结果

参数详解¶

命令	参数	说明
`phase`	`--output`	输出VCF文件路径
`phase`	`--reference`	参考基因组（CRAM输入时必须）
`phase`	`--sample`	指定要分型的样本名（多样本VCF）
`phase`	`--indels`	v2.8已默认开启Indel分型
`phase`	`--only-snvs`	仅分型SNV（跳过Indel）
`phase`	`-B` / `--block-cut-sensitivity`	分型块敏感度0-5（值越小块越长，默认4）
`haplotag`	`--output`	输出BAM文件
`haplotag`	`--ignore-read-groups`	忽略RG信息（单样本BAM时简化）
`stats`	`--tsv`	输出TSV格式统计表

实战案例¶

# 完整单倍型分型流程（长读ONT数据）

REF="hg38.fasta"
# 假设已用Clair3检测出变异
VCF="clair3_output/merge_output.vcf.gz"
BAM="ont_reads.bam"

# 1. 单倍型分型
whatshap phase \
    --output phased_variants.vcf.gz \
    --reference $REF \
    $VCF \
    $BAM

# 建索引
tabix -p vcf phased_variants.vcf.gz

# 2. 验证分型结果（查看GT字段是否有"|"分隔的单倍型）
zcat phased_variants.vcf.gz | grep -v "^#" | \
    awk '{print $9,$10}' | head -20
# 未分型：0/1  (斜杠表示未知哪条染色体)
# 已分型：0|1  (竖线表示已确定单倍型：染色体1为ref，染色体2为alt)

# 3. 统计分型覆盖率
whatshap stats phased_variants.vcf.gz
# 重点关注：
# - phased heterozygous SNPs: 多少杂合SNP被分型
# - phased fraction: 分型比例（长读通常>80%，短读可能<50%）
# - Longest block: 最长单倍型块长度

# 4. 给reads打上单倍型标签（用于后续分析如甲基化分析）
whatshap haplotag \
    --output haplotagged_reads.bam \
    --reference $REF \
    phased_variants.vcf.gz \
    $BAM

samtools index haplotagged_reads.bam

# 5. 查看标签后的reads（检查HP标签）
samtools view haplotagged_reads.bam | head -5 | \
    grep -oP "HP:i:\K[0-9]+"    # 提取HP标签值（1或2）

# 6. 基于单倍型提取特定reads（例如仅提取单倍型1的reads）
samtools view -d HP:1 haplotagged_reads.bam -o haplotype1_reads.bam
samtools index haplotype1_reads.bam

常见报错与解决¶

报错1：ERROR: No variants could be phased - 原因：VCF中无杂合变异（GT=0/0或1/1），或BAM覆盖度太低 - 解决：检查VCF杂合变异数bcftools view -g het phased.vcf.gz | grep -v "^#" | wc -l；确认测序深度

报错2：VCF file must be indexed - 原因：输入VCF缺少tabix索引 - 解决：tabix -p vcf variants.vcf.gz（VCF文件需先用bgzip压缩）

报错3：分型比例很低（<20%） - 原因：使用短读数据时，read长度不足以跨越两个杂合位点 - 解决：短读分型比例低是正常现象；使用长读（ONT/PacBio）可大幅提高分型率；或使用群体信息工具（SHAPEIT）补充

速查表¶

命令	说明
`whatshap phase --output out.vcf.gz vars.vcf.gz reads.bam`	基本单倍型分型
`whatshap haplotag --output tagged.bam phased.vcf.gz reads.bam`	给reads打单倍型标签
`whatshap stats phased.vcf.gz`	分型统计（覆盖率、块长度等）
`whatshap compare phased1.vcf.gz phased2.vcf.gz`	对比两种分型结果
`0\\|1` vs `0/1`	已分型（竖线）vs 未分型（斜杠）
HP:i:1 / HP:i:2	BAM read中的单倍型标签
`-B 0`	最激进的分型（最长块，可能错误多）
`-B 5`	最保守的分型（最短块，错误少）