The ability of a population of bacteria to survive the lethal effects of an antibiotic

Previous studies of large populations have established the ‘minimum inhibitory concentration’ (MIC; the lowest concentration of the drug that inhibits population growth) as the most critical parameter for characterizing the dynamics of a bacterial population under antibiotics (Regoes et al., 2004; Czock et al., 2009; Craig, 1998; Falagas et al., 2012). The dynamics of bacterial populations exposed to different concentrations of antibiotics have been examined and modeled deterministically in relation to the MIC, as follows. First, without drugs, the growth rate of cells, λ, is higher than the death rate, ϕ (i.e., λ > ϕ), and thus a bacterial population always grows. When drug concentration increases, as long as the concentration remains below the MIC (i.e., sub-MIC), growth rate is higher than death rate (λ > ϕ), and thus a population still grows, albeit at slower rates. When the drug concentration increases further and reaches the MIC, growth rate becomes equal to death rate (λ = ϕ), and the population size is maintained at a constant level. Only at drug concentrations above the MIC does a bacterial population decline. Extrapolating this deterministic knowledge to population extinction, studies often claimed that maintaining drug concentrations above the MIC was absolutely essential to eradicate bacterial populations. As will become evident later, our data challenge this deterministic framework.

As a first step to examine the dynamics of a small population in relation to the MIC, we used a plate assay and characterized how single E. coli cells grew and formed colonies at various antibiotic concentrations. Antibiotic-susceptible, wild-type (WT) E. coli cells were cultured in liquid LB medium without antibiotics and then spread on LB agar plates containing increasing concentrations of antibiotics. After 18 hr of incubation, the number of colony-forming units (NCFU) was determined. By normalizing the NCFU to that for an antibiotic-free plate (N0CFU), we then obtained the plating efficiency (=NCFU/N0CFU), which indicates the fraction of cells forming colonies. Following the definition of MIC as the lowest drug concentration that inhibits population growth, the lowest concentration yielding no visible colonies on the plates was defined as the MIC here. See Figure 1—figure supplement 1 for a detailed illustration of this procedure and Supplementary file 1 for the MIC values for all of the antibiotics examined.

When we performed this plate assay for various antibiotic drugs, we observed two strikingly distinct trends, which depended on whether the drug used was bacteriostatic (which suppresses cell growth) or bactericidal (which induces cell death). For bacteriostatic drugs, at increasing concentrations, the plating efficiency remained nearly constant and abruptly dropped to zero when the drug concentration reached the MIC (Figure 1A); the grey line was obtained from a linear regression analysis of the whole data set below 0.75 × MIC (see Figure 1A caption for details). This trend indicates that almost every single cell spread on the plate grew and formed colonies at a wide range of sub-MIC drug concentration, and no cells formed colonies at (and above) MIC. This observation, suggesting homogeneous population dynamics, agrees with the deterministic prediction discussed above. Additionally, we observed a decrease in colony size at increasing drug concentrations (Figure 1—figure supplement 2).

(A) When we performed the plate assay using bacteriostatic drugs, NCFU was generally maintained at increasing drug concentrations. See Figure 1—figure supplement 1 for a detailed illustration of how this plot was made. A linear regression analysis was performed for the whole data set below 0.75 × MIC, resulting in the grey line (slope = −0.05, intercept = 0.9942, R-squared = 0.99). Colony size decreased at increasing drug concentrations (Figure 1—figure supplement 2). (B) For bactericidal drugs, NCFU decreased at increasing drug concentrations, indicating heterogeneous population dynamics of bacteria. A linear regression analysis was performed for the whole data set, and the result was plotted as a grey line (slope = −1.07, intercept = 1.10, R-squared = 0.79). For both groups of drugs, the lowest concentration yielding no colonies was defined as the MIC. The Supplementary File 1 lists the MICs for all drugs examined in this study. We plotted plating efficiency for each antibiotic in separate panels in Figure 1—figure supplement 4. (C) For each bactericidal drug used in the experiment (for which the results were shown in Figure 1B), we picked a few colonies from the agar plate exhibiting a plating efficiency of ~0.5 (e.g., near 0.6 × MIC) and plated them immediately on fresh agar plates containing various concentrations of the same drug. The plating efficiency was similar or marginally lower on the second plating, possibly because exposure to the antibiotics on the first plate adversely affected the cells and rendered them more susceptible to the antibiotics on the second plating. See Figure 1—figure supplement 3 for similar results for other drugs. Therefore, the ability of bacteria to grow and form colonies on plates containing bactericidal drugs was not heritable. We performed at least two biological repeats for all the experiments and plotted the mean here. The error bars represent one standard deviation from the repeats.

For bactericidal drugs, at increasing concentrations, the plating efficiency decreased gradually from 1 to 0 (Figure 1B); the grey line was obtained from a linear regression analysis of the whole data set (see Figure 1B caption for details). This trend contrasts with our finding for bacteriostatic drugs (compare the grey lines in Figure 1A and B) and cannot be explained by the deterministic model. In the literature, other studies have reported a similar gradual decrease in the plating efficiency (Liu et al., 2011; Ernst et al., 2002; Dong et al., 2000). However, those studies primarily concerned how to better determine the MIC in the face of such a gradual decrease, and have not characterized population dynamics underlying the gradual decrease.

The plating efficiency between 0 and 1 indicates heterogeneous colony formation. In the plate assay above, we found that at ~0.6 × MIC, the plating efficiency was ~0.5, meaning that approximately 50% of the cells plated formed colonies and 50% did not. One possible explanation is that the colony-forming cells were intrinsically more resistant to the drugs than the cells that did not form colonies, subsequently giving rise to resistant daughter cells (i.e., heritable resistance). To examine this possibility, for each bactericidal drug used in the experiment (for which the results are shown in Figure 1B), we picked colonies from agar plates exhibiting a plating efficiency of ~0.5 (near 0.6 × MIC), suspended them in liquid medium, and immediately plated them on fresh agar plates containing various concentrations of the same drug. The results were plotted in Figure 1C and Figure 1—figure supplement 3. Contrary to our expectation, the plating efficiency of the second plating was about the same as that of the first plating, or in some cases, marginally lower (possibly because cells were challenged with drugs twice consecutively). This observation rules out heritable resistance as an explanation for heterogeneous colony formation at sub-MIC drug concentrations.

Our data above showing the absence of heritable resistance in the surviving populations suggest the possible involvement of stochasticity. That is, bactericidal drugs might induce stochastic fluctuations in the bacterial population size. To investigate this possibility, we performed a plate assay as above and followed the population dynamics of growing micro-colonies at single-cell resolution; we spread cells on agar plates, and rather than waiting for 18 hr and counting colonies visible to the naked eye, we examined how isolated single cells grew to form micro-colonies using time-lapse microscopy. First, as a control, we examined the dynamics in the absence of antibiotics (Video 1). We counted the number of growing cells in each colony and plotted the number versus time (Figure 2A). The colonies proliferated homogeneously, meaning that different colonies grew similarly. We then repeated the experiment using a bacteriostatic drug (chloramphenicol, thiolutin, or tetracycline) at a sub-MIC level. Visually inspecting the image sequences, we found that cells stably grew, albeit at lower rates (Video 2). The number of growing cells in each micro-colony increased homogeneously and similarly (Figure 2B and Figure 2—figure supplement 1). Previously, we developed a microfluidic chemostat for cell culture (Kim et al., 2012; Deris et al., 2013). When we repeated these experiments using this device, we again observed the same homogeneous population dynamics (Figure 2—figure supplement 2A).

We characterized how cells grow and form micro-colonies on LB agar plates containing different concentrations of antibiotics using time-lapse microscopy. At time zero, we transferred cells growing in antibiotic-free LB liquid medium to a LB agar plate, and confirmed that individual cells were spread out and isolated. Setting up a time-lapse imaging experiment after the transfer took 1 ~ 2 hr, which is why there is a gap in the data immediately after time zero (dashed lines). The experiments were independently repeated twice (biological repeats) and the data from one experiment are shown here. (A, B) We counted the number of growing cells in each micro-colony (represented by a different color). In the absence of antibiotics (panel A) or with a sub-MIC level of a bacteriostatic drug (panel B, 0.7 × MIC of chloramphenicol), the number increased similarly for different colonies, revealing homogeneous population growth. Example image sequences are shown in Videos 1 and 2. The data are replotted on a semi-log scale (insets). Such homogeneous population growth was observed for other bacteriostatic drugs (thiolutin and tetracycline) as well; see Figure 2—figure supplement 1. When we repeated the experiment using a microfluidic device, we again observed the same homogeneous population dynamics (Figure 2—figure supplement 2A). (C) In contrast, the population dynamics of bacteria exposed to a bactericidal drug were highly heterogeneous. An example image sequence was shown in Video 3. The number of growing cells within micro-colonies at 0.8 × MIC of cefsulodin is plotted here; the result shows stochastic population fluctuations. Such population fluctuations were again observed when experiments were repeated using other bactericidal drugs (ofloxacin, kanamycin, and 6-APA); see Figure 2—figure supplement 3A–C. The light red-shaded region indicates the number equal to zero (i.e., population extinction). When we repeated these experiments using a microfluidic device, we again observed population fluctuations (Figure 2—figure supplement 2B–D). (D) Our model predicts that the probability of population extinction increases linearly with death rate, with the slope being 1 / λ (the solid line). We experimentally characterized the extinction probability (Figure 2—figure supplement 4), and the death rate (Figure 2—figure supplements 5–6), at different concentrations of bactericidal drugs. Using these data, we obtained the relationship between the extinction probability and the death rate, and plotted it here (green circles: cefsulodin, pink squares: ofloxacin, and grey triangles: 6-APA). We found good agreement between the model prediction and experimental data. Note that at increasing drug concentrations, λ changed little (Figure 2—figure supplement 6B), and thus was taken as a constant in the analysis here.

This video cannot be played in place because your browser does support HTML5 video. You may still download the video for offline viewing.

We characterized how cells grew and formed micro-colonies on LB agar using time-lapse microscopy. 0.7 × MIC of chloramphenicol was used. An example image sequence is shown here. The time interval between each frame is 19 min.

Next, we characterized population dynamics for a bactericidal drug (cefsulodin, ofloxacin, kanamycin, or 6-APA), at a sub-MIC level. We found that the population dynamics were highly stochastic (Video 3). Visual inspection of such image sequences indicated that within a given population, some cells were killed stochastically, whereas other cells survived and divided. Such demographic stochasticity would lead to random fluctuations in the population size. Indeed, the number of growing cells in each colony fluctuated randomly over time (Figure 2C and Figure 2—figure supplement 3A–C). These fluctuations led to dramatically different dynamics for different colonies, even though they originated from genetically identical cells and were cultured under homogeneous antibiotic conditions. When we repeated these experiments using the microfluidic chemostat, we again observed significant population fluctuations (Figure 2—figure supplement 2B–D).

This video cannot be played in place because your browser does support HTML5 video. You may still download the video for offline viewing.

We characterized how cells grew and formed micro-colonies on LB agar using time-lapse microscopy. 0.8 × MIC of cefsulodin was used. An example image sequence is shown here. The time interval between each frame is 15 min.

Importantly, the fluctuations drove some colonies into extinction (the light red-shaded area in Figure 2C). We counted the number of colonies that went extinct and plotted the probability of colony extinction at various drug concentrations; here, we are interested in colony extinction because it is equivalent to bacterial clearance, which means treatment success. Figure 2—figure supplement 4 showed that the extinction probability increased with increasing drug concentrations. This increase in the extinction probability agrees with the trend of decreasing plating efficiency we found above (Figure 1B).

A population will undergo extinction if cells die more frequently than divide. Because bactericidal drugs induce cell death, an increase in extinction probability at higher drug concentrations is expected to be due to an increase in the rate of cell death, ϕ. Additionally, bactericidal drugs inflict damage on cells (Belenky et al., 2015; Lobritz et al., 2015). Thus, the rate of cell growth, λ, might decrease at higher drug concentrations, which could also contribute to colony extinction. We next sought to determine how bactericidal drugs affect ϕ and λ. Previous studies of population growth have shown that at a higher concentration of bactericidal drugs, the ‘net growth rate’, which is equal to λ− ϕ, decreases (Regoes et al., 2004). But, to separately resolve changes in ϕ and λ, the growth and death of cells must be examined at single-cell resolution. We analyzed the single-cell-level image sequences we obtained above (see Figure 2—figure supplement 5 for details of the analysis), and determined λ and ϕ at various concentrations of bactericidal drugs. We found that at increasing drug concentrations, ϕ increased (Figure 2—figure supplement 6A). But λ changed little, remaining nearly constant (Figure 2—figure supplement 6B). Therefore, bactericidal drugs have significant effects on cell death, but not on cell growth.

Our data above suggest that an increase in the probability of population extinction at higher drug concentrations (Figure 2—figure supplement 4) is likely due to an increase in the rate of cell death (Figure 2—figure supplement 6). To quantitatively understand the relationship between the extinction probability and death rate, we employed a stochastic model, known as the Markovian birth-and-death process, that has been widely used to study the basic features of stochastic population dynamics (Novozhilov et al., 2006; Pavel Krapivsky and Ben-Naim, 2010; Kendall, 1948). This model contains two parameters, the rate of cell growth and death, λ and ϕ, respectively. Each individual cell can divide or die stochastically with the probabilities determined by these parameters. Due to this demographic stochasticity, the number of cells within a population, n, fluctuates over time. Thus, n cannot be predicted deterministically but only probabilistically, and the probability is described as follows,

(1) Pn.=λ(n−1)Pn−1−(λ+ϕ)nPn+ϕ(n+1)Pn+1,

where Pn refers to the probability of n cells being present in a population. The key boundary condition in this model is that once n reaches 0, it cannot change afterward. Known as an ‘absorbing boundary’, this condition reflects the biological reality that once a population goes extinct, it cannot revive. Therefore, a solution of this model for n = 0 (i.e., Pn=0 or simply P0) describes the dynamics of population extinction. Please see Appendix 1 for the detail and solution.

We first considered the extinction probability P0 at low drug concentrations in which the death rate is lower than the growth rate (ϕ < λ). In this range, P0 is given by their ratio, P0 = ϕ / λ; see Equation A6. Thus, if the death rate is zero (ϕ = 0), the extinction probability P0 is 0, meaning n always increases (this makes intuitive sense). As the death rate increases (0 < ϕ < λ), P0 increases and becomes non-zero, meaning that n may stochastically reach the absorbing boundary, agreeing with our observation of stochastic population extinction at sub-MIC drug concentrations; Figure 2C shows that some populations reached n = 0 (marked by the light red-shaded area), while other populations thrived (also see Figure 2—figure supplement 3A–C). The solution P0 = ϕ / λ predicts that the extinction probability increases linearly with death rate, with the slope being 1 / λ (the solid line in Figure 2D). We sought to test this prediction quantitatively by comparing it with experimental data. Above, analyzing time-lapse microscope images, we obtained the probability of population extinction (Figure 2—figure supplement 4), and the death rate (Figure 2—figure supplement 6), at different concentrations of bactericidal drugs. Using these data, we obtained the relationship between the probability of population extinction and the death rate, and plotted it in Figure 2D. We found good agreement between the model prediction and experimental data (compare the solid line and symbols in Figure 2D).

Next, using the quantitative relationship we found above (P0 = ϕ / λ), we will specify the condition for the MIC. In our plate assay (Figure 1B), we observed that the plating efficiency decreases at higher drug concentrations (in the sub-MIC range) and reaches zero at the MIC. Also, the quantitative relationship we found above showed that extinction probability increases at higher drug concentrations (consistent with a decrease in the plating efficiency), reaching one when the death rate is equal to the growth rate (see the dashed line in Figure 2D); thus, P0 = 1 at ϕ = λ. The extinction probability equal to one (P0 = 1) means that all colonies go extinct, which corresponds to zero plating efficiency. The drug concentration at which the plating efficiency reaches zero is the MIC (Figure 1B). Taken together, at the MIC, the plating efficiency is zero because extinction probability is one (P0 = 1), and the extinction probability is one because the growth rate and death rate are equal to each other (ϕ = λ). In short, ϕ = λ at the MIC.

Next, we considered drug concentrations above the MIC, where the death rate is higher than the growth rate (ϕ > λ). In this range, the model predicts that all populations eventually go extinct (P0 = 1 in Equation A6); this makes intuitive sense. Importantly, due to population fluctuations, populations are expected to go extinct at various times (Equation A5), meaning that the number of live populations (the populations that have not undergone extinction yet) decreases gradually over time. The model predicts that this decrease can be approximated by an exponential decay in the long time limit (t≫1/|ϕ−λ|); see Equation A7. We tested these model predictions by repeating time-lapse microscope experiments at drug concentrations above the MIC. All the populations indeed went extinct at various times (Figure 3A and Figure 2—figure supplement 3D–F). When we counted the number of live colonies, this number decreased gradually over time (Figure 3B). In this semi-log plot, the decrease was linear (compare it with the dashed line), consistent with the model prediction of an exponential decay (Equation A7).

(A, B) In Figure 2C and Figure 2—figure supplement 3A–C, we showed the population dynamics of bacteria exposed to sub-MIC levels of bactericidal drugs. Here, we show the dynamics at drug concentrations above MIC (1.2 × MIC). The number of growing cells within micro-colonies for cefsulodin is plotted in panel A. The results for other bactericidal drugs were plotted in Figure 2—figure supplement 3D–F. We observed that all populations went extinct at various times. We then counted the number of live colonies (colonies that have not undergone extinction yet) at various times (~150 colonies monitored). The number decreased gradually over time (green circles: cefsulodin, pink squares: ofloxacin, and grey triangles: 6-APA); see panel B. The decrease was linear in a semi-log plot, consistent with the model prediction of an exponential decay (Equation A7). (C) We characterized the extinction probability of populations starting with ~640 cells. We prepared a large volume of a cell culture with the cell density of ~640 cells/ 200 µL, supplemented the culture with a low concentration of cefsulodin (23 µg/ml), and then distributed 200 µL of the cell culture equally to 12 isolated chambers in a microtiter plate. We repeated this procedure using higher cefsulodin concentrations (but maintaining the inoculum size). We then incubated the microtiter plate overnight. By counting the chambers that became turbid or clear, we calculated the extinction probability. We used a heat map to graphically represent the probability; for illustration purpose, we also added the values of the probability to the plot. At low cefsulodin concentrations (to the left of the dashed line), all chambers became turbid (P0 = 0). At high concentrations (to the right of the solid line), all chambers were clear (P0 = 1); thus, the solid line indicates MIC. At intermediate concentrations (between the dashed and solid lines), we observed heterogeneous results (only some chambers were clear, 0 < P0 < 1). (D) We repeated this experiment using different inoculum sizes and bactericidal drugs. Please note that although we prepared a large volume culture and distributed it equally to chambers, the number of cells in each chamber might vary. We found that the variation was ~10% or less. See Appendix 3 for details.

Our findings above indicate that the simple stochastic model can adequately capture the extinction dynamics of populations exposed to bactericidal drugs. What is striking in our findings is that, due to drug-induced population fluctuations, a bacterial population may undergo extinction even at sub-MIC concentrations, the concentrations that were previously deemed inefficacious to clear bacteria. We have established this result by examining the dynamics of colonies originated from single bacterial cells, the smallest possible inoculum size. Clinical studies have characterized the bacterial inoculum size that can produce infections (i.e., infectious dose). The infectious dose can be as low as one (Jones et al., 2006; Haas and Rose, 1994; Jones et al., 2005), but is generally 10–100 (Tuttle et al., 1999; DuPont et al., 1989; Hara-Kudo and Takatori, 2011; Kaiser et al., 1992) or larger (Kothary and BABU, 2001; Gama et al., 2012). It is expected that, with larger inoculum size, a population experiences less fluctuations (because demographic stochasticity gets averaged out). When we estimated the magnitudes of fluctuations by the coefficient of variation (CV), that is, the standard deviation divided by the mean, using our model, Equation A16 shows decreasing CV with increasing inoculum size, supporting the expectation above. Interestingly, Equation A16 also predicted that the magnitude of population fluctuations depends on rates of cell growth and death as well (Equation A16); CV increases as the death rate approaches the growth rate, meaning that population fluctuations become intensified when neither growth nor death is a dominant event. This prediction, together with our finding above that the death rate approaches the growth rate as the drug concentration increases to the MIC (Figure 2D), suggests that near the MIC, a population with relatively large inoculum size may still be prone to stochastic extinction.

To test this possibility, we experimentally characterized stochastic clearance of a bacterial population starting with different inoculum sizes. First, we prepared a large volume of a cell culture with the cell density of ~640 cells/ 200 µL, supplemented the culture with a low concentration of cefsulodin (23 µg/ml), and then distributed 200 µL of the cell culture equally to 12 isolated chambers in a microtiter plate. Here, an isolated cell culture in each chamber represents a separate population. We repeated this procedure using higher cefsulodin concentrations (but maintaining the inoculum size). We then incubated the microtiter plate overnight. We found that, at low cefsulodin concentrations (to the left of the dashed line in Figure 3C), all chambers became turbid, meaning all populations grew. Thus, the probability of population extinction was zero (P0 = 0). In Figure 3C, we used a heat map to graphically represent the probability (we also added the values of the probability in the graph for additional clarification). At high concentrations (to the right of the solid line), all chambers were clear (P0 = 1); thus, the solid line indicates the MIC. At intermediate concentrations (between the dashed and solid lines), we observed heterogeneous population growth; some chambers were clear while others were turbid (0 < P0 < 1). Subsequent plating of the clear cultures on drug-free LB agar plates yielded no colonies, indicating population extinction.

We then repeated this experiment by using different inoculum sizes and bactericidal drugs. The results were plotted as heat maps in Figure 3D; here, the Y axis represents different inoculum sizes. As above, a solid line represents the MIC, above which all chambers were clear (P0 = 1). MIC values were higher at higher cell densities. We note that the higher MIC at higher cell density (inoculum effect) has been observed previously and is being actively studied by others in the field (Artemova et al., 2015; Brook, 1989; Tan et al., 2012; Karslake et al., 2016). Thus, it is not the focus of our study; rather, we focus on stochastic clearance below the MIC. As above, a dashed line represents the concentrations below which all chambers were turbid (P0 = 0). The area between the dashed and solid lines indicates the range of drug concentrations and inoculum sizes that exhibited heterogeneous population growth (0 < P0 < 1), meaning stochastic clearance. Figure 3D shows that stochastic clearance occurs even for a population starting with as large as ~20,000 cells, inoculum size much larger than infectious doses for many infectious diseases.

Previously, antibiotic treatment at sub-MIC levels was not considered a viable option for bacterial eradication, because the deterministic model predicts that all bacterial populations should grow at sub-MIC levels (i.e., antibiotic treatment failure). However, our experimental results and stochastic model above indicate that at sub-MIC levels of bactericidal drugs, a population might undergo extinction stochastically. An increase in this probability while keeping the drug concentration low would be therapeutically useful; with an increased probability of extinction, sub-MIC ranges of drugs could be used to eradicate bacteria reliably. We therefore employed our model to explore how the extinction probability can be altered by means other than changing the bactericidal drug concentration.

Our model indicates that the extinction probability is determined by the ratio of the death and growth rates (P0 = ϕ /λ; Equation A6). Thus, based on the model, a reduction in growth rate (denominator) should lead to an increase in the extinction probability. Growth rate can be reduced by using poor growth media, or alternatively using bacteriostatic drugs. This means, for a sub-MIC concentration of a bactericidal drug (for which the extinction probability is less than 1), either a switch to poor growth media or addition of a sub-MIC level of a bacteriostatic drug would lead to an increase in the extinction probability. We note that the latter represents combination therapies, and other studies have characterized bacterial responses to combination therapies (Bollenbach, 2015). However, these studies primarily concerned deterministic changes in the MIC of a large population. Conversely, our study focuses on how combination therapies affect stochastic occurrence of population extinction. Another difference is our focus on sub-MIC drug ranges, an important point given previous research showing that the effects of drug combinations at the MIC might differ from those at sub-MIC levels (Ocampo et al., 2014).

To characterize the extinction probability, we introduced the plating inefficiency (= 1 – plating efficiency); the plating efficiency reflects the probability that a bacterial cell forms a population of a bacterial colony, and therefore, the plating inefficiency reflects the probability of population extinction. We first calculated the plating inefficiency using the plate assay results in Figure 1B and plotted it as white columns in Figure 4A; thus, white columns represent the extinction probability in rich growth media (LB) with bactericidal drugs alone. To test the effect of our growth reduction strategy, we repeated the plate assay, either by replacing LB with casamino acids (a poor nutrient source that leads to slower growth than LB; see Figure 4—figure supplement 1), or by adding a sub-MIC concentration of chloramphenicol or tetracycline. The plating inefficiency obtained with these treatments was plotted as solid columns in Figure 4A. The rise of solid columns above white columns indicates that growth reduction indeed led to an increase in the extinction probability, in agreement with our prediction.

(A, B) Previously, we plated WT (antibiotic-susceptible) E. coli cells on LB agar plates at various concentrations of cefsulodin and plotted the plating efficiency in Figure 1B. Here, we plotted the plating inefficiency (1- plating efficiency) as white columns in the top row; the plating inefficiency represents the extinction probability. We then repeated a plate assay, either by using casamino acids agar plates (left), or by adding a sub-MIC concentration of a bacteriostatic drug, chloramphenicol (0.7 × MIC, middle) or tetracycline (0.5 × MIC, right), to LB agar plates. Note that casamino acids lead to slower growth than LB; see Figure 4—figure supplement 1. We plotted the plating inefficiency obtained with these treatments as solid columns in panel A. We repeated this procedure for rifampicin (middle row) and vancomycin (bottom row), for kanamycin and 6-APA (Figure 4—figure supplement 3), and for ciprofloxacin and ofloxacin (Figure 4—figure supplement 4). Note that how generally growth reduction leads to an increase in the extinction probability depends on whether ϕ is dependent on λ or not; see the main text and Figure 4—figure supplement 2. We then used antibiotic-resistant strains and repeated these experiments (panel B). Please note that the MICs of these mutants were five to ten fold higher than those of the WT strain. The rise of solid columns above white columns indicates an increase in the extinction probability. We performed at least two biological repeats for all the experiments and plotted the mean. The error bars represent one standard deviation from the repeats.

We next examined how generally such an increase in the extinction probability might occur. The extinction probability depends on the growth rate (λ) and death rate (ϕ), which we assumed to be independent. Because the probability is equal to their ratio (ϕ/λ; Equation A6), this growth reduction strategy might not work when ϕ is not independent but decreases in response to a decrease in λ. Such coupling between ϕ and λ could occur for bactericidal drugs that kill only growing cells, possibly because these drugs target processes critical for cell growth. This means that, for bactericidal drugs that exhibit a killing rate of zero for non-growing cells (i.e., ϕ→ 0 when λ→ 0), neither the switch to poor growth medium nor addition of bacteriostatic drugs would increase the extinction probability. To test this possibility, we first identified such bactericidal drugs; we stopped cell growth in cultures by depriving the cells of nutrients, added bactericidal drugs at concentrations capable of eradicating growing cells, and then determined which drugs were no longer capable of killing the bacteria. As shown in Figure 4—figure supplement 2, killing was completely abolished for kanamycin and 6-APA, indicating that ϕ→ 0 when λ→ 0. As expected, for these drugs, neither the switch to poor growth medium nor addition of bacteriostatic drugs resulted in an increase in the extinction probability (Figure 4—figure supplement 3).

Furthermore, our model does not consider specific drug-drug interactions. For example, a previous study showed bacteriostatic translation-inhibiting drugs and bactericidal quinolone drugs affect gene expression in a way to negate their effects (Bollenbach et al., 2009). Consistent with this study, we failed to observe significant changes in the extinction probability when ofloxacin or ciprofloxacin (quinolone drugs) is used with and without bacteriostatic translation-inhibiting drugs (Figure 4—figure supplement 4).

Above, we tested the growth reduction strategy for a WT (antibiotic-susceptible) E. coli strain. Although the strategy did not work for some drugs due to their complex effects on cells, for those that worked, the strategy substantially increased the chance of bacterial clearance at sub-MIC drug concentrations. We wondered if this strategy could be applicable to antibiotic-resistant strains. Resistant strains are difficult to eradicate because their MICs are very high, to levels that are toxic to hosts. Therefore, the development of therapies that utilize sub-MIC doses would be highly useful. To test if our growth-reduction strategy would work for antibiotic-resistant bacteria, we repeated a plate assay using resistant strains. Laboratory evolution of rifampicin resistance has been frequently reported in the literature (Goldstein, 2014). By plating WT (antibiotic-susceptible) E. coli cells on a LB agar plate containing the rifampicin concentration equal to 2 × MIC, we isolated a rifampicin-resistant (Rifr) mutant; the MIC of this mutant was ~10 fold higher than that of the parent strain. We then plated the mutant on LB agar plates containing various rifampicin concentrations. Additionally, we acquired clinically-isolated E. coli strains that were resistant to either cefsulodin or vancomycin (see Methods), and plated them on LB agar plates with increasing concentrations of cefsulodin or vancomycin. These resistant strains exhibited non-zero plating inefficiency at sub-MIC concentrations of bactericidal drugs (white columns in Figure 4B), indicating stochastic clearance. We then repeated the experiments either by replacing LB with casamino acids or by adding a sub-MIC concentration of a bacteriostatic drug. As with the antibiotic-susceptible strain, these treatments led to an increase in the plating inefficiency (Figure 4B), indicating that our growth reduction strategy facilitated the clearance of antibiotic-resistant bacteria at sub-MIC concentrations.

Experiments were conducted using E. coli strain NCM3722 (Soupene et al., 2003; Lyons et al., 2011; Brown and Jun, 2015). Bacteria were grown in 5 mL of Lysogeny Broth (LB, Fisher Bioreagents) in 20 mL borosilicate glass culture tubes at 37° C with shaking (250 rpm). Our typical experimental procedure is as follows. Cells were first cultured in LB broth overnight (pre-culture). The next morning, the cells were sub-cultured into pre-warmed LB broth at the optical density (OD600) of ~0.001 (measured using a Genesys20 spectrophotometer, Thermo-Fisher) and allowed to grow exponentially. The culture at the OD600 of ~0.4 was used for a plate assay or for microscope experiments.

We evolved a rifampicin-resistant mutant from our WT strain (NCM3722). We plated ~108 WT cells on a LB agar plate containing 20 μg/ml (2 × MICWT). 18 colonies were formed next day. We chose one colony and purified it by re-streaking. The mutant (EMK32) had the MIC of 160 μg/ml. We obtained the cefsulodin-resistant strain (EMK35, MIC = 170 μg/ml) and vancomycin-resistant strain (EMK36, MIC = 250 μg/ml) from Georgia Emerging Infections Program MuGSI collection (their MuGSI strain numbers were Mu519 and Mu107, respectively.)

Cells were spread on LB agar plates containing different concentrations of various antibiotics (see below). Through serial dilutions (using 1.16% (w/v) NaCl solution), we ensured the number of colonies to be between 50 and 250 on a plate (100 × 15 mm Petri dish). The plates were then incubated at 37°C. After 18 hr of incubation, the number of visible colonies on plates was determined. As indicated in the main text, we also used casamino acids agar plates. We dissolved 2% casamino acids in N-C- minimal medium (Csonka et al., 1994) and filtered the medium. We separately autoclaved agar, and when agar was cooled and felt warm to the touch, the casamino acids medium was added.

Stock solutions of ciprofloxacin (1 mg/ml), kanamycin (50 mg/ml), 6-APA (2 mg/ml), streptomycin (25 mg/ml), ofloxacin (1 mg/ml), erythromycin (10 mg/ml), vancomycin (100 mg/ml), and cefsulodin (30 mg/ml) were prepared in sterilized water. Stock solutions of tetracycline (50 mg/ml) and chloramphenicol (10 mM) were prepared in methanol. Stock solutions of thiolutin (2 mg/ml) and rifampicin (50 mg/ml) were prepared in DMSO. Antibiotics were purchased from Biobasic Inc (Canada), Acros Organics, or Sigma-Aldrich (St. Louis, MO.).

Cells were cultured as described above without antibiotics first. At OD600 of ~0.4, 2 µL of aliquot from a culture was loaded into a pre-warmed 35 mm glass-bottom Petri dish (InVitro Scientific). A pre-warmed LB agarose pad containing antibiotics was placed over them (this procedure marked time zero). The dish was then moved into a pre-warmed (at 37°C) inverted microscope (Olympus IX83), and appropriate stage positions (typically ~50 positions) were selected. Selection of multiple state positions and setting up the software for automatic image acquisition normally took 1 ~ 2 hr. The microscope had an automated mechanical XY stage and auto-focus, and was controlled by the MetaMorph software (Molecular Devices). Also, it was housed by a microscope incubator (InVivo Scientific) which maintained the temperature of samples at 37°C during experiments. An oil immersion phase-contrast 60 × objective was used for imaging. Images were captured using a Neo 5.5 sCMOS camera (Andor). We also cultured cells in a microfluidic chemostat (using LB broth as growth medium). The detailed procedure for the microfluidic experiments was published in our previous articles (Kim et al., 2012; Deris et al., 2013). Images were analyzed using MicrobeJ, a plug-in for the ImageJ software (Ducret et al., 2016), and the analysis results were validated manually. In our experiments, all cells that stopped growing eventually underwent lysis. Although dormant cells that survive antibiotics (i.e., persisters) may complicate our analysis, their frequency is very low, ~10−5 (i.e., 1 in 105 cells) (Lewis, 2010). Because our study involves small populations (less than 100 cells), dormancy has little relevance to our study.

We first cultured cells in 5 mL LB medium to the OD600 of ~0.4, as described above. Previously, using a plate assay, we determined that OD600 of 0.4 contains 2 × 108 cells per 1 mL. Using this relation, we diluted the culture in 10 mL LB broth such that there were ~640 cells per 200 μL (we separately confirmed this density by plating cells on LB agar plates containing no antibiotics and counting colony-forming units). We supplemented this culture with various concentrations of antibiotics. We then equally distributed 200 μL of this culture to 12 chambers in a microtiter plate. We repeated this procedure by varying cell density. Next day, we measured the OD600 of each chamber.

If cells replicate at rate, λ, and die at rate, ϕ, a general deterministic model of population growth can be described by

where n is the number of cells in the population. Its solution is

Previous studies modeled the population dynamics of bacteria exposed to antibiotics using this solution; e.g., see (Regoes et al., 2004; Czock et al., 2009). At drug concentrations below the MIC (sub-MIC), growth rate is higher than death rate (λ > ϕ); thus, a bacterial population grows. When the drug concentration reaches the MIC, growth rate and death rate are equal (λ = ϕ); thus, the population size is maintained. At drugs concentrations above the MIC, growth rate is lower than death rate (λ < ϕ); thus, a bacterial population declines.

We found that bactericidal drugs induce demographic stochasticity (i.e., cells may replicate or die stochastically), and thus the population size fluctuates over time. To account for such fluctuations, we use the Markovian birth-and-death process model, which represents the most fundamental formulation of stochastic population growth (Novozhilov et al., 2006; Pavel Krapivsky and Ben-Naim, 2010; Kendall, 1948) and consider the probability of there being n cells at time t in a population, Pn(t). Temporal evolution of Pn(t) is given by the following master equation,

(A3) Pn.=λ(n−1)Pn−1−(λ+ϕ)nPn+ϕ(n+1)Pn+1.

This equation can be solved analytically; see (Novozhilov et al., 2006; Pavel Krapivsky and Ben-Naim, 2010; Kendall, 1948). For a population that originates from a single cell (e.g., a colony on agar originates from a single bacterium in our experiments), the initial condition is P1(0)=1. Then, the solution is

where A(t)=E(t)⋅(1−ϕ/λE(t)−ϕ/λ)2, a(t)=(E(t)−1E(t)−ϕ/λ), and E(t)=e(λ−ϕ)⋅t.

For n = 0,

(A5) P0(t)=ϕλ(E(t)−1E(t)−ϕ/λ),

which represents the probability that a population has gone extinct before time t.

In the long-time limit, P0(t) converges to

(A6) P0={1, if ϕ≥λϕλ, if ϕ<λ.

Because bactericidal drugs induce cell death, the death rate ϕ can be increased by using higher concentrations of drugs (Figure 2—figure supplement 6A). The solution (Equation A6) indicates that, if ϕ is equal to or greater than λ, a population will always go extinct, meaning that the lowest value of ϕ that can consistently achieve population extinction is equal to λ. This means that the lowest drug concentration that guarantees bacterial eradication, c1, is determined by ϕ (c1) = λ. In our plate assay (Figure 1B), c1 corresponds to the lowest antibiotic concentration at which no cells form colonies, that is, MIC; thus c1 = MIC. Alternatively put, at the MIC, ϕ = λ; see the main text for details.

Above the MIC (ϕ > λ), all populations eventually undergo extinction (Equation A6), but the extinction time cannot be predetermined. Equation A5 predicts that as time goes by, the probability of population extinction converges to P0(t)≈1+(λϕ−1)⋅e(λ−ϕ)⋅t. In our experiments, we have counted the number of live colonies, i.e., colonies that have not gone extinct yet (i.e., n≠ 0). The probability of finding live colonies at time t, Pn≠0(t), which is equal to 1−P0(t), is

(A7) Pn≠0(t)=1−P0(t)≈(1−λϕ)⋅e(λ−ϕ)⋅t.

Thus, at later times, the probability of finding live colonies can be approximated by an exponential decay.

Using Equation A3, moment dynamics can be derived (Novozhilov et al., 2006; Pavel Krapivsky and Ben-Naim, 2010; Kendall, 1948). The first moment, ∑n=0∞nP, is equal to the mean of n, ⟨n⟩, and the equation for the first moment is

In our experiments, the initial condition, n(0), is the same across the populations. For example, for all our microscope experiments, n(0)=1.

Its solution is

(A9) ⟨n(t)⟩=n(0)⋅e(λ−ϕ)⋅t.

These equation and solution are the same as the deterministic equation and solution we discussed above (Equation A1 and A2); therefore, the mean of the stochastic variable n obeys the deterministic model.

Next, we consider the mean number of dead cells ⟨ nD ⟩. Because live cells die at the rate of ϕ, we have

Solving Equation A10 together with Equation A9, we have

(A11) ⟨nD(t)⟩=ϕ⋅⟨n(0)⟩λ−ϕ⋅[e(λ−ϕ)⋅t−1].

or alternatively,

(A12) ⟨nD(t)⟩=ϕλ−ϕ⋅[⟨n(t)⟩−n(0)].

We can also write an equation for the second moment, ∑n=0n2Pn,

(A13) d⟨n2⟩dt=2(λ−ϕ)⟨n2⟩+(λ+ϕ)⟨n⟩.

The solution is

(A14) ⟨n2⟩=n(0)λ+ϕλ−ϕe(λ−ϕ)t[e(λ−ϕ)t−1]+n(0)2e2(λ−ϕ)t.

Using Equation A9 and S14, we can calculate the variance Var[n]=⟨ n2 ⟩−⟨ n ⟩2,

(A15) Var[n]=n(0)λ+ϕλ−ϕ⋅e(λ−ϕ)⋅t[e(λ−ϕ)⋅t−1].

Next, we quantified population fluctuations by calculating the coefficient of variation (CV), that is, standard deviation relative to the mean. The standard deviation is equal to Var[n] and can be calculated using Equation A15. The mean is given by Equation A9. Thus, we have

(A16) CV=Var[n]⟨n⟩≈λ+ϕλ−ϕ⋅n(0)−1/2.

Equation A16 holds for a sub-MIC region (λ > ϕ) and considers the long-time limit. It indicates that the magnitude of fluctuations in population size depends on inoculum size (n(0)), and the rates of growth and death.

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

We thank Jeff Gore, Véronique Perrot, and Bruce Levin for helpful comments on the manuscript. We also thank Emily Crispell and David Weiss for sharing strains with us. This work was funded by Research Corporation for Science Advancement (#24097) and the Human Frontier Science Program (RGY0072/2015).

1. Aleksandra M Walczak, Ecole Normale Superieure, France

© 2018, Coates 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.