Suppressing evolution through environmental switching

Ecology and evolution under changing environments are important in many subfields of biology with implications for medicine. Here, we explore an example: the consequences of fluctuating environments on the emergence of antibiotic resistance, which is an immense and growing problem. Typically, high doses of antibiotics are employed to eliminate the infection quickly and minimize the time under which resistance may emerge. However, this strategy may not be optimal. Since resistance typically has a reproductive cost to the microbes, there can be trade-offs between resistance and competitiveness depending on the environment. Here we show conditions under which such trade-offs can be exploited to prevent the emergence of resistance. We develop a stochastic Lotka-Volterra model of a microbial system with competing phenotypes: a wild strain susceptible to the antibiotic, and a mutant strain that is resistant. We investigate the impact of various pulsed applications of antibiotics on population suppression. Leveraging the trade-offs between susceptibility and resistance when the antibiotic is applied or is not, we show how a strategy of environmental switching can suppress the infection while avoiding resistant mutants. We discuss limitations of the procedure depending on the microbe and pharmacodynamics and methods to ameliorate them.


Introduction
Populations will face a variety of environmental fluctuations of both biotic and abiotic nature. Since phenotypes typically have different reproductive success in differing environments, the dynamics of these fluctuations can be crucial in determining phenotypic composition. Here, we consider the effects of varying environments on the emergence and maintenance of antibiotic resistance.
The rise of microbial resistance is a looming catastrophe, and prudential use of antimicrobials is a fundamental means to prevent it [30]. Such strategies to limit the chance of resistance can be made at all levels of disease dynamics, from population level protocols to individual patient therapies. Studies of antibiotic resistance in vivo, in hospitals, and in the community at large using simple mathematical models can help address the pharmacodynamics, pharmacokinetics, and epidemiology of resistance [4,9,14,19,31,36]. Such models have found use in effectively modelling real-world experimental data [8,35,44,50]. In particular, modelling has been used in identifying dosing regimens that suppress the emergence of resistance [48,49].
There are several mechanisms by which bacteria can be resistant to antibiotics [42], an example of which is overexpression of the efflux pump [10,47,54], which bacteria use to expel antibiotics. Typically, such resistant mechanisms have a fitness cost, which can result in trade-offs between resistance and growth [7,16,33,53]. Example costs of resistance include less energy for other cellular processes, and impaired motility. When the antibiotic is present, the resistant mechanism pays for its cost by providing a fitness advantage relative to the susceptible strain. However, when the antibiotic is not present, the resistant mechanism incurs a fitness disadvantage. Resistance, therefore, can be reversed under evolutionary forces by altering the environment [3]. As such, a pulsed protocol, where the antibiotic is periodically applied so that the environment switches from antibiotic to antibiotic-free regimes, may be able to eliminate the bacteria. However, there is a risk that resistant mutants evolve to reduce the fitness cost of resistance rather than lose the resistance mechanism, whereby they could be competitive whether the antibiotic is present or not [39]. Further, mutations in regulatory genes can produce phenotypes of irreversible resistance [51]. These risks can cripple pulsed protocols aimed at controlling the infection while preventing resistance. Preventing a sustained presence of resistance is therefore a high priority.
Here, we develop a mathematical model of pulsed protocols of antibiotic and antibiotic-free regimes, switching rapidly from one environment to the other, to control a bacterial population of Escherichia coli. We consider concentration-independent (i.e. time-dependent) bactericides such as -lactams (e.g. penicillins and cephalosporins), which require high maintained concentrations to be effective [2]. We assume that there is a maximum benefit to the concentration amount (due either to the pharmacodynamics of the bactericidal mechanism or tolerance of the patient to the antimicrobial). Therefore, we fix the dose concentration when applied, and find the proper periods for each regime that prevent the emergence of resistance while eliminating the infection.
Previous theoretical studies have shown that pulsed protocols of antibiotics can eliminate bacteria [1,12,13,29]. However, these studies feature only a "persistent phenotype", that neither grows nor dies under application of the antimicrobial agent [6,56]. Bacteria may transition between the persister type and wild type, depending on the environmental conditions. The number of persisters remain at low levels and act as a staging ground for the bacteria to repopulate after the antimicrobial is removed. Pulsed protocols of antibiotics, however, can disrupt this process and lead to the elimination of the bacteria. Experimental studies have shown that pulsed protocols can be effective in controlling such a system [46].
Although pulsed protocols can eliminate non-persister and persister colonies, they have more difficulty in eliminating colonies with resistant phenotypes that can grow when antibiotics are present. However, these protocols have been shown to be effective in containing an infection both theoretically and experimentally [5,22]. In such cases, antimicrobials can act as ecological disturbances and can be approximately as effective as a constant application of the antimicrobial in controlling the bacterial load while also diminishing the probability of the emergence of resistance [5]. With short durations of high concentrations of drugs, the period under which resistance is selected for can be minimized.
The above studies have explored pulsed protocols in different ways: controlling persisters, and controlling emergence of resistance. The resistant strain we consider here does grow in the presence of antibiotics, and thus are not persisters. Our scenario is thus more similar to and an extension of [5]. Our main contribution here is to show how leveraging competition can not only suppress the emergence of resistance as in [5], but also reduce the overall bacterial load. Additionally, we explore the impact of other important mechanisms on pulsed protocols including the evolvability of the bacteria and the lethality of the antibiotic.
The key mechanism of our models in suppressing the population is competition between the two phenotypes. Competition can be modelled in several ways such as through the competitive Lotka-Volterra equation or resource competition models. Competition, however, can be low when the total size of the population is small (e.g. the population is well below the carrying capacity or there is a high amount of resource relative to the number of bacteria). Yet, if competition is sufficiently strong even when the population is small, then the resistant type can be eliminated during the no-antibiotic phase of the pulse. We explore the scenarios where competition is and is not substantial at low population levels. Our models also feature stochasticity, which we develop in a stochastic kinetic framework [55]. We show that only when selection against the resistant type is high when the antibiotic is off can pulsed protocols effectively eliminate the population.

Methods
Stochastic birth-death processes are widely used in biological modelling [37], and, in particular, stochastic modelling of the Lotka-Volterra system [24]. Our stochastic model features microscopic processes of birth, death, competition, and mutation, as detailed in Box 1. These processes operate on two phenotypes and , which represent a wild-type strain, which is susceptible to the antibiotic, and a mutant strain, which is resistant, respectively.
Consider first the dynamics of the birth and death processes without the presence of antibiotics, i.e. in the antibiotic-free regime. Reaction set 1 represents these processes, where and are the birth and death rates, respectively, for the wild-type strain with > > 0. We assume that the death rate for the resistant strain is also . However, assuming a cost > 0 to resistance, the growth rate of the resistant strain is − (costs applied to birth rather than death rates have also been applied similarly in ecological games [23]).
In the presence of the antibiotic, the above processes still occur, but with an additional set of reactions involving the antibiotic. Since we are considering a concentration-independent or timedependent antibiotic, we will assume that the amount of antibiotic,¯ , remains unchanged while we are in the antibiotic regime. At the maximum dose¯ = 1, normalized. The antibiotic is bactericidal, and kills both types of bacteria. Reaction set 3 represents death from the antibiotic with rates and for the susceptible and resistant strains, respectively. Note that for to be susceptible and to be resistant, we must have > − > + .
The species also die from competition as represented in Reaction set 2. Let be the rate at which they come into contact and compete where one dies. We assume that the rate of death from competition is proportional to the difference between the mean growth rate (not including the effect of antibiotics) for the focal type vs. its competitor. Thus, for same type competition, we set the rate of death to 1. However, the rate of death of a mutant in competition with a wild-type is > 1 and the rate of death of a wild-type in competition with a mutant is 1/ .

Box 1: Stochastic Lotka-Volterra processes
Mutation: The bacteria can turn into the other type via mutation at rate , which occurs in both regimes (Reaction set 4). However, under the stress of the antibiotic, the susceptible type will mutate to at a higher rate > .
The environmental switching is controlled by a choice of the on and off durations of the drug. We assume 100% bioavailability of the drug at application, e.g. intravenous application. Thus, when the antibiotic is "turned on," its effects are immediate. Further, when it is "turned off," the dissipation of the antibiotic -i.e. the rate at which it breaks into ineffective material, is metabolized, etc. -is rapid, which is common for concentration-independent and time-dependent antibiotics. In the antibiotic regime, we apply the maximum effective dose. Thus, the set of pulsed protocols we consider are sequences of durations of the regimes, where¯ = 1 when the antibiotic is on and = 0 when it is off.
Parameter and default value Definition = 0.35 Growth rate. Stress induced mutation rate of susceptible to resistant bacteria. We conducted numerical simulations of the model to test the effects of various protocols. We average realizations for each parameter combination. We use the Gillespie algorithm [18] from Julia's DifferentialEquations and Catalyst packages [43]. Table 1 lists the default parameter values used with rates per 15 minutes. We vary these values to explore nearby parameter space. We assume that a new generation occurs after 1/( − ) ≈ 1 hour. We estimate that the relative fitness of the resistant strain in the antibiotic-free environment is ( − − )/( − ) ≈ 0.8 =⇒ = 0.05, which is within experimentally evaluated values [34]. E. Coli has an average rate of mutation per genome per generation on the order of 10 −5 [40]. Further, under stress from the antibiotic, the mutation rate can be even larger, up to ten times the non-stressed rate [28]. Thus, we consider = 10 . We assume that resistant bacteria die from the antibiotic at 1/10 th the rate susceptibles do [11], i.e.
= 0.1 . The initial condition is a population of 100% susceptible bacteria, 0 = 10 5 . We explore a variety of competition parameters . We consider fixed on/off durations, where we repeat switching until the population is extinct or = 14 days.

Results
We first consider the mean field (i.e. non-stochastic) behaviour of the system with average values and (for details of the derivation see Appendix A.) In either environment, both types will coexist at equilibrium due to mutations, though the less adapted type will remain at low frequency. When the population size is sufficiently low (i.e. below the zero isoclines), competition will also be low, which will allow both types to increase in abundance. However, the higher the competition term , the smaller this region is. A large will cause the mutant strain to be suppressed in the antibiotic-free environment, and in the stochastic scenario, this effect increases the chance that the mutant strain be eliminated. In the remainder of the results, we detail the effects of switching the drug on and off, competition, and stochasticity in the stochastic scenario. Figure 1 depicts a representative time series for a switching protocol vs. a constant application of the antibiotic. With a sufficient competitive disadvantage for resistance, i.e. high , we can effectively control the bacterial load and antibiotic resistance. However, due to stochastic effects, resistance can briefly rise as seen in Figure 1c before it is brought under control again. This phenomenon emerges from the interactions between the two types at different time scales. Oscillations between the types emerge at a time scale longer than the period of the protocol. When the mutants are at their peak, the wild-type strain is not, and vice verse. However, peaks of the wild-type strain are relatively low as can be seen in panel d. Between day 10 until sometime during day 12, the mutants are suppressed and the wild-type are at a relative high. This high, however, is not much higher than when they are suppressed. As such, the total bacterial load is relatively low. Regardless of these waves of resistant bacteria, the time average of the total bacterial load is still less than when the mutant becomes established under constant application. Figures 1b and 1d are specific instances where the mutant becomes established under the constant application protocol. However, the mutations and fluctuatios in abundances are stochastic, so it is also possible that we are fortunate and constant application drives the population extinct before resistance emerges. Thus, to better understand the effectiveness of therapies we must evaluate the statistics of the bacterial load as a function of system parameters.
Averaged over 50 realizations, we calculate the average total bacterial load over time for pulsed protocols with various on and off durations and compare these to the bacterial load for constant application of the antibiotics. We plot these results in heat maps, where the colour indicates the long-term bacterial load relative to the outcome from constant antibiotic application. Red indicates that the therapy is on average worse than constant application, yellow is on average equal, and blue indicates that it is on average better.
We observe that pulsed protocols along a diagonal do best. One reason for this is that the switching times explored here are much less than the time to reach carrying capacity in either regime. For example, even a day-period protocol will not reach carrying capacity (the expected time to reach carrying capacity is between one and two days). In such a case, the population can swing from predominantly one type to the other (see Appendix B for an example time series). However, this behaviour can still be beneficial, since each application of the antibiotic is another chance of eliminating the population, since switching environments drives the dominant type down potentially to extinction before the other type can become established. In addition to plotting heat maps, we plot the total bacterial load over time for constant and pulsed (2hrs on and off each) for various values of each parameter averaged over 50 realizations. Figure 2 shows that the higher the competition, the lower the diagonal (i.e. the best results come from protocols where the duration on is greater than the duration off). The increased competition suppresses the emergence of resistance even in the antibiotic environment, and thus the duration of application can be longer. However, we note that the effectiveness of the optimal pulsed protocols drastically falls if the duration on relative to off is too high. We do not see this sensitivity when reducing the time on. We also observe an intermediate level of competition is best for pulsed protocols relative to constant application. We can see this effect in Figure 2g. Increasing decreases the mean bacterial load for the constant application as we would expect. Since, high competition between the types will suppress the emergence of resistant mutants (which is true in both environments). However, increasing competition has an initially steeper effect upon pulsed protocols before it levels off. A sufficient amount of competition is required for pulsed protocols to work. As is increased, the difference between the outcomes of the two approaches decreases.
Another reason for angle of the optimal diagonal of successful protocols in Figure 2, is due to the relationships between the mean growth rates and the antibiotic kill rates. Figure 3 depicts the results for various values of . The higher the antibiotic kill rate, the shorter the duration on for the most successful protocols. Like in the case of , impacts the effectiveness of pulsed protocols nonlinearly. Figure 3g depicts the mean results for various . The higher the , the better constant application does. However, this is not true for pulsed protocols. An intermediate value is best. This result is due to the impact of on competitiveness. If is too high, then the wild-type is suppressed too much, and thus cannot be used to suppress the mutant strain through competition.
To explore how robust our results are to mutation rates, we considered various values of and . We can see the effects of various in the rows of Figure 4, which show that the pulsed protocols are more effective under a higher mutation rate. The first row depicts the case where there is a stress induced mutation rate from wild-type to mutant ( = 10 ). The second row depicts the results where stress does not increase the mutation rate ( = ). The stress induced mutation makes the antibiotic environment more conducive to generating resistance, and thus makes it harder to control the emergence of resistance. The impact is a small relative effect on the pulsed protocols vs. constant application protocols. Figure 4g shows that the mutation rate impacts the constant application more so than the pulsed protocol. This result matches intuition; the higher the mutation rate, the less likely a constant application can eliminate the colony before a mutant arises and becomes established. In summary, the more evolvable the system, the better switching environments works. In Figure 5 we explored the effect of varying the contact rate , and observed that switching was more effective for an intermediate . For high , the region where both types can grow is small, which magnifies the impact of stochastic effects leading to elimination of emerging mutants. Further, same-type competitive interactions are also more intense, and thus the population is driven to extinction more quickly. Therefore, a constant application is best. If is too low, however, there is not sufficient levels of cooperation for pulsed protocols to outperform constant ones. This observation has broader implications to theory such as the paradox of the plankton. In a switching Here what we mean is that switching environments which do poorly either lead to extinction or the resistant strain dominating (i.e less diversity than if both types are controlled at low levels). Figure  5d shows that the effectiveness of both constant application and pulsed protocols increases as increases (and thus the carrying capacity is lower). However, the gap between the two shrinks, and thus for a high system, the absolute bacterial load of the constant application is similarly to the low level of the pulsed protocol.

Discussion
Our model furthers the contribution of mathematical modelling to the fight against resistance [38]. We found a control strategy that can eliminate the infection while avoiding resistant mutants. Our study differs from [5] in that their study showed how resistance could be greatly reduced while keeping bacterial load steady. However, this only works when competition is high, else at low population levels, resistance can be maintained. Previous empirical and theoretical research has also found the importance of high competition in containing an infection [22]. Further, pulsed protocols and competition can be effective in containing an infection even in a well-mixed population [22]. Spatial effects, such as those found in biofilms, could heighten the degree of competition and thus the effectiveness of the pulsed protocols. Since, spatial heterogeneity due to clumping could keep the competition level between different types high even when the population is small relative to the carrying capacity. The competition effects we consider can be interpreted crudely as arising from such effects. Bactericides with significant post-antibiotic effects (PAE), such as fluoroquinolones, may hamper our control strategy. Such antimicrobials can impact the bacteria at sub-MIC levels long after they have been removed from the system [2], and as such can select for resistant bacteria after the antibiotic has been turned off. Additionally, sub MIC levels can lead to multidrug resistance through radical-induced mutagensis [26]. Therefore, we must have a rapid dissipation of the antibiotic once below the MIC to prevent selection for the mutant (the range between the MIC and the point at which the susceptible strain is selected for) [20]. PAE is frequently caused by antibiotics that impair DNA functioning. Hence, -lactams, which inhibit cell wall production, are a good choice for our therapies. Antibiotics with a short half life would also be effective with our protocols.
Intermediate drug concentration is a dimension we did not consider here. The effects of concentration on the spread of the disease due to within host and between host dynamics have been found to be an important factor [45]. It is not always beneficial to have a high concentration due to a U-shaped probability of resistance emerging vs. drug concentration [15]. Here we only considered a specific concentration for all applications.
We restrict ourselves to the case where the system cannot be well monitored. Clearly, the best strategy is to alter the durations dependent on the state of the system, which can both maximize the duration of the antibiotic regime and prevent the emergence of resistance. The more repeated applications of the treatment, the less likely it will work. However, such observations may not be feasible, especially when the bacterial load is small and heterogeneous.
Future models could incorporate other biological factors. For example, the setting and source of resistance (source-sink dynamics) are important factors in controlling antibiotic resistance [41]. The method by which resistance is spread is another important factor such as where plasmids confer resistance. In such a scenario, resistance can reemerge rapidly. Since, plasmids can remain in the population due to horizontal transfer even when the plasmid confers a cost [32]. Although, the transfer rate has been show to dramatically fall once the population is low [21]. Future models could also incorporate more of the complexity of interactions between the bacteria and the patients' natural flora [17,52].
Our work has applications to cancer [25, 27] and which families of chemotherapies will work best. More generally, it applies to species conflict in varying environments and the impact of such variations on extinction and coexistence. Though our aim here has been to lower the overall bacterial load, in some simulations (see Figure 1d), alternating environments prevents one species dominating thereby sustaining coexistence. This temporal heterogeneity in competitiveness can thus act as a stabilizing mechanism that promotes diversity (as measured by the relative proportions of each type over time).

B Further simulations results
Here we depicts other simulations results. In Figure 6 we depict an example time series for a one day on and off protocol where the bacteria is not eliminated. Longer durations on or off can lead to the population reaching the carrying capacity, after which they will be markedly less successful at suppressing the bacteria than the constant application therapy. = 4 and the remaining parameters are from Table 1.
In order to explore the effects of greater stochasticity in the birth and death rates on the outcomes of these therapies, we ran simulations with birth and death rates + and + , which keeps the mean growth rates unchanged. Figure 7 depicts the results. We see that if the higher stochasticity can reduce the relative effectiveness of pulsed protocols. However, this effect is not large. In Figure  7g, we can see that the gulf between the outcomes of constant applications and pulsed protocols decreases as increases.