Stochastic modelling, Bayesian inference, and new in vivo measurements elucidate the debated mtDNA bottleneck mechanism

  1. Iain G Johnston
  2. Joerg P Burgstaller
  3. Vitezslav Havlicek
  4. Thomas Kolbe
  5. Thomas Rülicke
  6. Gottfried Brem
  7. Jo Poulton
  8. Nick S Jones  Is a corresponding author
  1. Imperial College London, United Kingdom
  2. IFA Tulln, Austria
  3. University of Veterinary Medicine Vienna, Austria
  4. University of Veterinary Medicine, Austria
  5. University of Natural Resources and Life Sciences, Austria
  6. University of Oxford, United Kingdom
10 figures and 3 tables


The mitochondrial bottleneck, and elements of a general model for bottlenecking mechanisms.

(A) The mitochondrial DNA (mtDNA) bottleneck acts to produce a population of oocytes with varying heteroplasmies from a single initial oocyte with a specific heteroplasmy value. During development, mtDNA copy number per cell decreases (by a debated amount, which we address; see Main text) then recovers, suggesting a ‘bottleneck’ of cellular mtDNA populations. (B) Cellular mtDNA populations during the bottleneck are modelled as containing wildtype and mutant mtDNAs. MtDNAs can replicate and degrade within a cell cycle, with rates λ and ν respectively. (C) At cell divisions, the mtDNA population is partitioned between two daughter cells either deterministically, binomially, or through the binomial partitioning of mtDNA clusters. (D) Symbols used to represent quantities and model parameters used in the Main text, and their biological interpretations.
Different mechanisms for the mtDNA bottleneck.

(A) Trajectories of mean copy number E(m) and normalised heteroplasmy variance V(h) arising from the models described in the text, optimised with respect to data from experimental studies. Birth-death-partition (BDP) denotes the BDP model, encompassing Cree and Wai (A) mechanisms. Left plots show trajectories during development; right plots show behaviour in mature oocytes in the next generation. * denotes measurements in mature oocytes, modelled as 100 dpc (see ‘Materials and methods’). (B) Statistical support for different mechanisms from approximate Bayesian computation (ABC) model selection with thresholds ϵ1,2,3,4 = 75, 60, 50, 45. As the threshold decreases, forcing a stricter agreement with experiment (thinner, darker columns), support converges on the BDP model.
Parameterisation of the BDP model and inferred details of bottleneck mechanism.

Trajectories of (A) mean copy number E(m) and (B) normalised heteroplasmy variance V(h) resulting from BDP model parameterisations sampled using ABC with a threshold ϵ = 40. * denotes measurements in mature oocytes, modelled as 100 dpc (see ‘Materials and methods’). Note: the range in (B) does not correspond to a credibility interval on individual measurements, but rather on an expected underlying (population) variance, from which individual variance measurements are sampled. We thus expect to see, for example, several measurements lower than this range due to sampling limitations (see text). (C) Posterior distributions on mtDNA turnover ν with time. (D) Posterior distribution on min E(m), the minimum mtDNA copy number reached during development. (E) Posterior distribution on σ=i=36τiνi, a measure of the total amount of mtDNA turnover.
Predictions and experimental verification of the BDP model.

(A) New V(h) measurements from the HB mouse system, with optimised fits for the BDP, Wai (b) and Cao models. (B) Posterior probabilities of each model given this data under decreasing ABC threshold: ϵ = {50, 40, 30, 25}. (C) All V(h) measurements from the HB model (points) with inferred V(h) behaviour from ABC applied to the BDP model (red curves). As in Figure 3, this range does not correspond to a credibility interval on individual measurements, but rather on an expected underlying (population) variance, from which individual variance measurements are sampled. The inferred behaviour strongly overlaps with the inferred behaviour for the BALB/c system (blue curves), suggesting that the BDP model applies to a genetically diverse range of systems. (D) Heteroplasmy distributions. The transformation h=ln|(h11)E(h)/(1E(h))| (Burgstaller et al., 2014) is used to compare distributions with different mean heteroplasmy. Red jitter points are samples from sets used to parameterise the BDP model; red curves show the 95% range on transformed heteroplasmy with time inferred from these samples. Blue jitter points are samples withheld independent from this parameterisation; their distributuions fall within the independently inferred range. Insets show, in untransformed space, distributions of the withheld heteroplasmy measurements (blue) compared to parameterised predictions (red); no withheld datasets show significant support against the predicted distribution (Anderson-Darling test, p < 0.05).
Figure 4—source data 1

Individual heteroplasmy measurements in the HB mouse model contributing to the new heteroplasmy variance data used to test our theory.
Quantitative influences and clinical results from our bottlenecking model.

(AC) Trajectories of copy number E(m) and normalised heteroplasmy variance V(h) resulting from perturbing different physical parameters. Trajectory C labels the ‘control’ trajectory resulting from a fixed parameterisation; black dots show experimental data; * denotes measurements from primary oocytes, modelled at 100 dpc. (A) Increasing (T+) and decreasing (T) mtDNA turnover (both mtDNA replication and degradation) by 20%. (B) Increasing (M+) and decreasing (M) mtDNA degradation throughout development by a constant value (2 × 10−4, in units of day−1), while keeping replication constant. (C) Applying a positive (S+) and negative (S) selective pressure to mutant mtDNA by 5 × 10−6 day−1. (D) Probability of crossing different heteroplasmy thresholds h* with time, starting with initial heteroplasmy h0 = 0.3. (E) Probability distributions over embryonic heteroplasmy h given a measurement hm from preimplantation sampling (** hm = 0.1; *** hm = 0.4) at different times.
Model for the mtDNA bottleneck.

A summary of our findings. (A) There is most statistical support for a bottlenecking mechanism whereby mtDNA dynamics is stochastic within a cell cycle, involving random replication and degradation of mtDNA, and mtDNAs are binomially partitioned at cell divisions. (B) This mechanism results in heteroplasmy variance increasing both due to stochastic partitioning at divisions and due to random turnover. The absolute magnitude of the copy number bottleneck is not critical: a range of bottleneck sizes can give rise to observed dynamics. Random turnover of mtDNA increases heteroplasmy variance through folliculogenesis and germline development.
Appendix figure 1
Heteroplasmy variance in a model system under several different group-inheritance regimes.

V(h) over many cell divisions when the elements of inheritance are heteroplasmic or homoplasmic groups of different size. Groups may be quenched (Q; constituents remain the same across cell divisions) or unquenched (UQ; constituents are randomly resampled from the cellular population each cell cycle); for homoplasmic clusters, an unquenched protocol yields identical results to the quenched protocol. V(h) behaviour differing from binomial partitioning (c = 1) is only observed for homoplasmic groups with c ≥ 2. Points for heteroplasmic groups are slightly offset in the x-direction for clarity.
Appendix figure 2
Residual distributions at different ABC thresholds ϵ.

The distribution of squared residuals corresponding to individual experimental datapoints compared to an ensemble of simulated trajectories for (top) logE(m) (bottom) V(h). The V(h) residuals are scaled by A1 = 103 to ensure that the two sets of measurements are compared on a quantitatively equal footing. As ϵ is decreased (ϵ1,2,3,4 = 40, 50, 75, 100), distributions of residuals from accepted trajectories tighten around zero.
Appendix figure 3
Posterior distributions on model parameters.

The posterior distributions on individual model parameters, assuming the inferred BDP bottlenecking mechanism. Replication rates are presented as κ = λν, thus representing overall proliferation rates of mtDNA. Units are omitted for clarity. Pale, single-values distributions correspond to parameter values fixed within the model (κ6 = 0 to prevent mtDNA proliferation after development; τ1 = 7 hr, τ2 = 16 hr fixed by data on cell doubling times; h0 = 0.2 fixed for simplicity as heteroplasmy variances are normalised; δλ = 0 fixed to avoid varying selective pressure). The ‘turnover’ parameter, described in the text, is i=36τiνi, a measure of the total random turnover in the mtDNA population.
Appendix figure 4
Comparison of truncated Normal approximation with exact heteroplasmy distribution.

Representations of heteroplasmy distributions at a time t = 21 dpc, with various starting heteroplasmies, using (as an example) the maximum likelihood parameterisation emerging from the inference procedure in the Main text. Dark lines and bars show exact distributions from Equation 72; pale lines and bars show distributions arising from the truncated Normal distribution described in the text.


Table 1

Source data used in this study
Time/dpcE(m)NSource study
14.51.8e32615Cree, §
5.51.5e385Cao, §
6.51.6e343Cao, §
7.52.0e353Cao, §
7.752.0e342Cao, §
8.52.1e382Cao, §
9.52.2e393Cao, §
10.52.1e374Cao, §
11.52.0e367Cao, §
12.51.9e3124Cao, §
13.52.1e371Cao, §
Time/dpcV(h)NSource study
*0.14026Jenuth, **
*0.01724Jenuth, **
*0.14431Jenuth, **
*0.02149Jenuth, **
*0.03331Jenuth, **
8.50.01620Wai††, ‡‡
8.50.01820Wai††, ‡‡
9.50.01820Wai††, ‡‡
10.50.01720Wai††, ‡‡
10.50.02520Wai††, ‡‡
10.50.00820Wai††, ‡‡
12.50.01720Wai††, ‡‡
13.50.01420Wai††, ‡‡
13.50.01520Wai††, ‡‡
13.50.02020Wai††, ‡‡
13.50.02120Wai††, ‡‡
14.50.00420Wai††, ‡‡
16.50.01520Wai††, ‡‡
25.00.00220Wai††, ‡‡
25.00.00220Wai††, ‡‡
29.00.01320Wai††, ‡‡
29.00.02720Wai††, ‡‡
32.00.01620Wai††, ‡‡
32.00.02320Wai††, ‡‡
35.00.02620Wai††, ‡‡
35.00.02920Wai††, ‡‡
50.00.03920Wai††, ‡‡
51.10.04320Wai††, ‡‡
65.00.02620Wai††, ‡‡
  1. Data referenced by number of cells post–conception is assigned a time measurement assuming the 7 hr → 16 hr doubling times from Lawson and Hage (1994).

  2. Mean copy number taken directly from tabulated data.

  3. §

    (Weighted) average over germline cell classes presented at this time point.

  4. #

    Extracted from data in figures; n not explicitly available so estimated as n = 20 from accompanying histograms and discussion.

  5. Manually normalised from given data.

  6. **

    Data from mature oocytes in next generation: time in dpc not available.

  7. ††

    Extracted from data in figure in correspondence following study.

  8. ‡‡

    n not explicitly available so estimated as n = 20 from accompanying histograms and discussion in original paper.

Table 2

New heteroplasmy measurements from the HB model system
h × 10040.826.57.416.726.925.
  1. Heteroplasmy measurements and statistics from the HB model system. Ages are given in days after birth.

Table 3

Experiments for further elucidation of the mtDNA bottleneck
MtDNA copy number before and after cell divisions and/or variance of copy number between daughter cellsTo elucidate mechanism of mtDNA partitioning and whether this partitioning is deterministic or stochastic
Copy number trajectories with different mtDNA heteroplasmiesTo assess the modulation of copy number dynamics by mtDNA heteroplasmy via retrograde signalling
Measurement of mean heteroplasmy through development, with a variety of mtDNA type pairingsTo assess and quantify to what extent selection modulates mtDNA dynamics during germline development
Copy number measurements after upregulation of mitophagyTo assess the presence and strength of compensatory mechanisms that may act to preserve mtDNA copy number—and hence whether upregulating mitophagy will act to increase mtDNA turnover or simply lower copy number
Heteroplasmy variance after upregulation of mitophagyTo assess the efficacy of mitophagy for increasing the power of the bottleneck
Heteroplasmy distribution in cells after the bottleneck from sampled/known initial heteroplasmyTo confirm predictions for threshold crossing and statistics between generations
BrU incorporation in oocytes between 30 and 40 dpcTo confirm the random turnover mechanism: we expect a large proportion of BrU incorporation subset of mtDNAs to be observed in this time period (see section ‘Observation of a subset of replicating genomes’)
Mitochondrial ultrastructure and mtDNA localisation during developmentTo assess and characterise any potential modulation of the size of units of mitochondrial inheritance by mitochondrial dynamics through development, in particular, investigating whether there is time-varying modulation of cluster size at points of division

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Iain G Johnston
  2. Joerg P Burgstaller
  3. Vitezslav Havlicek
  4. Thomas Kolbe
  5. Thomas Rülicke
  6. Gottfried Brem
  7. Jo Poulton
  8. Nick S Jones
Stochastic modelling, Bayesian inference, and new in vivo measurements elucidate the debated mtDNA bottleneck mechanism
eLife 4:e07464.