Medaka Variant Calling: A Practical Guide for Genomic Studies

By Dr. Zubair Khalid, DVM, MS, PhD ·

Medaka Variant Calling: A Practical Guide for Genomic Studies

Introduction to Medaka Variant Calling

Medaka variant calling is the computational process of identifying single nucleotide polymorphisms (SNPs), insertions and deletions (indels), and structural variants from sequencing data of the Japanese rice fish Oryzias latipes. This small freshwater teleost has become a cornerstone of vertebrate genetics due to its short generation time, transparent embryos, and the availability of highly contiguous reference genomes. The term "medaka variant calling" encompasses the entire bioinformatic pipeline—from raw sequencing reads to annotated, filtered variant calls—that enables researchers to associate genetic variation with phenotype.

The medaka genome is approximately 800 megabases (Mb) distributed across 24 chromosomes. Its compact size, relative to zebrafish (1.4 Gb) or mammals (3 Gb), makes it tractable for whole-genome resequencing at moderate cost. The current reference genome assembly, HdrR (strain HNI-II), provides chromosome-level contiguity with an N50 exceeding 30 Mb, which is sufficient for most variant discovery applications. However, the medaka genome also contains regions of high GC content, segmental duplications, and a substantial transposable element fraction (~20%), all of which present specific challenges for read alignment and variant detection.

The core workflow for medaka variant calling proceeds through five stages: sequencing data generation, read alignment to a reference genome, post-alignment preprocessing, variant discovery, and variant filtering/annotation. Each stage introduces potential sources of error, and the choices made at each step propagate downstream. Understanding the mechanistic basis of these choices—rather than simply running default parameters—is essential for producing reliable variant calls.

Applications in Evolutionary and Disease Genomics

Medaka variant calling underpins several active research areas. In evolutionary genomics, medaka populations from Japan, Korea, and China show deep phylogeographic structure, and whole-genome resequencing of wild-caught individuals has revealed signatures of local adaptation, including genes involved in temperature tolerance and pigmentation. The olvas (vasa) gene promoter, commonly used for germline transgenesis, shows allelic variation that affects expression levels—a finding that emerged directly from population-level variant calling.

In disease genomics, medaka serves as a model for human congenital disorders, particularly those affecting craniofacial development and pigmentation. The slc45a2 gene, which encodes a membrane-associated transporter protein, harbors a well-characterized 4-base pair deletion in the i-3 albino mutant; identifying such causal variants requires precise indel calling. Additionally, medaka is used in chemical mutagenesis screens (ENU-based), where variant calling must distinguish induced point mutations from background polymorphisms. This application demands stringent filtering to avoid false positives that would otherwise waste months of phenotypic analysis.

Sequencing Technologies and Data Generation

The choice of sequencing platform fundamentally constrains what variants can be detected. Medaka variant calling is most commonly performed on short-read Illumina data, but long-read technologies are increasingly used to resolve complex regions.

Short-Read vs. Long-Read Sequencing

Illumina short-read sequencing (150 bp paired-end) remains the workhorse for medaka variant calling. The error rate of ~0.1% per base is low enough that, at 30× coverage, true heterozygous variants can be distinguished from sequencing errors with high confidence. Paired-end reads with an insert size of 300–500 bp provide sufficient fragment length to map uniquely across most of the medaka genome. However, short reads cannot span long repetitive elements or segmental duplications, leading to alignment ambiguities in these regions.

Oxford Nanopore Technologies (ONT) and Pacific Biosciences (PacBio) long-read sequencing produce reads of 10–100 kb, which can span entire transposable elements and resolve structural variants. The trade-off is a higher raw error rate (5–15% for ONT, 1–2% for PacBio HiFi). For medaka, where the genome contains many young transposable element insertions that differ between strains, long reads are particularly valuable for detecting presence/absence variants. Hybrid approaches—combining short-read accuracy with long-read contiguity—are emerging as the gold standard for de novo assembly and variant calling in non-model strains.

Coverage and Quality Considerations

Coverage depth directly determines variant calling sensitivity. For diploid medaka, 30× genome-wide coverage is generally sufficient to detect heterozygous SNPs with >99% sensitivity using GATK HaplotypeCaller. Lower coverage (10–15×) can still identify common variants but will miss a substantial fraction of heterozygotes. Coverage should be calculated as the number of sequenced bases divided by genome size, accounting for the fact that not all reads map uniquely.

Sequencing quality metrics—per-base quality scores (Phred scores), GC bias, and duplication rates—should be assessed before proceeding. The FastQC tool provides per-base quality plots, adapter content, and GC distribution. For Illumina data, a median Phred score above 30 (0.1% error rate) is expected. GC bias, where regions with extreme GC content (below 20% or above 80%) are under-represented, can be mitigated by library preparation methods such as the KAPA HyperPrep kit, which uses an optimized PCR amplification step. The medaka genome has a mean GC content of ~40%, but gene-dense regions on chromosome 1 can exceed 60% GC, making them prone to coverage drops.

Read Alignment and Preprocessing

Read alignment is the process of determining the genomic origin of each sequencing read. This step is computationally intensive and requires careful parameter selection to balance sensitivity against mapping errors.

Reference Genome Selection

The choice of reference genome is the single most important decision in medaka variant calling. The HdrR reference (assembly version ASM223467v1) is the standard for most studies, but its utility depends on the genetic distance between the reference strain and the samples being analyzed. For wild medaka populations from the Southern Japan clade, HdrR provides a close match, and variant calling will yield mostly true polymorphisms. For Northern Japanese or Chinese populations, which diverged from HdrR approximately 4 million years ago, the increased divergence leads to reference bias—reads carrying alternative alleles map less efficiently, causing false homozygosity calls.

If studying a strain that is highly divergent from HdrR, consider generating a strain-specific reference assembly using long-read sequencing, or at minimum using a reference-free variant calling approach such as de novo assembly of pooled reads followed by variant detection against the assembled contigs. For most applications, however, the HdrR reference is adequate, and the Variant Calling in Genomics overview provides context on how reference choice affects downstream analyses.

Alignment Tools (BWA, Minimap2)

For Illumina short reads, BWA-MEM (Burrows-Wheeler Aligner with Maximal Exact Matches) is the standard aligner. BWA-MEM uses a seed-and-extend strategy: it identifies short exact matches (seeds) between the read and the reference, then extends these seeds using a Smith-Waterman algorithm to produce a full alignment. The default parameters (minimum seed length of 19 bp, band width of 100) work well for medaka, but two adjustments improve performance:

  1. -M flag: Mark shorter split hits as secondary alignments, which is required for compatibility with Picard tools.
  2. -R flag: Add a read group header, which is essential for GATK downstream processing.

For long reads (ONT, PacBio), Minimap2 is the preferred aligner. Minimap2 uses minimizers—sampled k-mers that reduce the search space—to find candidate alignment regions, then performs base-level alignment. The -ax map-ont preset for ONT and -ax map-hifi for PacBio HiFi reads optimize the scoring parameters for the error profiles of these platforms. Minimap2 is substantially faster than BWA-MEM for long reads and handles the higher error rates gracefully.

Post-Alignment Processing (MarkDuplicates, BQSR)

After alignment, the BAM (Binary Alignment Map) file requires several preprocessing steps before variant calling:

  1. Sort and index: Use samtools sort to coordinate-sort the BAM file and samtools index to create a .bai index file. This enables random access to specific genomic regions during variant calling.
  1. Mark duplicates: PCR amplification during library preparation creates duplicate reads that originate from the same original fragment. These duplicates inflate coverage estimates and violate the independence assumption of variant calling algorithms. Picard's MarkDuplicates tool identifies read pairs with identical 5' coordinates and insert sizes, marking all but one as duplicates. For medaka, where input DNA is often limiting (e.g., from fin clips), duplicate rates can exceed 20%; marking rather than removing duplicates preserves information while preventing bias.
  1. Base Quality Score Recalibration (BQSR): This GATK tool empirically recalibrates base quality scores by comparing observed mismatches to known variant sites. The rationale is that sequencing machines systematically overestimate or underestimate base qualities in certain contexts (e.g., at the ends of reads, or following specific dinucleotide motifs). BQSR requires a set of known variants—for medaka, this can be generated from a pilot variant calling run or obtained from public databases such as the Medaka Variant Database. The recalibration produces a recalibration table that is applied to the BAM file, producing more accurate quality scores that improve variant calling specificity.

Variant Calling Algorithms and Tools

Variant calling is the statistical inference of genotypes at each genomic position from aligned reads. The core challenge is distinguishing true genetic variation from sequencing errors and alignment artifacts.

Germline vs. Somatic Calling

Germline variant calling assumes that all cells in the sample carry the same genotype—either homozygous reference, heterozygous, or homozygous alternate. The genotype likelihood is calculated from the base qualities of all reads covering a position, and the most likely genotype is reported with a quality score (Phred-scaled probability that the genotype is wrong).

Somatic variant calling, used for cancer samples or mosaic mutations, must account for the possibility that a variant is present in only a fraction of cells. This requires a different statistical framework that estimates variant allele fraction (VAF) and detects variants at low VAF (as low as 1–5%). For medaka, somatic calling is relevant when studying spontaneous mutations in aging studies or ENU-induced mosaicism in F0 founders. Tools such as Mutect2 (GATK) and VarScan2 implement somatic calling algorithms, but they require matched normal samples to distinguish true somatic mutations from sequencing artifacts.

Tool-Specific Parameters

GATK HaplotypeCaller is the most widely used tool for medaka germline variant calling. Its algorithm operates in two stages:

  1. Haplotype assembly: For each active region (determined by the presence of reads with mismatches), HaplotypeCaller reassembles the reads into a local de Bruijn graph. This graph represents all possible haplotypes in the region, allowing the tool to identify variants that span multiple polymorphisms (e.g., two SNPs in the same read).
  1. Genotype likelihood calculation: Each read is then aligned against all possible haplotypes, and the genotype likelihood is computed using a Bayesian model that incorporates base qualities and mapping qualities.

Key parameters for HaplotypeCaller include:

  • --minimum-mapping-quality (default 20): Reads with mapping quality below this threshold are ignored.
  • --standard-min-confidence-threshold-for-calling (default 30): The minimum Phred-scaled confidence for a variant to be emitted.
  • --ploidy (default 2): Set to 2 for diploid medaka.

FreeBayes is a Bayesian haplotype-based caller that does not require a reference genome for haplotype assembly. It uses a different statistical model that can handle pooled samples and polyploid genomes. For medaka, FreeBayes is often used as a secondary caller to validate GATK results, particularly for indels where the two tools show moderate discordance.

DeepVariant uses a convolutional neural network to classify candidate variants from read pileups. It transforms aligned reads into a multi-channel image (one channel per base, with intensity representing base quality) and trains a network to distinguish true variants from artifacts. DeepVariant consistently outperforms GATK and FreeBayes on human data, and it performs well on medaka when a suitable training set is available. However, it requires substantial computational resources (GPU recommended) and a high-quality truth set for fine-tuning.

The choice of tool depends on the application. For population-scale medaka studies where computational cost matters, GATK HaplotypeCaller with GVCF mode (which emits per-sample genotype likelihoods that can be combined across samples) is the standard. For a detailed comparison of variant calling strategies across species, the Variant Calling in Genomics resource provides additional context.

Variant Filtering and Quality Control

Raw variant calls contain a substantial fraction of false positives—typically 10–30% depending on coverage and tool. Filtering is therefore not optional; it is an integral part of the variant calling process.

Hard Filters vs. Machine Learning (VQSR)

Hard filters apply fixed thresholds to variant quality metrics. Common filters for medaka variant calling include:

MetricDescriptionTypical Threshold
QUALPhred-scaled variant confidence> 30
DPRead depth at variant site10–100 (not too low, not too high)
QDQUAL divided by depth (variant quality per read)> 2.0
FSFisher strand bias (probability of strand bias)< 60
MQRoot mean square mapping quality> 40
MQRankSumMapping quality rank sum test (z-score)> -12.5
ReadPosRankSumRead position rank sum test (z-score)> -8.0

The strand bias metrics (FS, MQ, ReadPosRankSum) detect systematic errors where variants appear predominantly on one strand or at read ends—hallmarks of alignment artifacts. The rank sum tests compare the distribution of mapping qualities or read positions between reference and alternate alleles; significant deviations suggest that the alternate allele is supported by low-quality reads.

Variant Quality Score Recalibration (VQSR) is a machine learning approach that trains a Gaussian mixture model on known true variants (from databases or validated sets) and known false variants (from sites that are monomorphic in many samples). VQSR is more powerful than hard filters because it captures correlations between metrics, but it requires a large truth set (at least 10,000 known variants) to be effective. For medaka, where public variant databases are less comprehensive than for humans, hard filters are often more practical.

Benchmarking with Known Variants

The most reliable way to assess variant calling accuracy is to compare against a set of known true variants. For medaka, this can be generated by:

  1. Sanger sequencing validation: Randomly select 50–100 variants and validate by PCR amplification and Sanger sequencing. The Sanger Sequencing Protocol provides the experimental details. This approach yields a precise estimate of the false positive rate.
  1. Cross-tool concordance: Call variants with two independent tools (e.g., GATK and FreeBayes) and retain only variants called by both. The concordant set has a higher precision, though it may miss some true variants (lower recall).
  1. Parent-offspring trios: If family samples are available, Mendelian inheritance checks can identify genotyping errors. Variants that violate Mendelian transmission are likely false positives or alignment artifacts.

Annotation and Functional Impact Prediction

Once a high-confidence variant set is obtained, the next step is to understand what these variants do. Annotation maps variants to genomic features (genes, regulatory regions, splice sites) and predicts their functional consequences.

Variant Annotation Databases

The medaka reference genome annotation (Ensembl release 101) contains approximately 20,000 protein-coding genes. Variant annotation tools such as SnpEff and ANNOVAR use this annotation to classify variants as:

  • Synonymous: Base change that does not alter the amino acid sequence.
  • Missense: Base change that alters one amino acid.
  • Nonsense: Base change that introduces a premature stop codon.
  • Frameshift: Insertion or deletion that shifts the reading frame.
  • Splice site: Variant within 2 bp of an exon-intron boundary.
  • 5' UTR / 3' UTR: Variant in untranslated regions.
  • Intergenic: Variant outside any annotated gene.

SnpEff requires a pre-built database for the medaka genome. If the database is not available, it can be built from the Ensembl GTF file using the java -jar snpEff.jar build -gtf22 -v medaka command. ANNOVAR uses a different format and requires downloading the medaka annotation files from the UCSC Genome Browser or Ensembl.

Predicting Deleteriousness (SIFT, PolyPhen)

For missense variants, predicting whether the amino acid change is deleterious is a key step in prioritizing candidates. SIFT (Sorting Intolerant From Tolerant) uses sequence conservation across related species to predict whether an amino acid substitution is likely to affect protein function. PolyPhen-2 uses a combination of sequence conservation and structural features (if a protein structure is available) to classify variants as benign, possibly damaging, or probably damaging.

Both tools require multiple sequence alignments of medaka proteins with orthologs from other species. For medaka, the Ensembl Compara database provides these alignments. The tools output a probability score (SIFT: 0 = deleterious, 1 = tolerated; PolyPhen: 0 = benign, 1 = damaging) that can be used to rank candidate variants.

It is important to note that these predictions are not definitive—they identify variants that may be deleterious, but functional validation (e.g., CRISPR knockout or transgenic rescue) is required for confirmation. Combining multiple prediction tools (SIFT, PolyPhen, MutationTaster, CADD) and requiring consensus across tools increases confidence.

Common Pitfalls and Troubleshooting

Even with careful pipeline design, medaka variant calling can fail in predictable ways. Understanding these failure modes is essential for troubleshooting.

Reference Bias and Mapping Errors

Reference bias occurs when reads carrying the alternate allele map less efficiently to the reference genome than reads carrying the reference allele. This is particularly problematic for divergent populations or for variants in repetitive regions. The consequence is an overestimation of homozygosity for the reference allele and a systematic loss of alternate alleles.

Symptoms: Excess of homozygous reference calls, low VAF for known heterozygotes, and a deficit of variants in regions with high divergence.

Solutions:

  • Use a reference genome from the same population as the samples.
  • Increase the --min-seed-length parameter in BWA-MEM to reduce spurious alignments.
  • Apply a mapping quality filter (MQ > 20) to remove reads that align ambiguously.
  • Consider using a variant caller that is less sensitive to reference bias, such as FreeBayes, which does not rely on a reference for haplotype assembly.

Handling Low-Complexity Regions

The medaka genome contains microsatellites (di- and trinucleotide repeats), homopolymer runs, and transposable element remnants. These regions are prone to alignment errors and PCR slippage during library preparation, leading to false indel calls.

Symptoms: Clusters of indels in repeat regions, high indel-to-SNP ratio, and poor concordance between replicates.

Solutions:

  • Mask low-complexity regions in the reference genome using RepeatMasker or Tandem Repeats Finder before alignment.
  • Apply a variant quality filter that is more stringent for indels than SNPs (e.g., QD > 5 for indels vs. QD > 2 for SNPs).
  • Use a caller that models repeat-induced errors, such as GATK HaplotypeCaller, which reassembles reads and can phase indels across repeats.

Strand Bias and Read-End Artifacts

Systematic errors often appear as strand bias—variants supported predominantly by reads on one strand. This can result from oxidative damage during DNA extraction (creating 8-oxo-guanine, which mispairs with adenine) or from sequencing chemistry artifacts.

Symptoms: High FS scores, variants that fail the strand bias filters, and enrichment of C>T and G>A transitions (the signature of oxidative damage).

Solutions:

  • Use a DNA extraction method that minimizes oxidation (e.g., QIAGEN DNeasy rather than phenol-chloroform).
  • Apply the FS filter (FS < 60) and the ReadPosRankSum filter (ReadPosRankSum > -8.0).
  • For FFPE-like damage (not typical for medaka, but relevant for archival samples), use the GATK BaseRecalibrator with a damage model.

Practical Workflow Summary

The following workflow represents a best-practice medaka variant calling pipeline. It assumes Illumina paired-end data at 30× coverage and the HdrR reference genome.

Example Pipeline Commands

Step 1: Quality control and trimming

fastqc sample_R1.fastq.gz sample_R2.fastq.gz
trimmomatic PE -phred33 sample_R1.fastq.gz sample_R2.fastq.gz \
  sample_R1_trim.fastq.gz sample_R1_unpaired.fastq.gz \
  sample_R2_trim.fastq.gz sample_R2_unpaired.fastq.gz \
  ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36

Step 2: Alignment and preprocessing

bwa mem -M -R "@RG\tID:sample\tSM:sample\tPL:ILLUMINA" \
  medaka_hdrr.fa sample_R1_trim.fastq.gz sample_R2_trim.fastq.gz | \
  samtools sort -o sample.sorted.bam
samtools index sample.sorted.bam
picard MarkDuplicates I=sample.sorted.bam O=sample.dedup.bam M=dup_metrics.txt
gatk BaseRecalibrator -I sample.dedup.bam -R medaka_hdrr.fa \
  --known-sites medaka_known_variants.vcf -O recal_table.txt
gatk ApplyBQSR -I sample.dedup.bam -R medaka_hdrr.fa \
  --bqsr-recal-file recal_table.txt -O sample.recal.bam

Step 3: Variant calling

gatk HaplotypeCaller -I sample.recal.bam -R medaka_hdrr.fa \
  -O sample.g.vcf.gz -ERC GVCF
gatk GenotypeGVCFs -R medaka_hdrr.fa -V sample.g.vcf.gz -O sample.vcf.gz

Step 4: Filtering

gatk VariantFiltration -V sample.vcf.gz -O sample.filtered.vcf.gz \
  --filter-expression "QD < 2.0 || FS > 60.0 || MQ < 40.0" \
  --filter-name "StandardFilter"

Step 5: Annotation

java -jar snpEff.jar eff -v medaka_hdrr sample.filtered.vcf.gz > sample.annotated.vcf

Key Parameters to Optimize

The following parameters have the largest impact on variant calling quality:

  1. Read trimming thresholds: Trimming to a minimum length of 36 bp and a sliding window quality of 15 (Phred) removes low-quality bases without discarding too much data.
  1. Mapping quality threshold: Filtering reads with MQ < 20 removes ambiguously mapped reads. For repetitive regions, increasing to MQ < 30 improves specificity at the cost of sensitivity.
  1. Variant quality threshold: A QUAL > 30 corresponds to a 1 in 1000 chance that the variant is a false positive. For low-coverage studies, lowering to QUAL > 20 may be necessary, but this increases the false positive rate.
  1. Allele balance filter: For heterozygous variants, the alternate allele fraction should be approximately 0.5. Filtering variants with allele balance < 0.2 or > 0.8 removes many false positives caused by alignment errors.

Frequently Asked Questions

What is the best tool for medaka variant calling?

There is no single "best" tool; the optimal choice depends on your data and goals. GATK HaplotypeCaller is the most robust and widely used for germline SNP and indel calling in medaka, particularly for population studies where GVCF mode enables joint calling across samples. DeepVariant achieves higher accuracy but requires GPU resources and a training set. FreeBayes is a good secondary caller for validation. For somatic variants (e.g., ENU-induced mutations), Mutect2 is the standard. We recommend using GATK HaplotypeCaller as the primary caller and validating a subset of variants by Sanger sequencing.

How do I choose a reference genome for medaka variant calling?

The HdrR reference (ASM223467v1) is the default for most studies. If your samples are from the Southern Japan clade (the origin of HdrR), this reference is appropriate. For Northern Japanese or Chinese populations, consider whether the divergence level (approximately 1–2% sequence divergence) will introduce reference bias. If you are studying a highly divergent strain, options include generating a strain-specific reference assembly or using a reference-free approach. Always check the assembly version and annotation quality—older assemblies may have misassembled regions that produce spurious variant calls.

What depth of coverage is needed for reliable medaka variant calling?

For germline SNP detection, 30× genome-wide coverage is the recommended minimum. This provides >99% sensitivity for heterozygous variants and allows accurate genotyping. For indel detection, higher coverage (40–50×) is advisable because indels are harder to call and require more supporting reads. For population studies where cost is a concern, 15–20× coverage can detect common variants (minor allele frequency > 5%) but will miss rare variants. For somatic variant detection, coverage must be high enough to detect low VAF variants—typically 100× or more.

How do I filter medaka variants to reduce false positives?

Apply a combination of hard filters based on variant quality metrics: QUAL > 30, QD > 2.0, FS < 60, MQ > 40, and ReadPosRankSum > -8.0. Additionally, filter by depth (DP > 10 and DP < 100) to remove sites with insufficient or excessive coverage. For indels, use more stringent thresholds (QD > 5.0). Validate a subset of variants by Sanger sequencing to estimate the false positive rate. If you have a large number of samples, consider using VQSR, but this requires a substantial truth set.

What is the difference between SNP and indel calling in medaka?

SNP calling identifies single base substitutions, while indel calling identifies insertions or deletions of one or more bases. Indel calling is more challenging because indels cause alignment shifts (frameshifts in the alignment) that can confuse read mapping. GATK HaplotypeCaller handles both by reassembling reads in active regions, which allows it to detect indels that would be missed by simple pileup-based callers. Indels in repetitive regions are particularly error-prone; these should be filtered more stringently and validated by Sanger sequencing or long-read sequencing.

Can I use RNA-seq data for medaka variant calling?

RNA-seq data can be used for variant calling, but with important caveats. RNA-seq reads originate from expressed genes, so variants are only detected in transcribed regions. Additionally, RNA editing (A-to-I editing) can create apparent variants that are not present in the genome. Splice junctions cause reads to map across exon boundaries, which complicates alignment. If you must use RNA-seq data, use a splice-aware aligner such as STAR or HISAT2, and be aware that the Process RNA-seq Data workflow differs from DNA-based variant calling. For most applications, DNA sequencing is strongly preferred for variant discovery.

How do I annotate medaka variants after calling?

Use SnpEff or ANNOVAR with the medaka reference annotation. SnpEff requires building a database from the Ensembl GTF file (java -jar snpEff.jar build -gtf22 -v medaka). ANNOVAR requires downloading the medaka annotation from the UCSC Genome Browser. Both tools classify variants by functional impact (synonymous, missense, nonsense, etc.) and provide gene names. For predicting deleteriousness of missense variants, run SIFT and PolyPhen-2, which require multiple sequence alignments of medaka proteins with orthologs. The Variant Calling in Genomics resource provides additional guidance on annotation strategies.

Key Takeaways

  • Medaka variant calling is a multi-stage pipeline—sequencing, alignment, preprocessing, variant discovery, filtering, and annotation—where each stage affects downstream accuracy.
  • The HdrR reference genome is suitable for most studies, but reference bias becomes problematic for highly divergent populations; consider strain-specific references when necessary.
  • 30× Illumina coverage is the recommended minimum for reliable germline SNP and indel detection; higher coverage is needed for somatic variants or complex regions.
  • GATK HaplotypeCaller is the standard tool for medaka variant calling, with FreeBayes and DeepVariant as alternatives; validation by Sanger sequencing is essential for high-stakes applications.
  • Hard filtering using quality metrics (QUAL, QD, FS, MQ) is practical for medaka, where VQSR truth sets are limited; strand bias and depth filters are particularly important.
  • Annotation with SnpEff or ANNOVAR, followed by functional prediction with SIFT and PolyPhen-2, is necessary to prioritize candidate variants for downstream validation.
  • Common pitfalls—reference bias, low-complexity region artifacts, and strand bias—can be mitigated through careful parameter selection and targeted filtering.

Related Clinical & Scientific Guides