GATK变异检测示例

# === GATK最佳实践流程 ===
# 1. 标记重复（MarkDuplicates）
gatk MarkDuplicates \
    -I sample.sorted.bam \
    -O sample.dedup.bam \
    -M sample.metrics.txt \
    --REMOVE_DUPLICATES false  # 标记但不删除，GATK推荐
samtools index sample.dedup.bam
# 2. 碱基质量校正（BQSR）
# 需要已知变异位点（如dbSNP）作为训练集
gatk BaseRecalibrator \
    -I sample.dedup.bam \
    -R reference.fa \
    --known-sites dbsnp.vcf.gz \
    -O sample.recal.table
gatk ApplyBQSR \
    -R reference.fa \
    -I sample.dedup.bam \
    --bqsr-recal-file sample.recal.table \
    -O sample.recal.bam
# 3. 变异检测（HaplotypeCaller）
gatk HaplotypeCaller \
    -R reference.fa \
    -I sample.recal.bam \
    -O sample.g.vcf.gz \
    -ERC GVCF  # 输出gVCF格式（适合多样本联合分析）
# 4. 多样本联合基因型（如果有多个样本）
# GenomicsDBImport导入gVCF
gatk GenomicsDBImport \
    -V sample1.g.vcf.gz \
    -V sample2.g.vcf.gz \
    -V sample3.g.vcf.gz \
    --genomicsdb-workspace-path cohort_db \
    -L intervals.list
# GenotypeGVCFs输出最终VCF
gatk GenotypeGVCFs \
    -R reference.fa \
    -V gendb://cohort_db \
    -O cohort.vcf.gz
# 5. 变异质控（VQSR，可选）
gatk VariantRecalibrator \
    -R reference.fa \
    -V cohort.vcf.gz \
    --resource:hapmap,known=false,training=true,truth=true,prior=15.0 hapmap.vcf.gz \
    --resource:omni,known=false,training=true,truth=false,prior=12.0 omni.vcf.gz \
    -an QD -an MQ -an MQRankSum -an ReadPosRankSum \
    -mode SNP \
    -O cohort.recal \
    --tranches-file cohort.tranches
gatk ApplyVQSR \
    -R reference.fa \
    -V cohort.vcf.gz \
    --recal-file cohort.recal \
    --tranches-file cohort.tranches \
    --truth-sensitivity-filter-level 99.5 \
    -mode SNP \
    -O cohort.vqsr.vcf.gz
# 6. 硬过滤（如果VQSR不适用，如小样本）
gatk VariantFiltration \
    -R reference.fa \
    -V cohort.vcf.gz \
    -O cohort.filtered.vcf.gz \
    --filter-expression "QD < 2.0 || MQ < 40.0 || FS > 60.0 || SOR > 3.0" \
    --filter-name "basic_filters"