A Step-by-Step Guide to Running FreeBayes for Germline Variant Calling: Parameters, Performance, and Pitfalls

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

A Step-by-Step Guide to Running FreeBayes for Germline Variant Calling: Parameters, Performance, and Pitfalls

Key Takeaways

  • FreeBayes germline variant calling necessitates precise input data preparation, including coordinate-sorted and indexed BAM files, and alignment to the GRCh38 reference genome to avoid liftover errors and capture novel loci.
  • Optimal germline calling requires accurate ploidy specification, particularly for sex chromosomes (e.g., using a --ploidy-map for males), and careful tuning of --min-alternate-fraction (default 0.2) and --min-alternate-count (default 2) to balance sensitivity and specificity.
  • Joint calling across multiple samples, achieved by listing multiple --bam files, enhances variant detection sensitivity by leveraging cohort-wide evidence, but requires unique read group identifiers for each sample.
  • Post-calling quality control is critical, involving hard filtering based on genotype quality (e.g., >20), read depth (e.g., 10-100), and alternate allele fraction (e.g., 0.2-0.8), followed by variant annotation for functional interpretation.
  • Common pitfalls include low variant yield due to overly stringent parameters and excessive false positives from insufficient quality filtering or unhandled duplicate reads, necessitating careful parameter selection and rigorous QC.
  • Reproducibility is paramount, requiring meticulous documentation of FreeBayes version, reference genome, and all parameter settings, ideally within a containerized environment (e.g., Docker) or version-controlled scripts.

Germline variant calling with FreeBayes requires a deliberate workflow that begins with aligned sequencing data and ends with a filtered variant call format (VCF) file suitable for downstream analysis. This guide walks through the complete process, from input preparation through parameter selection, execution, and quality control, with attention to the common failure points that compromise variant detection accuracy. The intended reader is a bioinformatics practitioner who has aligned reads to a reference genome and needs practical direction on running FreeBayes effectively for germline samples.

Understanding FreeBayes in the Germline Calling Context

FreeBayes is a Bayesian genetic variant detector designed to identify single nucleotide polymorphisms (SNPs), insertions and deletions (indels), and multi-nucleotide polymorphisms from aligned sequencing data. For germline calling, the tool analyzes each genomic position independently, using the observed allele frequencies in the sample to estimate the probability that a variant is present. This approach differs from somatic calling, where tumor and normal samples are compared to distinguish acquired mutations from inherited variation.

The distinction between germline and somatic variant calling matters for parameter selection. Germline variants are expected to be present in either a heterozygous state, with approximately 50 percent alternate allele fraction, or a homozygous state, with near 100 percent alternate allele fraction. Somatic variants, by contrast, often appear at lower allele fractions because they are present in only a subset of cells within a tumor sample. FreeBayes can be configured for both applications, but the default parameters and filtering thresholds differ based on the expected allele frequency spectrum.

The input requirements for FreeBayes are straightforward. The tool accepts BAM files that have been aligned to a reference genome and sorted by genomic coordinate. Duplicate reads should be marked or removed prior to variant calling, as PCR duplicates artificially inflate read depth at specific positions and can bias allele frequency estimates. Base quality score recalibration is recommended but not strictly required, particularly for datasets generated on newer sequencing platforms with improved base calling accuracy.

Reference genome choice influences variant calling results. The GRCh38 assembly represents the current standard for human genomics, and variant calling directly on this assembly avoids the complications of coordinate liftover from older builds. A 2019 effort that generated variant calls on GRCh38 from 2,548 samples across 26 populations demonstrated that direct calling on the current assembly captures novel, medically relevant loci that are not represented in older reference builds (Wellcome Open Research). This work highlighted the practical advantages of using the current reference instead of relying on coordinate conversions from GRCh37.

At a Glance: FreeBayes Germline Calling Decision Table

Decision PointRecommended SettingRationaleCommon Pitfall
Reference genomeGRCh38 primary assemblyDirect calling avoids liftover errors and captures novel lociUsing GRCh37-aligned reads with GRCh38 reference
Ploidy2 for autosomes, 1 for sex chromosomes via ploidy mapMatches expected germline allele fractionsTreating male X chromosome as diploid
Minimum alternate fraction0.2 defaultAccommodates heterozygous imbalance while excluding noiseSetting too high and missing true heterozygotes
Minimum alternate count2 readsPrevents single-read artifactsLowering to 1 in low-depth samples
Minimum coverage10 for WGS, higher for targeted panelsEnsures confident genotype callsSetting too high and excluding valid regions
Duplicate handlingMark duplicates before callingPrevents inflated depth and false variantsSkipping duplicate marking entirely
Joint callingMultiple BAMs in one commandImproves sensitivity across cohortMixing samples without unique read groups

Preparing Input Data for FreeBayes

Alignment Files and Reference Requirements

FreeBayes requires a reference genome in FASTA format and one or more aligned BAM files. The reference must be indexed with both samtools faidx and bwa index or an equivalent indexing tool, because FreeBayes accesses the reference through these index files during variant calling. The BAM files must be sorted by coordinate and indexed with samtools index before FreeBayes can process them efficiently.

The alignment process itself affects downstream variant calling quality. Reads aligned with BWA-MEM to the GRCh38 reference produce BAM files that are suitable for FreeBayes input. The alignment should use the same reference version that will be used for variant calling. Mixing references between alignment and variant calling steps introduces systematic errors because read mapping positions and base qualities are reference-dependent.

For germline calling, the sample BAM file should represent a single individual. Pooled samples or mixed populations complicate the allele frequency interpretation because FreeBayes assumes that observed allele frequencies reflect the genotype of one individual. If the input contains reads from multiple individuals, the tool will call variants that represent the mixture, which is not appropriate for germline analysis.

Read Group and Sample Metadata

FreeBayes uses read group information from the BAM file header to assign reads to samples. Each sample must have a unique read group identifier, and the sample name in the read group must match the sample name used in downstream analysis. When processing multiple samples simultaneously, FreeBayes can call variants jointly across all samples, which improves sensitivity for low-frequency alleles that might be present in one sample but supported by reads from another.

The BAM file header should contain complete read group information, including the sample name, library, and platform. This metadata becomes part of the VCF output and is essential for downstream filtering and interpretation. Missing or inconsistent read group information is a common source of errors in variant calling pipelines, particularly when samples are merged or processed in batches.

Coverage and Depth Considerations

Sequencing depth directly influences variant calling confidence. For germline calling, the recommended depth depends on the application. Whole genome sequencing at 30x coverage provides sufficient depth for most germline variant detection, while targeted panels may require higher depth to achieve confidence in low-frequency variant detection. FreeBayes uses the observed read depth at each position to calculate genotype likelihoods, so positions with inadequate depth will have uncertain genotype calls.

The distribution of coverage across the genome matters as much as the average depth. Regions with extreme GC content, repetitive sequences, or structural variation often have reduced mappability and lower coverage. These regions will have fewer reads supporting variant calls, and the genotype quality scores will reflect this uncertainty. Practitioners should examine coverage statistics before variant calling to identify regions that may produce unreliable calls.

Core FreeBayes Parameters for Germline Calling

Ploidy and Population Model Settings

The --ploidy parameter specifies the expected ploidy of the sample. For human germline calling, the default ploidy of 2 is appropriate for autosomal chromosomes. Sex chromosomes require special consideration. Males have a single X chromosome and a single Y chromosome, so the ploidy for these chromosomes should be set to 1. FreeBayes allows per-chromosome ploidy specification through the --ploidy-map option, which enables accurate calling on sex chromosomes without treating the single X chromosome as a homozygous diploid region.

The --pvar parameter sets the prior probability that a site is variant. The default value of 0.001 reflects the expectation that most positions in the genome are invariant. This prior influences the posterior probability calculation and affects the tradeoff between sensitivity and specificity. Lowering the prior increases specificity by requiring stronger evidence for variant calls, while raising the prior increases sensitivity at the cost of more false positives.

Allele Frequency and Genotype Parameters

The --min-alternate-fraction parameter sets the minimum alternate allele fraction required to consider a variant. For germline calling, the default value of 0.2 is appropriate because heterozygous variants are expected at approximately 0.5 allele fraction, and the threshold provides tolerance for sequencing error and allelic imbalance. Lowering this threshold increases sensitivity for variants with skewed allele fractions but also increases false positive calls from sequencing artifacts.

The --min-alternate-count parameter specifies the minimum number of reads supporting the alternate allele. The default value of 2 requires at least two reads to support a variant call. This threshold prevents single-read artifacts from being called as variants. For low-depth samples, this parameter may need adjustment, but lowering it below 2 risks calling sequencing errors as genuine variants.

The --min-coverage parameter sets the minimum depth required to make a genotype call at a position. Positions with coverage below this threshold are reported as missing genotypes. The default value of 0 allows calls at any depth, but practitioners often set this to a higher value to avoid low-confidence calls in regions with sparse coverage.

Quality and Filtering Parameters

The --min-base-quality parameter sets the minimum base quality score required for a read to contribute to variant calling. The default value of 13 corresponds to a base call accuracy of approximately 95 percent. Raising this threshold excludes lower quality bases from the analysis, which can reduce false positives in regions with systematic base quality issues.

The --min-mapping-quality parameter sets the minimum mapping quality required for a read to be considered. The default value of 1 is permissive, allowing reads with low mapping quality to contribute to variant calls. For germline calling, raising this threshold to 20 or 30 excludes reads that map ambiguously to multiple genomic locations, reducing false positive calls in repetitive regions.

The --genotype-qualities parameter enables output of genotype quality scores in the VCF. These scores represent the probability that the called genotype is correct given the observed data. Genotype quality scores are essential for downstream filtering and for identifying low-confidence calls that should be validated by an orthogonal method.

Running FreeBayes: Command Structure and Execution

Basic Command for Single Sample Germline Calling

A basic FreeBayes command for a single germline sample follows this structure:

freebayes \
  --fasta-reference reference.fasta \
  --bam sample.bam \
  --ploidy 2 \
  --min-alternate-fraction 0.2 \
  --min-alternate-count 2 \
  --min-coverage 10 \
  --genotype-qualities \
  --output-vcf sample.vcf

This command calls variants on a diploid sample with a minimum coverage of 10 reads, requiring at least 2 reads supporting the alternate allele at a fraction of at least 20 percent. The genotype quality scores are included in the output VCF for downstream filtering.

Joint Calling for Multiple Samples

Joint calling across multiple samples improves variant detection sensitivity and produces a unified VCF that simplifies downstream analysis. The command structure for joint calling uses multiple --bam arguments:

freebayes \
  --fasta-reference reference.fasta \
  --bam sample1.bam \
  --bam sample2.bam \
  --bam sample3.bam \
  --ploidy 2 \
  --min-alternate-fraction 0.2 \
  --min-alternate-count 2 \
  --genotype-qualities \
  --output-vcf cohort.vcf

In joint calling mode, FreeBayes evaluates each genomic position across all samples simultaneously. A variant is reported if the combined evidence across samples supports it, even if no single sample has strong evidence. This approach is particularly valuable for family-based studies where a variant may be present in a parent at low frequency but inherited by a child at higher frequency.

Parallelization and Performance Optimization

FreeBayes supports parallel execution through region-based splitting. The --region parameter restricts variant calling to a specific genomic interval, allowing multiple FreeBayes processes to run concurrently on different regions of the genome. The output VCF files from each region are then concatenated to produce the complete variant set.

The --region parameter accepts coordinates in the format chromosome:start-end. For whole genome analysis, the genome is typically split into intervals of 10 to 50 megabases, depending on available computational resources. Each interval can be processed independently, and the results are merged after all processes complete.

The --threads parameter enables multithreading within a single FreeBayes process. This parameter controls the number of threads used for parallel processing of genomic regions. The optimal thread count depends on the available CPU resources and the size of the input BAM files. Performance scaling is not linear, and practitioners should benchmark their specific dataset to determine the optimal thread configuration.

Variant Filtering and Quality Control

Understanding VCF Output Fields

The FreeBayes VCF output contains essential information for variant filtering. The QUAL field represents the probability that the site is variant, expressed as a Phred-scaled score. Higher QUAL values indicate stronger evidence for a variant. The INFO field contains allele frequency, read depth, and other site-level statistics. The FORMAT field contains per-sample genotype information, including genotype calls, genotype quality, and read depth.

The genotype quality score in the FORMAT field represents the confidence in the specific genotype call. A genotype quality of 20 corresponds to a 1 percent error rate, while a quality of 30 corresponds to a 0.1 percent error rate. These thresholds are commonly used for filtering low-confidence genotype calls.

Hard Filtering Approaches

Hard filtering applies fixed thresholds to VCF fields to remove low-quality variant calls. Common hard filters for germline FreeBayes output include:

  • Genotype quality greater than 20
  • Read depth between 10 and 100
  • Alternate allele fraction between 0.2 and 0.8 for heterozygous calls
  • Mapping quality greater than 30
  • Base quality greater than 20

These thresholds are applied using tools such as vcffilter from vcflib or bcftools filter. The specific thresholds should be adjusted based on the sequencing platform, coverage, and application. A practical framework for variant interrogation in tumor samples emphasized that filtering and validation are critical phases in variant analysis, requiring a systematic approach to prioritize meaningful variants while excluding artifacts (PLOS Computational Biology).

Variant Annotation and Interpretation

After filtering, variants should be annotated with functional information to support interpretation. Annotation tools add information about gene location, amino acid changes, and population frequency. This annotation step is essential for identifying clinically relevant variants and for prioritizing candidates for validation.

Population frequency databases provide context for interpreting variant significance. Variants that are common in the general population are less likely to be pathogenic, while rare variants require additional evidence for clinical interpretation. The 1000 Genomes Project data, which includes variants from 2,548 samples across 26 populations, provides a valuable reference for population frequency (Wellcome Open Research). This resource is particularly useful for filtering common polymorphisms from rare disease-associated variants.

Common Failure Patterns and Troubleshooting

Low Variant Yield

A common problem in FreeBayes germline calling is unexpectedly low variant yield. This issue often stems from overly stringent filtering parameters. The --min-alternate-fraction parameter set too high will exclude genuine heterozygous variants with skewed allele fractions. The --min-coverage parameter set too high will exclude variants in regions with below-average coverage.

Diagnosing low variant yield requires examining the distribution of allele fractions and coverage in the output VCF. If the allele fraction distribution shows a gap around 0.5, the filtering parameters may be excluding genuine heterozygous calls. If coverage is consistently below the minimum threshold across the genome, the sequencing depth may be inadequate for the intended application.

Excessive False Positive Calls

High false positive rates in FreeBayes output often result from insufficient quality filtering. Reads with low mapping quality in repetitive regions produce spurious variant calls. Base quality errors at read ends contribute false positive SNPs. The --min-mapping-quality and --min-base-quality parameters should be raised to exclude these artifacts.

Another source of false positives is the inclusion of duplicate reads. PCR duplicates inflate the read depth at specific positions and can create the appearance of a variant where none exists. Marking duplicates with picard MarkDuplicates or samtools markdup before variant calling reduces this artifact.

Strand Bias and Read Position Bias

Systematic sequencing errors often show strand bias, where the alternate allele appears predominantly on one strand. FreeBayes does not apply strand bias filtering by default, so these artifacts can pass through the initial variant calling step. Downstream filtering with tools that assess strand bias, such as bcftools filter with the SB annotation, can remove these artifacts.

Read position bias occurs when the alternate allele appears predominantly at read ends, where base calling errors are more common. This artifact is particularly problematic for indels, which are difficult to align at read ends. Filtering based on read position bias requires tools that annotate the position of variant-supporting reads within their full read length.

Reference Bias in Variant Representation

FreeBayes reports variants relative to the reference genome, and the representation of indels can vary depending on the local sequence context. The same biological variant can be represented as different alleles in different samples if the alignment differs. This representation issue complicates comparison across samples and requires normalization of variant representation before downstream analysis.

Variant normalization tools, such as bcftools norm, standardize the representation of variants by left-aligning indels and splitting multi-nucleotide variants. This normalization step is essential for comparing variants across samples and for merging variant calls from different callers.

Performance Considerations and Resource Management

Memory and Storage Requirements

FreeBayes memory usage scales with the size of the input BAM files and the number of samples being processed jointly. For a single whole genome sample, FreeBayes typically requires 4 to 8 gigabytes of RAM. Joint calling across multiple samples increases memory requirements proportionally. The --region parameter can reduce memory usage by restricting analysis to a single genomic interval.

Storage requirements for variant calling output are substantial. A whole genome VCF file can exceed 100 gigabytes for a single sample, and joint calling across many samples produces even larger files. Compressed VCF files using bgzip reduce storage requirements by approximately 80 percent. The compressed VCF can be indexed with tabix for efficient random access.

Runtime Expectations

FreeBayes runtime depends on the genome size, coverage, and available computational resources. A whole human genome at 30x coverage typically requires 10 to 20 hours of single-threaded runtime. Parallelization across multiple regions and threads reduces wall clock time substantially. The --region parameter enables embarrassingly parallel execution across compute nodes, making FreeBayes suitable for cluster environments.

The runtime also depends on the variant density in the sample. Samples with high heterozygosity produce more candidate variants and require more computation. The --pvar parameter influences runtime because lower priors reduce the number of sites that pass the initial variant probability threshold.

Scaling to Large Cohorts

For large cohort studies, joint calling across hundreds or thousands of samples requires careful resource planning. The computational cost of joint calling scales superlinearly with sample number because each position must be evaluated across all samples simultaneously. An alternative approach is to call variants per sample and then merge the results, but this approach loses the sensitivity benefits of joint calling.

Community-developed workflow frameworks provide structured approaches to large-scale variant calling. The nf-core documentation maintains standardized pipelines for genomic analysis that include variant calling steps. These pipelines handle the complexity of parallelization, resource allocation, and reproducibility, allowing researchers to focus on the biological interpretation of results.

Reproducibility and Documentation

Version Control and Environment Management

Reproducible variant calling requires documentation of the software environment. FreeBayes version, reference genome version, and parameter settings all influence the output. Container technologies such as Docker and Singularity provide a mechanism for capturing the complete software environment, ensuring that the same analysis can be reproduced at a later time.

The Carpentries lessons on version control with Git provide foundational training for managing analysis code. Tracking the exact commands and parameters used for variant calling in a version-controlled repository enables reproducibility and facilitates collaboration. The analysis code should include the FreeBayes version, the reference genome version, and all parameter settings.

Documentation of Analysis Steps

The variant calling workflow should be documented in sufficient detail that another researcher can reproduce the analysis. This documentation should include the exact FreeBayes command, the reference genome version, the alignment and preprocessing steps, and the filtering thresholds applied. The documentation should also note any deviations from standard protocols and the rationale for those deviations.

The Galaxy Training Network provides tutorials on variant analysis workflows that emphasize reproducibility. These tutorials demonstrate how to construct analysis workflows that can be shared and rerun by other researchers. The workflow approach ensures that all analysis steps are documented and that the analysis can be reproduced with the same parameters.

Data Management and Storage

Variant calling produces large intermediate files that require management. The raw sequencing data, aligned BAM files, and variant call files should be stored according to a data management plan. The NCBI data resources provide repositories for sequencing data, including the Sequence Read Archive for raw data and dbSNP for variant data. Depositing data in these repositories ensures long-term access and enables data sharing with the research community.

The storage requirements for genomic data are substantial, and the EMBL-EBI training resources provide guidance on data management practices for bioinformatics. These resources cover data formats, storage options, and data sharing considerations. Proper data management ensures that the results of variant calling analyses remain accessible and interpretable over time.

Limitations and Interpretation Boundaries

Technical Limitations of Short-Read Variant Calling

FreeBayes operates on short-read sequencing data, which has inherent limitations for variant detection. Structural variants, including large insertions, deletions, and rearrangements, are poorly detected by short-read approaches. Repetitive regions of the genome, including centromeres and telomeres, have limited mappability and produce unreliable variant calls.

Long-read sequencing technologies address some of these limitations by generating reads that span repetitive regions and structural variant breakpoints. A 2019 study using long-read sequencing to investigate the IGH-DUX4 translocation in B-cell acute lymphoblastic leukemia demonstrated the value of long reads for resolving complex genomic rearrangements (Nature Communications). These approaches complement short-read variant calling but require different analysis tools and workflows.

Population Diversity and Reference Bias

Variant calling accuracy depends on the reference genome and the population being studied. Reference genomes are built from a limited number of individuals, and populations that are genetically distant from the reference population have more variants relative to the reference. This reference bias affects variant calling sensitivity and accuracy.

A 2019 review of whole genome data interpretation highlighted the historical lack of population diversity in genotype-phenotype studies, which limits the generalizability of findings (Human Genetics). The All of Us research program aims to address this gap by ascertaining a diverse cohort of at least one million participants. For variant calling, this diversity means that reference genomes and variant databases must represent the populations being studied to avoid systematic bias.

Validation Requirements for Clinical Applications

Variant calls from FreeBayes require validation before clinical use. The filtering and validation phase of variant analysis is critical for distinguishing true variants from artifacts (PLOS Computational Biology). Validation typically involves orthogonal methods such as Sanger sequencing or an independent variant calling approach. The validation results should be documented and reported alongside the variant calls.

The clinical relevance of a variant depends on its functional impact and population frequency. Variant annotation tools provide information about the predicted effect of a variant on protein function, but these predictions require experimental validation. The interpretation of variant significance should follow established guidelines and should be performed by qualified personnel.

Professional Escalation Criteria

When to Seek Expert Consultation

Certain variant calling scenarios warrant consultation with a bioinformatics specialist or clinical geneticist. These include:

  • Variants in regions with complex structure or low mappability
  • Discrepancies between variant callers or between sequencing platforms
  • Variants with potential clinical significance that require validation
  • Unexpected patterns in variant calls that suggest systematic errors

The complexity of variant analysis pipelines and tool selection remains a barrier for researchers new to the field. A practical framework that guides researchers through the critical steps of variant interrogation, from planning through dissemination, can help navigate this complexity (PLOS Computational Biology). When the analysis exceeds the practitioner's expertise, consultation with a specialist is appropriate.

Documentation for Escalation

When escalating a variant calling issue, the practitioner should provide complete documentation of the analysis. This documentation should include the input data, the exact commands used, the parameter settings, and the output files. The documentation should also include any quality metrics that were calculated and any anomalies observed during the analysis.

The escalation documentation should be organized to facilitate efficient review by the specialist. The VCF file should be filtered to the variants of interest, and the supporting evidence for those variants should be summarized. The documentation should clearly state the question or concern that prompted the escalation.

Practical Implementation Steps

Step 1: Verify Input Data Integrity

Before running FreeBayes, confirm that the BAM files are sorted, indexed, and contain complete read group information. Check that the reference FASTA file has both .fai and .bwt index files. Validate that the reference version used for alignment matches the reference that will be used for variant calling. Run samtools quickcheck on each BAM file to identify corrupted files before starting the analysis.

Step 2: Assess Coverage Distribution

Calculate per-base coverage statistics across the genome using samtools depth or mosdepth. Examine the distribution to identify regions with unusually low or high coverage. Low coverage regions will produce unreliable genotype calls, while extreme high coverage may indicate duplicate reads or mapping artifacts. Document the median coverage and the proportion of the genome below the minimum coverage threshold.

Step 3: Configure Parameters Based on Application

Select parameters based on the specific germline calling application. For whole genome sequencing, use the default minimum alternate fraction of 0.2 and set minimum coverage to 10. For targeted panels with higher depth, increase the minimum coverage threshold to match the expected depth. Set ploidy to 2 for autosomes and configure a ploidy map for sex chromosomes when analyzing male samples.

Step 4: Execute FreeBayes with Region-Based Parallelization

Split the genome into intervals of 10 to 50 megabases and run FreeBayes on each interval in parallel. Use the --region parameter to specify each interval and the --threads parameter to enable multithreading within each process. Monitor runtime and memory usage to identify intervals that require more computational resources.

Step 5: Merge and Normalize Variant Calls

Concatenate the per-region VCF files using bcftools concat and sort the merged file by genomic coordinate. Normalize the variant representation using bcftools norm to left-align indels and split multi-nucleotide variants. This step ensures consistent variant representation across samples and enables accurate comparison with reference databases.

Step 6: Apply Hard Filters and Annotate

Apply hard filters to remove low-quality variant calls based on genotype quality, read depth, and allele fraction. Annotate the filtered variants with functional information using tools such as SnpEff or VEP. Add population frequency information from public databases to support variant interpretation.

Step 7: Document and Archive the Analysis

Record the exact FreeBayes version, reference genome version, parameter settings, and filtering thresholds in a version-controlled analysis script. Store the final VCF file in a compressed and indexed format. Deposit the raw data and final variant calls in appropriate repositories to ensure long-term access and reproducibility.

Records and Measurements

Essential Records for Variant Calling Runs

Maintain a laboratory notebook or electronic record that captures the following information for each FreeBayes run:

  • FreeBayes software version and installation method
  • Reference genome version and source
  • Input BAM file names and alignment parameters
  • Complete command line with all parameter settings
  • Runtime and peak memory usage
  • Number of variants called before and after filtering
  • Transition to transversion ratio as a quality metric
  • Genotype concordance with known control samples when available

These records enable troubleshooting when results are unexpected and provide the documentation needed for publication and clinical reporting.

Quality Metrics to Track

Track the following metrics for each variant calling run to assess data quality:

  • Total number of SNPs and indels called
  • Transition to transversion ratio, typically between 2.0 and 2.2 for human germline data
  • Heterozygous to homozygous variant ratio
  • Number of variants in repetitive or low complexity regions
  • Genotype quality distribution across all calls
  • Proportion of variants with population frequency above 1 percent

Deviations from expected values for these metrics indicate potential issues with the input data or parameter settings.

Common Failure Patterns and Troubleshooting

Low Variant Yield

A common problem in FreeBayes germline calling is unexpectedly low variant yield. This issue often stems from overly stringent filtering parameters. The --min-alternate-fraction parameter set too high will exclude genuine heterozygous variants with skewed allele fractions. The --min-coverage parameter set too high will exclude variants in regions with below-average coverage.

Diagnosing low variant yield requires examining the distribution of allele fractions and coverage in the output VCF. If the allele fraction distribution shows a gap around 0.5, the filtering parameters may be excluding genuine heterozygous calls. If coverage is consistently below the minimum threshold across the genome, the sequencing depth may be inadequate for the intended application.

Excessive False Positive Calls

High false positive rates in FreeBayes output often result from insufficient quality filtering. Reads with low mapping quality in repetitive regions produce spurious variant calls. Base quality errors at read ends contribute false positive SNPs. The --min-mapping-quality and --min-base-quality parameters should be raised to exclude these artifacts.

Another source of false positives is the inclusion of duplicate reads. PCR duplicates inflate the read depth at specific positions and can create the appearance of a variant where none exists. Marking duplicates with picard MarkDuplicates or samtools markdup before variant calling reduces this artifact.

Strand Bias and Read Position Bias

Systematic sequencing errors often show strand bias, where the alternate allele appears predominantly on one strand. FreeBayes does not apply strand bias filtering by default, so these artifacts can pass through the initial variant calling step. Downstream filtering with tools that assess strand bias, such as bcftools filter with the SB annotation, can remove these artifacts.

Read position bias occurs when the alternate allele appears predominantly at read ends, where base calling errors are more common. This artifact is particularly problematic for indels, which are difficult to align at read ends. Filtering based on read position bias requires tools that annotate the position of variant-supporting reads within their full read length.

Reference Bias in Variant Representation

FreeBayes reports variants relative to the reference genome, and the representation of indels can vary depending on the local sequence context. The same biological variant can be represented as different alleles in different samples if the alignment differs. This representation issue complicates comparison across samples and requires normalization of variant representation before downstream analysis.

Variant normalization tools, such as bcftools norm, standardize the representation of variants by left-aligning indels and splitting multi-nucleotide variants. This normalization step is essential for comparing variants across samples and for merging variant calls from different callers.

Frequently Asked Questions

What is the difference between germline and somatic variant calling with FreeBayes?

Germline variant calling identifies variants that are inherited and present in all cells of an individual. These variants are expected at approximately 50 percent allele fraction for heterozygous calls and 100 percent for homozygous calls. Somatic variant calling identifies acquired mutations that are present in only a subset of cells, typically in tumor samples. Somatic variants often have lower allele fractions, and somatic calling requires comparison between tumor and normal samples to distinguish acquired mutations from inherited variation. FreeBayes can be used for both applications, but the parameters and interpretation differ based on the expected allele frequency spectrum.

How should I set the minimum alternate allele fraction for germline calling?

The default minimum alternate allele fraction of 0.2 is appropriate for most germline calling applications. This threshold accommodates heterozygous variants that may show allelic imbalance due to sequencing error or mapping bias. Lowering the threshold to 0.1 increases sensitivity for variants with skewed allele fractions but also increases false positive calls. Raising the threshold to 0.3 reduces false positives but may miss genuine heterozygous variants with significant allelic imbalance. The optimal threshold depends on the sequencing platform, coverage, and the tolerance for false positives in the specific application.

Why is my FreeBayes output missing known variants?

Missing known variants in FreeBayes output typically results from overly stringent filtering parameters or inadequate sequencing coverage. The --min-coverage parameter set too high will exclude variants in regions with below-average coverage. The --min-alternate-fraction parameter set too high will exclude heterozygous variants with skewed allele fractions. The --min-alternate-count parameter set too high will exclude variants supported by few reads. Examine the distribution of coverage and allele fractions in your data to identify which parameter is excluding the expected variants.

How do I run FreeBayes on multiple samples jointly?

Joint calling across multiple samples uses multiple --bam arguments in a single FreeBayes command. Each BAM file must have unique read group information with the sample name specified in the read group header. Joint calling evaluates each genomic position across all samples simultaneously, which improves sensitivity for variants present in one sample but supported by reads from another. The output VCF contains genotype information for all samples, enabling direct comparison of variant genotypes across the cohort.

What is the recommended approach for filtering FreeBayes output?

The recommended filtering approach combines hard filters with annotation-based prioritization. Hard filters remove variants with low genotype quality, extreme read depth, or low mapping quality. Annotation adds functional information about gene location, amino acid changes, and population frequency. The filtered and annotated variant set is then prioritized based on the application, with rare, protein-altering variants receiving the highest priority for validation. The specific filtering thresholds should be adjusted based on the sequencing platform and the tolerance for false positives in the application.

How does reference genome choice affect FreeBayes results?

The reference genome version directly affects variant calling results because variants are reported relative to the reference sequence. Using the current GRCh38 assembly avoids the complications of coordinate liftover from older builds and captures novel, medically relevant loci that are not represented in older references (Wellcome Open Research). The reference used for alignment must match the reference used for variant calling. Mixing references between alignment and variant calling introduces systematic errors because read mapping positions and base qualities are reference-dependent.

What are the main causes of false positive variant calls in FreeBayes?

False positive variant calls in FreeBayes output result from sequencing errors, mapping errors, and systematic artifacts. Sequencing errors at read ends produce spurious SNPs. Reads with low mapping quality in repetitive regions produce false variants. PCR duplicates inflate read depth and create the appearance of variants where none exist. Strand bias and read position bias indicate systematic errors that should be filtered. Raising the minimum base quality and mapping quality thresholds and marking duplicates before variant calling reduce false positive calls.

When should I use an alternative variant caller instead of FreeBayes?

Alternative variant callers may be appropriate when FreeBayes does not meet the needs of the specific application. Callers that use haplotype-based assembly, such as GATK HaplotypeCaller, may perform better in regions with complex variation. Callers designed for specific data types, such as long-read sequencing, are appropriate when using those platforms. The choice of variant caller should be based on the specific application, the data type, and the performance characteristics of each tool. Comparing results from multiple callers can identify variants that are consistently called across approaches, which increases confidence in the variant set.

Related Bioinformatics Guides

Related Clinical & Scientific Guides

References and Further Reading

This article is educational and does not replace validated analysis plans, institutional policy, clinical interpretation, or specialist review.