Understanding the impact of numerical solvers on inference for differential equation models

Most ordinary differential equation (ODE) models used to describe biological or physical systems must be solved approximately using numerical methods. Perniciously, even those solvers that seem sufficiently accurate for the forward problem, i.e. for obtaining an accurate simulation, might not be sufficiently accurate for the inverse problem, i.e. for inferring the model parameters from data. We show that for both fixed step and adaptive step ODE solvers, solving the forward problem with insufficient accuracy can distort likelihood surfaces, which might become jagged, causing inference algorithms to get stuck in local ‘phantom’ optima. We demonstrate that biases in inference arising from numerical approximation of ODEs are potentially most severe in systems involving low noise and rapid nonlinear dynamics. We reanalyse an ODE change point model previously fit to the COVID-19 outbreak in Germany and show the effect of the step size on simulation and inference results. We then fit a more complicated rainfall run-off model to hydrological data and illustrate the importance of tuning solver tolerances to avoid distorted likelihood surfaces. Our results indicate that, when performing inference for ODE model parameters, adaptive step size solver tolerances must be set cautiously and likelihood surfaces should be inspected for characteristic signs of numerical issues

A Monte Carlo method to estimate cell population heterogeneity from cell snapshot data

Variation is characteristic of all living systems. Laboratory techniques such as flow cytometry can probe individual cells, and, after decades of experimentation, it is clear that even members of genetically identical cell populations can exhibit differences. To understand whether variation is biologically meaningful, it is essential to discern its source. Mathematical models of biological systems are tools that can be used to investigate causes of cell-to-cell variation. From mathematical analysis and simulation of these models, biological hypotheses can be posed and investigated, then parameter inference can determine which of these is compatible with experimental data. Data from laboratory experiments often consist of “snapshots” representing distributions of cellular properties at different points in time, rather than individual cell trajectories. These data are not straightforward to fit using hierarchical Bayesian methods, which require the number of cell population clusters to be chosen a priori. Nor are they amenable to standard nonlinear mixed effect methods, since a single observation per cell is typically too few to estimate parameter variability. Here, we introduce a computational sampling method named “Contour Monte Carlo” (CMC) for estimating mathematical model parameters from snapshot distributions, which is straightforward to implement and does not require that cells be assigned to predefined categories. The CMC algorithm fits to snapshot probability distributions rather than raw data, which means its computational burden does not, like existing approaches, increase with the number of cells observed. Our method is appropriate for underdetermined systems, where there are fewer distinct types of observations than parameters to be determined, and where observed variation is mostly due to variability in cellular processes rather than experimental measurement error. This may be the case for many systems due to continued improvements in resolution of laboratory techniques. In this paper, we apply our method to quantify cellular variation for three biological systems of interest and provide Julia code enabling others to use this method

The power of evolutionary rescue is constrained by genetic load

International audienceThe risk of extinction faced by small isolated populations in changing environments can be reduced by rapid adaptation and subsequent growth to larger, less vulnerable sizes. Whether this process, called evolutionary rescue, is able to reduce extinction risk and sustain population growth over multiple generations is largely unknown. To understand the consequences of adaptive evolution as well as maladaptive processes in small isolated populations, we subjected experimental Tribolium castaneum populations founded with 10 or 40 individuals to novel environments, one more favorable, and one resource poor, and either allowed evolution, or constrained it by replacing individuals one-for-one each generation with those from a large population maintained in the natal environment. Replacement individuals spent one generation in the target novel environment before use to standardize effects due to the parental environment. After eight generations we mixed a subset of surviving populations to facilitate admixture, allowing us to estimate drift load by comparing performance of mixed to unmixed groups. Evolving populations had reduced extinction rates, and increased population sizes in the first four to five generations compared to populations where evolution was constrained. Performance of evolving populations subsequently declined. Admixture restored their performance, indicating high drift load that may have overwhelmed the beneficial effects of adaptation in evolving populations. Our results indicate that evolution may quickly reduce extinction risk and increase population sizes, but suggest that relying solely on adaptation from standing genetic variation may not provide long-term benefits to small isolated populations of diploid sexual species, and that active management facilitating gene flow may be necessary for longer term persistence

Data from: The power of evolutionary rescue is constrained by genetic load

Extinction risk of small isolated populations in changing environments can be reduced by rapid adaptation and subsequent growth to larger, less vulnerable sizes. Whether this process, called evolutionary rescue, is able to reduce extinction risk and sustain population growth over multiple generations is largely unknown. To understand the consequences of adaptive evolution as well as maladaptive processes in small isolated populations, we subjected experimental Tribolium castaneum populations founded with 10 or 40 individuals to novel environments, one more favorable, and one resource poor, and either allowed evolution, or constrained it by replacing individuals one-for-one each generation from a non-adapting large population to minimize both adaptive and non-adaptive evolutionary processes. Replacement individuals spent one generation in the target novel environment before use to standardize effects due to the parental environment. After 8 generations we mixed a subset of surviving populations to facilitate admixture, allowing us to estimate drift load by comparing performance of mixed to unmixed groups. Evolving populations had reduced extinction rates, and increased population sizes in the first four to five generations compared to populations where evolution was constrained. Performance of evolving populations subsequently declined. Admixture restored their performance, indicating high drift load that may have overwhelmed the beneficial effects of adaptation in evolving populations. Our results indicate that evolution may quickly reduce extinction risk and increase population sizes, but suggest that relying solely on adaptation from standing genetic variation may not provide long-term benefits to small isolated populations of diploid sexual species, and that active management facilitating gene flow may be necessary for longer-term persistence