Variant Calling Workflow: Design Decisions Before You Start
A variant calling workflow maps raw sequencing reads to a reference genome, identifies differences, and produces a list of candidate variants. The design decisions you make before running any tool determine whether those candidates reflect true biology or technical artifacts. This guide is for bioinformatics trainees and lab researchers who need to build a defensible germline or somatic variant calling pipeline. We focus on the choices that separate technical confidence from biological confidence: alignment strategy, base quality calibration, joint versus per sample calling, filtration thresholds, and validation approaches.
At a Glance
| Decision Point | Options | Key Consideration | Impact on Confidence |
|---|---|---|---|
| Read alignment | BWA MEM, Bowtie2, minimap2 | Accuracy for indels and structural variants, computational cost | Technical: poor alignment inflates false positives |
| Base quality recalibration | GATK BQSR, no recalibration | Systematic instrument errors | Technical: recalibration reduces false SNPs in high coverage regions |
| Calling approach | Single sample GVCF, joint genotyping, cohort or trio calling | Population priors, batch effects, rare variant discovery | Biological: joint calling boosts sensitivity for low frequency alleles |
| Variant filtration | Hard filters (QD, FS, SOR) vs VQSR | Training set availability, cohort size | Technical: overfiltering removes real variants, underfiltering adds noise |
| Validation | PCR Sanger, orthogonal sequencing, in silico replicas | Cost, throughput, variant type | Biological: validation confirms true positives and reveals systematic errors |
Decision Criteria
Study Type Determines the Core Pipeline
Germline variant calling expects heterozygous and homozygous sites with allele frequencies near 0.5 or 1.0. Somatic calling must detect subclonal variants at allele fractions as low as 1% to 5% and must cope with tumor heterogeneity and sample contamination. The EMBL EBI training materials emphasize using a specialized somatic caller such as Mutect2 for cancer samples and HaplotypeCaller for germline projects. Mixing these workflows is a common source of misinterpreted results.
Sample Size and Read Depth
Joint calling across many samples improves specificity for rare variants because population level information helps distinguish sequencing errors from true polymorphisms. In small cohorts (fewer than 30 samples) the gain is minimal and the extra computational overhead may not be justified. Deeply sequenced targeted panels (over 500x) require different filtering thresholds than whole exome data (50x to 100x). The Galaxy Training Network provides hands on examples of adjusting hard filter parameters for coverage.
Reference Genome Selection
Using the most current reference build (GRCh38 or T2T CHM13) reduces the number of apparent variants that are actually reference errors. If you must use an older build for compatibility with clinical databases, document the lift over strategy and expect higher false positive rates in regions that were refined in newer releases. The NCBI Bookshelf contains authoritative guidance on reference genome choices.
Practical Workflow or Implementation Sequence
Step 1: Read Alignment and Preprocessing
Map paired end FASTQ files with BWA MEM for Illumina data or minimap2 for long reads. Sort the output BAM files by coordinate and mark duplicate reads. Duplicate marking is not optional for accurate variant allele fraction estimation. According to a study of a custom myeloid assay, deduplication reduced false positive calls by more than 40% in high coverage amplicon data [6].
Step 2: Base Quality Score Recalibration
Run GATK BaseRecalibrator using known variant sites from dbSNP and a set of known indels from Mills and 1000 Genomes. Recalibration corrects for systematic errors in base quality scores that vary by machine cycle, read group, and dinucleotide context. The Bioconductor project provides R packages such as Rsamtools that can help inspect recalibration tables. Skipping this step leads to inflated variant quality scores and more false positives in the final VCF.
Step 3: Variant Calling and GVCF Generation
For germline analysis, run HaplotypeCaller in GVCF mode on each sample. For somatic analysis, run Mutect2 on tumor normal pairs. GVCF mode produces a record for every genomic position, not just variant sites. This is essential for joint genotyping. Do not use legacy caller modes that emit records only at variant positions if you plan to call in a cohort.
Step 4: Joint Genotyping
Combine sample GVCFs with CombineGVCFs, then run GenotypeGVCFs. Joint genotyping uses the allele counts across all samples to refine genotype likelihoods. In a study of whole exome sequencing for orofacial clefts in a Sub Saharan African cohort [8], joint calling improved detection of rare variants that would have been missed in single sample calls. For somatic projects, joint calling is typically not performed because tumor normal pairs are processed independently.
Step 5: Variant Filtration
Apply hard filters when you lack a large training set. The GATK recommended hard filter parameters for SNPs are QD less than 2.0, FS greater than 60.0, MQ less than 40.0, and SOR greater than 3.0. For indels use QD less than 2.0, FS greater than 200.0, and ReadPosRankSum less than negative 20.0. These thresholds work well for high quality whole exome data but should be adjusted for targeted panels. For variant quality score recalibration (VQSR), you need at least two truth sets (e.g., HapMap and Omni) and a large cohort of 30 or more samples. VQSR often outperforms hard filters on whole genome data but can overfilter if the training resources are mismatched to the study population.
Step 6: Annotation and Validation
Annotate variants with SnpEff or VEP to predict functional impact. Validate a subset of clinically significant calls with orthogonal methods such as Sanger sequencing or digital droplet PCR. In a diagnostic setting, the segmental copy number variant detection demonstrated that amplicon based NGS panels require specific validation for CNVs because read depth based methods have high false positive rates in GC rich regions [11]. For research, you can use technical replicates to estimate false discovery rates without orthogonal assays.
Common Mistakes
Mistake 1: Skipping base recalibration. Many users assume that modern sequencers produce accurate base qualities. In practice, systematic errors persist and recalibration markedly improves variant quality scores. A review of AI driven neoantigen identification [9] noted that unrecalibrated data introduced artifactual neopeptides, wasting downstream validation effort.
Mistake 2: Using default filters without inspection. The GATK best practices filters are designed for 30x whole exome data. Running them unchanged on 500x targeted panel data removes many real variants. You must calculate the distribution of filter parameters on your own data and set thresholds accordingly.
Mistake 3: Confusing technical confidence with biological confidence. A variant that passes all quality filters may still be a germline polymorphism in a somatic study, or an artifact from FFPE deamination. Technical confidence means the sequencing reads consistently support the call. Biological confidence requires allele frequency, population frequency, and functional evidence to align with the hypothesis.
Mistake 4: Ignoring copy number variation. Standard single nucleotide variant pipelines miss large events. In a study of developmental and epileptic encephalopathy [10], copy number variants accounted for 12% of diagnostic findings. If your question includes structural variation, you must integrate a CNV caller such as GATK gCNV or Canvas.
Limits and Uncertainty
No variant caller achieves 100% sensitivity and specificity. False negatives occur in low complexity regions, homopolymer runs, and areas with poor mappability. The NCBI Sequence Read Archive contains billions of reads that can be used to benchmark callers, but benchmark results from one sequencing platform or library preparation method do not transfer directly to another. Somatic variant calling has higher uncertainty because low allele frequency variants are difficult to distinguish from sequencing noise. A study of antibiotic resistance gene tracking [7] showed that even with deep sequencing, certain mobile elements were missed because they reside in repetitive regions that are collapsed during alignment. Always report the regions of the genome that were excluded from analysis and provide confidence intervals for variant allele fractions.
Frequently Asked Questions
What is the difference between joint calling and per sample calling?
Joint calling combines genotype likelihoods across multiple samples to infer variants using population information. Per sample calling treats each sample independently. Joint calling improves sensitivity for rare variants but requires GVCF files and more compute time. It is the recommended approach for germline cohort studies with more than 10 samples.
How do I choose between hard filters and VQSR?
Use hard filters when you have fewer than 30 samples or when you lack reliable truth sets for your specific platform. VQSR works well for whole genome data with many samples and known variant resources. For targeted panels, hard filters with empirically adjusted thresholds often give better results.
What is the role of base quality score recalibration?
Recalibration corrects systematic errors in base quality scores that are introduced by the sequencer. Without recalibration, variant quality scores are systematically inflated, leading to more false positive calls. The correction is essential for achieving the specificity reported in GATK best practices.
How do I distinguish technical artifacts from true variants?
Compare the variant to a panel of normals to subtract sequencing artifacts. Check strand bias, read position bias, and mapping quality. Validate a random subset with an orthogonal method. For research projects, require at least two independent read pairs supporting the call and ensure the variant is present in both forward and reverse strand reads.
References and Further Reading
- NCBI Bookshelf: Genome Analysis Toolkit Best Practices (authoritative reference for GATK workflows)
- EMBL EBI Training: Variant Calling and Annotation (hands on courses with practical examples)
- Galaxy Training Network: Variant Analysis (open workflows and tutorials)
- Bioconductor: VariantAnnotation Package (R tools for manipulating and annotating VCF files)
- NCBI Sequence Read Archive (public repository for benchmarking data)
- Performance Evaluation of a Custom Myeloid Assay (technical validation of an amplicon panel) [6]
- AI Driven Neoantigen Identification Review (discusses impact of variant calling quality on downstream applications) [9]
- Segmental Copy Number Variant Detection Using Amplicon NGS (validation strategies for CNV in targeted panels) [11]
- Whole Exome Sequencing in Sub Saharan African Cohort (example of joint calling in a real world study) [8]
- Developmental and Epileptic Encephalopathy CNV Study (importance of including CNV analysis) [10]
Related Articles
- VCF Annotation: How to Read Variant Consequence Reports
- Genome Browsers for Researchers: A Guide to Inspecting Genomic Evidence
- Sequence Alignment: Choosing the Right Method for Your Biological Question
- Multiple Sequence Alignment: Common Pitfalls and Quality Checks
- Phylogenetic Tree Workflow: From Aligned Sequences to a Defensible Figure