The B-value calculator: expected diversity under background selection

This paper introduces Bvalcalc, a Python-based command-line tool that efficiently calculates genome-wide expected diversity under background selection (B-values) at single base-pair resolution by integrating analytical theory with recombination maps and various biological parameters, and validates its accuracy through simulations and applications to human, fruit fly, and Arabidopsis genomes.

The Gist

Imagine your genome (your complete set of DNA instructions) as a massive, bustling city. In this city, most of the buildings are just empty lots or simple houses (neutral sites), but there are also critical skyscrapers and power plants (conserved sites) that are absolutely vital for the city to function.

The Problem: The "Bad Neighbor" Effect
In this city, the "power plants" are constantly under attack by vandals (deleterious mutations). The city's security force (natural selection) has to work overtime to catch these vandals and remove them.

Here's the catch: When the security force is busy chasing a vandal near a power plant, they are less available to patrol the empty lots nearby. Because the security is distracted, the empty lots near the power plants end up looking a bit messier and less diverse than they would if the security force were free.

In genetics, this is called Background Selection (BGS). It's the idea that the struggle to keep important parts of the genome clean accidentally "drags down" the diversity of the neutral parts nearby. Scientists call the measure of this effect the B-value.

• B = 1.0: Everything is perfect; no one is distracted.
• B = 0.5: The area is only half as diverse as it should be because the "security" is too busy.

The Old Way vs. The New Tool
For a long time, calculating exactly how much diversity is lost in every single neighborhood of this genetic city was like trying to count every grain of sand on a beach by hand. It was slow, required a PhD in math, and the tools were clunky.

Enter Bvalcalc (The B-value Calculator). Think of this as a high-tech drone that flies over the entire genetic city in seconds.

What Bvalcalc Does (The Analogy)
Instead of counting sand grains, Bvalcalc uses a sophisticated map to predict exactly how "messy" any specific spot in the genome will be, based on how close it is to the "power plants" (genes) and how busy the "security" is.

It accounts for several real-world complications that the old models missed:

1. Gene Conversion: Imagine if two neighbors swapped their mailboxes. Sometimes, this helps clean up the mess faster. Bvalcalc knows how to factor this in.
2. Self-Fertilization: Some organisms (like certain plants) are their own neighbors. This changes how the "security" patrols. Bvalcalc adjusts the math for this.
3. Population Boom/Bust: If the city suddenly doubles in size or shrinks by half, the "security" force changes its strategy. Bvalcalc can calculate how this affects the messiness of the genome.
4. The "Unlinked" Effect: Sometimes, a problem on one side of the city (one chromosome) affects the whole city's security budget, even if you are on the other side. Bvalcalc calculates this "global distraction."

Why Do We Need This?
Scientists use this tool to solve two main mysteries:

1. Finding the "Bad Guys" (Selective Sweeps): Sometimes, a new, beneficial mutation spreads through the city like wildfire, wiping out diversity. This looks a lot like the "messiness" caused by the bad neighbors (BGS). If you don't know where the "bad neighbors" are, you might think you found a beneficial mutation when you actually just found a messy neighborhood. Bvalcalc draws a map of the "expected messiness," so scientists can subtract it and see the real beneficial mutations.
2. Reading the City's History (Demography): Scientists try to guess how the city's population grew or shrank over time by looking at how diverse the empty lots are. But if they don't account for the "bad neighbors," they might think the city shrank when it actually just had a lot of vandals. Bvalcalc helps them get the history right.

The Results
The authors tested their new drone (Bvalcalc) by simulating millions of genetic cities on a computer. The drone's predictions matched the computer simulations almost perfectly.

They then used it to map three real species:

• Humans: A huge city with many districts. The "global distraction" (unlinked effects) is huge here.
• Fruit Flies: A smaller city where the "bad neighbor" effect is very strong and obvious.
• Thale Cress (a plant): A city where everyone lives in one big family (selfing), making the "security" patrol very different.

The Bottom Line
Before Bvalcalc, understanding how the "struggle for survival" at important genes affects the rest of the genome was like trying to navigate a dark room with a flashlight. Bvalcalc turns on the lights. It gives researchers an easy, free, and accurate way to see the invisible forces shaping our DNA, helping them understand evolution, disease, and history with much greater clarity.

Technical Summary

1. Problem Statement

Background Selection (BGS) is a fundamental evolutionary force where purifying selection at conserved sites reduces genetic diversity at linked neutral sites. This reduction in effective population size ($N_e$) and distortion of the Site Frequency Spectrum (SFS) confounds demographic inference and the detection of selective sweeps.

• The Challenge: While theoretical models for BGS exist (e.g., Charlesworth et al., Hudson and Kaplan), integrating these into unified, accessible tools for real-world genomes is difficult. Existing software (e.g., calc_bkgd) is often technically demanding and lacks modern extensions like gene conversion, self-fertilization, or complex demographic histories.
• The Gap: There is a need for an accessible, analytical tool that can calculate expected neutral diversity (the B-value, defined as $B = \pi / \pi_0$) at single-base-pair resolution across entire genomes, accounting for modern population genetic complexities without requiring computationally intensive forward simulations.

2. Methodology

The authors developed Bvalcalc, a Python-based command-line interface (CLI) and API that analytically calculates B-values. The methodology integrates several theoretical extensions to classic BGS models:

• Core Calculation: The total B-value ($B_{tot}$) at a focal site is calculated as the product of linked effects from conserved elements on the same chromosome and unlinked effects from other chromosomes:
  $$B_{tot} = \left( \prod_{i=1}^{E} B_i \right) \times B_{unlinked}$$
• Distribution of Fitness Effects (DFE): The tool integrates over a DFE modeled as four non-overlapping uniform distributions (effectively neutral, weakly, moderately, and strongly deleterious mutations) to calculate the cumulative effect of purifying selection.
• Key Extensions Implemented:
  • Gene Conversion: New analytical expressions account for gene conversion rates and tract lengths, modifying the crossover/gene conversion parameters in the integration.
  • Self-Fertilization (Selfing): The model adjusts population size, crossover rates, gene conversion rates, and dominance coefficients using Wright's inbreeding coefficient ($F$).
  • Unlinked Effects: A derived analytical expression (Equation 4) calculates the reduction in diversity caused by conserved sites on other chromosomes, correcting previous assumptions that failed for strongly deleterious mutations.
  • Demography: A two-epoch model allows for instantaneous population size changes (bottlenecks or expansions), scaling the B-value based on the time since the change.
  • Hill-Robertson Interference (HRI): While the default model assumes no HRI, the authors provide an experimental API (calculateB_hri) for non-recombining regions using coarse-graining methods to estimate $B'$ (B-values with interference).
• Validation: The analytical results were validated against forward-in-time simulations using SLiM 4.0.1. Parameters included varying DFEs, gene conversion, population size changes, and selfing rates.

3. Key Contributions

• Bvalcalc Software: A user-friendly, open-source Python package that enables the analytical calculation of B-maps (genome-wide B-values) at single-base resolution for both model and non-model species.
• Theoretical Unification: The paper unifies classic BGS theory with modern extensions (gene conversion, selfing, unlinked effects, and demography) into a single coherent framework.
• Pre-built Templates: The authors provide pre-configured parameter templates for Homo sapiens, Drosophila melanogaster, and Arabidopsis thaliana, lowering the barrier to entry for empirical researchers.
• New Analytical Derivations: The authors derived new expressions for unlinked BGS effects that remain accurate for strongly deleterious mutations, correcting limitations in previous literature (e.g., Charlesworth 2012).

4. Results

• Validation against Simulations:
  • Bvalcalc accurately recovered mean B-values from SLiM simulations across 10,000 replicates for basic models, gene conversion, and population expansions.
  • Demographic Nuance: The tool performed well for population expansions but showed deviations (underestimation of $B$) for severe population contractions, particularly when the contraction occurred long ago ($1N_{anc}$ generations), suggesting limits in the current two-epoch approximation for complex bottlenecks.
  • Unlinked Effects: The new Equation 4 accurately predicted unlinked BGS effects across various DFEs, whereas previous methods significantly underestimated diversity for strongly deleterious mutations.
  • HRI: The HRI extension ($B'$) was accurate for single-class DFEs but underestimated diversity when multiple mutation classes were present in non-recombining regions.
• B-Maps for Model Organisms:
  • Humans: Showed the highest mean $B$ (0.89 for CDS), but unlinked effects contributed 87% to the total reduction in diversity due to large genome size.
  • D. melanogaster: Showed strong BGS effects ($B \approx 0.64$ when including UTRs/phastCons). The B-map showed a high correlation ($R^2 = 0.76$ at 1Mb) with empirical nucleotide diversity, outperforming previous maps at finer scales (100kb/10kb).
  • A. thaliana: Exhibited the lowest mean $B$ (0.65) due to high selfing rates ($F \approx 0.97$), which intensifies BGS. This is the first published B-map for this species.
• Correlation with Empirical Data: The correlation between predicted $B$ and observed diversity was strongest in D. melanogaster and weakest in humans. The authors attribute the lower human correlation to the dominance of unlinked effects (reducing local variance), potential inaccuracies in mutation rate maps, and temporal shifts in recombination hotspots.

5. Significance

• Improved Inference: By providing accurate B-maps, Bvalcalc allows researchers to develop better null models for demographic inference (e.g., correcting PSMC+ or SFS-based methods) and selective sweep detection (e.g., SweepFinder2), reducing false positives/negatives caused by BGS.
• Accessibility: The tool democratizes complex population genetic modeling, making it accessible to non-technical biologists through comprehensive documentation and pre-built templates.
• Scalability: It offers a computationally efficient alternative to whole-genome forward simulations, which are often infeasible for large genomes or complex parameter spaces.
• Future Directions: The paper highlights that while HRI and complex demography remain challenges, Bvalcalc provides a robust foundation. It enables the study of non-model species where B-maps were previously unavailable, facilitating a deeper understanding of how selection shapes genomic variation across the tree of life.

Availability: The software and documentation are available at johrilab.github.io/Bvalcalc/, and pre-computed B-maps are hosted on Zenodo.