5. Somatic & Cancer Variant Calling

GATK Mutect2 calls somatic SNVs and indels by assembling active regions in tumor and matched-normal samples. It uses a Panel of Normals (PoN) to remove recurrent germline and sequencing artifacts, and applies orientation-bias and contamination models.

# 1. Mutect2 Somatic Calling with Matched Normal & PoN
gatk Mutect2 \
  -R GRCh38.fa \
  -I tumor.dedup.bam -tumor TUMOR_SAMPLE \
  -I normal.dedup.bam -normal NORMAL_SAMPLE \
  --panel-of-normals pon.vcf.gz \
  --germline-resource af-only-gnomad.vcf.gz \
  -O somatic_raw.vcf.gz \
  --f1r2-tar-gz f1r2.tar.gz

# 2. Learn Read Orientation Model (filters FFPE oxidation artifacts)
gatk LearnReadOrientationModel -I f1r2.tar.gz -O read_orientation_model.tar.gz

# 3. Filter Somatic Calls with Contamination and Orientation Priors
gatk FilterMutectCalls \
  -R GRCh38.fa \
  -V somatic_raw.vcf.gz \
  --ob-priors read_orientation_model.tar.gz \
  -O somatic_filtered_pass.vcf.gz

Strelka2 is an ultra-fast and accurate somatic and germline variant caller developed by Illumina. It uses tiered haplotype estimation and Empirical Variant Scoring (EVS) to call somatic SNVs and indels in minutes.

# 1. Configure Somatic Tumor-Normal Workflow
configureStrelkaSomaticWorkflow.py \
  --tumorBam=tumor.dedup.bam \
  --normalBam=normal.dedup.bam \
  --referenceFasta=GRCh38.fa \
  --runDir=strelka_somatic_run \
  --exome

# 2. Execute Strelka2 with multi-core parallelism
python strelka_somatic_run/runWorkflow.py -m local -j 16

# Outputs available at:
# strelka_somatic_run/results/variants/somatic.snvs.vcf.gz
# strelka_somatic_run/results/variants/somatic.indels.vcf.gz

VarScan 2 uses heuristic thresholds and Fisher’s exact test to detect somatic SNVs, indels, and copy number alterations/loss of heterozygosity (LOH) from simultaneous samtools mpileup streams of tumor and normal tissues.

# 1. Generate synchronized mpileup stream for normal and tumor
samtools mpileup -B -q 20 -f GRCh38.fa normal.dedup.bam tumor.dedup.bam > sync.mpileup

# 2. Run VarScan somatic mode
java -jar VarScan.jar somatic sync.mpileup varscan_out \
  --mpileup 1 \
  --min-coverage-normal 10 \
  --min-coverage-tumor 15 \
  --min-var-freq 0.05 \
  --somatic-p-value 0.05 \
  --output-vcf 1

# 3. Process somatic, germline, and LOH calls
java -jar VarScan.jar processSomatic varscan_out.snp.vcf

VarDict is an ultra-sensitive variant caller capable of detecting SNVs, complex indels, structural variations, and loss of heterozygosity down to 0.5% variant allele frequency in cell-free DNA (cfDNA) and liquid biopsy panels.

# 1. Run VarDictJava in Tumor-Normal paired mode
vardict-java -G GRCh38.fa -N "TUMOR_cfDNA" \
  -b "tumor_cfdna.bam|normal.bam" \
  -z 1 -c 1 -S 2 -E 3 -g 4 \
  targeted_panel_regions.bed \
  -f 0.005 \
  | teststrandbias.R \
  | var2vcf_paired.pl -N "TUMOR_cfDNA|NORMAL" -f 0.005 > vardict_paired.vcf

# 2. Compress and index output
bgzip vardict_paired.vcf && tabix -p vcf vardict_paired.vcf.gz

LoFreq is a model-free variant caller that models sequencing error rates to accurately call somatic SNVs and indels at very low allele frequencies (<1%) in deep sequencing (>1000x), viral quasispecies, and subclonal cancer biopsies.

# 1. Insert indel qualities (required prior to LoFreq somatic calling)
lofreq indelqual --dindel -f GRCh38.fa -o tumor.lofreq.bam tumor.dedup.bam
lofreq indelqual --dindel -f GRCh38.fa -o normal.lofreq.bam normal.dedup.bam
samtools index tumor.lofreq.bam && samtools index normal.lofreq.bam

# 2. Execute two-sample somatic calling
lofreq somatic -n normal.lofreq.bam -t tumor.lofreq.bam \
  -f GRCh38.fa \
  -o lofreq_somatic_ \
  --threads 12 \
  --dbsnp dbsnp_146.hg38.vcf.gz