GATK Reference

HaplotypeCaller is designed for germline variant calling, detecting inherited SNPs and Indels present in the sample's genome. Mutect2 is designed for somatic variant calling, detecting mutations acquired in tumor cells that are absent from matched normal tissue. HaplotypeCaller supports GVCF mode for scalable joint genotyping across many samples, while Mutect2 supports tumor-normal pair analysis with Panel of Normals and contamination estimation.

Base Quality Score Recalibration (BQSR) corrects systematic errors in the quality scores assigned by the sequencing machine. BaseRecalibrator builds a recalibration model using known variant sites (dbSNP, Mills indels, 1000G SNPs) to learn error patterns, then ApplyBQSR writes a new BAM with corrected quality scores. This improves variant calling accuracy by ensuring quality scores accurately reflect the true probability of sequencing errors.

VQSR (Variant Quality Score Recalibration) is preferred for whole-genome and large whole-exome datasets because it uses machine learning with known variant training resources (HapMap, Omni, 1000G, dbSNP) to learn the profile of true variants. Hard filtering with VariantFiltration is recommended when you have too few variants for VQSR to build a reliable model, such as small targeted panels or single-gene studies. The reference includes both SNP and Indel hard filter thresholds.

First, HaplotypeCaller runs on each sample individually in GVCF mode (-ERC GVCF), producing per-sample GVCF files that include reference confidence blocks. Then GenomicsDBImport consolidates all GVCFs into a GenomicsDB workspace. Finally, GenotypeGVCFs performs joint genotyping across all samples simultaneously. This approach is scalable because new samples can be added to GenomicsDB without reprocessing existing ones.

Key INFO fields include QD (quality by depth), FS (Fisher strand bias), SOR (strand odds ratio), MQ (mapping quality), MQRankSum and ReadPosRankSum (rank sum tests). Key FORMAT fields include GT (genotype: 0/0, 0/1, 1/1), AD (allelic depth for ref and alt), DP (read depth), GQ (genotype quality), and PL (Phred-scaled genotype likelihoods). These annotations are used by both hard filtering and VQSR.

The RNA-seq variant calling pipeline differs from DNA-seq in several preprocessing steps: 1) Use STAR 2-pass alignment instead of BWA-MEM2, 2) Run MarkDuplicates, 3) Run SplitNCigarReads to handle splice junction reads with N-in-CIGAR, 4) Perform BQSR, 5) Run HaplotypeCaller. SplitNCigarReads is the critical RNA-seq-specific step that splits reads spanning introns into separate supplementary alignments.

A Panel of Normals (PoN) is a VCF created from Mutect2 calls on multiple normal (non-tumor) samples. It captures recurrent technical artifacts and germline variants. During somatic calling, Mutect2 uses the PoN (--panel-of-normals) to filter out these artifacts, improving the specificity of true somatic variant detection. CreateSomaticPanelOfNormals merges the individual normal VCFs into the final PoN.

Use --java-options "-Xmx16g" to set Java heap memory (adjust based on your data size and available RAM). For threading, HaplotypeCaller supports --native-pair-hmm-threads for PairHMM parallelization. MarkDuplicatesSpark uses Spark-based parallelism with --conf spark.executor.cores. The --tmp-dir flag specifies a temporary directory for intermediate files, useful when the default /tmp has limited space.