🧬 Gold-Standard Germline & Somatic Variant Discovery
The Broad Institute GATK (Genome Analysis Toolkit) Best Practices pipeline is the internationally accepted standard for discovering single nucleotide polymorphisms (SNPs) and small indels from high-throughput whole genome (WGS) and whole exome (WES) sequencing data.
🔢 Phred Quality Score & Error Probability Converter
Interactive ToolPhred Quality Score Input
Computed Statistical Metrics
⚖️ GATK Hard Filtering Command Generator
Filter GeneratorVariant Class
Generated GATK VariantFiltration Expression
0. GATK Best Practices Architecture Overview
The standard pipeline progresses through three phases: 1) Pre-processing (alignment, duplicate marking, BQSR); 2) Variant discovery (HaplotypeCaller per-sample GVCFs); 3) Joint cohort genotyping and statistical filtration.
1. Reference Genome Pre-processing
Generate FASTA index (.fai) and sequence dictionary (.dict) required by all GATK tools.
# Index reference FASTA samtools faidx hg38.fasta # Create Sequence Dictionary gatk CreateSequenceDictionary -R hg38.fasta -O hg38.dict
2. BWA-MEM Alignment with Read Group Tags (@RG)
Align reads with BWA-MEM while injecting mandatory @RG tags identifying sample (SM), library (LB), and platform (PL).
bwa mem -t 8 \ -R "@RG\tID:FLOWCELL_L1\tSM:NA12878\tLB:LIB_1\tPL:ILLUMINA" \ hg38.fasta sample_R1.fq.gz sample_R2.fq.gz | \ samtools sort -@ 8 -m 2G -o sample.sorted.bam -
3. MarkDuplicatesSpark
Identify PCR and optical duplicates with multithreaded Spark acceleration.
gatk MarkDuplicatesSpark \ -I sample.sorted.bam \ -O sample.dedup.bam \ -M sample.dup_metrics.txt \ --spark-master local[8]
4. Base Quality Score Recalibration (BQSR)
Correct systematic platform errors in sequencer Phred scores by training a machine learning model on known polymorphic sites (dbSNP, 1000 Genomes).
# Step 1: Build Recalibration Model gatk BaseRecalibrator \ -R hg38.fasta \ -I sample.dedup.bam \ --known-sites dbsnp_146.hg38.vcf.gz \ --known-sites Mills_and_1000G_gold_standard.indels.hg38.vcf.gz \ -O recal_data.table # Step 2: Apply Recalibration Table gatk ApplyBQSR \ -R hg38.fasta \ -I sample.dedup.bam \ --bqsr-recal-file recal_data.table \ -O sample.bqsr.bam
5. Germline HaplotypeCaller in GVCF Mode
Perform de novo local haplotype assembly around active regions and output comprehensive Genomic VCF records (-ERC GVCF).
gatk HaplotypeCaller \ -R hg38.fasta \ -I sample.bqsr.bam \ -O sample.g.vcf.gz \ -ERC GVCF
6. Multi-Sample Joint Genotyping (GenomicsDBImport)
Import GVCFs into a high-performance TileDB database and joint-call genotypes across the cohort to distinguish homozygous reference from missing data.
# Import into GenomicsDB gatk GenomicsDBImport \ --genomicsdb-workspace-path cohort_genomics_db \ -L chr20 \ -V sample1.g.vcf.gz -V sample2.g.vcf.gz # Joint Genotyping gatk GenotypeGVCFs \ -R hg38.fasta \ -V gendb://cohort_genomics_db \ -O cohort_raw_variants.vcf.gz
7. Somatic Variant Calling (Mutect2)
Detect low-frequency somatic mutations in matched Tumor-Normal samples using a Panel of Normals (PoN) to suppress recurring sequencing artifacts.
gatk Mutect2 \ -R hg38.fasta \ -I tumor.bam -tumor TUMOR_SAMPLE \ -I normal.bam -normal NORMAL_SAMPLE \ --panel-of-normals pon.vcf.gz \ -O somatic_raw.vcf.gz
8. Variant Filtration (Hard Filtering vs VQSR)
Filter raw candidate variants to remove false positives caused by strand bias, low coverage depth, or misaligned reads.
# Apply SNV Hard Filters gatk VariantFiltration \ -R hg38.fasta \ -V cohort_raw_snps.vcf.gz \ -O cohort_filtered_snps.vcf.gz \ --filter-name "GATK_SNP_HardFilter" \ --filter-expression "QD < 2.0 || FS > 60.0 || MQ < 40.0 || MQRankSum < -12.5 || ReadPosRankSum < -8.0 || SOR > 3.0"
9. Functional Annotation with Ensembl VEP & SnpEff
Annotate functional consequences (missense, nonsense, frameshift, splice donor/acceptor) and population frequencies (gnomAD).
vep -i cohort_filtered_snps.vcf.gz \ --assembly GRCh38 --cache --offline \ --everything --vcf -o cohort_annotated.vcf.gz
10. VCF Querying & Extraction with BCFtools
Extract passing high-impact mutations for clinical interpretation or population genetic analysis.
# Query passing variants with high impact bcftools view -f PASS cohort_annotated.vcf.gz | \ bcftools query -f '%CHROM\t%POS\t%REF\t%ALT\t%QUAL\t%INFO/CSQ\n' > impact_variants.tsv
🔀 Decision Architecture: Variant Filtration Strategy
Decision TreeFor cohorts with under 30 samples, VQSR lacks sufficient statistical power to train machine learning Gaussian mixture models. Apply standard hard filter thresholds (QD, FS, MQ, SOR).