
Mathematical biology has long relied on mechanistic models, including ordinary and partial differential equations, stochastic systems, and agent-based models, to study biological processes across scales. These approaches remain central because they provide structure, interpretability, and biological insight. However, modern biological and biomedical data, including longitudinal clinical records, medical imaging, and multi-omics measurements, are often high-dimensional, noisy, heterogeneous, and incomplete. These features make model calibration, simulation, and uncertainty quantification increasingly difficult. As a result, artificial intelligence (AI) and machine learning (ML) are playing a growing role in mathematical biology, not as replacements for mechanistic modeling, but as complementary tools that can support prediction, hybrid modeling, inference, and control. In this review, we organize these uses of AI from a mathematical biology perspective and examine how data-driven and mechanistic approaches interact across different modeling tasks. The selected examples span biomedical, epidemiological, ecological, evolutionary, biochemical, and population-level systems, while the review is intended to be representative rather than exhaustive. We emphasize recurring challenges that strongly affect biological credibility and practical usefulness, including interpretability, identifiability, generalization, and uncertainty quantification. Our goal is to clarify the roles AI can play in mathematical biology and to highlight the opportunities and limitations that arise when flexible learning methods are integrated with biologically grounded modeling.
This paper addresses the increasing need for comprehensive mathematical descriptions of cell organization by examining the algebraic structure of mitochondrial network dynamics. Mitochondria are cellular structures involved in metabolism that take the form of a network of membrane-based tubes that undergo continuous re-arrangement by a set of morphological processes, including fission and fusion, carried out by protein-based machinery. Because of their network structure, mitochondria can be represented as graphs, and the morphological operations that take place in the cell, referred to as mitochondrial dynamics, can be represented by changes to the graphs. Prior studies have classified mitochondrial graphs based on graph-theoretic features, but an alternative approach is to focus not on the graphs themselves but on the set of morphological operations inducing mitochondrial dynamics, since this may provide a simpler representation. Moreover, the operations are what determine the graphs that will be generated in a biological system. Here we show that mitochondrial dynamics give rise to a category in which the objects are equivalence classes of graphs defined by one of the morphological operations and morphisms are mappings between these equivalence classes defined by the remaining morphological operations. For mitochondria consisting of a single component this gives rise to a particularly simple representation. Using these formalisms we define a distance metric for similarity between mitochondrial structures based on an edit distance, and demonstrate how this representation can be used for visualization and statistical analysis of biological data. In the course of defining these structures we provide a mathematical motivation for new experimental questions regarding mitochondrial fusion, the impacts of cell division on mitochondrial morphology, and the presence of a single giant component in some cell types. This work points to a general strategy for formulating a cell structure state-space, based not on the shapes of cellular structures, but on relations between the dynamic operations that produce them.
Tumour recurrence after oncolytic virotherapy is a critical clinical challenge, as solid cancers tend to persist and spread after such therapies alone. Often, mathematical models consider tumour growth without taking possible cooperative cell behaviours into account. Instead, the well-known Allee effect has the ability to drive low-density tumour dynamics and affect therapeutic outcomes. For this reason, we here explore the contribution of Allee effects into two of the most common cancer growth paradigms, e.g. the logistic and the Gompertz frameworks. Our analysis reveals that density-dependent cooperation fundamentally alters virotherapy efficacy. Weak effects can enable pseudo-extinction states where tumour populations collapse to undetectable levels, whilst strong effects can interestingly create non trivial, bistable outcomes. Initial tumour burden seems to be determinant: bifurcation analysis shows scenarios where modest increase in infection or decrease in clearance rates can shift outcomes dramatically. Notably, the Gompertz model exhibits multistability under marginal Allee thresholds, pointing at a possible explanation for spontaneous remission that have been clinically observed. Overall, these findings suggest that Allee effects may be an important factor in the future to adjust dosing schedules, reduce viral loads, lower recurrence risk and shorten therapeutic times. Our framework may also offer quantitative guidance for patient-specific regimes for virotherapy and combination therapies.
The mechanistic target of rapamycin complex 1 (mTORC1) has been implicated in coronavirus pathogenesis, yet its precise role in shaping antiviral defenses in pneumocytes remains unresolved. This study combines ex vivo human lung tissue assays with a literature-curated, logic-based (Boolean) signaling model to study how mTORC1 influences SARS-CoV-2 replication and type-I interferon (IFN) responses. Our in-vitro data showed that pharmacologic inhibition of mTORC1 with sirolimus was associated with reduced viral replication and increased IFN- β expression across donor samples. To interpret these data, we constructed a network integrating inflammation, stress, apoptosis, and RIG-I–IFN signaling, and evaluated four mechanistic hypotheses for explaining our experiments in which mTORC1 either promotes replication, inhibits IFN, both, or neither. The model was constructed according to the best evidence in the literature and evaluated against 56 independent protein/phosphoprotein readouts from 26 studies not used for model construction; the model achieved high qualitative accuracy overall (F1 ≳ 0.75 ) and on live-virus experiments (F1 ≳ 0.9 ). We found that the hypotheses where mTORC1 promotes viral replication, inhibits IFN, or both were able to qualitatively reproduce the empirical data, i.e., mTORC1 inhibition by sirolimus leads to lower viral replication and higher IFN expression. We queried the robustness of those predictions through a systematic edge knock-in/knockout screen generating 6,328 network variants per scenario (25,312 total). We found that the scenario where mTORC1 enhances viral replication is the more plausible, due to parsimony (this mechanism is also observed in other viral infections), higher robustness (more variants still reproduce the result), and stronger inhibition that is more likely to be observed in biological experiments with low statistical power. This scenario supports the idea that the mTORC1 effects on IFN expression are indirectly mediated by viral replication itself.
Cassava mosaic disease (CMD), caused by the cassava mosaic virus (CMV) and transmitted primarily by the whitefly Bemisia tabaci, reflects the most severe and widespread viral threat to cassava cultivation. The emergence of CMD with a nonlinear saturated incidence of Holling type II form is stochastically modeled and analyzed in this paper utilizing the continuous-time Markov chain (CTMC) modeling approach. The most significant distinction between deterministic and stochastic models is that the deterministic model predicts disease persistence when the basic reproduction number R 0 > 1 , whereas the stochastic model, characterized by the branching-process threshold ρ ( M ) , indicates that disease extinction can still occur in finite time due to random fluctuations, even when R 0 > 1 . The likelihood of disease extinction is determined using the Galton-Watson branching process (GWbp) approximation and compared to the estimated probability obtained from the stochastic model's 10,000 sample paths. This comparison reveals a substantial link between these probabilities. It is found that the disease has a higher probability of extinction when transmission occurs solely through infected vectors, as opposed to transmission through infected plants or both infected plants and vectors. In the stochastic settings, the implicit equation for the mean first passage time is derived to assess the typical time until the first state transition. Additionally, we evaluate both the quasi-stationary distribution of infected individuals and the probability distribution of the epidemic's ultimate size.
Climate variability plays a fundamental role in desert locust population dynamics, yet most existing modelling and optimal control frameworks rely on simplified representations of seasonal forcing that inadequately capture short-term climatic variability. To address this limitation, we develop a climate-informed predictive optimal control framework that integrates a mechanistic stage- and phase-structured desert locust model with machine-learning-based climate forecasting. Long short-term memory networks and gradient-boosting regression are employed to predict temperature and rainfall, respectively, providing data-driven climatic inputs to a non-autonomous system describing locust development, reproduction, mortality, vegetation dynamics, and phase transitions. The mathematical analysis establishes the well-posedness of the model through the existence, uniqueness, positivity, and boundedness of solutions, while an autonomous reduction yields a closed-form basic offspring number that provides analytical insight into invasion potential under representative climatic conditions. A predictive optimal control problem is formulated to evaluate stage-specific physical, biological, and chemical interventions targeting juvenile and adult locust populations. Numerical simulations demonstrate that machine-learning-derived climate forcing more accurately reproduces observed climatic variability than conventional harmonic forcing, thereby providing improved inputs for predictive population modelling. Integrated juvenile–adult intervention strategies consistently outperform single-stage controls by reducing locust abundance, preserving vegetation, and lowering the composite management index. Among the intervention strategies considered, chemical control provides the greatest short-term suppression, biological control offers a more environmentally sustainable alternative, and physical control is most effective as a complementary measure during the early stages of population growth. Robustness analyses under deterministic climate perturbations and stochastic forecast errors show that the comparative ranking of intervention strategies remains stable under realistic climate forecast uncertainty. These findings demonstrate that integrating machine-learning climate prediction with predictive optimal control provides a mathematically rigorous and flexible framework for climate-informed decision support in desert locust management and establishes a foundation for future developments incorporating spatial dynamics, probabilistic forecasting, and real-time surveillance.
We develop a novel, comprehensive, and rigorously validated mathematical framework to investigate the kinetics of amyloid- β (A β ) aggregation in the presence of biologically relevant metal ions, chelating agents, and inhibitor drugs. Building upon and extending existing aggregation models, our approach integrates metal-assisted aggregation, A β self-assembly, and therapeutic interventions within a unified and mechanistically consistent formulation. The model captures the microscopic reaction pathways governing A β dynamics and explicitly incorporates the catalytic roles of copper, zinc, and iron ions-key contributors to neurotoxic plaque formation in Alzheimer's disease. Distinctively, the framework combines dual therapeutic strategies: (i) metal chelation therapy, which sequesters free metal ions, and (ii) direct inhibition of A β aggregation. Numerical simulations across multiple kinetic regimes reveal how these interventions modulate aggregation pathways, both independently and synergistically. To further validate the model, we perform a quantitative comparison with experimental data by reconstructing aggregate morphology distributions and benchmarking them against reported AFM measurements. The model successfully captures key experimental features, including peak structure and metal-dependent heterogeneity, thereby demonstrating its predictive capability. Overall, this work provides an extended and unified modeling platform that advances the quantitative understanding of metal-mediated amyloid aggregation and offers a predictive tool for evaluating and optimizing therapeutic strategies for Alzheimer's disease.
Medication non-adherence is recognized as a critical determinant of therapeutic efficacy and safety. However, its specific impact on the concentration fluctuations of chiral drugs remains poorly understood. By incorporating dual stochasticity in dosing intervals and dosages, we developed a stochastic pharmacokinetic model for chiral drugs comprising bio-active and bio-inactive enantiomers under multiple intravenous bolus administrations. Utilizing characteristic functions and second-type Volterra integral equations, we derived explicit expressions for the expectation and variance of both the active moiety and total measured concentrations, as well as their discrepancy. Moreover, leveraging Kesten theorem and ergodic theory, the existence and the geometric ergodicity of the stationary probability distribution were mathematically demonstrated. Taking ibuprofen as a case study, we simulated various medication non-adherence scenarios to quantify the enantiomeric randomness. The results demonstrated that stochastic concentrations not only cause significant deviation from the ideal levels under perfect adherence, but also induce non-negligible discrepancy between the effective and total measured concentrations. These findings of the work provide a rigorous theoretical guidance for the characterization of two enantiomers variability of chiral drugs induced by medication non-adherence, offering new insights for enhancing therapeutic efficacy and safety.
The aim of this study was to assess the potential of FTIR spectroscopy for monitoring biochemical changes in serum samples of individuals with carotid atherosclerosis following surgical intervention. Principal Component Analysis (PCA) of FTIR spectra from serum samples reveals distinct biochemical patterns at different time points: pre-surgery, 24 h post-surgery, and 48 h post-surgery. Two spectral ranges, 800–1800 cm−1 and 2800–3000 cm−1, were analyzed. PCA demonstrated that pre-surgery samples can be clearly differentiated from those taken 24 and 48 h post-surgery. However, no significant distinction was found between the 24-hour and 48-hour post-surgery samples. For the 800–1800 cm−1 range, the first principal component (PC1) explained 77.49% of the variance, highlighting the molecular vibrations of lipids, proteins, and carbohydrates. In the 2800–3000 cm−1 range, PC1 accounted for 94.89% of the variance, primarily reflecting lipid-related vibrations. These findings indicate a clear separation between pre-surgery and post-surgery samples, with the most significant variance explained by PC1. Additionally, the Boruta algorithm identified a key spectral range between 1506 cm−1 and 1673 cm−1, critical for distinguishing the samples. Classification models, including k-Nearest Neighbors, Gradient Boosting, Support Vector Machine, and Neural Network, demonstrated excellent performance in differentiating pre-surgery and post-surgery samples. However, the models struggled to distinguish between the 24-hour and 48-hour post-surgery time points. This suggests that FTIR spectroscopy may be useful for monitoring post-surgery recovery in carotid artery atherosclerosis, although subtle changes in the biochemical profile are challenging to detect between 24 and 48 h post-surgery.
Attenuated total reflectance-Fourier transform infrared spectroscopy (ATR-FTIR) provides information on the molecular composition and structure of samples. The use of ATR-FTIR was evaluated for biochemical analysis and taxonomic differentiation of entomopathogenic nematodes (EPNs). Spectra were obtained from a small sample (pellet) of a nematode population recovered from commercial EPN packages, which was placed directly on the ATR plate. Differences in signal intensity at multiple peaks associated with biomolecules critical to the survival of EPN (trehalose, glycogen, and triglyceride) were measured and visualized using Non-Metric Multidimensional Scaling (nMDS) and Principal Component Analysis (PCA). Statistically significant differences in peak signal intensity were observed between EPN species for each biochemical parameter, providing a basis for assessing the likelihood of their performance success in the field conditions. The present study also evaluated FTIR analysis of EPN for taxonomic differentiation. Results demonstrate that FTIR can be used to identify and differentiate Steinernema and Heterorhabditis genera/species, offering a potentially faster, less expensive alternative to molecular identification techniques. Ultimately, this study demonstrates the efficacy of ATR-FTIR as a reliable method for assessing the biochemical suitability of EPN products for field applications and differentiating between EPNs.
Endometrial cancer (EC) is increasingly prevalent worldwide, highlighting the need for non-invasive blood-based diagnostic triage tools. ATR-FTIR spectroscopy enables rapid, label-free biochemical profiling of plasma or serum for experimental cancer detection. To date, no systematic review or meta-analysis has evaluated the experimental performance of infrared spectroscopy for discriminating EC from non-cancer in blood-based samples. This study synthesizes available evidence to characterize the strength, consistency, and heterogeneity of the underlying spectroscopic signal across preclinical and proof-of-concept studies. MEDLINE, Web of Science, EMBASE, Scopus, Google Scholar, and CENTRAL were searched without language restrictions. Eligible studies evaluated ATR-FTIR spectroscopy of plasma or serum using histopathology as the reference standard. Pooled sensitivity, specificity, likelihood ratios, and diagnostic odds ratios were estimated using a bivariate random-effects model, with assessment of heterogeneity, threshold effects, and publication bias. Five case–control studies comprising 1376 participants were included. For plasma-based analyses, pooled sensitivity was 0.61 (95% CI: 0.59–0.68) and specificity was 0.73 (95% CI: 0.69–0.76), with a diagnostic odds ratio of 4.23 (95% CI: 3.33–5.37). For serum-based analyses, pooled sensitivity and specificity were both 0.62 (95% CI: 0.59–0.65), with a diagnostic odds ratio of 2.65 (95% CI: 2.16–3.25). Substantial heterogeneity and significant threshold effects were observed. Current evidence supports reproducible spectroscopic differences between EC and non-cancer blood samples under experimental conditions. However, methodological heterogeneity and retrospective case–control study designs limit clinical interpretability. These findings provide a benchmark for future prospective validation rather than immediate clinical application.
Vaccine-related communication can be harnessed to curb infection; yet, in practice, it often undermines vaccine uptake and sustains transmission even when effective vaccines are available. Our goal is to understand how vaccine information dynamics may shape infection spread within a coupled information–infection modeling framework. We present a coupled information-infection modeling framework that links a standard SVIRS infection-spread model with an integrate-and-fire-inspired information-dissemination model. Vaccine-positive and vaccine-critical active groups disseminate competing information that shapes non-active/hesitant individuals’ vaccine attitudes through direct peer influence, threshold-based acceptance, and persistence of engagement, while infection prevalence can feed back by amplifying caution. The evolving vaccine attitudes modulate vaccination uptake, and the model tracks the joint evolution of information dynamics and epidemic trajectories. Our analysis shows that, under the assumed coupling, information dynamics can shift the system among qualitatively distinct infection outcomes. In particular, reducing resistance to vaccine-positive information and sustaining vaccine-positive engagement can move the system toward lower-endemic regimes more reliably than changes focused only on weakening vaccine-critical engagement, for the parameter ranges considered here. These findings highlight the potential importance of information dynamics in epidemic modeling and suggest that sustained vaccine-positive engagement can be an important qualitative mechanism for reducing long-term infection burden.
Ecological interactions shape the dynamics of natural populations in the wild. Density-dependent processes are widespread and may change the respective effects of populations on one another, for instance by shifting interactions from mutualistic to parasitic relationships. Here, we develop a general deterministic model of two interacting populations, assuming density-dependent costs and benefits for the interacting individuals within and between species. This framework aims at generalizing pre-existing population dynamics models involving competition, predation, mutualism and parasitism, by allowing ecological interactions to transition when the respective densities of interacting species change. Through ordinary differential equations and phase portrait analysis, we derive general principles governing these systems, identifying constraints on the organization of equilibria and sufficient conditions for the emergence of certain dynamic behaviors. In particular, we show that equilibrium indices alternate along isoclines under broad geometric assumptions, and that limit cycles can arise when interactions include mutualistic and parasitic phases, while they cannot be generated locally in strictly mutualistic regions where the relevant interaction signs remain fixed. This framework provides a general approach for characterizing the population dynamics of interacting species and highlights the effect of the density dependent transitions in ecological interactions.
Life and evolution require precise yet imperfect transmission of genetic information across generations. Accurate transmission maintains the genetic blueprint for phenotypes that succeed in current conditions, while some transmission error is necessary to generate heritable diversity that can further increase fitness and hedge against future environmental change. This background noise, however, risks excessive mutation accumulation that degrades fitness (“Muller’s ratchet”) or even destroys genetic information beyond recovery (“error catastrophe”). Across the history of life these competing pressures resolve into two regimes: during environmental fluctuation, a higher mutation rate aids survival; during prolonged stability, populations converge toward an evolutionarily stable state (ESS) in which mutations can only reduce fitness. We model this tension with a Fisher Information framework in which genetic inheritance is a noisy two-state channel. The nontrivial eigenvalue of the channel’s symmetric circulant transition operator governs the decay of allelic contrast across generations; selection modifies this eigenvalue by suppressing the channel’s switching (mutation) probability, acting as an anti-depolarizing force. We show this compensation is bounded: an instability boundary at p(1-p)=1/9 marks the point beyond which selection can no longer maintain informational stability — a boundary that coincides with the loss of quantum coherence in the formally equivalent depolarizing channel and with thermodynamic fine-graining failure, suggesting that genetic stability, regulatory fidelity, and structural order share a common information-geometric limit. We apply this framework to cancer, framing carcinogenesis as an informational transition: a shift from host-regulated, high-fidelity transmission that maintains tissue homeostasis to a regime in which cells retain only the genetic and epigenetic changes that maximize their own proliferation — including loss of functions that serve the host and gain of functions that improve competition for space and nutrients and evasion of the predator-like immune response. Because cancer cells operate farther from this information-geometric limit than the normal cells around them, we propose that they may pursue an informational “niche construction” strategy: producing a mutagenic, acidic, hypoxic microenvironment that pushes neighboring host cells past the same collapse boundary, driving loss of function (e.g., “T cell exhaustion”) in infiltrating immune cells. We present empirical support for this framework together with specific, experimentally testable predictions.
Polyploidy occurs in plants and animals, and is an important force in speciation and genome evolution. The main focus of this paper is the following fundamental question that was recently posed by Huber and Maher: Given the ploidy numbers of a collection of extant species, or their ploidy profile, what is the smallest number of hybridizations needed in any evolutionary history for these species to completely represent these numbers? In this paper, we shall show that this question can be rephrased in terms of addition chains and the closely related addition sequences, which have been studied for over a century in mathematics and computer science. These are sequences of natural numbers that start with 1, so that each number in the sequence larger than 1 is the sum of two other numbers arising earlier in the sequence. In our first main result, we show that finding the smallest number of hybridization events to explain a ploidy profile, or the hybrid number, is equivalent to solving the so-called addition sequence problem. This immediately implies that computing the hybridization number is computationally intractable. Even so, it also leads to new connections to representing polyploid evolution using networks. More specifically, in our second main result we show that ploidy profiles representable by tree-child networks are exactly the addition chains, implying a polynomial-time algorithm for identifying these profiles. We then consider beaded tree-child networks, which permit the representation of autopolyploidy events, and in our third main result we provide a greedy polynomial-time algorithm to decide whether a given profile can be realized by such a network. We expect that our results can be leveraged in future work through, for example, making use of known algorithms for computing short addition sequences to give bounds for the hybrid number, and in guiding network reconstruction for polyploid species.
Cancer mortality remains high in part because tumor-immune interactions can be unpredictable and may exhibit multistability, making malignant progression difficult to anticipate. We develop and analyze a tumor-immune model incorporating chemotherapy-induced toxicity and a Norton-Simon type tumor response, and show via deterministic analysis and continuation that increasing toxicity can generate hysteresis and bistability separating tumor dormancy from an uncontrolled full-growth state, so that small perturbations may precipitate abrupt progression. To capture uncertainty, we study environmental and demographic stochasticity and observe noise-induced switching in both cases; however, demographic noise sustains higher resilience of the tumor-dominant state under bistability, whereas stronger environmental noise tends to suppress uncontrolled tumor growth. We further assess early-warning signals using single and composite rolling-window indicators and find that they can provide advance warning of transitions from dormancy to full growth, with composite measures offering greater robustness across noise levels. Finally, we formulate an optimal control problem for combination immunotherapy and radiotherapy and demonstrate that appropriately timed treatment can substantially reduce tumor burden when initiated from either dormancy or full-growth conditions, highlighting how stochasticity-aware monitoring and optimized interventions may help prevent catastrophic tumor progression.
Cancer diagnostic methods based on Raman spectroscopy are being actively investigated, and there is a strong need for simple approaches to amplify the intensity of Raman scattered light from biological fluids such as serum and urine, which contain only trace amounts of nucleic acids, proteins, amino acids, and other analytes together with highly autofluorescent background components. We evaluated two measurement methods. One was the needle method (NM), in which a laser irradiates a droplet of liquid sample held at the tip of a fine-diameter stainless-steel needle. The other was the quartz glass fiber sheet method (QSM), in which a quartz glass fiber sheet is imbued with a liquid sample, allowed to dry, and then irradiated at the sheet surface. Raman spectra of sodium benzoate, sodium sulfate, human serum, and human urine were recorded. For the model compounds, spectra obtained by QSM reproduced the Raman shifts of the solid state, whereas spectra of aqueous solutions measured by NM showed clear peak shifts, and the scattered-light intensity increased monotonically with the number of drops on the sheet. Based on these findings, we infer that the samples crystallize and become concentrated within the quartz glass fiber sheet, enabling acquisition of spectra with high scattered-light intensity even from low-concentration solutions. For human serum and urine, QSM increased the intensity of characteristic bands by up to about seven-fold compared with NM while preserving the spectral fingerprints. Our results indicate that a quartz glass fiber sheet is a practical low-background substrate for obtaining FT-Raman spectra of liquid biological samples whose components are present at low concentrations.
Deterministic dynamical models are widely used in infectious disease modelling, but they often become systematically biased when key drivers such as seasonality and random environmental variation are simplified or omitted. This paper proposes a practical way to account for structural bias while preserving the underlying mechanistic model. We model the observed time series as the sum of (i) a deterministic transmission component given by a reduced Ross malaria model and (ii) a latent stochastic seasonal component that captures unresolved seasonal forcing and other unmodelled variability. The seasonal component is defined as a mean-reverting stochastic differential equation with periodic forcing, which can be interpreted as a seasonally forced Ornstein–Uhlenbeck (OU) process and, at the same time, as a dynamically constrained model-discrepancy term. The combined model forms an additive Bayesian state-space system. A key challenge in additive decompositions is identifiability: many combinations of the deterministic and stochastic components can explain the same observations. We therefore use informative priors to stabilise this decomposition, together with a non-centred parameterisation (a reparameterisation that improves MCMC efficiency by sampling standardised noise terms instead of states directly) of the latent stochastic differential equation (SDE) states that enables efficient joint inference with Hamiltonian Monte Carlo (NUTS) in Stan. We validate the approach in three steps: (1) an OU example that contrasts parameter-only inference with latent state-space inference, (2) synthetic experiments that have a good agreement with the observed signal while highlighting the expected negative posterior dependence between components, and (3) an application to five years of monthly malaria case reports from Delta State, Nigeria. In the real-data analysis, the proposed additive ODE–SDE model produces a close fit with coherent uncertainty quantification and a flexible seasonal reconstruction, while keeping the mechanistic transmission model interpretable. Overall, the framework provides a flexible and transferable method for accounting for structural model discrepancy in misspecified dynamical models using a structured stochastic discrepancy.
Mathematical models, e.g., differential equations and stochastic processes, have gained considerable attention for understanding evolution of antibiotic resistance. However, most existing models assume standing genetic variation and do not consider the possibility of random or drug-induced mutation of reference bacterial strains. Therefore, we propose a pharmacokinetics/pharmacodynamics (PK/PD)-based continuous-time Markov chain considering the competition and mutation between sensitive and resistant bacterial within an infected host during treatment. The proposed model is approximated as a generalized birth–death process with immigration, allowing for explicit derivation of the probability resistant population establishes during treatment. Besides capturing the stochasticity of de novo emergence of a resistant bacterial strain, we explore the effects of different antibiotic modes of action, horizontal gene transfer, nutrient availability and drug pharmacokinetics on antibiotic resistance. We find that replication-targeting (biostatic) drugs suppress resistance more than death-targeting (biocidal) drugs. Like prior works, we obtain maximized resistance at intermediate drug concentrations, however the consideration of de novo mutation magnifies the superiority of higher doses in preventing resistance emergence.