Microbiomepathogen interactions drive epidemiological dynamics of antibiotic resistance: A modeling study applied to nosocomial pathogen control
Abstract
The human microbiome can protect against colonization with pathogenic antibioticresistant bacteria (ARB), but its impacts on the spread of antibiotic resistance are poorly understood. We propose a mathematical modeling framework for ARB epidemiology formalizing withinhost ARBmicrobiome competition, and impacts of antibiotic consumption on microbiome function. Applied to the healthcare setting, we demonstrate a tradeoff whereby antibiotics simultaneously clear bacterial pathogens and increase host susceptibility to their colonization, and compare this framework with a traditional strainbased approach. At the population level, microbiome interactions drive ARB incidence, but not resistance rates, reflecting distinct epidemiological relevance of different forces of competition. Simulating a range of public health interventions (contact precautions, antibiotic stewardship, microbiome recovery therapy) and pathogens (Clostridioides difficile, methicillinresistant Staphylococcus aureus, multidrugresistant Enterobacteriaceae) highlights how speciesspecific withinhost ecological interactions drive intervention efficacy. We find limited impact of contact precautions for Enterobacteriaceae prevention, and a promising role for microbiometargeted interventions to limit ARB spread.
Introduction
Bacteria are fundamental drivers of human health and disease. On one hand, bacterial pathogens are leading causes of global infectious disease burden, with antibioticresistant and healthcareassociated infections posing significant risks to patient safety and public health (Cassini et al., 2019; Cassini et al., 2016; Friedman et al., 2016; Naylor et al., 2018). On the other, the bacterial microbiome – the trillions of individual bacteria that collectively inhabit the human body – provides support to development and homeostasis, facilitates core physiological processes like digestion, and protects against diseases ranging from colitis to cancer (Bäckhed et al., 2005; Kamada et al., 2013b; Lynch and Pedersen, 2016; Round and Mazmanian, 2009; Roy and Trinchieri, 2017). The microbiome can also protect against colonization with infectious bacterial pathogens, a phenomenon known as colonization resistance, limiting their capacities to establish colonies, grow, persist, and transmit (Bäumler and Sperandio, 2016; Buffie and Pamer, 2013). This is perhaps best exemplified by the canonical Clostridioides difficile: growth within the intestine is inhibited by secondary metabolites of commensal gut bacteria, protecting colonized hosts from disease and limiting propagation of the infectious spores that drive betweenhost transmission (Kamada et al., 2013a; Pamer, 2016). The human microbiome is purported to play a defensive role against a range of bacterial pathogens, but the mechanistic nature and epidemiological consequences of withinhost microbiomepathogen interactions are poorly understood.
This relationship between microbiome ecology and bacterial epidemiology is important in the context of widespread antibiotic use and global dissemination of highrisk pathogenic antibioticresistant bacteria (ARB). When prescribed appropriately, antibiotics target particular bacterial pathogens, but cocolonizing microbiota are also exposed (Tedijanto et al., 2018). This can unintentionally destabilize healthy microbial communities, resulting in dysbiosis, a state of population dynamic disequilibrium (Bhalodi et al., 2019; Coyte et al., 2015; Dethlefsen and Relman, 2011). Microbiome dysbiosis is associated with reduced abundance and diversity of commensal bacteria, impaired host immune responses, and loss of colonization resistance, altogether increasing host susceptibility to ARB colonization (Kim et al., 2017; Sorbara and Pamer, 2019; Zhang et al., 2015). Antibioticinduced dysbiosis may further result in elevated expression of antibiotic resistance genes, increased rates of horizontal transfer of such genes, and ecological release, whereby subdominant ARB are released from competition via clearance of drugsensitive bacteria, growing out into dominant colonies (Doan et al., 2019; Letten et al., 2021; Ruppé et al., 2019; Stecher et al., 2013). These associations may be particularly relevant for healthcare settings, where antibiotic use is widespread, and where microbiome dysbiosis is increasingly recognized as a key driver of ARB colonization and infection (Baggs et al., 2018; Prescott et al., 2015; Ravi et al., 2019). From a clinical perspective, this motivates a need for public health interventions that minimize or reverse harm to patient microbiota, from antibiotic stewardship, to fecal microbiota transplantation, to microbiome protective therapies. (de Gunzburg et al., 2018; Relman and Lipsitch, 2018).
In the face of uncertainty, mathematical modeling is a useful tool to analyze epidemiological dynamics of antibiotic resistance and evaluate control measures (Heesterbeek et al., 2015; Opatowski et al., 2011). However, an incomplete understanding of key ecoevolutionary principles is highlighted as a limitation to using modeling to predict future trends and inform decisionmaking (Birkegård et al., 2018; Knight et al., 2019). Disruption of the host microbiome is a longstanding theory explaining how antibiotics select for the spread of resistance at both the individual and population levels (Lipsitch and Samore, 2002), but most epidemiological models consider just one species of bacteria at a time, under the traditional assumption that antibiotic selection for resistance results from intraspecific competition between cocirculating strains (Blanquart, 2019; Ramsay et al., 2018; Spicknall et al., 2013). This simple framework has been particularly useful for bacteria like Streptococcus pneumoniae and Staphylococcus aureus, in which different strains – sometimes conceptualized as drugsensitive vs. drugresistant, or communityassociated vs. healthcareassociated – are believed to be in close ecological competition (Blanquart, 2019; Domenech de Cellès et al., 2011; KardaśSłoma et al., 2011; Pressley et al., 2010; van Kleef et al., 2013).
Accounting for other forms of complexity in epidemiological models – from treatment intensity, to ageassortative contact behavior, to hospital referral networks, to animalhuman interactions, to genetic linkage between resistance and nonresistance genes – has helped to unravel the many, disparate forces that contribute to drive the spread of resistance (Blanquart et al., 2018; Cobey et al., 2017; Colijn and Cohen, 2015; Donker et al., 2017; Lehtinen et al., 2017; van Bunnik and Woolhouse, 2017). Nevertheless, withinhost bacterial competition remains a key mechanism of selection for antibiotic resistance dissemination, and an active area of research at the forefront of resistance modeling (Lipsitch and Samore, 2002; Mulberry et al., 2020; Spicknall et al., 2013). For instance, the ‘mixedcarriage’ model by Davies et al. demonstrates how intraspecific competition results in negative frequencydependent selection for either of two competing strains, and provides a satisfying mechanistic explanation for widespread strain coexistence at the population level (Davies et al., 2019). However, contemporary work has stopped short of evaluating consequences of betweenspecies competition on resistance epidemiology. Yet for many ARB, including emerging highpriority multidrugresistant bacteria like extendedspectrum betalactamase (ESBL) producing Enterobacteriaceae, interactions with the host microbiome may be important mediators of nosocomial colonization dynamics (Kim et al., 2017; Lerminiaux and Cameron, 2019; Pilmis et al., 2020).
Here, we use mathematical modeling to evaluate how microbiome ecology and antibiotic consumption combine to drive the spread and control of antibioticresistant bacteria in the healthcare setting. This is presented in two parts. First, we propose a modeling framework for ARB colonization dynamics, accounting for different withinhost ecological interactions – including intraspecific pathogen strain competition, interspecific microbiomepathogen competition, and horizontal gene transfer (HGT) – in the context of antibiotic treatment. Synthesizing these into a final model, we show how different combinations of ecological interactions drive antibiotic selection for the spread of resistance, with heterogeneous impacts on classic epidemiological indicators. Second, using parameter estimates from the literature, we apply this framework to simulate the nosocomial epidemiology of four highrisk ARB: C. difficile, methicillinresistant S. aureus (MRSA), ESBLproducing Escherichia coli (ESBLEC) and carbapenemaseproducing Klebsiella pneumoniae (CPKP). Expert elicitation interviews were conducted to characterize the clinical relevance of microbiome dysbiosis for each species, and to qualify and quantify interaction coefficients with uncertainty. By simulating a range of different public health interventions, we demonstrate the theoretical importance of microbiomepathogen interactions as mediators of ARB epidemiology, and determining factors in the control of resistance dissemination.
Model and results
Part 1: A modeling framework for antibiotic resistance epidemiology in healthcare settings
We propose a series of five models describing colonization dynamics of an antibioticresistant bacterial pathogen, denoted P^{R}, among hospital inpatients in an acute care setting. Each model accounts for different withinhost ecological interactions between P^{R} and other bacteria. Models are described using systems of ordinary differential equations (ODEs), and are evaluated deterministically using numerical integration. Across models, three primary epidemiological outcomes are calculated at steadystate equilibrium: P^{R} prevalence (the proportion of patients colonized), P^{R} incidence (the daily rate of colonization acquisition within the hospital), and the pathogen resistance rate (the proportion of patients colonized with the focal antibioticresistant strain P^{R} relative to a competing drugsensitive strain P^{S}). We also derive and evaluate the basic reproduction number R_{0} for P^{R}, an indicator of pathogen epidemicity representing the average number of patients expected to acquire a novel pathogen from an initial index patient. See Materials and methods for technical details, and the supplementary appendix for the complete modeling framework and assumptions (Appendix Equation A1).
A simple transmission model for bacterial colonization
We start with a SusceptibleColonized transmission model (Figure 1A) representing a population of N hospital patients as either susceptible to colonization (S) or colonized (C^{R}) by P^{R}, the focal strain or species:
This model is adapted from classic colonization models of antibioticresistant bacteria (Austin et al., 1997), includes no ecological interactions with nonfocal bacteria, and reflects a suite of common assumptions relevant to the healthcare setting, including: (i) a symmetric rate of patient admission and discharge μ, holding N constant; (ii) a proportion of patients colonized upon admission f, reflecting pathogen prevalence in the community; (iii) a dynamic rate of colonization acquisition λ_{R}=β×C^{R}/N, for hosttohost transmission; (iv) a static rate of acquisition α_{R}, for endogenous routes of acquisition; (v) a rate of natural clearance γ_{R}; and (vi) a rate of effective antibiotic treatment σ_{R}. Note that all state variables are functions of time t, though this is omitted from ODEs for brevity.
Efficacy of antibiotic treatment is assumed to depend both on the distribution of antibiotics consumed in the hospital and on the intrinsic antibiotic resistance profile of P^{R}. We express this as
where $\mathrm{a}$ is the hospital population’s antibiotic exposure prevalence (the proportion of patients exposed to antibiotics at any t), θ_{C} is the antibioticinduced clearance rate (the rate at which effective antibiotics clear pathogen colonization), and r_{R} is the antibiotic resistance level (the proportion of antibiotics that are ineffective against P^{R}). Modeling the latter as a continuous proportion reflects that bacteria are not necessarily fully drugsensitive (r_{R} = 0) nor resistant (r_{R} = 1), but can range in their sensitivity to different antibiotics (0 ≤ r_{R} ≤ 1). The resistance level r_{R} is thus a model input interpreted as an overall measure of the pathogen’s innate degree of resistance to the particular antibiotics to which it is exposed.
Following models build upon these assumptions, representing the same pathogen P^{R}, but altering its ecological interactions with other bacteria from one model to the next. Models were evaluated over the same generic parameter space, to isolate impacts of model structure on epidemiological outcomes in the context of antibiotic use (see parameters in Appendix 1—table 1).
Antibiotic selection for resistance: the role of strain competition
Antibiotic consumption selects for the epidemiological spread of antibioticresistant bacteria (Chatterjee et al., 2018). To explain this mechanistically, a classic modeling assumption is that selection results from intraspecific competition between two or more drugsensitive strains P^{S} and drugresistant strains P^{R} (Spicknall et al., 2013). The reasoning goes: strains of the same species occupy the same ecological niche, so colonization with one strain inhibits colonization with another. In turn, antibiotics that preferentially clear P^{S} render the withinhost niche available to potential colonization with cocirculating P^{R}, indirectly favouring P^{R} spread through the host population. A simple twostrain ‘exclusive colonization’ model (Figure 1B) can be written as:
where patients can be colonized (C^{S}, C^{R}) by either strain (P^{S}, P^{R}). Subscripts S and R denote strainspecific rates, accounting for ecological differences between strains (e.g. antibiotic resistance levels, fitness costs of resistance; see Appendix Equations 3–7). Strains are labeled as sensitive or resistant, but here this is interpreted as relative (r_{S} < r_{R}). For simplicity, patient demography (admission and discharge) is given as Δ_{j} for each compartment j.
Unlike the singlestrain SusceptibleColonized model (Figure 1D and G), the strain competition model can explain how antibiotic consumption selects for the spread of resistance (Figures 1E,H). When P^{R} is resistant to all antibiotics (r_{R} = 1), its prevalence increases monotonically with increasing antibiotic exposure (Appendix 1—figure 1). When resistance is partial – when P^{R} is still cleared by antibiotics, but at a lower rate than P^{S} (r_{S} < r_{R} < 1) – selection for resistance can peak at intermediate antibiotic exposure, owing to a tradeoff in how antibiotics both clear and facilitate P^{R}. As such, R_{0} for P^{R} tends to increase with antibiotic use when P^{R} is highly resistant to antibiotics (high r_{R}), but decrease when still largely sensitive (low r_{R}).
Most contemporary antibiotic resistance models are variants of this typology (Blanquart, 2019; Spicknall et al., 2013), but intraspecific strain competition is not the only mechanism by which antibiotic consumption can drive ARB spread (Lipsitch and Samore, 2002), and may have limited relevance for certain species, settings, timescales and epidemiological indicators.
Antibiotic selection for resistance: the role of microbiome competition
Antibioticinduced microbiome dysbiosis is an alternative mechanism by which antibiotics may select for the spread of resistance. We propose a model in which (i) bacterial pathogens and commensal microbiota compete ecologically within the host, and (ii) antibioticinduced dysbiosis disrupts these interactions, predisposing hosts to pathogen colonization. We consider three competitive withinhost microbiomepathogen interactions and conceptualize how they affect pathogen epidemiology (Figure 2; Buffie et al., 2015; Hooper et al., 2003; Kamada et al., 2013a). First, a stable microbiome can act as a barrier preventing introduced pathogens from establishing colonies. We modeled this colonization resistance (ε) as a reduced rate of pathogen transmission (β) to hosts with stable microbiota, β_{ε} = (1 – ε) × β. Second, cocolonizing bacteria can compete for space, nutrients and other limited resources within the host. In this case, antibiotics that reduce the density of potential competitors may favour pathogen persistence. We modeled this resource competition (η) as a reduced rate of pathogen clearance (γ) among hosts undergoing dysbiosis, γ_{η} = (1 – η) × γ. Lastly, microbiome dysbiosis can favor the emergence or outgrowth of subdominant pathogen colonies, and this ecological release (ϕ) was modeled as an increased rate of endogenous pathogen acquisition (α) upon dysbiosis, α_{ϕ} = ϕ × α.
We integrate the microbiome, these three interactions, and antibioticinduced microbiome dysbiosis into Equation 1. The resulting ‘microbiome competition model’ (Figure 1C) is given by
describing colonization (C^{R}) with a single pathogen strain (P^{R}) across two host types: patients with a microbiome at dynamic equilibrium (subscript e) and those undergoing dysbiosis (subscript d). Antibiotics induce dysbiosis at a rate σ_{m} given by
such that microbiome dysbiosis depends on both antibiotic exposure prevalence ($\mathrm{a}$) and the rate at which antibiotic exposure causes dysbiosis (θ_{m}). Accordingly, the same level of antibiotic exposure ($\mathrm{a}$) can have asymmetric effects on microbiome stability (via θ_{m}) and P^{R} colonization (via (1 – r_{R}) × θ_{c}, as in Equation 2). Dysbiosis can result in longterm changes to microbiome composition, but ecological function and population dynamic stability tend to recover in the days or weeks following antibiotic therapy (Lozupone et al., 2012), represented here by microbiome recovery rate δ.
Included microbiomepathogen interactions predict in different ways how antibiotic consumption selects for the epidemiological spread of antibioticresistant bacterial pathogens (Figure 1F). Like strain competition, microbiome competition underlies a tradeoff in antibiotic selection: P^{R} prevalence and R_{0} increase monotonically with antibiotic exposure when P^{R} is completely antibioticresistant (r_{R} = 1) (Appendix 1—figure 1), but can peak at intermediate antibiotic exposure when resistance is partial (0 < r_{R} < 1) because antibiotics simultaneously clear pathogen colonization and induce greater host susceptibility through dysbiosis (Figure 1I). This tradeoff is more pronounced when multiple interactions are combined: when microbiota simultaneously limit multiple colonization processes (transmission, emergence, persistence), patients with stable microbiota are more protected from P^{R} colonization, but antibiotic use is predicted to have greater epidemiological costs, selecting more strongly for P^{R} spread (Appendix 1—figure 2).
Antibiotic selection for resistance: combining strain and microbiome competition
Strain competition and microbiome competition are not mutually exclusive: when combined in a twostrain ‘microbiomestrain competition’ model, different strains of the same pathogen species compete for hosts whose microbiota and history of antibiotic consumption influence susceptibility to colonization. We assume that microbiomepathogen interactions are species and not strainspecific, that is that ε, η and ϕ apply equally to P^{S} and P^{R}. This is given by:
Introducing a drugsensitive strain P^{S} to the microbiome model dampens R_{0} of the focal strain P^{R}, because fewer patients are susceptible to colonization when the competing strain cocirculates in the population (Appendix 1—figure 3). However, antibiotic use makes way for P^{R} not only through preferential clearance of P^{S}, but remaining P^{R} also benefit from increased host susceptibility to colonization when antibiotics cause dysbiosis. Accordingly, antibiotics that both disrupt host microbiota and clear drugsensitive pathogen strains tend to select more strongly for the spread of resistant strains than antibiotics that only target one or the other (Figure 3). Overall, antibiotics with stronger impacts on pathogen clearance (higher θ_{C}) lead to increased resistance rates, while consequences for P^{R} colonization prevalence are modest and depend on potential interactions with the microbiome. Conversely, antibiotics with stronger impacts on microbiome stability (higher θ_{m}) lead to higher prevalence, with only modest effects on resistance rates. Further, impacts of antibiotic treatment on resistance rates are greater when P^{R} resists a greater share of antibiotics (higher r_{R}), while impacts on prevalence are greater when microbiomepathogen interactions are stronger (higher ε, η, ϕ) (Appendix 1—figure 4). Different microbiomepathogen interactions also underlie distinct dynamic responses to theoretical public health interventions. For the same generic pathogen P^{R}, antibiotic stewardship interventions generally, but do not always prevent colonization, with predictions depending on the ecological interactions in effect (e.g. strain competition, microbiome competition), the impact of the intervention (e.g. reduced microbiome disruption, reduced overall prescribing), and the epidemiological outcome considered (e.g. colonization prevalence, resistance rate) (Appendix 1—figure 5).
Antibiotic selection for resistance: introducing interspecific horizontal gene transfer
Horizontal transfer of resistanceencoding genes may also contribute to hospital resistance dynamics (Evans et al., 2020; Lerminiaux and Cameron, 2019). We extend our model to include interspecific HGT, conceptualized as twoway withinhost transfer of a focal resistance gene, either from a resistant pathogen strain P^{R} to cocolonizing microbiota, or from resistancebearing microbiota to a cocolonizing drugsensitive pathogen strain P^{S}. (See final model and assumptions, Appendix Equations A1–A11). Overall, an increasing rate of withinhost HGT (χ) is found to drive increasing P^{R} prevalence (C^{R}), regardless of other microbiomepathogen interactions (Figure 4). Under our modeling assumptions, HGTdriven gains in C^{R} result from symmetric declines in C^{S}, such that HGT’s potential impact depends on the presence of both sufficient resistance donors and sufficient recipients. This results in nonlinear feedbacks between HGT and other processes that drive strain prevalence (Appendix 1—figure 7). For instance, impacts of HGT on C^{R} are greatest at intermediate antibiotic exposure ($\mathrm{a}$). This is linked to another antibiotic selection tradeoff: higher $\mathrm{a}$ affords a selective advantage to P^{R} relative to P^{S}, increasing the potential impact of HGT on C^{R}; but higher $\mathrm{a}$ also reduces the pool of recipient C^{S}, ultimately limiting HGT’s ability to drive C^{R}.
Part 2: Model application to highrisk nosocomial pathogens
We applied this modeling framework to simulate colonization dynamics of four nosocomial pathogens in the hospital setting: C. difficile, methicillinresistant S. aureus (MRSA), extendedspectrum betalactamaseproducing E. coli (ESBLEC), and carbapenemaseproducing K. pneumoniae (CPKP). Data from the literature were used to parameterize the model to each pathogen (Appendix 1—figure 5, Appendix 1—table 2–6). Literature estimates for microbiomepathogen interaction coefficients are scarce, and the speciesspecific relevance of different withinhost interactions (intraspecific strain competition, microbiome competition, HGT) in the hospital environment are not welldefined. To characterize ecological interactions for each species and inform model structure, we conducted interviews with a panel of subjectmatter experts in medical microbiology and antibiotic resistance epidemiology (details in Materials and methods). Based on their beliefs, all pathogens were assumed to compete with microbiota; MRSA, ESBLEC and CPKP were further assumed to compete intraspecifically with nonfocal strains, for simplicity characterized as methicillinsensitive S. aureus (MSSA), E. coli (EC), and K. pneumoniae (KP); and both ESBL resistance and carbapenem resistance were assumed to be borne by plasmids capable of horizontal transfer between patient microbiota and, respectively, EC and KP (Figure 5A). To quantify speciesspecific strengths of microbiomepathogen interactions, withinhost ecological parameters were translated into clinical parameters (Appendix 1—table 7), and experts were asked to quantify these using standardized expert elicitation methodology (Figure 5B, Appendix 1—figures 8–10).
Colonization dynamics for each ARB were modeled against a common backdrop of antibiotic consumption. This was parameterized at the level of antibiotic class, using national data from French hospitals in 2016 (Appendix 1—table 8; Agence nationale de sécurité du médicament et des produits de santé, 2017). Using data from the literature, antibiotic classes were assumed to vary in their impact on microbiome dysbiosis (very low, low, medium, or high rate of inducing dysbiosis) and on each pathogen strain (classified as sensitive, intermediate, or resistant) (Figure 5C and D; Baggs et al., 2018; Brown et al., 2013; McCormack and Lalji, 2015). For each of 10,000 Monte Carlo simulations, in which parameter values were sampled randomly from their respective probability distributions, epidemiological outcomes were evaluated at population dynamic equilibrium. For all outcomes, we compare results from simulations that include microbiomepathogen interactions and dysbiosis (‘microbiome simulations’) and those that exclude them (‘singlespecies simulations’). Multivariate sensitivity analyses were also conducted, using partial rank correlation coefficients (PRCCs) to evaluate impacts of parameter uncertainty on model outcomes (details in Methods).
Speciesspecific hospital colonization dynamics
Across simulations, C. difficile was the most prevalent pathogen (Figure 6A), MRSA had the highest resistance rate (Figure 6B), and ESBLEC had the highest rate of incidence within the hospital (Figure 6C). CPKP had the lowest prevalence, resistance rate and incidence, but was the pathogen most favoured by the hospital environment, its prevalence increasing by approximately 5.4fold (95% uncertainty interval: 2.1–10.9) among hospital patients relative to baseline prevalence in the community (Appendix 1—figure 12A). For pathogens subject to intraspecific strain competition, resistance rates in the hospital also tended to exceed rates in the community (Appendix 1—figure 12B). Patienttopatient transmission was the primary route of MRSA acquisition, while endogenous acquisition was the primary route for the enteric pathogens C. difficile, ESBLEC and CPKP (Appendix 1—figure 13). HGT played a potentially important but highly uncertain role for ESBLEC and CPKP, accounting for 8.7% (<0.01–49.7%) and 2.1% (<0.01–22.8%) of acquisition events, respectively. In multivariate sensitivity analysis, community prevalence (f_{C} for C. difficile, f_{R} for others) and rates of endogenous acquisition (α_{R}) had overall the strongest positive impacts (highest PRCCs) on hospital colonization prevalence across ARB, while rates of hospital admission/discharge (μ) and microbiome recovery (δ) had the strongest negative impacts (lowest PRCCs; Appendix 1—figure 14A). For resistance rates, parameters with the strongest positive impacts (highest PRCCs) were community prevalence (f_{R}), the rate of endogenous acquisition (α_{R}) and the rate of antibioticinduced pathogen clearance (θ_{C}). Parameters with the strongest negative impacts (lowest PRCCs) were rates of hospital admission/discharge (μ) and endogenous acquisition of competing drugsensitive strains (α_{S}) (Appendix 1—figure 14B). Across ARB, prevalence estimates, but not resistance rates, were generally sensitive to microbiome parameters.
Microbiome ecology underlies epidemiological responses to public health interventions
We simulated three types of public health intervention, each across three levels of intervention compliance: (i) contact precautions, which reduced transmission rates of all strains by 20%, 35%, or 50%; (ii) antibiotic stewardship interventions, which reduced or modified hospital antibiotic consumption patterns by 20%, 35%, or 50%; and (iii) a theoretical ‘microbiome recovery therapy’ intervention that facilitates recovery from dysbiosis. We assumed a mean 2 day delay to microbiome recovery, and compliance levels corresponding to use among 10%, 30%, or 50% of antibioticexposed patients (see Materials and methods). Intervention efficacy was evaluated using: reduction in colonization incidence, 1 – IRR (where IRR is the incidence rate ratio of postintervention to preintervention incidence); and reduction in the resistance rate, 1 – RRR (where RRR is the resistance rate ratio, the ratio of the postintervention to preintervention resistance rate).
Efficacies of simulated interventions varied considerably by pathogen and type of intervention (Figure 7). Contact precautions were highly effective for reducing MRSA incidence, of intermediate efficacy for C. difficile, and minimally effective for ESBLEC and CPKP. Contact precautions had comparatively little impact on resistance rates, with a median 2–4% reduction across simulations and compliance levels for MRSA, 0–2% for CPKP and negligible impact for ESBLEC (Appendix 1—figure 15). These interventions were overall less effective in microbiome simulations, which tended to limit the role of betweenhost transmission (via colonization resistance) and favor the role of endogenous acquisition (via ecological release) in the hospital environment when compared to singlespecies simulations (Appendix 1—figure 13).
Antibiotic stewardship interventions led to substantial reductions in nosocomial incidence for all pathogens, but only when the microbiome was taken into account (Figure 7). Unlike contact precautions, which only reduced incidence via transmission, stewardship reduced incidence through all acquisition routes, including HGT (Appendix 1—figure 16). Overall efficacy estimates and speciesspecific responses were similar across three types of stewardship considered (Appendix 1—figure 17). Pooling these together under intermediate compliance, colonization incidence was reduced by a median 20% for CPKP, 18% for C. difficile, 15% for ESBLEC, and 10% for MRSA. Singlespecies simulations excluding microbiome competition predicted negligible efficacy of all stewardship interventions for reducing incidence, and nonefficacy for C. difficile. Stewardship interventions also had a substantial impact on resistance rates, with overall greater reductions for CPKP than MRSA and ESBLEC, and similar outcomes across microbiome and singlespecies simulations (Appendix 1—figure 15).
Lastly, microbiome recovery therapy was potentially highly effective for limiting pathogen incidence, but efficacy varied greatly across different levels of intervention compliance (Figure 7). This intervention was most effective against C. difficile, of similar efficacy against ESBLEC and CPKP, and comparatively least effective against MRSA. Across pathogens, microbiome recovery therapy reduced incidence through all acquisition routes,(Appendix 1—figure 16) but had no clear impact on resistance rates (Appendix 1—figure 15).
Discussion
Antibiotics are essential medicines for the treatment and prevention of bacterial infections, but their use selects for the spread of pathogenic antibioticresistant bacteria (ARB), and can inadvertently disrupt the host microbiome and its associated immune function (Buffie and Pamer, 2013; Chatterjee et al., 2018; Kim et al., 2017). Withinhost ecological interactions between cocolonizing bacteria can have important consequences for their colonization dynamics, which likely extend to influence the epidemiology of human pathogens in clinical settings. Yet microbiome ecology remains largely absent from the epidemiological modeling of antibiotic resistance (Assab et al., 2017; Birkegård et al., 2018; Blanquart, 2019; Niewiadomska et al., 2019), suggesting a need to better understand withinhost competition between ARB and the host microbiome, its potential epidemiological consequences, and more broadly how antibiotics exert selection pressure on resistant bacteria (Knight et al., 2019). We present a modeling framework that includes withinhost ecological costs of antibiotic use in the form of microbiome dysbiosis, incorporating a leading hypothesis for antibiotic selection into classical models of resistance epidemiology (Austin et al., 1997; Lipsitch and Samore, 2002). We formalize three examples of microbiomepathogen competition potentially affected by dysbiosis, and show how they, either separately or in combination with other forces of selection, help explain how antibiotic use drives the spread of ARB in healthcare settings. We use probabilistic simulations to apply this framework across a panel of characteristic nosocomial pathogens and compare results to a traditional strainbased framework, demonstrating utility of a microbiomeoriented approach for modeling antibiotic resistance epidemiology and associated interventions.
MRSA, C. difficile and ESBLproducing Enterobacteriaceae are leading causes of antibioticresistant and healthcareassociated infection, while carbapenemaseproducing Enterobacteriaceae represent emerging threats of particular concern due to limited therapeutic options for effective treatment of invasive infection (Cassini et al., 2019; Jernigan et al., 2020; Miller et al., 2011; RodríguezBaño et al., 2018). Antibiotic stewardship is a core component of public health efforts to limit the emergence and spread of these ARB in clinical settings (Baur et al., 2017), and an important focus of antibiotic resistance modeling (Niewiadomska et al., 2019). Different antibiotics vary in which organisms they intrinsically target, as well as pharmacodynamic factors like route of administration (e.g. oral vs. systemic), clearance mechanism (e.g. biliary vs. renal excretion) and site of absorption (e.g. small vs. large intestine). These differences affect the degree to which bacteria in different host niches are exposed to and cleared by different antibiotics, and are increasingly well described (Bhalodi et al., 2019; Kim et al., 2017; Zhang et al., 2015). Our findings suggest that asymmetric impacts of different antibiotics on competing commensal and pathogenic bacteria can drive antibioticdriven selection for these highrisk ARB, with important consequences for resistance dynamics and stewardship efficacy.
Findings also suggest promise for interventions that effectively restore microbiome stability and associated colonization resistance as a means to control ARB spread. Fecal microbiota transplantation is already used to treat recurrent C. difficile infection, and is under investigation for multidrugresistant Enterobacteriaceae decolonization (Davido et al., 2019a; Kassam et al., 2013; Saha et al., 2019). However, its appropriateness for dysbiosis recovery in the absence of other clinical indications is unclear. Transplantation requires rigorous donor screening and close longitudinal followup, and cases of donor stool contaminated with toxicogenic and multidrugresistant bacteria highlight nonnegligible risks (Gupta et al., 2021; Zellmer et al., 2021). Alternative microbiome protective therapies now exist, like DAV132, a novel activatedcharcoal product currently undergoing clinical trials. When coadministered with antibiotics by the oral route, DAV132 has been shown to absorb antibiotic residues in the colon and preserve the richness and composition of intestinal microbiota, while maintaining systemic antibiotic exposure (de Gunzburg et al., 2018; de Gunzburg et al., 2015; Pinquier et al., 2021). Modeling has been used previously to evaluate impacts of such microbiomeoriented interventions at the withinhost level (Guittar et al., 2021; Guk et al., 2021), but knockon impacts on ARB transmission dynamics and epidemiological burden are unknown. Our simulations were limited to a select few ARB, but our framework and findings likely have relevance for other bacteria known to interact with the microbiome, including vancomycinresistant Enterococci and other multidrugresistant Enterobacteriaceae (Davido et al., 2019b; Stecher et al., 2013). and could be further extended to explore impacts of the microbiome on resistance dynamics and intervention efficacy beyond healthcare settings.
Our simulations predicted colonization dynamics broadly consistent with previous findings from the literature. Input from the community was the main driver of hospital prevalence, and both prevalence and resistance rates tended to increase in the hospital relative to the community, as estimated in previous modeling studies (Knight et al., 2018; MacFadden et al., 2019). Findings reflected an important role for betweenhost nosocomial transmission for MRSA, as observed clinically (Khader et al., 2019; Nadimpalli et al., 2020), and a comparatively important role for endogenous acquisition for enteric ARB, as estimated elsewhere (Bootsma et al., 2007; Gurieva et al., 2018). Our estimates of intervention efficacy were more consistent with previous findings from the literature when microbiome interactions were taken into account. Contact precautions were effective for reducing incidence of MRSA and to a lesser extent C. difficile, but had limited impact against Enterobacteriaceae, consistent with clinical trials and modeling estimates. (Khader et al., 2021; Luangasanatip et al., 2015; Maechler et al., 2020). By contrast, simulated antibiotic stewardship interventions were broadly effective for reducing incidence across all included ARB. This is consistent with findings from a metaanalysis of clinical trials evaluating the efficacy of hospital antibiotic stewardship interventions for reducing incidence of ARB colonization and infection (Baur et al., 2017), which we updated to exclude studies coimplementing stewardship with alternative interventions (Appendix 1—figure 18). In comparison to our findings, estimates from the metaanalysis were associated with greater uncertainty across more heterogeneous interventions, but predicted the same rank order of efficacy across included ARB, and similar mean efficacy for C. difficile (19% from n=seven studies, vs. 18% under intermediate compliance in our simulations), ESBLEC (18% from n=four studies, vs. 15%) and MRSA (12% from n=10 studies, vs. 10%), but higher efficacy for CPKP (54% from n=one study, vs. 20%).
When excluding microbiome interactions, simulations predicted negligible efficacy of antibiotic stewardship interventions for controlling ARB incidence. Previous models have predicted efficacy using predominantly strainbased approaches, but often focus on resistance rates as the primary outcome, and in many cases assume that patienttopatient transmission is the only route of colonization acquisition (Niewiadomska et al., 2019). Here, microbiome competition was found to have a large impact on incidence but comparatively little impact on resistance rates – both for the theoretical pathogen evaluated in Part 1 (Figures 3 and 4) and for the four ARB simulated in Part 2 (Appendix 1—figures 12 and 14) – underlying why stewardship interventions were of similar efficacy for reducing resistance rates across singlespecies and microbiome simulations (Appendix 1—figure 15). This reflects the importance of different forces of antibiotic selection for different epidemiological outcomes: while strainbased competition explains relative ecological dynamics of cocirculating strains, microbiome interactions may be better suited to explain how increasing antibiotic use favours ARB incidence, and how antibiotic stewardship and microbiome recovery interventions can help prevent colonization acquisition.
Under strain competition, pathogen colonization is limited by closely related strains sharing an ecological niche, with the same approximate epidemiological profile and transmission characteristics. Competition against ARB is assumed to depend on the epidemic spread of competing strains, and removal of drugsensitive strains releases antibiotic selection for resistance. By contrast, microbiome population structure is inherently stable, depending less on the epidemiological transmission of particular taxa, and more on host factors like diet, maternal inheritance, genetics, and antibiotic exposure (Brito and Alm, 2016). Despite great interindividual diversity in microbiome composition, there is functional redundancy from one host to the next, such that colonization resistance and other forms of microbiome competition are shared across healthy hosts colonized with different taxa (Bashan et al., 2016; Kim et al., 2017). For these reasons, microbiome stability is modeled as a host trait reflecting the functional ecology of the microbiome in different population dynamic states, as opposed to a more traditional bacterial colonization process governed by rates of acquisition and clearance. This is clearly an oversimplification of real microbiome dynamics and complexity (Hooks and O’Malley, 2017), but is a useful approximation for the needs of epidemiological modeling, particularly in the absence of data, and reflects the universality of both human microbiome function and of the ecological impacts that antibiotics have on microbiome stability (Bashan et al., 2016).
Previous studies have modeled antibiotic exposure as a risk factor for ARB colonization, which can be interpreted as indirectly accounting for microbiome dysbiosis. In transmission models of C. difficile, ESBLEC and Pseudomonas aeruginosa, among others, patients undergoing antibiotic therapy are assumed to be at greater risk of colonization and/or infection (Gingras et al., 2016; Hurford et al., 2012; MacFadden et al., 2019). An alternative approach has been to use antibiotic exposure as a coefficient on epidemiological parameters (e.g. transmission, endogenous acquisition), allowing ARB colonization rates to scale with antibiotic use (Knight et al., 2018). These strategies reflect widespread recognition that antibiotic use favours ARB acquisition, through erosion of colonization resistance or other supposed mechanisms, and independent of potential competition with other strains. The present work formalizes examples of the microbiomepathogen interactions that underlie these assumptions, demonstrating their relevance to various epidemiological outcomes, distinguishing them from strainbased selection, and providing a framework for their application.
This study has a number of limitations. First, hospitals and healthcare settings are heterogeneous environments with nonrandom contact patterns and relatively small population sizes. Stochastic, individualbased models accounting for these factors reproduce more realistic nosocomial transmission dynamics than deterministic ODE simulations, allowing for local extinction events, superspreaders, and other inherently random epidemiological phenomena. Nonetheless, our goal was to study how ecological mechanisms impact average epidemiological outcomes in the context of different model assumptions and parameter uncertainty, and in this context, ODE modeling was the more appropriate tool, particularly for widely endemic ARB like C. difficile, MRSA and ESBLEC. Still, further insights could certainly be gained by accounting for additional complexity and stochastic heterogeneity in future work, from withinhost spatial organization (Estrela et al., 2015), to patient and staff contact behavior (Duval et al., 2019), to interinstitutional or interward metapopulation dynamics (Shapiro et al., 2020). These distinctions may be particularly important for rare or nonendemic ARB (e.g. CPKP in some regions). Second, our evaluation of strain competition was limited to exclusive colonization, a widely used approach (an estimated 12% of published strain competition models allow cocolonization or infection) (Niewiadomska et al., 2019). Yet alternative models predict unique impacts on resistance dynamics (Spicknall et al., 2013), and explicit consideration of higher resolution withinhost population dynamics has been shown to better reproduce empirical findings in previous work (Davies et al., 2019). Further, our exclusive colonization approach precluded assessment of intraspecific HGT, which may have different impacts on resistance dynamics than interspecific HGT.
Third, Monte Carlo simulations were limited by the availability of speciesspecific model parameters from the literature, in some instances necessitating use of previous modeling results, approximations, or estimates from small studies in specific locations, making the generalizability of results unclear. For instance, Khader et al. estimated a fourfold difference in MRSA transmission rates between hospitals and nursing homes (Khader et al., 2019). Such differences could have a substantial impact on dynamics and estimated intervention efficacy, with higher transmission rates favouring use of contact precautions, and higher rates of endogenous acquisition favouring antibiotic stewardship (in the context of a high ecological release coefficient). Uncertainty in endogenous acquisition rates may be particularly important: in multivariate sensitivity analyses, this parameter emerged as a key driver of both colonization prevalence and resistance rates across ARB (Appendix 1—figure 14). Further, data used to estimate classspecific rates of microbiome dysbiosis are specific to the gut; classspecific data for dysbiosis of the skin, the preferred niche of S. aureus, were not available. This may overestimate impacts of dysbiosis on MRSA colonization dynamics.
Finally, the nature of microbiomepathogen interactions and their epidemiological consequences remain poorly understood and largely unquantified. We show in theory why these interactions matter at the population level, but empirical data were unavailable to inform their parameterization. Instead, we translated microbiomepathogen competition coefficients into clinical parameters, designed a structured expert elicitation protocol, and conducted interviews allowing subjectmatter experts to quantify their beliefs. Although these estimates are subject to substantial bias and uncertainty, they facilitated speciesspecific characterization of the epidemiological impact of microbiome dysbiosis, and represent useful proxy measures in the absence of clinical data. More broadly, our characterizations of microbiomepathogen interactions are conceptual, and were mapped mechanistically to particular colonization processes (transmission, clearance, endogenous acquisition), but we note that in other contexts terms like colonization resistance, resource competition and ecological release may map to specific biochemical processes that could affect epidemiological parameters in different ways.
Despite data limitations, epidemiological conclusions from Monte Carlo simulations were largely consistent with empirical findings (discussed above), suggesting that final parameter distributions were reasonable approximations. Uncertainty in parameter inputs translated to uncertainty in model outputs, reflecting the knowledge gaps underlying our simulations (Appendix 1—figure 14). This is exemplified by HGT and its highly uncertain role in driving colonization incidence; to date, HGT modeling has largely been limited to withinhost dynamics, and impacts on epidemiological dynamics are only just beginning to come to light (Evans et al., 2020; Leclerc et al., 2019; Lerminiaux and Cameron, 2019). Increasing availability and synthesis of highquality withinhost microbiological data will help to further characterize epidemiological impacts of microbiomepathogen interactions. Studies are needed that describe ecological impacts of antibiotic exposure on microbiome population structure across control and treatment groups, with longitudinal followup evaluating subsequent nosocomial ARB colonization risk. In the absence of clinical data, insights from experiments and withinhost models nonetheless suggest that antibiotic disruption of microbiomepathogen competition is a key driver of selection for resistance (Baumgartner et al., 2020; Estrela and Brown, 2018; O’Brien et al., 2021; Shaw et al., 2019; Stein et al., 2013; Tepekule et al., 2019). The present work highlights the importance of extending these withinhost concepts to the population level. Links between within and betweenhost dynamics have been widely studied in various contexts, including the theory underlying the evolution of parasite life history and antimicrobial resistance (Day et al., 2011; Greischar et al., 2019; zur Wiesch et al., 2011), but remain largely absent from clinical models of antibiotic resistance (Blanquart, 2019). A clear extension of the present work is to explicitly account for simultaneous within and betweenhost processes using nested models (Birkegård et al., 2018; Niewiadomska et al., 2019). However, such models do not necessarily provide more epidemiological clarity, especially when data are lacking and simple heuristic parameters can capture epidemiological consequences of withinhost processes (Gog et al., 2015; Mideo et al., 2008).
In conclusion, we have proposed a mathematical modeling framework for the epidemiology of antibioticresistant bacteria that accounts for their potential interactions with the host microbiome. This model simplifies into accessible epidemiological parameters what are in reality highly complex ecological systems, comprising a staggering diversity of microbes and interactions among them. We demonstrate that accounting for at least some of this ecological complexity may help to explain how antibiotics select for the epidemiological spread of resistance, how antibiotic stewardship works to reduce pathogen colonization incidence, and how interventions favouring healthy microbiome function may help to mitigate the epidemiological burden of antibiotic resistance.
Materials and methods
Mathematical models of bacterial colonization
Request a detailed protocolWe evaluated ODE systems describing colonization dynamics of bacterial pathogens in the healthcare setting. The final model, comprising pathogen colonization (Equation 1), intraspecific strain competition (Equation 3), microbiomepathogen competition (Equation 4), and horizontal gene transfer (HGT), is given alongside all assumptions in the supplementary appendix (Appendix Equation A1). ODEs were integrated numerically to calculate steadystate epidemiological outcomes for nosocomial P^{R} colonization: colonization prevalence (the sum of all compartments C^{R}), colonization incidence (the daily rate of C^{R} acquisition), and the resistance rate (C^{R}/(C^{S} + C^{R})). (See Appendix 1 for technical details, and R and Mathematica files available online at https://github.com/drmsmith/microbiomeR [Smith, 2021; copy archived at swh:1:rev:a3682a24970d79e4f748952ecc49fcdb16adf48f].) For each model, outcomes were evaluated over the same parameter space representing a generic pathogen P^{R} (parameters in Appendix 1—table 1), while varying specific parameters through univariate and bivariate analysis to assess their impacts on dynamic equilibria in the context of different modeling assumptions. We focused on impacts of the patient population’s antibiotic exposure prevalence ($\mathrm{a}$), rates of antibioticinduced pathogen clearance (θ_{c}) and microbiome dysbiosis (θ_{m}), the focal pathogen’s intrinsic antibiotic resistance level (r_{R}), and mediating impacts of microbiomepathogen competition (ε, η, ϕ) and horizontal gene transfer (χ_{e}, χ_{d}).
Monte Carlo simulations over parameter distributions
Request a detailed protocolWe applied the final model to simulate epidemiological dynamics of four highrisk nosocomial pathogens, varying the withinhost ecological interactions in effect for each (Figure 5A). For MRSA, ESBLEC and CPKP, the focal pathogen was taken as the ‘drugresistant’ strain P^{R}, while the ‘drugsensitive’ strain P^{S} was taken to represent all other cocirculating strains of the same species. We simulated colonization dynamics using these model characterizations and accounted for parameter uncertainty using Monte Carlo methods. Specifically, 10,000 unique parameter vectors Ω were created by drawing random values for each parameter from its respective probability distribution (parameters in Appendix 1—table 2–6). For each Ω, epidemiological outcomes (prevalence, incidence, resistance rate) were calculated as above. Final outcome distributions were reported as the median and 95% uncertainty interval, that is the 50th (2.5th–97.5th) percentiles across all simulations.
Parameterizing models for application to specific nosocomial pathogens
Request a detailed protocolModels were parameterized using estimates from the literature, prioritizing clinical studies from the French hospital setting where available. We used expert elicitation to inform model structure and quantify parameter values for microbiomepathogen interactions. This involved development of a protocol and questionnaire (provided separately), designed using established expert elicitation methodologies for quantitative estimation of unknown parameter values (Johnson et al., 2010). To facilitate more reliable parameter interpretation, microbiomepathogen coefficients were translated into clinical parameters, for example relative risks of pathogen colonization among hospital patients undergoing dysbiosis relative to those with stable microbiota (Appendix 1—table 7). Expert estimates were quantified with uncertainty using the MATCH Uncertainty Elicitation Tool (Morris et al., 2014), and final parameter distributions were generated by pooling their individual distributions while maintaining estimated species rank order (Appendix 1—figure 8–10).
Characterizing ecological impacts of antibiotic consumption
Request a detailed protocolAntibiotic exposure prevalence was quantified at the level of antibiotic class using nationally representative antibiotic consumption data from French hospitals in 2016 (Agence nationale de sécurité du médicament et des produits de santé, 2017). A study from American hospitals in 2006–2012 was used to further stratify antibiotics with ATC code J01F into J01FA and J01FF, J01X into J01XA and J01XD, and J01DD+DE into J01DD and J01DE (Baggs et al., 2016). Final antibiotic consumption data are presented in Appendix 1—table 8. To simulate classspecific rates of microbiome dysbiosis, we used a fourpoint loglinear scale of intestinal microbiome disruption from Brown et al., 2013. For interpretation, we present classes as inducing dysbiosis at a high, medium, low, or very low rate. This scale was supplemented with data for additional antibiotic classes from Baggs et al., 2018. Applied to our model, the rate that antibiotic treatment induces dysbiosis (σ_{m}) is given by
where a_{k} is the proportion of antibiotics consumed of each group k, and where the most ecologically disruptive group (k=0) causes dysbiosis a mean 12 hr after antibiotic exposure (θ_{m} = 2 day^{−1}), with classes in less disruptive groups (k=1,2,3) causing dysbiosis at successively slower rates (Figure 5C).
To characterize classspecific effects on pathogen clearance, characteristic antibiograms for all strains were adapted from an online compendium from the Therapeutics Education Collaboration (Figure 5D; McCormack and Lalji, 2015). Each strain i was classified as sensitive (r_{i,j} = 0), of intermediate sensitivity (0 < r_{i,j} < 1), or resistant (r_{i,j} = 1) to each antibiotic class j. Overall, the rate that antibiotic treatment clears each strain is given by:
across the included classes. Under these assumptions, MRSA was resistant to the greatest proportion of antibiotics consumed in hospital (median r_{R} = 94.5%), followed by C. difficile (94.3%), CPKP (91.7%) and ESBLEC (84.8%); competing ‘drugsensitive’ strains bore considerably less resistance median r_{S} = 33.0% for MSSA and 23.0% for both E. coli and K. pneumoniae.
Sensitivity analyses
Request a detailed protocolTwo distinct sensitivity analyses were conducted. First, to evaluate the impact of microbiome competition on model outcomes, a second lot of ‘singlespecies simulations’ was run after removing microbiomepathogen interactions from all Ω (ε = 0, η = 0, ϕ = 1, χ = 0). Second, the impact of parameter uncertainty on model outcomes was evaluated. For each pathogen, model parameter values were resampled from their distributions (Appendix 1—table 2–6) using Latin Hypercube Sampling over 10,000 iterations, epidemiological outcomes were recalculated at population dynamic equilibrium, and partial rank correlation coefficients were calculated between each parameter and the pathogen’s (i) colonization prevalence and (ii) resistance rate (using the R package pse) (Chalom and Prado, 2013; Marino et al., 2008).
Public health interventions
Request a detailed protocolThree public health interventions were incorporated into the final model (see Appendix Equations A12–A15). First, contact precautions were assumed to represent physical or behavioral barriers that block opportunities for transmission, reducing transmission rates by the same fraction τ_{ipc} across all pathogens relative to baseline. Second, antibiotic stewardship programmes were assumed to alter antibiotic consumption patterns in the hospital. Two main types were considered: antibiotic reduction, which limits overall antibiotic prescribing by a fraction τ_{asp}, and antibiotic restriction, which adjusts the distribution of antibiotic classes consumed in the hospital by the same fraction. We considered two types of restriction, the first favouring classically narrowspectrum antibiotics (e.g. macrolides) over broadspectrum (e.g. quinolones), and the second favouring antibiotics that cause dysbiosis at very low or low rates (k={3,2}) over those causing dysbiosis at medium or high rates (k={1,0}) (see Appendix 1—table 8 for antibiotic classification). Third, microbiome recovery therapy was assumed to trigger recovery at rate 0.5 day^{−1} and was apportioned to the fraction τ_{pbt} of patients. Overall, different values assumed for intervention parameters (τ_{ipc}, τ_{asp}, τ_{pbt}) are interpreted as different levels of compliance to the respective interventions. Intervention efficacy was evaluated using the IRR and RRR (defined in Model and Results). Outcomes were matched across Monte Carlo simulations, such that IRRs and RRRs were calculated for each intervention and compliance level for each Ω. The distribution of outcomes is expressed as the median and 95% uncertainty interval.
Appendix 1
Model and assumptions
ODEs
The final ODE system, which incorporates horizontal gene transfer (HGT) of a transmissible resistance gene R into the mixed microbiomestrain competition model (Equation 6 in the main text), is given by:
This system of equations describes epidemiological colonization dynamics of a bacterial pathogen among host classes H_{j,k}, corresponding to patients of pathogen colonization status H (S, C^{S}, C^{R}), microbiome status j (e, equilibrium; or d, dysbiosis) and microbiome resistance profile k (s, not bearing R; r, bearing R). Equations can be found in an R file and Mathematica notebook available online (https://github.com/drmsmith/microbiomeR), and modeling assumptions are listed below.
Patient demography
We assume constant population size (N=1), balanced by a daily rate of hospital admission and discharge μ. The number of patients in each compartment is expressed as proportions of the total population. The proportion of patients entering each compartment j is given by admission fractions f_{j}. These are interpreted as representing prevalence of different host types in the community, which are assumed to be stable over time and unaffected by hospital dynamics. For simplicity, demography is expressed as Δ_{j} for each class j, which expand to:
where f_{C} is the proportion of patients colonized with the pathogen upon admission, f_{R} is the proportion of colonized patients bearing the focal resistant strain (i.e. the community resistance rate), f_{d} is the proportion of patients undergoing dysbiosis upon hospital admission, and f_{ω} is the proportion of admitted patients with microbiota bearing R.
Pathogen epidemiology
Pathogen colonization can be acquired through dynamic patienttopatient transmission via the dynamic force of infection λ, which is defined as the product of the pathogen’s intrinsic transmission rate (β) and the prevalence of patients colonized at a particular time t. For each strain P^{j}, with j in {S, R}, this is given by
In hospital environments, pathogen colonization incidence also results from endogenous acquisition (α). This subsumes alternative routes of acquisition not resulting from direct patienttopatient transmission, and can include processes like translocation and withinhost outgrowth of a subdominant/nondetectable/nontransmissible colony (into a dominant/detectable/transmissible colony) (Archambaud et al., 2019; Bootsma et al., 2007; Duval et al., 2019; Gurieva et al., 2018). This reflects the tenet that ‘everything is everywhere, but the environment selects’ (de Wit and Bouvier, 2006). We assume that pathogen colonization is not acquired endogenously among patients undergoing antibiotic treatment capable of clearing that pathogen, such that for each strain P^{j},
where α’ is the baseline rate in untreated patients. For the hospital simulation study (below), available literature estimates for endogenous acquisition rates were not strainspecific. These were generated by scaling speciesspecific rates by the baseline prevalence of each strain in the community, or
Pathogen colonization is naturally cleared by the host at rate γ. To reflect potential metabolic costs associated with bearing antibiotic resistance genes (Melnyk et al., 2015), we assumed that C^{R} (colonization with the resistant strain) is naturally cleared at a faster rate than C^{S} (colonization with the sensitive strain), such that γ_{S} < γ_{R}. This is given by
where c is interpreted as the fitness cost of bearing R.
Antibiotic resistance
Pathogen colonization is cleared by antibiotics according to strainspecific rates of effective antibiotic treatment, given by
Antibiotic exposure prevalence $\mathrm{a}$ is applied independently of colonization status to reflect high estimated rates of bystander selection for ARB, which predominantly colonize patients asymptomatically, are rarely detected and only opportunistically cause disease (Tedijanto et al., 2018). Among patients exposed to effective antibiotic therapy, pathogen colonization is cleared at a constant rate θ_{C}. The proportion of antibiotics that are ineffective is given by the antibiotic resistance level r_{j}, which is a proportion reflecting the pathogen’s innate resistance to the antibiotics to which it is exposed, ranging from complete antibiotic sensitivity (r_{j} = 0) to complete resistance (r_{j} = 1), For model application, effective antibiotic treatment is considered at the level of antibiotic class, with rates for each class varying across strains and species (see main text Equation 8).
Microbiome dynamics
Patient microbiota are considered in terms of their population dynamic stability, and are assumed either to be at stable equilibrium or undergoing temporary dysbiosis (population dynamic disequilibrium). Microbiome stability is disrupted via antibiotic treatment according to
For model application, rates of antibioticinduced microbiome dysbiosis are also considered at the level of antibiotic class (see main text Equation 7). Microbiome stability recovers at rate δ, but only among patients not actively undergoing antibiotic treatment, such that
where δ ‘is the baseline rate in untreated hosts. Further, during model application we assume that microbiota can still recover among patients exposed to aminoglycosides and tetracyclines, the antibiotic classes inducing microbiome dysbiosis at the lowest rate (k=3, see Materials and methods), such that
Microbiomepathogen competition
We assume that the host microbiome and the focal pathogen compete ecologically within the host, such that host microbiota limit pathogen transmission, persistence and endogenous acquisition. We propose three examples of microbiomepathogen competition – colonization resistance (ε), resource competition (η), and ecological release (ϕ) – and conceptualize how they affect pathogen colonization processes, and by extension how microbiome dysbiosis predisposes patients to pathogen colonization (see Figure 2 in the main text). Interactions are given by
where, respectively, ε reduces the transmission rate to patients with stable microbiota, η reduces the clearance rate among patients undergoing dysbiosis, and ϕ favours endogenous acquisition among patients undergoing dysbiosis. We assume that microbiomepathogen interactions are species and not strainspecific, that is that ε, η and ϕ apply equally to P^{S} and P^{R}.
Horizontal gene transfer
HGT is conceptualized as twoway withinhost transfer of the focal resistance gene R, either from an Rbearing pathogen strain P^{R} to cocolonizing microbiota, or from Rbearing microbiota to a cocolonizing drugsensitive pathogen strain P^{S}. For simplicity, we assumed (i) a symmetric rate of HGT (χ) from each donor, (ii) no loss of resistance upon donation, (iii) no accumulation of resistance (e.g. plasmid copy number dependence), (iv) that all patients bear microbiota capable of receiving and transferring R, (v) that dysbiosis can accelerate the rate of transfer (χ_{d} ≥ χ_{e}), (vi) no impact of R on rates of microbiome dysbiosis and recovery, (vii) that a proportion f_{ω} of patients are admitted to hospital with microbiota bearing R, and lastly (viii) that microbiota of a proportion ω of patients can spontaneously acquire R subsequent to dysbiosis. The latter assumption was made to reflect increased expression of antibiotic resistance genes among host microbiota following antibiotic therapy, and can be interpreted as a corollary to endogenous pathogen acquisition (Ruppé et al., 2019).
Public health interventions
Three public health interventions τ were incorporated into Appendix Equation A1 for subsequent model application. First, contact precautions were assumed to represent physical or behavioral barriers that block opportunities for transmission, reducing transmission rates by the same fraction τ_{ipc} across all pathogens relative to baseline. This is given by
Second, antibiotic stewardship programmes were assumed to alter antibiotic consumption patterns in the hospital. Two main types were considered: antibiotic reduction, which limits overall antibiotic prescribing by a fraction τ_{asp}, given by
and antibiotic restriction, which adjusts the distribution of antibiotic classes consumed in the hospital. We considered two types of restriction, the first favouring classically narrowspectrum antibiotics (e.g. macrolides) over broadspectrum (e.g. quinolones), and the second favouring antibiotics that cause dysbiosis at very low or low rates (k={3,2}) over those causing dysbiosis at medium or high rates (k={1,0}). For both, the proportion of antibiotics in the restricted group (p_{restrict}) was reduced by the same fraction τ_{asp} without adjusting the distribution of antibiotics consumed within that group, given by
and the proportion of nonrestricted antibiotics was increased symmetrically. With this structure, τ_{asp} alters antibiotic consumption for the same proportion of hospital patients across all three stewardship interventions. In all simulations, p_{restrict} > τ_{asp}.
Third, microbiome recovery therapy was assumed to trigger recovery at rate 0.5 day^{−1} and was apportioned to the fraction τ_{pbt} of patients, such that the overall rate of microbiome recovery when including these interventions (δ_{pbt}) is given by
For simplicity, different values assumed for intervention parameters (τ_{ipc}, τ_{asp}, τ_{pbt}) are interpreted as different levels of compliance to the respective interventions.
Part 1: Model evaluation and parameterization
R_{0} expressions
R_{0} was calculated for P^{R} across Equations 1, 3, 4 and 6 from the main text (see Mathematica notebook), and was interpreted as the expected number of secondary patients colonized by an index patient in a fully susceptible (uncolonized) hospital population (i.e. at diseasefree equilibrium, DFE, as indicated by the equilibrium symbol ^ above state variables). To reflect this epidemiological context, we made two additional assumptions specifically for R_{0} calculations: no C^{R} input from the community (f_{R}=0), and that the rate of endogenous acquisition α_{R} scales linearly with current colonization prevalence C^{R}. R_{0} expressions were derived following van den Driessche, 2017.
For the susceptiblecolonized model,
where $\hat{S}=1$, such that pathogens with higher rates of transmission (β) and endogenous acquisition (α_{R}) have higher R_{0}, while those with higher rates of natural clearance (γ_{R}) have lower R_{0}. Higher rates of effective antibiotic treatment (σ_{R}) and patient admission and discharge (μ) also reduce R_{0}.
For the strain competition model, we derive R_{0} for P^{R} assuming that P^{S} is at endemic equilibrium with input from the community (f_{C} > 0, f_{R} = 0). In this context, the same R_{0} expression was found as for the susceptiblecolonized model (Appendix Equation A16), except R_{0} is reduced because of a lower equilibrium prevalence of susceptible hosts when P^{S} is endemic ($\hat{S}<1$). DFE were not found analytically for models with strain competition and were solved numerically (details below).
For the microbiome model, R_{0} was calculated using nextgeneration theory as the spectral radius of the nextgeneration matrix (NGM),
where NGM is the product of the transmission matrix F (describing the rates at which existing colonies cause colonization in new hosts) and the inverse transition matrix V (describing the rates at which colonized hosts shift between colonized classes or are removed). These are given by
and
which give
where the two terms in the numerator correspond to pathogen acquisition, respectively, in hosts with a stable microbiome and in those undergoing dysbiosis. Antibioticinduced dysbiosis (σ_{m}) and treatment (σ_{R}) are found in both numerator and denominator.
When assuming no microbiomepathogen interactions (ε = 0, η = 0, ϕ = 1), DFE is
and R_{0} reduces to the same expression as given by the susceptiblecolonized model
With two strains, as in the microbiomestrain competition model, the same R_{0} expression is found for P^{R} as the singlestrain model, but when P^{S} is endemic,
so strain competition dampens R_{0} relative to the singlestrain microbiome competition model. R_{0} expressions are evaluated numerically in Figure 1 and Appendix 1—figure 2 and 3.
Epidemiological outcomes
Colonization prevalence
Colonization prevalence was defined as the proportion of hospital patients colonized (C^{R}) by the focal drugresistant strain (P^{R}) at endemic equilibrium:
Colonization prevalence of the drugsensitive strain P^{S} was also calculated:
These were found by substituting a corresponding vector of numerical parameter values Ω and solving ODE systems through numerical integration. Solutions were found using the lsoda method of the function ode from the package deSolve in R (version 3.6.0). For evaluation of the theoretical pathogen (part 1), solutions were corroborated using the function NSolve in Mathematica.
Resistance rate
The resistance rate was defined as the proportion of patients colonized with the focal resistant strain P^{R} relative to the sensitive strain P^{S}, calculated using equilibrium prevalence values as C^{R}/(C^{S} + C^{R}).
Colonization incidence
Nosocomial colonization incidence was defined as the daily rate of colonization acquisition within the hospital, calculated separately for each route of acquisition. For P^{R}, incidence rates corresponding to transmission (inc_{β}), endogenous acquisition (inc_{α}) and HGT (inc_{χ}) were calculated as
We evaluated incidence at dynamic equilibrium by solving ODE systems numerically for each Ω, then using resulting equilibria as initial values for subsequent numerical integration to calculate incidence from time t to t +1.
Parameter space for P^{R} colonization
Part 1: Supplementary results
Here, we expand on results for Part one from the main text. For both strain competition and microbiome competition models, C^{R} increases monotonically with antibiotic use under the assumption of complete antibiotic resistance (r_{R} = 1) (Appendix 1—figure 1). In Appendix 1—figure 2, we evaluate R_{0} in the context of each microbiomepathogen interaction separately, and in concert, over a range of assumed values. This demonstrates how colonization resistance dampens pathogen R_{0}, while resource competition and ecological release augment it. Taken together, strong microbiomepathogen interactions can protect patients from P^{R} colonization in some conditions (e.g. at low $\mathrm{a}$ and r_{R}), while predisposing them to P^{R} colonization in others (high $\mathrm{a}$ and r_{R}). In Appendix 1—figure 3, we demonstrate numerically that introducing a drugsensitive strain P^{S} to the microbiome model reduces R_{0} of P^{R}. In Appendix 1—figure 4, we again demonstrate how the strength of a pathogen’s microbiome interactions (ε, η, ϕ) and the level of antibiotic resistance for P^{R} (r_{R}) drive epidemiological consequences of antibiotic use, here disentangling the role of antibioticinduced pathogen clearance (θ_{c}) from antibioticinduced microbiome dysbiosis (θ_{m}). In Appendix 1—figure 5, we simulate dynamic responses of P^{R} to public health interventions. First, dynamic equilibria were found, from which ODEs were integrated for an additional 365 days, evaluating C^{R} and the resistance rate over time in the context of each microbiomepathogen interaction (separately and in concert). We introduced theoretical interventions at days 90 (halving r_{R} from 0.8 to 0.4), 180 (halving θ_{m} from 1.0 to 0.5), and 270 (halving $\mathrm{a}$ from 0.2 to 0.1), demonstrating that an otherwise identical pathogen can experience diverse, and sometimes opposing responses to public health interventions in the context of different microbiomepathogen interactions and epidemiological outcomes.
In Appendix 1—figure 6, we show how HGT impacts P^{R} colonization incidence, in contrast to prevalence and resistance rate as presented in Figure 4 in the main text. In Appendix 1—figure 7, we demonstrate effects of HGT on pathogen colonization outcomes. First, we vary key epidemiological parameters over broad ranges to show how the absolute impact of HGT depends on how other parameters mediate colonization dynamics (Appendix 1—figure 7A). Second, we show how increasing the rate of HGT in hosts undergoing dysbiosis (χ_{d}), while holding the rate constant in hosts with stable microbiota (χ_{e}), augments HGT’s impact but does not substantively shift the level of antibiotic exposure that maximizes C^{R} (Appendix 1—figure 7B). Third, we demonstrate how effects of HGT depend on characteristics of the resistance gene being transmitted (Appendix 1—figure 7C). For instance, a metabolically costly but highly drugresistant gene R (c = 2, r_{R} = 0.8) is potentially disadvantageous for the pathogen species as a whole at low antibiotic use (reducing total pathogen prevalence) but advantageous at high antibiotic use (increasing total prevalence).
Part 2: Model application, species characterization, and parameterization
Host and pathogen parameterization
Expert elicitation
Protocol
Expert elicitation is a scientific consensus methodology involving the estimation of unknown parameter values from subjectmatter experts. We developed an expert elicitation protocol, using published recommendations for elicitation of clinical parameters,(Johnson et al., 2010) to characterize the role of microbiome dysbiosis in driving hospital colonization dynamics across the ARB included in this study (provided separately). First, the scientific context of the study was introduced, including the potential epidemiological relevance of microbiomepathogen interactions, and potential roles for intraspecific strain competition and microbiome dysbiosis in explaining how antibiotics select for resistance. Second, the MATCH uncertainty elicitation tool and ‘chips and bins’ method were introduced, allowing experts to build distributions visually in the form of histograms to quantify parameter values and associated uncertainty. (Morris et al., 2014) Third, potential cognitive biases were reviewed and tips for reliable parameter estimation were provided. Finally, experts were asked to respond to questions and quantify their beliefs and uncertainty about how microbiome dysbiosis among hospital patients affects pathogen colonization processes (acquisition, clearance, growth and horizontal gene transfer). It was stressed that the goal of the elicitation was not to characterize precise numerical estimates of particular mechanistic interactions, but rather to quantify the relative impact of microbiome dysbiosis on different colonization outcomes in the hospital setting. Experts were initially identified through a review of authors of relevant publications in the field. At the end of each interview, experts were further asked to refer us to additional experts in the form of snowball sampling. Elicitation results are reported anonymously, but all experts agreed to their acknowledgment in this article.
Results and interpretation
Of 23 invited experts, ten ultimately participated in elicitation interviews (response rate 43.5%). Their specialties include bacteriology, internal medicine, intensive care, infectious disease epidemiology and clinical microbiology. Model parameters, the clinical parameters that experts were asked to estimate, and assumed links between them are provided in Appendix 1—table 7. First, experts were asked about the relevance of intraspecific strain competition for the nosocomial epidemiology of each pathogen. Then, for key colonization parameters (β, α, γ, χ), experts were asked whether hospital patients experiencing antibioticinduced dysbiosis are more or less likely to experience the corresponding colonization outcome than patients not experiencing dysbiosis (Appendix 1—figure 8). If yes, they were asked to quantify the relative impact of dysbiosis on that outcome i for each ARB j as a histogram (estimated using MATCH) given by a vector x_{i,j,k} for each expert k. Experts made estimates for each ARB consecutively, such that parameters for each ARB were estimated relative to other included ARB.
Raw data for expert distributions are provided in Appendix 1—figure 9. Experts differed substantially in estimated ranges of different parameters, but Friedman’s tests using medians of each parameter distribution suggested that the species rank order was conserved across experts, with C. difficile generally having the strongest estimated microbiomepathogen interaction coefficients, and MRSA the weakest. To conserve species rank order when combining expert estimates x_{i,j,k} to form pooled histograms, distributions were recentred by the mean relative distance z_{i,j} between each j and the reference ARB (taken as MRSA) over all experts k (Appendix 1—figure 10). For C. difficile and MRSA, HGT was excluded from simulations and pooled distributions were not recentred. Each expert distribution was weighted equally, contributing to 10% of the final pooled distribution for each species and parameter. Pooled histograms were fit to six candidate distributions (normal, lognormal, Weibull, Cauchy, exponential and gamma), and the final distribution was selected as the fitted distribution giving the lowest Akaike Information Criterion (using the function fitdist from the R package fitdistrplus). Final fitted distributions for each parameter and ARB are given in Appendix 1—table 3–6, and for select parameters in Figure 5.
Among other observations, experts highlighted (i) uncertainty about an influence of strain competition on hospital colonization dynamics for all ARB, with the majority of experts believing strain competition was irrelevant for C. difficile but at least partially relevant for other ARB, (ii) unanimous certainty about a role for antibioticinduced microbiome dysbiosis as a driver of colonization for C. difficile, ESBLEC and CPKP, with somewhat less certainty for MRSA, and (iii) certainty about a role for horizontal gene transfer (HGT) for ESBLEC and CPKP (Appendix 1—figure 8).
Antibiotic exposure data
Part 2: Supplementary results
Baseline colonization dynamics
Intervention evaluation
Metaanalysis of hospital antibiotic stewardship interventions
To assess consistency of intervention outcomes with the literature, we compared antibiotic stewardship results to a systematic review and metaanalysis of hospital antibiotic stewardship interventions initially conducted by Baur et al., 2017. They included interventional studies published worldwide from Jan 1960 to May 2016 describing incidence of bacterial colonization or infection among hospital inpatients. Stewardship interventions were heterogeneous and included audits, antibiotic restriction, antibiotic cycling, antibiotic mixing, feedback, guideline implementation and education. The authors found no evidence of publication bias or small study effects.
Their analysis included 32 studies over 9 million patientdays, though stewardship interventions were coimplemented with other infection control measures (e.g. hand hygiene, patient isolation) in 31% of studies. When limited to stewardship, patients were still significantly less likely to experience colonization and infection (IRR=0.81, 0.67–0.97), but speciesspecific estimates were not made. When excluding studies that included nonstewardship interventions, ten studies remained for MRSA, (Arda et al., 2007; Frank et al., 1997; Mach et al., 2007; Marra et al., 2009; Meyer et al., 2010; Miyawaki et al., 2010; Peto et al., 2008; Smith et al., 2008; Yeo et al., 2012; Zou et al., 2015) seven for C. difficile, (Borde et al., 2015; Dubrovskaya et al., 2012; Frank et al., 1997; Leung et al., 2011; Lübbert et al., 2014; Malani et al., 2013; Schön et al., 2011) two for ESBLE. coli, (Arda et al., 2007; Chong et al., 2013) and one for CPK. pneumoniae (Arda et al., 2007) For ESBLE. coli, we further included studies of ESBLproducing Enterobacteriaceae and Gramnegative bacteria (Grohs et al., 2014; Takesue et al., 2010) Using their published data and methodology, we calculated incidence risk ratios (IRR) using the standard inversevariance method, and pooled IRRs across species using mixedeffects metaanalysis models (R package metafor). Results for each study and pooled results for each pathogen are presented in Appendix 1—figure 18.
Data availability
Model equations and parameter values are provided in the manuscript, as well as in supporting R files and a Mathematica notebook available online at https://github.com/drmsmith/microbiomeR (copy archived athttps://archive.softwareheritage.org/swh:1:rev:a3682a24970d79e4f748952ecc49fcdb16adf48f).
References

SoftwareLa Consommation D’antibiotiques en France en 2016 – Rapport De l’ANSMLa Consommation D’antibiotiques en France en 2016 – Rapport De l’ANSM.

Antibiotic use and good practice in 314 french hospitals: the 2010 SPA2 prevalence studyMédecine Et Maladies Infectieuses 45:475–480.https://doi.org/10.1016/j.medmal.2015.10.001

Mathematical models of infection transmission in healthcare settings: recent advances from the use of network structured dataCurrent Opinion in Infectious Diseases 30:410–418.https://doi.org/10.1097/QCO.0000000000000390

The transmission dynamics of antibioticresistant Bacteria: the relationship between resistance in commensal organisms and antibiotic consumptionProceedings of the Royal Society of London. Series B: Biological Sciences 264:1629–1638.https://doi.org/10.1098/rspb.1997.0227

Estimating national trends in inpatient antibiotic use among US hospitals from 2006 to 2012JAMA Internal Medicine 176:1639–1648.https://doi.org/10.1001/jamainternmed.2016.5651

Risk of subsequent Sepsis within 90 days after a hospital stay by type of antibiotic exposureClinical Infectious Diseases 66:1004–1012.https://doi.org/10.1093/cid/cix947

Natural history and decolonization strategies for ESBL/carbapenemresistant Enterobacteriaceae carriage: systematic review and metaanalysisJournal of Antimicrobial Chemotherapy 71:2729–2739.https://doi.org/10.1093/jac/dkw221

Prevalence and pathogenicity of Clostridium difficile in hospitalized patientsArchives of Internal Medicine 156:1449.https://doi.org/10.1001/archinte.1996.00440120107012

Outpatient antibiotic use in France between 2000 and 2010: after the nationwide campaign, it is time to focus on the elderlyAntimicrobial Agents and Chemotherapy 58:71–77.https://doi.org/10.1128/AAC.0181313

Impact of antimicrobial therapy on the gut microbiomeJournal of Antimicrobial Chemotherapy 74:i6–i15.https://doi.org/10.1093/jac/dky530

Send more data: a systematic review of mathematical models of antimicrobial resistanceAntimicrobial Resistance & Infection Control 7:117.https://doi.org/10.1186/s1375601804061

The evolution of antibiotic resistance in a structured host populationJournal of the Royal Society Interface 15:20180040.https://doi.org/10.1098/rsif.2018.0040

Evolutionary epidemiology models to predict the dynamics of antibiotic resistanceEvolutionary Applications 12:365–383.https://doi.org/10.1111/eva.12753

An algorithm to estimate the importance of bacterial acquisition routes in hospital settingsAmerican Journal of Epidemiology 166:841–851.https://doi.org/10.1093/aje/kwm149

Tracking strains in the microbiome: insights from metagenomics and modelsFrontiers in Microbiology 7:712.https://doi.org/10.3389/fmicb.2016.00712

MetaAnalysis of antibiotics and the risk of CommunityAssociated Clostridium difficile infectionAntimicrobial Agents and Chemotherapy 57:2326–2332.https://doi.org/10.1128/AAC.0217612

Microbiotamediated colonization resistance against intestinal pathogensNature Reviews Immunology 13:790–801.https://doi.org/10.1038/nri3535

Impact of antibiotic gut exposure on the temporal changes in microbiome diversityAntimicrobial Agents and Chemotherapy 63:19.https://doi.org/10.1128/AAC.0082019

SoftwarePse: Parameter Space Exploration with Latin HypercubesPse: Parameter Space Exploration with Latin Hypercubes.

Quantifying drivers of antibiotic resistance in humans: a systematic reviewThe Lancet Infectious Diseases 18:e368–e378.https://doi.org/10.1016/S14733099(18)302962

Host population structure and treatment frequency maintain balancing selection on drug resistanceJournal of the Royal Society Interface 14:20170295.https://doi.org/10.1098/rsif.2017.0295

Fifty shades of graft: how to improve the efficacy of faecal Microbiota transplantation for decolonization of antibioticresistant BacteriaInternational Journal of Antimicrobial Agents 53:553–556.https://doi.org/10.1016/j.ijantimicag.2019.03.008

Fecal Microbiota transplantation to eradicate vancomycinresistant enterococci colonization in case of an outbreakMédecine Et Maladies Infectieuses 49:214–218.https://doi.org/10.1016/j.medmal.2018.11.002

Withinhost dynamics shape antibiotic resistance in commensal BacteriaNature Ecology & Evolution 3:440–449.https://doi.org/10.1038/s415590180786x

Targeted adsorption of molecules in the Colon with the novel adsorbentbased medicinal product, DAV132: a proof of concept study in healthy subjectsThe Journal of Clinical Pharmacology 55:10–16.https://doi.org/10.1002/jcph.359

Protection of the human gut microbiome from antibioticsThe Journal of Infectious Diseases 217:628–636.https://doi.org/10.1093/infdis/jix604

Intrinsic epidemicity of Streptococcus pneumoniae depends on strain serotype and antibiotic susceptibility patternAntimicrobial Agents and Chemotherapy 55:5255–5261.https://doi.org/10.1128/AAC.0024911

Antibiotic stewardship for intraabdominal infections: early impact on antimicrobial use and patient outcomesInfection Control & Hospital Epidemiology 33:427–429.https://doi.org/10.1086/664765

Quantifying Transmission of Clostridium difficile within and outside Healthcare SettingsEmerging Infectious Diseases 22:608–616.https://doi.org/10.3201/eid2204.150455

The demographic determinants of human microbiome healthTrends in Microbiology 23:134–141.https://doi.org/10.1016/j.tim.2014.11.005

Community interactions and spatial structure shape selection on antibiotic resistant lineagesPLOS Computational Biology 14:e1006179.https://doi.org/10.1371/journal.pcbi.1006179

Decrease in expenditures and selected nosocomial infections following implementation of an antimicrobialprescribing improvement programClinical Performance and Quality Health Care 5:180–188.

The negative impact of antibiotic resistanceClinical Microbiology and Infection 22:416–422.https://doi.org/10.1016/j.cmi.2015.12.002

Fighting the spread of AmpChyperproducing Enterobacteriaceae: beneficial effect of replacing ceftriaxone with cefotaximeJournal of Antimicrobial Chemotherapy 69:786–789.https://doi.org/10.1093/jac/dkt403

Resource competition and host feedbacks underlie regime shifts in gut MicrobiotaThe American Naturalist 198:1–12.https://doi.org/10.1086/714527

Modeling the effect of DAV132, a novel ColonTargeted adsorbent, on fecal concentrations of moxifloxacin and gut Microbiota diversity in healthy volunteersClinical Pharmacology & Therapeutics 109:1045–1054.https://doi.org/10.1002/cpt.1977

Fecal Microbiota transplantation: the evolving risk landscapeAmerican Journal of Gastroenterology 116:647–656.https://doi.org/10.14309/ajg.0000000000001075

The transmissibility of AntibioticResistant Enterobacteriaceae in intensive care unitsClinical Infectious Diseases 66:489–493.https://doi.org/10.1093/cid/cix825

Angiogenins: a new class of microbicidal proteins involved in innate immunityNature Immunology 4:269–273.https://doi.org/10.1038/ni888

MultidrugResistant bacterial infections in U.S. hospitalized patients, 20122017New England Journal of Medicine 382:1309–1319.https://doi.org/10.1056/NEJMoa1914433

A valid and reliable belief elicitation method for bayesian priorsJournal of Clinical Epidemiology 63:370–383.https://doi.org/10.1016/j.jclinepi.2009.08.005

Control of pathogens and pathobionts by the gut MicrobiotaNature Immunology 14:685–690.https://doi.org/10.1038/ni.2608

Role of the gut Microbiota in immunity and inflammatory diseaseNature Reviews Immunology 13:321–335.https://doi.org/10.1038/nri3430

Impact of antibiotic exposure patterns on selection of CommunityAssociated MethicillinResistant Staphylococcus aureus in hospital settingsAntimicrobial Agents and Chemotherapy 55:4888–4895.https://doi.org/10.1128/AAC.0162610

Fecal Microbiota transplantation for Clostridium difficile infection: systematic review and metaanalysisAmerican Journal of Gastroenterology 108:500–508.https://doi.org/10.1038/ajg.2013.59

The intestinal Microbiota: antibiotics, colonization resistance, and enteric pathogensImmunological Reviews 279:90–105.https://doi.org/10.1111/imr.12563

Mathematical modelling for antibiotic resistance control policy: do we know enough?BMC Infectious Diseases 19:1011.https://doi.org/10.1186/s128790194630y

Fitness and competitive growth advantage of new gentamicinsusceptible MRSA clones spreading in french hospitalsJournal of Antimicrobial Chemotherapy 47:277–283.https://doi.org/10.1093/jac/47.3.277

Mathematical modelling to study the horizontal transfer of antimicrobial resistance genes in Bacteria: current state of the field and recommendationsJournal of the Royal Society Interface 16:20190260.https://doi.org/10.1098/rsif.2019.0260

Horizontal transfer of antibiotic resistance genes in clinical environmentsCanadian Journal of Microbiology 65:34–44.https://doi.org/10.1139/cjm20180275

Growing a "positive culture" of antimicrobial stewardship in a community hospitalThe Canadian Journal of Hospital Pharmacy 64:314–320.https://doi.org/10.4212/cjhp.v64i5.1065

Antimicrobial use and antimicrobial resistance: a population perspectiveEmerging Infectious Diseases 8:347–354.https://doi.org/10.3201/eid0804.010312

Can the antibiotic prescription practice in a hospital be influenced by inhouse guidelines? an interventional study at the university hospital halle GermanyDeutsche Medizinische Wochenschrift 139:2578–2584.https://doi.org/10.1055/s00341387220

The human intestinal microbiome in health and diseaseNew England Journal of Medicine 375:2369–2379.https://doi.org/10.1056/NEJMra1600266

Clinical and economic outcomes from a community hospital's antimicrobial stewardship programAmerican Journal of Infection Control 41:145–148.https://doi.org/10.1016/j.ajic.2012.02.021

A methodology for performing global uncertainty and sensitivity analysis in systems biologyJournal of Theoretical Biology 254:178–196.https://doi.org/10.1016/j.jtbi.2008.04.011

The effect of limiting antimicrobial therapy duration on antimicrobial resistance in the critical care settingAmerican Journal of Infection Control 37:204–209.https://doi.org/10.1016/j.ajic.2008.06.008

The fitness costs of antibiotic resistance mutationsEvolutionary Applications 8:273–283.https://doi.org/10.1111/eva.12196

Linking within and betweenhost dynamics in the evolutionary epidemiology of infectious diseasesTrends in Ecology & Evolution 23:511–517.https://doi.org/10.1016/j.tree.2008.05.009

Comparison of the Burdens of HospitalOnset, Healthcare FacilityAssociated Clostridium difficile Infection and of HealthcareAssociated Infection due to MethicillinResistant Staphylococcus aureus in Community HospitalsInfection Control & Hospital Epidemiology 32:387–390.https://doi.org/10.1086/659156

A webbased tool for eliciting probability distributions from expertsEnvironmental Modelling & Software 52:1–4.https://doi.org/10.1016/j.envsoft.2013.10.010

Systematic comparison of coexistence in models of drugsensitive and drugresistant pathogen strainsTheoretical Population Biology 133:150–158.https://doi.org/10.1016/j.tpb.2019.12.001

Patient to healthcare personnel transmission of MRSA in the nonintensive care unit settingInfection Control & Hospital Epidemiology 41:601–603.https://doi.org/10.1017/ice.2020.10

Estimating the burden of antimicrobial resistance: a systematic literature reviewAntimicrobial Resistance & Infection Control 7:58.https://doi.org/10.1186/s137560180336y

Contribution of mathematical modeling to the fight against bacterial antibiotic resistanceCurrent Opinion in Infectious Diseases 24:279–287.https://doi.org/10.1097/QCO.0b013e3283462362

Species interactions drive the spread of ampicillin resistance in humanassociated gut MicrobiotaEvolution, Medicine, and Public Health 9:256–266.https://doi.org/10.1093/emph/eoab020

Faecal carriage of carbapenemaseproducing Gramnegative bacilli in hospital settings in southern franceEuropean Journal of Clinical Microbiology & Infectious Diseases 34:899–904.https://doi.org/10.1007/s1009601422981

Carriage of ESBLproducing Enterobacteriaceae in French hospitals: the PORTABLSE studyJournal of Hospital Infection 98:247–252.https://doi.org/10.1016/j.jhin.2017.11.022

A ColonTargeted adsorbent (DAV132) does not affect the pharmacokinetics of warfarin or clonazepam in healthy subjectsClinical Pharmacology in Drug Development 10:908–917.https://doi.org/10.1002/cpdd.901

Hospitalization type and subsequent severe sepsisAmerican Journal of Respiratory and Critical Care Medicine 192:581–588.https://doi.org/10.1164/rccm.2015030483OC

Effects of treatment with antimicrobial agents on the human colonic microfloraTherapeutics and Clinical Risk Management 4:1343–1358.https://doi.org/10.2147/tcrm.s4328

Application of dynamic modelling techniques to the problem of antibacterial use and resistance: a scoping reviewEpidemiology and Infection 146:2014–2027.https://doi.org/10.1017/S0950268818002091

The gut Microbiota shapes intestinal immune responses during health and diseaseNature Reviews Immunology 9:313–323.https://doi.org/10.1038/nri2515

Microbiota: a key orchestrator of Cancer therapyNature Reviews Cancer 17:271–285.https://doi.org/10.1038/nrc.2017.13

Faecal Microbiota transplantation for eradicating carriage of multidrugresistant organisms: a systematic reviewClinical Microbiology and Infection 25:958–963.https://doi.org/10.1016/j.cmi.2019.04.006

A comparative study of three methods to evaluate an intervention to improve empirical antibiotic therapy for acute bacterial infections in hospitalized patientsScandinavian Journal of Infectious Diseases 43:251–257.https://doi.org/10.3109/00365548.2010.544326

Infection due to Clostridium difficile among elderly residents of a longtermcare facilityClinical Infectious Diseases 17:672–678.https://doi.org/10.1093/clinids/17.4.672

A modeling framework for the evolution and spread of antibiotic resistance: literature review and model categorizationAmerican Journal of Epidemiology 178:508–520.https://doi.org/10.1093/aje/kwt017

'Blooming' in the gut: how dysbiosis might contribute to pathogen evolutionNature Reviews Microbiology 11:277–284.https://doi.org/10.1038/nrmicro2989

Variability in antibiotic use across Ontario acute care hospitalsJournal of Antimicrobial Chemotherapy 72:554–563.https://doi.org/10.1093/jac/dkw454

A payer perspective of the hospital inpatient additional care costs of antimicrobial resistance in France: a matched CaseControl studyApplied Health Economics and Health Policy 17:381–389.https://doi.org/10.1007/s4025801804511

Modelling the impact of curtailing antibiotic usage in food animals on antibiotic resistance in humansRoyal Society Open Science 4:161067.https://doi.org/10.1098/rsos.161067

Reproduction numbers of infectious disease modelsInfectious Disease Modelling 2:288–303.https://doi.org/10.1016/j.idm.2017.06.002

Faecal carriage of multidrugresistant Gramnegative bacilli during a nonoutbreak situation in a french university hospitalJournal of Antimicrobial Chemotherapy 65:2455–2458.https://doi.org/10.1093/jac/dkq333

Prospective audit and feedback on antibiotic prescription in an adult hematologyoncology unit in SingaporeEuropean Journal of Clinical Microbiology & Infectious Diseases 31:583–590.https://doi.org/10.1007/s1009601113516

Shiga ToxinProducing Escherichia coli transmission via fecal Microbiota transplantClinical Infectious Diseases : An Official Publication of the Infectious Diseases Society of America 72:e876–e880.https://doi.org/10.1093/cid/ciaa1486

Trends and correlation of antibacterial usage and bacterial resistance: time series analysis for antibacterial stewardship in a chinese teaching hospital (20092013)European Journal of Clinical Microbiology & Infectious Diseases 34:795–803.https://doi.org/10.1007/s1009601422936

Population biological principles of drugresistance evolution in infectious diseasesThe Lancet Infectious Diseases 11:236–247.https://doi.org/10.1016/S14733099(10)702644
Decision letter

Gwenan M KnightReviewing Editor; London School of Hygiene and Tropical Medicine, United Kingdom

George H PerrySenior Editor; Pennsylvania State University, United States

Esther van KleefReviewer; Institute of Tropical Medicine, Belgium

Erik S WrightReviewer; University of Pittsburgh, United States
In the interests of transparency, eLife publishes the most substantive revision requests and the accompanying author responses.
Acceptance summary:
This paper provides a mathematical modelling framework for a stepwise incorporation of ecological complexity up to the effects of host microbiome on the spread of antibioticresistance. This provides a key first step in our understanding of the heterogeneous impact of the host microbiome on the spread of resistance, as demonstrated by the author's application of the model to four key pathogens.
Decision letter after peer review:
Thank you for submitting your article "Microbiomepathogen interactions drive epidemiological dynamics of antibiotic resistance: modelling insights for infection control" for consideration by eLife. Your article has been reviewed by 3 peer reviewers, and the evaluation has been overseen by a Reviewing Editor and George Perry as the Senior Editor. The following individuals involved in review of your submission have agreed to reveal their identity: Esther van Kleef (Reviewer #1); Erik S Wright (Reviewer #2).
The reviewers have discussed their reviews with one another, and the Reviewing Editor has drafted this to help you prepare a revised submission.
Summary:
The reviewers agreed that the approach is thoughtful, and the analysis performed is well complemented by simulations using relevant parameters in healthcare settings in hospitals. The big problem with this area, which is stressed by the authors, is the lack of information to support parameter choice. Despite being stressed, it was often not completely clear how and why parameters were chosen. Further clarity is required to support publication of this work which does offer a rational framework for linking the epidemiology and the microbiome, with an important hypothesis for antibiotic resistance.
Essential revisions:
The main revisions required are
– A restructuring of the methods to provide greater clarity of model structures and justification of parameter choices (is there no data? why was this chosen if so?). The driving force behind this should be the reproducibility of the work and the exploration of the results based on the parameter choices.
– A rewriting of the introduction and discussion to include more context of prior modelling work as well as to provide justification for the focus of the work
– In line with (1), all reviewers were in agreement that this work needs a clear single model description: a formal description of the model framework in one place. This could be in the supplementary, but it needs to be a well referenced place.
1) In general, the model the authors use is highly parameterized. There is a worry about circular reasoning when modelling. Parameters are derived from observations, then a model is constructed to recapitulate this observation. Then similar observations are used to show the model's validity. How did the authors ensure this type of bias was not incorporated into their modelling? This can become a bigger issue with more highly parameterized models, and that is why simpler models are preferable. Another approach is to use crossvalidation, although we are not sure how that would work here.
Abstract
2) The Abstract is vague and somewhat convoluted as written. It did not do the best job of selling the paper and does not convey any specific information.
Introduction
3) The authors need to tailor this to their specific question: where is the knowledge gap? We recommend that the authors think carefully about the statements they are making and consider revising the Introduction to include a more thorough analysis of prior modelling work, which is largely absent from the Introduction.
4) The study focuses on four particular pathogens, however, no motivation of why they are chosen is provided, nor the relevance of their analysis in a more general overview is discussed. The authors could expand on how their results could be extended to other contexts/species in which the role of microbiome on resistance is also relevant. This could also be included in the discussion.
5) The authors conclude that disruption of the microbiome on an individual level (due to antibiotics) can result in selection for AMR on a population level. In a way, the work is, at least in part, a modelbased implementation of the theoretical perspective of Lipsitch and Samore 2002 EID doi: 10.3201/eid0804.010312 (which describes how antibiotic effects on an individual, withinhost level, can affect population level transmission dynamics). Could the comparison to Lipsitch and Samore and the novelty of the improved datadriven approach (although still limited) be further discussed?
6) Antibiotic stewardship is treated a single entity, whereas in reality there are countless types of stewardship interventions and a concomitant variety of outcomes. Comparing the model's results directly to an analysis of previous stewardship interventions and attempting to show alignment does not make sense to me. How were these stewardship interventions reflective of the model? Why would we expect agreement or disagreement a priori? The fact that there was agreement with the model was somewhat concerning given the complexity of many stewardship interventions.
Methods
7) In attempting to replicate this study, it became intractably difficult to recreate the ODEs because the equations and variables are interspersed throughout the main text and supplementary appendix. For example, in Figure 1e the authors mention differences in resistance levels (γ) and not transmission rate (λ, related to β as in Equation (3)), the value of which could not be found anywhere in the text. Similarly, where is patient demography (Δ) incorporated into equations 3 or 4? There are certainly answers to these questions, but it would have been helpful to have the model more clearly presented in a central location. Finding the parameters and making sense of them was extremely challenging, if not infeasible. Therefore, we were unable to replicate the results to any extent.
8) It is unclear how the authors arrived at their parameter values throughout the manuscript. For example, colonization resistance (epsilon) as it is formulated is not compelling, and where did the value of this parameter come from? It is listed in Table S5 as being drawn from a Cauchy distribution. Why? Where did these numbers come from? Some of the parameters seem justified by literature, but the rest are perhaps debatable. This might be inevitable with a complex model, but the authors should minimally provide a robustness analysis to each of their lesssupported parameter choices. What would be the outcome if the sampled distributions had been wider (e.g., if the true value had been different by an order of magnitude)?
9) The inclusion of the 'r' parameter, allowing for partial resistance in the resistant strain, is interesting (and innovative for modelling frameworks to our knowledge). Nonetheless, it is not clear why the authors use a baseline of r = 0.8 to illustrate the tradeoff between antibiotic induced clearance and selection for Cr (Figure 1F), whereas in figure 2, the authors use r=0.4, i.e. higher levels of sensitivity in Cr. I suppose with a lower r, strain coexistence is more likely, and thus higher likelihood of Horizontal Gene Transfer (HGT) (?), but it would be, for consistency and comparison of the different withinhost dynamics, more transparent to use the same baseline parameter values, unless the authors can provide a good reason why different values of r are assumed at baseline between Figure 1 and Figure 2B?
10) The model presented in Eq. (1) needs a better introduction. Is this a new model or based on other models previously studied? This is explained later in the text citing references (30) and (31), but we recommend to introduce them earlier. It would be also good to point the reader earlier in the main text to the first sections of the supplement for a better understanding of the model assumptions and parameters employed.
11) Also, should the parameters be more consistently labelled across frameworks? For example, we recommend writing λ(N,C), or something similar, in all the equations to explicitly state it depends on these parameters. The way it is written gives the impression that λ is a constant parameter.
12) The r vs resistance rate (Cr/Cr+Cs) caused some confusion. This as a resistance rate of Cr = 0.8 but an r = 0.4 could actually mean that 0.8*0.4 = 0.08 of pathogens carried are fully resistant against the antibiotic, while this is 0.64 when r=0.8. Can the same baseline values for r be used in Figure 1 and Figure 2 or the differences justified?
13) Furthermore, in figure 1F: this is representing a one strain model. Is the Ce+Cd strain similar to the Cr? Can this be clarified?
14) In line with this. In Figure S3, R0 seems to represent the R0 of Cr. However, for Figure 1H and 1I, how should R0 be interpreted respectively? And for Figure S2? Please clarify how similar Cr (Figure 1H) and the one strain modelled (which could be similar to the Cr strain of Figure 1H) for figure 1I are
15) Figure 2A, the minimum resistance rate is 0.35 (see legend, and this appeared the case when no antibiotic induced clearance nor antibiotic induced microbiome disruption). This seems rather high. In particular as in Figure 2B, minimum levels of less than 0.1 are shown. Could the authors explain where this discrepancy is coming from and to what extend these high baseline resistance rates are representative?
16) For parameterisation of the models for the different pathogens, estimates from literature, notably existing modelling studies are chosen. However, these estimates, at least for C. difficile and MRSA, are coming from models that don't incorporate endogenous acquisition explicitly (for the Enterobacterieacia, model estimates from Gurieva are used, which do use a modelling framework incorporating an exo and endogenous acquisition). Therefore, the acquisition rates may be overestimated for these grampositives, in particular for C. diff the endogenous acquisition, which is listed as the main acquisition route (pp 14 line 330).
This may affect the estimated intervention effectiveness (in particular for antibiotic stewardship (less effective) and contact precautions (more effective)) under assumptions of the microbiome model. Could the authors use more realistic estimates (not sure models with explicit exo and endogenous acquisition exist for MRSA and C. diff), or at least, reflect on how different values of α and β affect the model results?
17) It would be helpful to explicitly say how the simulations were performed, i.e., mention that the ODEs where integrated. My first impression was that the authors performed stochastic simulations using the processes illustrated in panels A, B, C from Figure 1. This is probably an issue that only mathematical modellers would have, but could the authors please add this clarification.
Discussion
18) The model is described by a set of differential equations. This means the model ignores stochasticity of the population size dynamics of each strain. The authors could expand a bit more on the assumptions of the model and why fluctuations are ignored.
19) Also, spatial structure is particularly important when taking into account interactions between microbiome and pathogens, but it is also neglected in the model. It would be useful if the authors could discuss this as well.
20) The results shown often focus on steadystate quantities (such as colonisation prevalence), but the temporal evolution of the model is not discussed. It would be interesting if the authors could investigate the timescale of the different interactions taken into account, and how relevant they may be in the healthcare settings they investigate. Is the temporal evolution of the model any relevant for their conclusions?
21) Can the authors please add examples of the hypothetical microbiome recovery interventions they have in mind? Are the authors thinking of things like faecal transplantation? Please add, also to provide more practical interpretation of the work.
Reviewer #1:
The authors use a theoretical dynamic transmission modelling framework, incorporating both between and withinhost dynamics of antimicrobial resistant (AMR) pathogens, to stress the importance of the microbiome in the transmission dynamics of bacterial species. The work considered five different colonisation models, comprising different variations of 'traditional' epidemiological models, which, according to the authors, often do not include withinhost microbiomepathogen interactions, and add such interactions. With regards to the latter, three different withinhost microbiomepathogen interactions, as well as horizontal gene transfer are considered, and their effect on population transmission as well as interventions, evaluated.
The authors have done an extensive amount of work, and have done a great job in documenting all the model assumptions, data used and methods employed. Also, visualisations and illustrations are nice and informative.
However, the main difficulty with incorporating withinhost dynamics of AMR pathogens in a humantohuman transmission modelling framework is the lack of data to inform key parameters. The authors stress this limitation is also part of their work. The authors emphasise that they aimed to show in theory the importance of microbiomepathogen interactions, and have conducted an expert elicitation to partly inform parameters representing these dynamics.
The authors conclude that disruption of the microbiome on an individual level (due to antibiotics) can result in selection for AMR on a population level. In a way, the work is, at least in part, a modelbased implementation of the theoretical perspective of Lipsitch and Samore 2002 EID doi: 10.3201/eid0804.010312 (which similar to here, describes how antibiotic effects on an individual, withinhost level, can affect population level transmission dynamics). Therefore, I am somewhere inclined to think, are the main points made by the authors new? As this work here does, similar to this earlier work, not provide a datadriven approach (although tries to some extend). On the other hand, the work does, in contrast to Lipsitch and Samore, provide a framework for how to model these interactions, and illustrate to some extend which parameter values require what data, as well as which dynamics are most likely at play for which pathogens. Moreover, what is new, is that the authors try to illustrate how such interactions may affect the effectiveness of interventions, which is an interesting and, to the best of my knowledge, a novel component. My comments are largely requests for clarification.
Reviewer #2:
Here the authors tackle the important challenge of incorporating the microbiome into a traditional susceptiblecolonized transmission model. The microbiome acts as a third party in infections, with roles in supplying or suppressing antibiotic resistance. The authors found that microbiome interventions hold the promise of avoiding antibiotic resistance dissemination. This theoretical work encourages continued research into how the microbiome could be used to mitigate antibiotic resistance. However, I am skeptical of the translatability of the model's predictions to real infections. The model is heavily parameterized, without clear reasoning for some parameter choices. The work's greatest contribution, in my opinion, is that it offers a rational framework for how the microbiome could be incorporated into epidemiology, and provides an intriguing hypothesis that the microbiome could play a substantial role in combating antibiotic resistance.
I would like to say the authors clearly did a good job of conceptualizing infections and resistance dissemination. This was a substantial body of work and the manuscript was interesting to read. My biggest concern is that some of the results intuitively reflect parameters of the model, suggesting they might be artifacts of the approach rather than independent results.
Reviewer #3:
This manuscript by Smith et al. proposes a novel mathematical framework that contributes to the understanding of the role of microbiome on the resistance colonisation of bacterial pathogen populations exposed to antibiotics. The model proposed incorporates microbiome competition and takes into account dysbiosis effects on the population dynamics. Extended versions of this model that incorporates strainmicrobiome competition and horizontal gene transfer contributions are also presented. The model is simple to understand and takes into account important considerations in healthcare settings usually ignored in the literature. The simulations performed use parameters of four pathogens obtained carefully from a panel of experts. The fact they used these parameters makes the implications of their analysis quite relevant in clinical infections. Different assumptions on the interactions of each of the pathogens are incorporated in the model, and different public health interventions are analysed. The study performed is quite complete and supports the authors objectives. Their work demonstrate how relevant is the accounting for microbiome interactions for the understanding of antibiotics incidence on antibiotic resistance. This is particularly relevant for future work on mathematical modelling of resistance.
Even though the work presented is well supported by the analysis performed, there are certain aspects regarding the model that need to expanded. In particular:
1) The study focuses on four particular pathogens, however, no motivation of why they are chosen is provided, nor the relevance of their analysis in a more general overview is discussed. The authors could expand on how their results could be extended to other contexts/species in which the role of microbiome on resistance is also relevant.
2) The model is described by a set of differential equations. This means the model ignores stochasticity of the population size dynamics of each strain. The authors could expand a bit more on the assumptions of the model and why fluctuations are ignored. Also, spatial structure is particularly important when taking into account interactions between microbiome and pathogens, but it is also neglected in the model. It would be useful if the authors could discuss this as well.
3) The results shown often focus on steadystate quantities (such as colonisation prevalence), but the temporal evolution of the model is not discussed. It would be interesting if the authors could investigate the timescale of the different interactions taken into account, and how relevant they may be in the healthcare settings they investigate. Is the temporal evolution of the model any relevant for their conclusions?
https://doi.org/10.7554/eLife.68764.sa1Author response
Essential revisions:
The main revisions required are
– A restructuring of the methods to provide greater clarity of model structures and justification of parameter choices (is there no data? why was this chosen if so?). The driving force behind this should be the reproducibility of the work and the exploration of the results based on the parameter choices.
We have significantly restructured the paper to provide greater clarity of model structures and parameter choices. Changes include:
– Explicitly separating the work into two parts: Part 1, Development of a novel ODE modelling framework, evaluated using a generic nosocomial pathogen denoted P^{R}; and Part 2, Model application to highrisk nosocomial pathogens with speciesspecific parameterization. This distinction is made clear from the end of the introduction.
– Part 1 now begins with a leading paragraph introducing the methods used and outcomes evaluated, with immediate reference to the Methods for more technical detail, and to the Supplement for the final model framework.
– The Methods section now begins by clearly explaining that the final model combines the 4 models described in the main text (bacterial colonization; strain competition; microbiome competition; strainmicrobiome competition), and directs readers to the final model 5 incorporating HGT in supplementary equation S1.
– Model parameterization, including our use of expert elicitation to approximate microbiomepathogen interaction coefficients, is expanded and more clearly presented from the end of the introduction, and when introducing Part 2 of the Model and Results section, with additional detail provided in the third paragraph of the Methods section, and with full context and methodological detail in the supplement (Appendix section 4.2, page 21), including raw elicitation data (Appendix figures 89) and final recentred data (Appendix figure 10).
To support reproducibility of the work, in addition to the Mathematica notebook provided previously, R code has now been made available online at https://github.com/drmsmith/microbiomeR. Code includes all ODEs, functions used to calculate model results, parameter values for the generic pathogen P^{R}, and reproduces figures from Part 1 of this work.
– A rewriting of the introduction and discussion to include more context of prior modelling work as well as to provide justification for the focus of the work
In the previous introduction paragraph describing existing modelling work, we now clearly state that antibioticinduced microbiome dysbiosis is a longstanding hypothesis for the spread of resistance (and not our own), dating to at least Lipsitch and Samore 2002 EID.
“Disruption of the host microbiome is a longstanding theory explaining how antibiotics select for the spread of resistance at both the individual and population levels,(34) but most epidemiological models consider just one species of bacteria at a time, under the traditional assumption that antibiotic selection for resistance results from intraspecific competition between cocirculating strains.” (35–37)
We have included a new introductory paragraph describing contemporary modelling studies in the field of antibiotic resistance epidemiology. We stress the importance of withinhost ecological competition among cocolonizing bacteria explored in contemporary work, e.g; Davies et al. 2019 Nature Ecol Evol.
“Accounting for other forms of complexity in epidemiological models – from treatment intensity, to ageassortative contact behaviour, to hospital referral networks, to animalhuman interactions, to genetic linkage between resistance and nonresistance genes – has helped to unravel the many, disparate forces that contribute to drive the spread of resistance.(42–47) Nevertheless, withinhost competition between cocolonizing bacteria is a key mechanism of selection for antibiotic resistance dissemination, and an active area of research at the forefront of resistance modelling.”
“For instance, the ‘mixedcarriage’ model by Davies et al. demonstrates how intraspecific competition results in negative frequencydependent selection for either of two competing strains, and provides a satisfying mechanistic explanation for widespread strain coexistence at the population level.(49) However, contemporary work has stopped short of evaluating consequences of betweenspecies competition on resistance epidemiology. Yet for many ARB, including emerging highpriority multidrugresistant bacteria like extendedspectrum βlactamase (ESBL) producing Enterobacteriaceae, interactions with the host microbiome may be important mediators of nosocomial colonization dynamics.”
Additional detail has also been provided in the Discussion to better contextualize our work, including:
– A sentence in the first paragraph broadly contextualizing our work:
“We present a modelling framework that includes withinhost ecological costs of antibiotic use in the form of microbiome dysbiosis, incorporating a leading hypothesis for antibiotic selection into classical models of resistance epidemiology.”
– The new paragraph on microbiome recovery interventions now references relevant modelling work in this context: two withinhost modelling papers newly released since initial submission of our work (Guittar et al. Am Nat 2021 doi: 10.1086/714527; Guk et al. Clin Pharmacol Ther 2021 doi: 10.1002/cpt.1977).
“Modelling has been used previously to evaluate impacts of such microbiomeoriented interventions at the withinhost level,(74,75) but knockon impacts on ARB transmission dynamics and epidemiological burden have not been evaluated previously.”
– The paragraph comparing modelling results to the literature now makes clear reference to which results come from modelling vs. clinical studies.
– The paragraph discussing alternative means of modeling antibiotic selection in previous work now makes a clear link to our work:
“These strategies reflect widespread recognition that antibiotic use favours ARB acquisition, through erosion of colonization resistance or other supposed mechanisms, and independent of potential competition with other strains. The present work formalizes examples of the microbiomepathogen interactions that underlie these assumptions, demonstrating their relevance to various epidemiological outcomes, distinguishing them from strainbased selection, and providing a framework for their application.”
– In the first limitations paragraph, we provide examples of alternative modelling methods commonly used in the literature, and reference alternative factors not considered here (e.g. space, behaviour), but which have been found to drive antibiotic resistance dynamics in other modelling work:
“This study has a number of limitations. First, hospitals and healthcare settings are heterogeneous environments with nonrandom contact patterns and relatively small population sizes. Stochastic, individualbased models accounting for these factors reproduce more realistic nosocomial transmission dynamics than deterministic ODE simulations, allowing for local extinction events, superspreaders, and other inherently random epidemiological phenomena. Nonetheless, our goal was to study how ecological mechanisms impact average epidemiological outcomes in the context of different model assumptions and parameter uncertainty, and in this context ODE modelling was the more appropriate tool, particularly for widely endemic ARB like C. difficile, MRSA and ESBLEC. Still, further insights could certainly be gained by accounting for additional complexity and stochastic heterogeneity in future work, from withinhost spatial organization,(91) to patient and staff contact behaviour,(92) to interinstitutional or interward metapopulation dynamics.(93) These distinctions may be particularly important for rare or nonendemic ARB (e.g. CPKP in some regions).”
To provide further justification for the focus of the work, the introduction has been expanded:
– In intro paragraph 2, we provide a more clear description of relationships between antibiotic consumption, microbiome dysbiosis, and consequences for susceptibility to colonization (foreshadowing model content):
“Microbiome dysbiosis is associated with reduced abundance and diversity of commensal bacteria, impaired host immune responses, and loss of colonization resistance, altogether increasing host susceptibility to ARB colonization.(18–20) Antibioticinduced dysbiosis may further result in elevated expression of antibiotic resistance genes, increased rates of horizontal transfer of such genes, and ecological release, whereby subdominant ARB are released from competition via clearance of drugsensitive bacteria, growing out into dominant colonies.(21–24)”
– In intro paragraph 4, we more clearly define the research gap we have aimed to fill:
“However, contemporary work has stopped short of evaluating consequences of betweenspecies competition on resistance epidemiology. Yet for many ARB, including emerging highpriority multidrugresistant bacteria like extendedspectrum βlactamase (ESBL) producing Enterobacteriaceae, interactions with the host microbiome may be important mediators of nosocomial colonization dynamics.(18,50,51)”
– In line with (1), all reviewers were in agreement that this work needs a clear single model description: a formal description of the model framework in one place. This could be in the supplementary, but it needs to be a well referenced place.
– A single model description is now provided from the beginning of the supplementary appendix (Appendix page 2), and is referenced throughout the main text, including at the end of the first paragraph of the Model and Results section; upon introduction of HGT at the end of part 1; and in first paragraph of the Methods, where we state:
“The final model, comprising pathogen colonization (equation 1), intraspecific strain competition (equation 3), microbiomepathogen competition (equation 4), and horizontal gene transfer (HGT), is given alongside all assumptions in the supplementary appendix (equation S1).”
– This entailed a complete reformulation of the supplementary appendix, which is now presented in a simplified format (see Table of Contents on Appendix page 1) as follows:
– Model and assumptions (including final 12compartment ODE system and describing all underlying modelling choices and assumptions).
– Part 1: Model evaluation and parameterization (including calculations for all outcomes evaluated).
– Part 1: Supplementary results.
– Part 2: Model application, species characterization and parameterization (including expert elicitation and parameter estimates from the literature).
– Part 2: Supplementary results.
1) In general, the model the authors use is highly parameterized. There is a worry about circular reasoning when modelling. Parameters are derived from observations, then a model is constructed to recapitulate this observation. Then similar observations are used to show the model's validity. How did the authors ensure this type of bias was not incorporated into their modelling? This can become a bigger issue with more highly parameterized models, and that is why simpler models are preferable. Another approach is to use crossvalidation, although we are not sure how that would work here.
Abstract
2) The Abstract is vague and somewhat convoluted as written. It did not do the best job of selling the paper and does not convey any specific information.
We agree that simple models are preferable, and that our final model (equation S1) combines a large number of assumptions into a single ODE system (although still fewer than a typical agentbased approach). This was our primary motivation for decomposing the final 12compartment model into individual submodels and presenting these gradually through Part 1, in order to evaluate each component ecological interaction separately and understand its impacts on resistance epidemiology.
– By parameterizing all models identically to the same generic pathogen P^{R}, by conducting identical univariate and bivariate analyses across models, and by evaluating a range of epidemiological outcomes, we believe that we have isolated and demonstrated epidemiological impacts of different ecological mechanisms that we have considered, and justified their inclusion in our final modelling framework (Figures 1, 2, 3)
By “Parameters are derived from observations, then a model is constructed to recapitulate this observation”, we presume that the reviewers refer to parameters like pathogen prevalence and resistance rates, which are simultaneously model input parameters (community burden) and epidemiological outputs (simulated hospital dynamics). We don’t disagree with the interpretation that there is inherently “circularity” here where model inputs drive outputs. However, we stress that such sink dynamics are a truism for healthcare institutions (see e.g. Blanquart et al. 2019 Evol App https://doi.org/10.1111/eva.12753): baseline prevalence and resistance rates for asymptomatic colonizers like MRSA and Enterobacteria inevitably reflect burden in the surrounding community.
– This is also why colonization incidence was evaluated as the primary outcome for evaluation of public health interventions, as this outcome reflects model dynamics but does not have a corresponding input parameter.
– We now better distinguish between interpretation of these input and output parameters in the caption for Figure 5:
“Dashed lines (model inputs) represent assumed community prevalence, i.e. the proportion of patients already colonized upon hospital admission (see Tables S3S6). Solid lines represent simulated prevalence within the hospital, as resulting from both importation from the community and withinhospital epidemiology.”
For part 2, we stress our use of multivariate sensitivity analysis (PRCC) to understand the impacts of parameter uncertainty on model outcomes, and our analysis of “singlespecies” vs. “microbiome” simulations to understand impacts of microbiome parameters on model outcomes
Title [editor request, included here by authors]:
Could you please edit your title so the second part gives clear indication of the study design under investigation, if appropriate, and avoiding abbreviations where possible. The remainder of the title should read as a single, concise sentence containing no dashes or questionanswer, and avoiding abbreviations where possible.
We propose the following title:
“Microbiomepathogen interactions drive epidemiological dynamics of antibiotic resistance: a modelling study applied to nosocomial pathogen control”
Introduction
3) The authors need to tailor this to their specific question: where is the knowledge gap? We recommend that the authors think carefully about the statements they are making and consider revising the Introduction to include a more thorough analysis of prior modelling work, which is largely absent from the Introduction.
A more thorough analysis of previous modelling work is now provided as the fourth paragraph of the introduction. In the context of this paragraph, we present the following knowledge gap.
“Accounting for other forms of complexity in epidemiological models – from treatment intensity, to ageassortative contact behaviour, to hospital referral networks, to animalhuman interactions, to genetic linkage between resistance and nonresistance genes – has helped to unravel the many, disparate forces that contribute to drive the spread of resistance.(42–47) Nevertheless, withinhost competition between cocolonizing bacteria is a key mechanism of selection for antibiotic resistance dissemination, and an active area of research at the forefront of resistance modelling.(34,36,48) For instance, the ‘mixedcarriage’ model by Davies et al. demonstrates how intraspecific competition results in negative frequencydependent selection for either of two competing strains, and provides a satisfying mechanistic explanation for widespread strain coexistence at the population level.(49) However, contemporary work has stopped short of evaluating consequences of betweenspecies competition on resistance epidemiology. Yet for many ARB, including emerging highpriority multidrugresistant bacteria like extendedspectrum βlactamase (ESBL) producing Enterobacteriaceae, interactions with the host microbiome may be important mediators of nosocomial colonization dynamics.”
4) The study focuses on four particular pathogens, however, no motivation of why they are chosen is provided, nor the relevance of their analysis in a more general overview is discussed. The authors could expand on how their results could be extended to other contexts/species in which the role of microbiome on resistance is also relevant. This could also be included in the discussion.
We thank the reviewer for this helpful comment which we have taken into account in several instances:
In the introduction, in addition to previous reference to C. difficile, we now reference the emerging threat of ESBLEnterobacteriaceae when introducing the knowledge gap motivating this work (see previous comment).
In the discussion, second paragraph:
– Additional context for these species as key drivers of nosocomial infection:
“MRSA, C. difficile and ESBLproducing Enterobacteriaceae are leading causes of antibioticresistant and healthcareassociated infection, while carbapenemaseproducing Enterobacteriaceae represent emerging threats of particular concern due to limited therapeutic options for effective treatment of invasive infection.” (4,63–65)
“Antibiotic stewardship is a core component of public health efforts to limit the emergence and spread of these ARB in clinical settings,(66) and an important focus of antibiotic resistance modelling.”
– Relevance for other species and settings:
“Our simulations were limited to a select few ARB, but our framework and findings likely have relevance for other bacteria known to interact with the microbiome, including vancomycinresistant Enterococci and other multidrugresistant Enterobacteriaceae,(23,76) and could be further extended to explore impacts of the microbiome on resistance dynamics and intervention efficacy beyond healthcare settings.”
5) The authors conclude that disruption of the microbiome on an individual level (due to antibiotics) can result in selection for AMR on a population level. In a way, the work is, at least in part, a modelbased implementation of the theoretical perspective of Lipsitch and Samore 2002 EID doi: 10.3201/eid0804.010312 (which describes how antibiotic effects on an individual, withinhost level, can affect population level transmission dynamics). Could the comparison to Lipsitch and Samore and the novelty of the improved datadriven approach (although still limited) be further discussed?
We agree that our work implements a longstanding hypothesis previously described in the perspective by written by Lipsitch and Samore 2002 EID. This is now emphasized first in the introduction:
“Disruption of the host microbiome is a longstanding theory explaining how antibiotics select for the spread of resistance at both the individual and population levels,(34) but most epidemiological models consider just one species of bacteria at a time, under the traditional assumption that antibiotic selection for resistance results from intraspecific competition between cocirculating strains.”
In the Model and Results, as previously:
“….intraspecific strain competition is not the only mechanism by which antibiotic consumption can drive ARB spread,(34) and may have limited relevance for certain species, settings, timescales and epidemiological indicators.”
And anew in the first paragraph of the discussion
“We present a modelling framework that includes withinhost ecological costs of antibiotic use in the form of microbiome dysbiosis, incorporating a leading hypothesis for antibiotic selection into classical models of resistance epidemiology.”
We also highlight the following passage from Lipsitch and Samore, which supports the present work:
“The use of mathematical models, and more generally the attempt to predict the relative merits of different interventions, will depend on an improved understanding of the mechanisms of antibiotic selection in particular organisms.”
The novelty of our approach is further emphasized in the discussion, when comparing our approach to previous models assuming increased susceptibility to colonization/infection among antibioticexposed individuals:
“These strategies reflect widespread recognition that antibiotic use favours ARB acquisition, through erosion of colonization resistance or other supposed mechanisms, and independent of potential competition with other strains. The present work formalizes examples of the microbiomepathogen interactions that underlie these assumptions, demonstrating their relevance to various epidemiological outcomes, distinguishing them from strainbased selection, and providing a framework for their application.”
In the second limitations paragraph of the discussion, we now highlight that, despite clear data limitations, our expert elicitationdriven approach provided useful proxy measures that allowed us to characterize the speciesspecific ecology of different ARB:
“Finally, the nature of microbiomepathogen interactions and their epidemiological consequences remain poorly understood and largely unquantified. We show in theory why these interactions matter at the population level, but empirical data were unavailable to inform their parameterization. Instead, we translated microbiomepathogen competition coefficients into clinical parameters, designed a structured expert elicitation protocol and conducted interviews allowing subjectmatter experts to quantify their beliefs. Although these estimates are subject to substantial bias and uncertainty, they facilitated speciesspecific characterization of the epidemiological impact of microbiome dysbiosis, and represent useful proxy measures in the absence of clinical data.”
And in the penultimate paragraph we further discuss our approach and findings in the context of limited data, including perspectives for future empirical work:
“Despite data limitations, epidemiological conclusions from Monte Carlo simulations were largely consistent with empirical findings (discussed above), suggesting that final parameter distributions were reasonable approximations. Uncertainty in parameter inputs translated to uncertainty in model outputs and reflects the knowledge gaps underlying our simulations (Figure S14). This is exemplified by HGT and its highly uncertain role in driving colonization incidence; to date, HGT modelling has largely been limited to withinhost dynamics, and impacts on epidemiological dynamics are only just beginning to come to light.(51,57,94) Increasing availability and synthesis of highquality withinhost microbiological data will help to further characterize epidemiological impacts of microbiomepathogen interactions. Studies are needed that describe ecological impacts of antibiotic exposure on microbiome population structure across control and treatment groups, with longitudinal followup evaluating subsequent nosocomial ARB colonization risk. In the absence of clinical data, insights from experiments and withinhost models nonetheless suggest that antibiotic disruption of microbiomepathogen competition is a key driver of selection for resistance.” (95–100).
6) Antibiotic stewardship is treated a single entity, whereas in reality there are countless types of stewardship interventions and a concomitant variety of outcomes. Comparing the model's results directly to an analysis of previous stewardship interventions and attempting to show alignment does not make sense to me. How were these stewardship interventions reflective of the model? Why would we expect agreement or disagreement a priori? The fact that there was agreement with the model was somewhat concerning given the complexity of many stewardship interventions.
– We totally agree that stewardship interventions are highly heterogeneous and thank the reviewer for raising this point. Heterogeneity likely drives high uncertainty in efficacy estimates from the metaanalysis by Baur et al., which groups interventions together as a single entity. Despite this high heterogeneity, stewardship interventions – whether reducing rates of antibiotic consumption, favouring lowerspectrum antibiotics, decreasing treatment duration, or otherwise limiting patient exposure to unnecessary antibiotics – should in theory all reduce populationlevel antibioticinduced dysbiosis (and hence reduce microbiomeinduced selection for resistance)
– Owing to our ODEbased approach, we were unable to reproduce highly granular stewardship interventions that reflect realworld heterogeneity, but we nonetheless simulated 3 types of stewardship, and found similar efficacy estimates (1 – IRR) for each (Appendix figure 17), which were also concordant with pooled estimates from the metaanalysis
– We agree that visualizing validation of model results to published results from Baur et al. was perhaps not appropriate, and have removed those estimates from Figure 6. Instead, we make comparison to their findings in the discussion:
“By contrast, simulated antibiotic stewardship interventions were broadly effective for reducing incidence across all included ARB. This is consistent with findings from a metaanalysis of clinical trials evaluating the efficacy of hospital antibiotic stewardship interventions for reducing incidence of ARB colonization and infection,(66) which we updated to exclude studies coimplementing stewardship with alternative interventions (Figure S18). In comparison to our findings, estimates from the metaanalysis were associated with greater uncertainty across more heterogeneous interventions, but predicted the same rank order of efficacy across included ARB, and similar mean efficacy for C. difficile (19% from n=7 studies, vs. 18% under intermediate compliance in our simulations), ESBLEC (18% from n=4 studies, vs. 15%) and MRSA (12% from n=10 studies, vs. 10%), but higher efficacy for CPKP (54% from n=1 study, vs. 20%).”
Methods7 In attempting to replicate this study, it became intractably difficult to recreate the ODEs because the equations and variables are interspersed throughout the main text and supplementary appendix. For example, in Figure 1e the authors mention differences in resistance levels (γ) and not transmission rate (λ, related to β as in Equation (3)), the value of which could not be found anywhere in the text. Similarly, where is patient demography (Δ) incorporated into equations 3 or 4? There are certainly answers to these questions, but it would have been helpful to have the model more clearly presented in a central location. Finding the parameters and making sense of them was extremely challenging, if not infeasible. Therefore, we were unable to replicate the results to any extent.
We appreciate that some aspects of our model were unclear to reviewers and regret any confusion caused. We have taken a number of steps to improve the clarity of our modelling assumptions and to improve reproducibility of our work:
– We provide equations for each of the analyzed submodels in the manuscript (the strainmicrobiome competition was previously in the supplement), more explicitly reference them in the text (equations 1, 3, 4 and 6, respectively), and include δ (demography) terms throughout
The final complete model and all of its underlying assumptions – the ‘framework’ – is now provided at the beginning of the supplement, delineating all assumptions, composite parameters and strainspecific differences. This final framework is referenced on multiple occasions in the manuscript, including:
– In the first paragraph of Methods and Results:
“See Methods for technical details, and the supplementary appendix for the complete modelling framework and assumptions (equation S1).”
– After introduction of HGT, to clarify where the final model including HGT can be found:
“(See final model and assumptions, Appendix equations 1 – 11.)”
– At the beginning of the methods, to clarify our approach of building subsequent submodels together to produce our final framework:
“We evaluated ODE systems describing colonization dynamics of bacterial pathogens in the healthcare setting. The final model, comprising pathogen colonization (equation 1), intraspecific strain competition (equation 3), microbiomepathogen competition (equation 4), and horizontal gene transfer (HGT), is given alongside all assumptions in the supplementary appendix (equation S1).”
In addition to the previously available Mathematica code, R code has been made available for the 5 ODE systems at https://github.com/drmsmith/microbiomeR, as referenced in the Methods:
“(See the supplementary appendix for technical details, and R and Mathematica files available online at https://github.com/drmsmith/microbiomeR.)”
In the caption for Figure 1E we have defined strainspecific differences in force of infection, as suggested:
“For E, P^{S} and P^{R} circulate simultaneously, assuming strainspecific differences in antibiotic resistance (r_{S} = 0, r_{R} = 0.8), natural clearance (γ_{S} = 0.03 day^{1}, γ_{R} = 0.06 day^{1}) and transmission (λ_{S} = β × C^{S}/N, λ_{R} = β × C^{R}/N).”
8) It is unclear how the authors arrived at their parameter values throughout the manuscript. For example, colonization resistance (epsilon) as it is formulated is not compelling, and where did the value of this parameter come from? It is listed in Table S5 as being drawn from a Cauchy distribution. Why? Where did these numbers come from? Some of the parameters seem justified by literature, but the rest are perhaps debatable. This might be inevitable with a complex model, but the authors should minimally provide a robustness analysis to each of their lesssupported parameter choices. What would be the outcome if the sampled distributions had been wider (e.g., if the true value had been different by an order of magnitude)?
We thank the reviewer for raising this important point and have taken a number of steps to address it. First, by explicitly splitting the work into 2 parts, we hope that it is now more clear that in Part 1 we evaluate impacts of different withinhost ecological interactions on a generic nosocomial pathogen P^{R} – hence with arbitrary parameter values, but using univariate and bivariate analysis to explore the epidemiology of this theoretical pathogen across wide parameter space – and that in Part 2 we tailored our model to particular ARB using speciesspecific parameter estimates from the literature and expert opinion
For consistency, we have introduced concepts of colonization resistance and ecological release – both wellestablished in cited literature – earlier in the manuscript:
“The microbiome can also protect against colonization with infectious bacterial pathogens, a phenomenon known as colonization resistance, limiting their capacities to establish colonies, grow, persist and transmit.”
“Antibioticinduced dysbiosis may further result in elevated expression of antibiotic resistance genes, increased rates of horizontal transfer of resistance genes, and ecological release, whereby subdominant ARB are released from competition via clearance of drugsensitive bacteria, growing out into dominant colonies.”
We now also discuss limitations of our conceptualization of these interactions in our modelling framework:
“More broadly, our characterizations of microbiomepathogen interactions are conceptual, and were mapped mechanistically to particular colonization processes (transmission, clearance, endogenous acquisition), but we note that in other contexts terms like colonization resistance, resource competition and ecological release may map to specific biochemical processes that could affect epidemiological parameters in different ways.”
We have taken a number of steps to more clearly describe parameterization of microbiomepathogen interaction coefficients and our use of expert elicitation:
– The protocol document for our expert elicitation exercise is provided as a separate supplementary file, an exact copy of which was used by experts during interviews.
– Use of expert elicitation to inform model parameters is now referenced as early as the introduction:
“Expert elicitation interviews were conducted to characterize the clinical relevance of microbiome dysbiosis for each species, and to qualify and quantify interaction coefficients with uncertainty.”
– In the Model and Results section of the manuscript, we more clearly reference the expert elicitation and what it was used for, with reference to Methods and relevant sections of the supplement:
“Data from the literature were used to parameterize the model to each pathogen (Figure 4, Tables S2 – S6). Literature estimates for microbiomepathogen interaction coefficients are scarce, and the speciesspecific relevance of different withinhost interactions (intraspecific strain competition, microbiome competition, HGT) in the hospital environment are not welldefined. To characterize ecological interactions for each species and inform model structure, we conducted interviews with a panel of subjectmatter experts in medical microbiology and antibiotic resistance epidemiology (details in Methods). Based on their beliefs, all pathogens were assumed to compete with microbiota; MRSA, ESBLEC and CPKP were further assumed to compete intraspecifically with nonfocal strains, for simplicity characterized as methicillinsensitive S. aureus (MSSA), E. coli (EC) and K. pneumoniae (KP); and both ESBL resistance and carbapenem resistance were assumed to be borne by plasmids capable of horizontal transfer between patient microbiota and, respectively, EC and KP (Figure 4A). To quantify speciesspecific strengths of microbiomepathogen interactions, withinhost ecological parameters were translated into clinical parameters (Table S7), and experts were asked to quantify these using standardized expert elicitation methodology (Figures 4B, S8 – S10).”
– In the Methods, description of the expert elicitation has been condensed into one succinct paragraph for clarity, while full methodological details have been provided in the corresponding section of the supplement (Appendix section 4.2, page 21), where we:
– Describe expert elicitation and our protocol.
– Describe elicitation results and their interpretation. Here we also detail how final pooled distributions were found (including the Cauchy distribution listed in Table S5 and referenced by the reviewer):
“Raw data for expert distributions are provided in Figure S9. Experts differed substantially in estimated ranges of different parameters, but Friedman’s tests using medians of each parameter distribution suggested that the species rank order was conserved across experts, with C. difficile generally having the strongest estimated microbiomepathogen interaction coefficients, and MRSA the weakest. To conserve species rank order when combining expert estimates x_{i,j,k} to form pooled histograms, distributions were recentred by the mean relative distance z_{i,j} between each j and the reference ARB (taken as MRSA) over all experts k (Figure S10). For C. difficile and MRSA, HGT was excluded from simulations and pooled distributions were not recentred. Each expert distribution was weighted equally, contributing to 10% of the final pooled distribution for each species and parameter. Pooled histograms were fit to six candidate distributions (normal, lognormal, Weibull, Cauchy, exponential and γ), and the final distribution was selected as the fitted distribution giving the lowest Akaike Information Criterion (using the function fitdist from the R package fitdistrplus). Final fitted distributions for each parameter and ARB are given in Tables S3S6, and for select parameters in Figure 4.”
– Provide Table S7, as previously, which translates the clinical parameters experts estimated into microbiomepathogen interaction coefficients.
– Provide Figure S8, as previously, showing expert opinion about the clinical relevance of different clinical parameters for different species.
– Provide Figure S9, the raw data resulting from the expert elicitation, i.e. the distributions created by each expert for each microbiomepathogen interaction for each species.
– Provide Figure S10, final recentred data used to generate pooled distributions while maintaining species rank order.
We agree that alternative model parameterization for each species would change modelling results and conclusions, but stress that the goals of Part 2 – evaluation of speciesspecific colonization dynamics and intervention efficacy – necessitates speciesspecific parameterization, which we estimated to the best of our ability given available data, while:
– accounting for parameter uncertainty;
– conducting multivariate sensitivity analysis with PRCC to quantify impacts of parameter uncertainty on model outcomes;
– conducting simulations without microbiome interactions, to show their impact on outcomes.
A key result of these analyses is that alternative speciesspecific model parameterization leads to different modelling results and estimates for intervention efficacy, and that these are further impacted by including microbiomepathogen interactions. We believe that additional analysis arbitrarily adjusting parameter distributions would not be appropriate:
– Widening distributions for all parameters should produce approximately the same median results with inflated uncertainty across outcomes.
– Univariate widening of parameter distributions is unnecessary given PRCC results (Appendix figure 14, page 30), which already convey which parameters have greatest impact on outcomes, and hence which parameters should most increase or decrease resistance if adjusting the value.
– Arbitrarily adjusting certain parameters goes against the stated goal of speciesspecific characterization. For example, it would be possible to run simulations with a high transmission rate for E. coli or low microbiomepathogen interaction coefficients for C. difficile, and outcomes would certainly change substantially. But these simulations would no longer represent these pathogens, and it is not clear to us that worthwhile insights would be gained that are not already demonstrated by PRCC results, by comparisons between the four included ARB, and by ‘singlespecies’ vs. ‘microbiome’ simulations
We nonetheless agree that imperfect parameterization is a key limitation to this study, and dedicate a new paragraph to discuss this, highlighting limited availability of several parameters in the literature, uncertain generalizability of results across settings, sensitivity of model conclusions to assumed parameters, the conceptual nature of our characterization of microbiomepathogen interactions, and the gutspecific nature of our microbiome dysbiosis scale (which could potentially bias results for predominantly skinresiding S. aureus).
Finally, for completeness, we reran our multivariate sensitivity analysis conducted previously (Latin Hypercube Sampling + PRCC) for an additional outcome: the pathogen resistance rate (Appendix figure 14B, page 30). Findings for both outcomes are now presented together in the Model and Results section:
“In multivariate sensitivity analysis, community prevalence (f_{C} for C. difficile, f_{R} for others) and rates of endogenous acquisition (α_{R}) had overall the strongest positive impacts (highest PRCCs) on hospital prevalence across ARB, while rates of hospital admission/discharge (μ) and microbiome recovery (δ) had the strongest negative impacts (lowest PRCCs; Figure S14A). For resistance rates, parameters with the strongest positive impacts (highest PRCCs) were community prevalence (f_{R}), the rate of endogenous acquisition (α_{R}) and the rate of antibioticinduced pathogen clearance (σ_{C}). Parameters with the strongest negative impacts (lowest PRCCs) were rates of hospital admission/discharge (μ) and endogenous acquisition of competing drugsensitive strains (α_{S}) (Figure S14B). Across ARB, prevalence estimates, but not resistance rates, were generally sensitive to microbiome parameters.”
9) The inclusion of the 'r' parameter, allowing for partial resistance in the resistant strain, is interesting (and innovative for modelling frameworks to our knowledge). Nonetheless, it is not clear why the authors use a baseline of r = 0.8 to illustrate the tradeoff between antibiotic induced clearance and selection for Cr (Figure 1F), whereas in figure 2, the authors use r=0.4, i.e. higher levels of sensitivity in Cr. I suppose with a lower r, strain coexistence is more likely, and thus higher likelihood of Horizontal Gene Transfer (HGT) (?), but it would be, for consistency and comparison of the different withinhost dynamics, more transparent to use the same baseline parameter values, unless the authors can provide a good reason why different values of r are assumed at baseline between Figure 1 and Figure 2B?
We thank you for raising this and agree that identical parameter values across plots for Part 1 makes for more consistent and transparent results. We have now reproduced all Figures using the same default parameter set (Appendix table 1 on page 11, results also provided in R code), such that all results for Part 1 now represent the same generic pathogen P^{R} (adjusting only the parameters explicitly made to vary for univariate or bivariate analysis in each figure). In particular, the previous Figure 2B has been expanded into a new Figure 3, which includes the default value of r_{R}=0.8, as well as results for lesser resistance levels (r_{R} = 0.2, 0.5) to demonstrate how different epidemiological outcomes (pathogen prevalence and resistance rate) vary with this important parameter.
10) The model presented in Eq. (1) needs a better introduction. Is this a new model or based on other models previously studied? This is explained later in the text citing references (30) and (31), but we recommend to introduce them earlier. It would be also good to point the reader earlier in the main text to the first sections of the supplement for a better understanding of the model assumptions and parameters employed.
This model is adapted from classic colonization models of antibioticresistant bacteria – the first example to our knowledge being Austin and Anderson Proc Biol Sci 1997 – with the key difference being our novel approach to modelling antibiotic treatment. We now introduce this model in this context, including the following:
“We start with a SusceptibleColonized transmission model (Figure 1A) representing a population of N hospitalized patients as either susceptible to colonization (S) or colonized (C^{R}) by P^{R}, the focal strain or species: [equations] This model is adapted from classic colonization models of antibioticresistant bacteria,(51) includes no ecological interactions with nonfocal bacteria, and reflects a suite of common assumptions relevant to the healthcare setting, including: …”
And subsequently highlighting how we model antibiotic exposure and its impacts on pathogen clearance:
“Modelling the latter as a continuous proportion reflects that bacteria are not necessarily fully drugsensitive (r_{R} = 0) nor resistant (r_{R} = 1), but can range in their sensitivity to different antibiotics (0 ≤ r_{R} ≤ 1). The resistance level r_{R} is thus a model input interpreted as an overall measure of the pathogen’s innate degree of resistance to the particular antibiotics to which it is exposed.”
Finally, we place this first model in the context of subsequent models
“Following models build upon these assumptions, representing the same pathogen P^{R}, but altering its ecological interactions with other bacteria from one model to the next. Models were evaluated over the same generic parameter space, to isolate impacts of model structure on epidemiological outcomes in the context of antibiotic use (see parameters in Table S1).”
11) Also, should the parameters be more consistently labelled across frameworks? For example, we recommend writing λ(N,C), or something similar, in all the equations to explicitly state it depends on these parameters. The way it is written gives the impression that λ is a constant parameter.
13) Furthermore, in figure 1F: this is representing a one strain model. Is the Ce+Cd strain similar to the Cr? Can this be clarified?
14) In line with this. In Figure S3, R0 seems to represent the R0 of Cr. However, for Figure 1H and 1I, how should R0 be interpreted respectively? And for Figure S2? Please clarify how similar Cr (Figure 1H) and the one strain modelled (which could be similar to the Cr strain of Figure 1H) for figure 1I are
We respond to these 3 comments simultaneously: yes, the focal antibioticresistant pathogen (previously denoted C in model 1, C^{R} in model 2 and C_{e}+C_{d} in model 3) is identical, and we understand the previous confusion. We thank the reviewer for raising these points and we have consequently made changes to clarify notation
We have rewritten the manuscript to explicitly distinguish between the pathogen in question and patient colonization with that pathogen: drugsensitive and resistant pathogen strains are denoted P^{S} and P^{R}, respectively, while patient colonization with those strains is denoted C^{S} and C^{R}
All ODE systems are rewritten to describe dynamics of C^{R}, i.e. all systems describe colonization with the same focal antibioticresistant pathogen P^{R}. Epidemiological parameters for P^{R} are always denoted by subscript R (e.g. α_{R}). Where strain competition is also included, ODE systems include C^{S} (representing patient colonization with the strain P^{S}), with epidemiological parameters for P^{S} denoted by subscript S (e.g. α_{S})
We now make the reader aware of these notations several times early in the paper. For instance we:
– Introduce P^{R} from the first sentence of the Model and Results:
“We propose a series of five models describing colonization dynamics of a focal antibioticresistant bacterial pathogen, denoted P^{R}, among hospital inpatients in an acute care setting. Each model accounts for different withinhost ecological interactions between P^{R} and other bacteria.”
– Distinguish between P^{R} and C^{R} upon introduction of the first model:
“We start with a SusceptibleColonized transmission model (Figure 1A) representing a population of N hospitalized patients as either susceptible to colonization (S) or colonized (C^{R}) by P^{R}, the focal strain or species:”
– Distinguish between P^{S} and P^{R} upon introduction of strain competition:
“… selection results from intraspecific competition between two or more drugsensitive strains P^{S} and drugresistant strains P^{R}.”
– And again upon introduction of the strain competition model:
“A simple twostrain ‘exclusive colonization’ model (Figure 1B) can be written as: [equations] where patients can be colonized (C^{S}, C^{R}) by either strain (P^{S}, P^{R}).”
– And so on throughout the manuscript.
We further make explicit that we are always evaluating the same P^{R} across models using the same parameter set, changing only its withinhost interactions with other bacteria across models:
“Following models build upon these assumptions, representing the same pathogen P^{R}, but altering its ecological interactions with other bacteria from one model to the next. Models were evaluated over the same generic parameter space, to isolate impacts of model structure on epidemiological outcomes in the context of antibiotic use (see parameters in Table S1).”
We have updated model diagrams in Figure 1 to represent these changes
We understand that not all readers may be familiar with the classic force of infection term λ, but in the final model λ already has up to 2 subscripts (λ_{S,ε}). We believe that adding more (e.g. λ_{S,ε,C,N}) could unnecessarily complexify model interpretation. Instead, we have proposed clear definitions of λ in 3 instances:
– as previously, upon introduction of model 1;
“(iii) a dynamic rate of colonization acquisition λ_{R}=β⨉C^{R}/N, for hosttohost transmission.”
– as suggested by reviewers, in the caption for Figure 1 for the strain competition model:
“For E, P^{S} and P^{R} circulate simultaneously, assuming strainspecific differences in antibiotic resistance (r_{S} = 0, r_{R} = 0.8), natural clearance (γ_{S} = 0.03 day^{1}, γ_{R} = 0.06 day^{1}) and transmission (λ_{S} = β × C^{S}/N, λ_{R} = β × C^{R}/N).”
– in the supplement, as the final generalized expression for any strain j. (Appendix equation 3 page 3)
12) The r vs resistance rate (Cr/Cr+Cs) caused some confusion. This as a resistance rate of Cr = 0.8 but an r = 0.4 could actually mean that 0.8*0.4 = 0.08 of pathogens carried are fully resistant against the antibiotic, while this is 0.64 when r=0.8. Can the same baseline values for r be used in Figure 1 and Figure 2 or the differences justified?
We regret the confusion. The same baseline value r_{R}=0.8 is now used throughout.
To clarify interpretation of the parameter r_{R} (and r_{S}), the following text has been added upon its introduction in the first model:
“…and r_{R} is the antibiotic resistance level (the proportion of antibiotics that are ineffective against P^{R}). Modelling the latter as a continuous proportion reflects that bacteria are not necessarily fully drugsensitive (r_{R} = 0) nor resistant (r_{R} = 1), but can range in their sensitivity to different antibiotics (0 ≤ r_{R} ≤ 1). The resistance level r_{R} is thus a model input interpreted as an overall measure of the pathogen’s innate degree of resistance to the particular antibiotics to which it is exposed.”
And we have also more clearly defined the resistance rate in the first paragraph of Model and Results:
“Across models, three primary epidemiological outcomes are calculated at steadystate equilibrium: P^{R} prevalence (the proportion of patients colonized), P^{R} incidence (the daily rate of colonization acquisition within the hospital), and the pathogen resistance rate (the proportion of patients colonized with the focal antibioticresistant strain P^{R} relative to a competing drugsensitive strain P^{S}).”
We have kept the initial name used for the r_{R} parameter: the antibiotic resistance level. This name was chosen to avoid confusion with alternative possible names like the antibiotic resistance rate/proportion/frequency/share (which are used elsewhere to describe the share of strains of a species – or the share of clinical isolates – that are drugresistant, as above). However, we understand that interpretation of this novel parameter may not be obvious and we propose the following alternatives to the reviewers and editors: antibiotic resistance coverage, antibiotic resistibility, antibiotic resistance degree, innate antibiotic resistance level/degree.
– It is not clear to us that any of these alternatives improve understanding, but our goal is to maximize clarity and interpretability of this important parameter, and we welcome any feedback and are happy to update it, if so desired by the editorial team.
15) Figure 2A, the minimum resistance rate is 0.35 (see legend, and this appeared the case when no antibiotic induced clearance nor antibiotic induced microbiome disruption). This seems rather high. In particular as in Figure 2B, minimum levels of less than 0.1 are shown. Could the authors explain where this discrepancy is coming from and to what extend these high baseline resistance rates are representative?
For ease of interpretation for the generic pathogen P^{R} in part 1, in our baseline parameter set (Table S1) we now assume that 50% of colonized individuals admitted from the community are colonized with the resistant strain (f_{R} = 0.5): when the steadystate resistance rate exceeds 50%, it can easily be understood that the hospital has provided a competitive advantage to the resistant strain despite assumed costs of resistance. We agree that this is a higher resistance rate than many realworld pathogens, examples of which are provided in Part 2 (see for e.g. Figure 5B).
In the previous Figure 2B, minimum levels of C^{R} prevalence were <0.1 (not the resistance rate, which was not shown), but we acknowledge that this may have been confusing as the adjacent Figure 2A presented both resistance rate and prevalence. We have now separated these figures: Figure 2A is presented as Figure 2; and Figure 2B has been expanded into a new Figure 3, which includes both prevalence and resistance rates to avoid potential misinterpretation, and to demonstrate asymmetry in these outcomes. We further demonstrate how outcomes related to HGT depend on the antibiotic resistance level r_{R} (columns). We have also added translucent grey bars to this figure representing assumed community P^{R} colonization prevalence (f_{C} × f_{R}) and resistance rate (f_{R}), to help the reader to contextualize results against baseline assumptions.
16) For parameterisation of the models for the different pathogens, estimates from literature, notably existing modelling studies are chosen. However, these estimates, at least for C. difficile and MRSA, are coming from models that don't incorporate endogenous acquisition explicitly (for the Enterobacterieacia, model estimates from Gurieva are used, which do use a modelling framework incorporating an exo and endogenous acquisition). Therefore, the acquisition rates may be overestimated for these grampositives, in particular for C. diff the endogenous acquisition, which is listed as the main acquisition route (pp 14 line 330).
This may affect the estimated intervention effectiveness (in particular for antibiotic stewardship (less effective) and contact precautions (more effective)) under assumptions of the microbiome model. Could the authors use more realistic estimates (not sure models with explicit exo and endogenous acquisition exist for MRSA and C. diff), or at least, reflect on how different values of α and β affect the model results?
Unfortunately, we are also unaware of models with explicit exoendo estimates for the grampositives, but we now dedicate substantial discussion to this limitation and discuss relevance to PRCC results from our sensitivity analysis:
“Third, Monte Carlo simulations were limited by the availability of speciesspecific model parameters from the literature, in some instances necessitating use of previous modelling results, approximations, or estimates from small studies in specific locations, making the generalizability of results unclear. For instance, Khader et al. estimated a fourfold difference in MRSA transmission rates between hospitals and nursing homes.(79) Such differences could have a substantial impact on dynamics and estimated intervention efficacy, with higher transmission rates favouring use of contact precautions, and higher rates of endogenous acquisition favouring antibiotic stewardship (in the context of a high ecological release coefficient). Uncertainty in endogenous acquisition rates may be particularly important: in multivariate sensitivity analyses, this parameter emerged as a key driver of both colonization prevalence and resistance rates across ARB (Figure S14).”
17) It would be helpful to explicitly say how the simulations were performed, i.e., mention that the ODEs where integrated. My first impression was that the authors performed stochastic simulations using the processes illustrated in panels A, B, C from Figure 1. This is probably an issue that only mathematical modellers would have, but could the authors please add this clarification.
We agree that more details on the modelling and simulations were missing and thank the reviewer for raising that point. This is now made explicit in the new introductory paragraph to the Methods and Results section:
“Models are described using systems of ordinary differential equations (ODEs) and are evaluated deterministically by numerical integration. Across models, three primary epidemiological outcomes are calculated at steadystate equilibrium: P^{R} prevalence (the proportion of patients colonized), P^{R} incidence (the daily rate of colonization acquisition within the hospital), and the pathogen resistance rate (the proportion of patients colonized with the antibioticresistant strain P^{R} relative to a drugsensitive strain P^{S}). We also derive and evaluate the basic reproduction number R_{0}, an indicator of pathogen epidemicity representing the average number of patients expected to acquire a novel pathogen from an initial index patient. See Methods for technical details, and the supplementary appendix for the complete modelling framework and assumptions (equation S1).”
In the caption for Figure 1:
“For all models, ODEs are integrated numerically using the same parameter values representing a generic nosocomial pathogen P^{R} (see Table S1).”
And again, in the methods:
“ODEs were integrated numerically to calculate steadystate epidemiological outcomes for nosocomial P^{R} colonization: colonization prevalence (the sum of all compartments C^{R}), colonization incidence (the daily rate of C^{R} acquisition), and the resistance rate (C^{R}/(C^{S} + C^{R})). (See the supplementary appendix for technical details, and R and Mathematica files available online at https://github.com/drmsmith/microbiomeR.) For each model, outcomes were evaluated over the same parameter space representing a generic pathogen P^{R} (parameters in Table S1), while varying specific parameters through univariate and bivariate analysis to assess their impacts on dynamic equilibria in the context of different modelling assumptions.”
Discussion
18) The model is described by a set of differential equations. This means the model ignores stochasticity of the population size dynamics of each strain. The authors could expand a bit more on the assumptions of the model and why fluctuations are ignored.
We have expanded upon these limitations as suggested, nuancing our statements with the intended goals of this study:
“This study has a number of limitations. First, hospitals and healthcare settings are heterogeneous environments with nonrandom contact patterns and relatively small population sizes. Stochastic, individualbased models accounting for these factors reproduce more realistic nosocomial transmission dynamics than deterministic ODE simulations, allowing for local extinction events, superspreaders, and other inherently random epidemiological phenomena. Nonetheless, our goal was to study how ecological mechanisms impact average epidemiological outcomes in the context of different model assumptions and parameter uncertainty, and in this context ODE modelling was the more appropriate tool, particularly for widely endemic ARB like C. difficile, MRSA and ESBLEC. Still, further insights could certainly be gained by accounting for additional complexity and stochastic heterogeneity in future work, from withinhost spatial organization,(91) to patient and staff contact behaviour,(92) to interinstitutional or interward metapopulation dynamics.(93) These distinctions may be particularly important for rare or nonendemic ARB (e.g. CPKP in some regions).”
19) Also, spatial structure is particularly important when taking into account interactions between microbiome and pathogens, but it is also neglected in the model. It would be useful if the authors could discuss this as well.
Importance of spatial structure to resistance dynamics is now referenced in the introduction:
“Accounting for additional forms of complexity in epidemiological models – from treatment intensity, to ageassortative contact behaviour, to hospital referral networks, to animalhuman interactions, to genetic linkage between resistance and nonresistance genes – has helped to unravel the many, disparate forces that contribute to drive the spread of resistance.”
and in the discussion:
“Still, further insights could certainly be gained by accounting for additional complexity and stochastic heterogeneity in future work, from withinhost spatial organization,(91) to patient and staff contact behaviour,(92) to interinstitutional or interward metapopulation dynamics.(93) These distinctions may be particularly important for rare or nonendemic ARB (e.g. CPKP in some regions).”
20) The results shown often focus on steadystate quantities (such as colonisation prevalence), but the temporal evolution of the model is not discussed. It would be interesting if the authors could investigate the timescale of the different interactions taken into account, and how relevant they may be in the healthcare settings they investigate. Is the temporal evolution of the model any relevant for their conclusions?
We include a new supplementary figure which evaluates P^{R} prevalence and resistance rate over time for the generic pathogen P^{R} evaluated in part 1, separating and then combining each microbiomepathogen interaction to assess their impacts on dynamics (Appendix Figure 5). We included three public health interventions, implemented at days 90, 180 and 270, to explore how pathogen responses to interventions vary with time. Associated code is provided in the R file. We include the following text in the main text with regards to these findings:
“Different microbiomepathogen interactions also underlie distinct dynamic responses to theoretical public health interventions. For the same generic pathogen P^{R}, antibiotic stewardship interventions generally, but do not always prevent colonization, with predictions depending on the ecological interactions in effect (e.g. strain competition, microbiome competition), the impact of the intervention (e.g. reduced microbiome disruption, reduced overall prescribing), and the epidemiological outcome considered (e.g. colonization prevalence, resistance rate) (Figure S5).”
As well as the following text in the supplement:
“In Figure S5, we simulate dynamic responses of P^{R} to public health interventions. First, dynamic equilibria were found, from which ODEs were integrated for an additional 365 days, evaluating C^{R} and the resistance rate over time in the context of each microbiomepathogen interaction (separately and in concert). We introduced theoretical interventions at days 90 (halving r_{R} from 0.8 to 0.4), 180 (halving θ_{m} from 1.0 to 0.5), and 270 (halving a from 0.2 to 0.1), demonstrating that an otherwise identical pathogen can experience diverse, and sometimes opposing responses to public health interventions in the context of different microbiomepathogen interactions and epidemiological outcomes.”
Systems generally returned to equilibria after 1 to 3 months. However, we avoid overinterpretation of timescales and lags in the manuscript, owing to limitations to using deterministic ODE modelling for predictive modelling in the hospital setting (see comments 18 and 19).
21) Can the authors please add examples of the hypothetical microbiome recovery interventions they have in mind? Are the authors thinking of things like faecal transplantation? Please add, also to provide more practical interpretation of the work.
We now include a statement in the introduction more clearly referencing these interventions:
– From a clinical perspective, this motivates a need for public health interventions that minimize or reverse harm to patient microbiota, from antibiotic stewardship, to faecal microbiota transplantation, to microbiome protective therapies. (28,29).
And devote a new discussion paragraph to discuss these interventions and relevance for our modelling framework and findings.
– Findings also suggest promise for interventions that effectively restore microbiome stability and associated colonization resistance as a means to control ARB spread. Fecal microbiota transplantation is already used to treat recurrent C. difficile infection, and is under investigation for multidrugresistant Enterobacteriaceae decolonization.(67–69) However, its appropriateness for dysbiosis recovery in the absence of other clinical indications is unclear. Transplantation requires rigorous donor screening and close longitudinal followup, and cases of donor stool contaminated with toxicogenic and multidrugresistant bacteria highlight nonnegligible risks.(70,71) Alternative microbiome protective therapies now exist, like DAV132, a novel activatedcharcoal product currently undergoing clinical trials. When coadministered with antibiotics by the oral route, DAV132 has been shown to absorb antibiotic residues in the colon and preserve the richness and composition of intestinal microbiota, while maintaining systemic antibiotic exposure.(28,72,73) Modelling has been used previously to evaluate impacts of such microbiomeoriented interventions at the withinhost level,(74,75) but knockon impacts on ARB transmission dynamics and epidemiological burden have not been evaluated previously. Our simulations were limited to a select few ARB, but our framework and findings likely have relevance for other bacteria known to interact with the microbiome, including vancomycinresistant Enterococci and other multidrugresistant Enterobacteriaceae,(23,76) and could be further extended to explore impacts of the microbiome on resistance dynamics and intervention efficacy beyond healthcare settings.
https://doi.org/10.7554/eLife.68764.sa2Article and author information
Author details
Funding
Agence Nationale de la Recherche (SPHINX17CE36000801)
 David RM Smith
 Laura Temime
 Lulla Opatowski
Canadian Institutes of Health Research (164263)
 David RM Smith
Agence Nationale de la Recherche (ANR10LABX62IBEID)
 Lulla Opatowski
The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.
Acknowledgements
We express our gratitude to the subjectmatter experts who participated in our expert elicitation exercise: Antoine Andremont, Christian BrunBuisson, Aurélien Dinh, Stephan Harbarth, JeanLouis Herrmann, Solen Kernéis, Alban Le Monnier, Benoît Pilmis, Etienne Ruppé and PaulLouis Woerther. The work was supported directly by internal resources from the French National Institute for Health and Medical Research, the Institut Pasteur, the Conservatoire National des Arts et Métiers, and the University of Versailles–SaintQuentinenYvelines / University of ParisSaclay. This study received funding from the French Government’s ‘Investissement d’Avenir’ program, Laboratoire d’Excellence ‘Integrative Biology of Emerging Infectious Diseases’ (Grant ANR10LABX62 IBEID). DS is supported by a Canadian Institutes of Health Research Doctoral Foreign Study Award (Funding Reference Number 164263) and all authors are supported by the French government through its National Research Agency project SPHINX17CE36000801.
Senior Editor
 George H Perry, Pennsylvania State University, United States
Reviewing Editor
 Gwenan M Knight, London School of Hygiene and Tropical Medicine, United Kingdom
Reviewers
 Esther van Kleef, Institute of Tropical Medicine, Belgium
 Erik S Wright, University of Pittsburgh, United States
Publication history
 Received: March 26, 2021
 Accepted: August 31, 2021
 Accepted Manuscript published: September 14, 2021 (version 1)
 Version of Record published: November 1, 2021 (version 2)
Copyright
© 2021, Smith et al.
This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.
Metrics

 1,355
 Page views

 245
 Downloads

 2
 Citations
Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.