Growth-dependent drug susceptibility can prevent or enhance spatial expansion of a bacterial population

As a population wave expands, organisms at the tip typically experience plentiful nutrients while those behind the front become nutrient-depleted. If the environment also contains a gradient of some inhibitor (e.g. a toxic drug), a tradeoff exists: the nutrient-rich tip is more exposed to the inhibitor, while the nutrient-starved region behind the front is less exposed. Here we show that this can lead to complex dynamics when the organism’s response to the inhibitory substance is coupled to nutrient availability. We model a bacterial population which expands in a spatial gradient of antibiotic, under conditions where either fast-growing bacteria at the wave’s tip, or slow-growing, resource-limited bacteria behind the front are more susceptible to the antibiotic. We find that growth-rate dependent susceptibility can have strong effects on the dynamics of the expanding population. If slow-growing bacteria are more susceptible, the population wave advances far into the inhibitory zone, leaving a trail of dead bacteria in its wake. In contrast, if fast-growing bacteria are more susceptible, the wave is blocked at a much lower concentration of antibiotic, but a large population of live bacteria remains behind the front. Our results may contribute to understanding the efficacy of different antimicrobials for spatially structured microbial populations such as biofilms, as well as the dynamics of ecological population expansions more generally.

Spatial gradients of environmental parameters play a central role in the ecology and evolution of both macro-and microorganisms [1][2][3][4][5][6][7][8][9][10]. A particularly important scenario is one in which a population wave is halted in its spatial expansion by an increasingly challenging environment, to which it must evolve resistance in order to advance further. Population genetic theory has shown that the sharpness of the boundary of the region colonised (the species' range) [11] depends on a complex interplay between flow of genetic variants from the main body of the population to its edges [12,13], evolution of genetic variance within the population [14], genetic drift [15][16][17] and the potential for interspecific competition [18]. In the specific context of a bacterial population which expands in a spatial gradient of an antimicrobial chemical (e.g. an antibiotic), experiments [19,20] and theory [21][22][23][24] show that the antimicrobial gradient can speed up the evolution of resistance; further theoretical work suggests that this depends both on the nature of the mutational pathway to resistance [23] and on the balance between the bacterial migration and mutation rates [22,24,25]. Similar spatial effects have been observed in the emergence of resistance to cancer therapy [26,27].
As a population wave expands in an environment that contains a gradient of an inhibitory substance, individuals within the population experience changes in both inhibitor and nutrient abundance. At the tip of the wave, organisms are exposed to high levels of the inhibitor, but nutrients are typically abundant because the population density is low. In contrast, behind the wave front, the environment is less inhibitory but organisms can become nutrient-limited because the population density is high. This scenario is especially interesting if a coupling exists between the organism's susceptibility to the inhibitor and the nutrient availability. For example, for bacteria exposed to antibiotics, it is known that nutrient conditions can strongly influence susceptibility to antibiotics. For most antibiotics, fast-growing Escherichia coli (on rich nutrient media) are more susceptible than slow-growing E. coli 2 P Sinclair et al (on poor nutrient media) [28][29][30][31][32]. For some antibiotics, however, the opposite is true: slow-growing bacteria are more susceptible [29][30][31]33]. Recent work by Greulich et al suggested that antibiotics within the same target class can be more effective on either fastor slow-growing bacteria depending on molecular parameters such as target-binding and transport rates [33,34]. For Pseudomonas aeruginosa growing in surface-associated biofilms, it is known that antibiotics such as ciprofloxacin, tetracycline and tobramycin are more active against the fast-growing subpopulation at the biofilm surface [35], while membrane-targeting antimicrobials such as the peptide colistin, ethylenediaminetetraacetic acid (EDTA) and sodium dodecyl sulfate (SDS) are more effective against the nongrowing subpopulation within the biofilm [35][36][37]. It remains unclear whether this classification of antimicrobials according to their growth-rate dependence is linked to the more traditional classification into bactericidal drugs, which cause bacterial death, and bacteriostatic drugs, which inhibit bacterial growth without killing [38].
Within a bacterial population wave, growthdependent susceptibility implies that fast-growing bacteria at the wave tip may show different susceptibility to antibiotics than slow-growing bacteria behind the wave front-in other words, there may be a coupling between the ability of the population to invade the antibiotic gradient and nutrient limitation within the population wave. Here, we investigate the implications of such a coupling, by modelling the spread of a bacterial population wave in a gradient of both bacteriostatic and bactericidal antibiotics, under conditions where either fast-or slow-growing bacteria are more susceptible to the antibiotic. Our simulations suggest significant differences in outcome in these two cases. If the fast-growing bacteria are more susceptible, the population wave stops at a low concentration of antibiotic, but, for a bactericidal antibiotic, many bacteria remain alive behind the front. If the slow-growing bacteria are more susceptible, the wave spreads much further, but may trail a wake of dead bacteria behind it. Our work reveals basic principles that may have implications for treating bacterial infections growing in spatially structured environments, such as biofilms.

Model
To model the expansion of a bacterial population in a spatial antibiotic gradient, we follow the work of Greulich et al [10,23] by discretising space as a chain of L microhabitats (figure 1), each of which can contain a different concentration of antibiotic. We assume that bacteria can proliferate within each microhabitat, consuming nutrient, can migrate to neighbouring microhabitats and can die (if the antibiotic is bactericidal). The variables in our model are the number N i of bacteria and the concentration s i of nutrient, for each microhabitat i ∈ [0, L]. The antibiotic concentration c i remains fixed for each microhabitat, as does the migration rate m (which is assumed to be the same for all microhabitats).
For a given microhabitat, the local bacterial growth rate b i depends on both the local nutrient concentration s i and the local antibiotic concentration c i : (1) In equation (1), g(s i ) is a Monod function describing the nutrient-concentration dependence of the growth rate [39]: K being the Monod constant. The function g increases with nutrient concentration s i and saturates for (1) describes bacterial inhibition by the antibiotic, with a minimum inhibitory concentration (MIC) β i . Importantly, growthdependent susceptibility means that the MIC is nutrient-dependent: β i depends on the local nutrientdependent growth function g(s i ) ( figure 1(b)). In this work, we consider two contrasting cases: 'fast-growth targeting antibiotics' (FGTA), for which fast-growing bacteria are more susceptible, and 'slow-growth targeting antibiotics' (SGTA), for which slow-growing bacteria are more susceptible. In both cases (and inspired by [33]), we assume a linear relation between MIC and growth rate: for the FGTA, the MIC decreases with growth rate: while for the SGTA the MIC increases with growth rate: (see figure 1 as well as the note on units and parameters below). We also model the situation where the MIC is independent of growth rate, which we term 'growthindependent antibiotic', or GIA; here we set β GIA = 5.
To model a bacteriostatic antibiotic, we assume a quadratic form for the pharmacodynamic function φ(c i , β i ) [23]: For a bactericidal antibiotic, we use the following pharmacodynamic function: corresponding to the general function proposed by Regoes et al [40] with parameters κ = 2, ψ max = 1 and ψ min = −5. If the local antibiotic concentration c i < β i , equation (6) is used as input to equation (1) to compute the local growth rate. For higher concentrations of the bactericidal antibiotic, c i > β i , φ becomes negative and bacteria do not grow but instead die, at rate We simulate the model using the stochastic Monte Carlo algorithm described by Greulich et al [23] (see Methods). All our simulations are initiated with a fixed nutrient concentration s max in all microhabitats, and with N 0 bacteria in the first microhabitat (and zero elsewhere). All model parameters are given in the Methods section. A note on units and parameters-our model aims to be generic, rather than to represent particular antibiotic-bacteria combinations. Nevertheless, it is important to consider the units of time, space, nutrient and antibiotic concentration in the model. The units of time in our model can be considered to be hours (this is consistent with the choice of g max = 1 h −1 , a typical growth rate for E. coli on minimal lab media [41]). We suppose that the volume of a microhabitat is ∼1 µl, or equivalently its size is ∼1 mm. For the migration rate, our chosen value of m = 0.1 then corresponds to an effective diffusion constant of m 2 , in units of mm 2 h -1 which is approximately 3 µm 2 s −1 [42]. We define the units of nutrient concentration in terms of the yield, or amount of nutrient needed to create a single bacterial cell. Thus the chosen value s max = 500 implies that the maximal bacterial density that can be reached in the model is 500 bacteria per microhabitat, or ∼5

For bacteriostatic antibiotics, growthdependent susceptibility affects population expansion Expansion in a uniform antibiotic concentration
We first study how a bacterial population invades an environment with a spatially uniform concentration of a bacteriostatic antibiotic (c i = c for i ∈ [0, L]). In the absence of growth-dependent susceptibility, we expect the model, in the continuum limit, to map onto the Fisher-KPP equation [23,[43][44][45][46]: where x is a spatial coordinate, n is the population density, D is an effective diffusion constant arising from migration between microhabitats, b max is the maximal growth rate and n max is a carrying capacity. Equation (7) has travelling wave solutions with wave speed 2 √ Db max n max [43,44]. In our model, the maximal growth rate scales with the pharmacodynamic function φ(c, β), which for the bacteriostatic antibiotic, obeys equation (5). Therefore we expect the speed v of spread of the bacterial population in our model to scale The concentration of antibiotic can be different in each microhabitat (and is fixed throughout the simulation). We simulate the dynamics of bacteria and nutrient within each microhabitat. As the bacteria multiply, the nutrient is consumed, which alters the growth rate of the bacteria. Growth-dependent susceptibility means that this in turn affects the minimal inhibitory concentration (MIC). Within the population, bacteria at the edge of the expanding population wave experience plentiful nutrients and grow fast, whereas those behind the wave are nutrient-depleted and grow slowly. (b) Growth-dependent susceptibility as implemented in our model. The MIC is plotted versus the Monod growth function g(s) (equation (2)), which depends on the nutrient concentration s, for the three types of antibiotics studied: fast-growth targeting antibiotics (FGTA), growth-independent antibiotics (GIA) and slow-growth targeting antibiotics (SGTA).
For concentrations above the MIC, the bacteria do not grow and so there is no wave propagation. Figure 2 shows that, for the GIA, the population indeed expands as a travelling wave for c < β GIA = 5 (figures 2(b), (e) and (h)), and that the speed v scales as expected (figures 2(j) and (k)).
Growth-dependent susceptibility changes the results quantitatively but not qualitatively. For both the fast-growth targeting antibiotic (FGTA) and the slow-growth targeting antibiotic (SGTA), our simulations show that the population expands as a travelling wave (figure 2), but the growth-dependent susceptibility alters the range of antibiotic concentrations over which invasion happens. For the FGTA, population expansion happens only for c < 1 (figures 2(a), (d) and (g)), whereas for the SGTA we observe expansion even for c = 5 (figures 2(c), (f) and (i)). Thus, for a uniform antibiotic concentration, the fast-growth targeting antibiotic is much more effective at preventing bacterial invasion.
This observation can be rationalised by noting that the propagation of Fisher-KPP waves is governed by the dynamics at the edge of the wave front [46,47]; i.e. these are 'pulled waves' [46][47][48]. The edge of the wave front is also where the population density is lowest, hence nutrients are abundant, and the local growth rate is highest. For the FGTA, the MIC, β FGTA i , decreases with the nutrient concentration, implying that bacteria at the edge of the wave . For clarity, every fourth snapshot is shown in black. Panels (j) and (k) show the average speed of expansion of the population wave as a function of antibiotic concentration, for the three antibiotic types. The dependence of wave speed on concentration is quadratic (inset to panel (k)), consistent with our expectation from Fisher-KPP wave theory. The wave travels faster in the SGTA case and slower in the FGTA case, but all the plots collapse onto one line when the antibiotic concentration is scaled by the value of the MIC for maximal nutrient concentration, β max ≡ β(s max ) (main plots in panels (j) and (k)). are maximally inhibited by the antibiotic. In contrast, for the SGTA, the MIC, β SGTA i , increases with the nutrient concentration, implying that bacteria at the edge of the wave are minimally inhibited by the antibiotic. This suggests that the invasion dynamics with growth-dependent susceptibility should map onto Fisher-KPP theory, but with a wave speed controlled by the MIC value at the wave edge. Figures 2(j) and (k) show that this is indeed the case; scaling the antibiotic concentration by β FGTA (s max ) = 1 for the FGTA and by β SGTA (s max ) = 10 for the SGTA results in a collapse of all of our simulation data onto a universal curve. Thus, in our model, the significant quantitative difference in the invasion of a uniform environment in the presence of FGTA and SGTA antibiotics is caused by the contrasting antibiotic susceptibilities of bacteria at the nutrient-rich edge of the travelling wave.

Bacterial expansion in an antibiotic gradient
When a bacterial population wave expands into a region with a gradient of antibiotic, we expect the wave to slow down and eventually stop when the population reaches an antibiotic concentration equal to its MIC [23]. Here, following Greulich et al [23], we assume an exponentially increasing concentration of antibiotic: This antibiotic profile is plotted in figure 3(a) for α = 0.0049 (purple line). This corresponds to a rather long decay length of ∼1/α ∼2 m (similar to that of the giant gradient plate of Baym et al [20]), but allows us to see clearly what is going on in our simulations (since the antibiotic concentration for the last microhabitat (i = 500) is equal to c 500 = 10.6, which is only just greater than our maximal MIC value). Figure 3(a) also shows the concentration profile for α = 0.02 (pink line; decay length ∼5 cm); we show in the supplementary material (stacks.iop.org/PhysBio/16/046001/mmedia) that we obtain qualitatively the same results for this steeper antibiotic gradient. Figure 3(b) (orange line) shows that, without growth-dependent susceptibility (i.e. the GIA), the population wave is indeed halted by the antibiotic gradient (the long timescales here arise from the shallowness of the gradient; in the supplementary material, figure S1, we show that the wave halts much earlier for a steeper gradient). Figure 3(d) shows the wave profiles for the GIA, plotted at different times; once the wave has stopped, further spread of the population occurs simply by diffusion. Figure 3(g) further illustrates the underlying cause of the halt in expansion: the growth rate b i decreases in front of the wave due to the increasing antibiotic concentration (note that the growth rate also decreases behind the wave due to nutrient depletion).
We expect growth-dependent susceptibility to alter this picture. For the FGTA, antibiotic susceptibility is greater at the advancing tip of the wave, due to the abundant nutrients. This leads to a lower growth rate at the wave tip (figure 3(f)), and therefore the wave advance is inhibited at a lower antibiotic concentration than for the GIA ( figure 3(c)). This can also be seen when we plot front position as a function of time ( figure 3(b)). We obtain contrasting results for the SGTA. Here, antibiotic susceptibility is reduced at the front of the wave, where nutrients are abundant, and increased behind the wave front, where nutrients are scarce (figure 3(h)). Figures 3(b) and (e) show that for the SGTA, the wave spreads further before being inhibited by the antibiotic.
In summary, our simulations show that growthdependent susceptibility can have important effects on bacterial population expansion in a drug gradient. If the antibiotic is more effective for fast-growing bacteria, we expect strong inhibition at the tip of the population wave, greatly limiting the population's ability to colonise the gradient. However, if the antibiotic is more effective for slow-growing bacteria, the antibiotic gradient will fail to inhibit bacteria at the wave tip, and the population will advance further into the gradient before being stopped.
Our conclusions are independent of the shape of the pharmacodynamic function Our choice of pharmacodynamic function φ(c, β) is somewhat arbitrary. The threshold-like form of equations (5) and (6) is realistic for some antibiotics, such as ciprofloxacin, rifampicin and streptomycin at sub-MIC concentrations [33,40]. However, other antibiotics such as tetracycline and chloramphenicol show smoother, Langmuir-like pharmacodynamic functions, for which the general form φ(c, β) = 1 1 + c β (9) has been found to provide a good fit to data [33].
To test the effect of the pharmacodynamic function, we repeated our simulations using the smoother form of φ(c, β) as in equation (9), but retaining the same relationships between the MIC and the nutrient concentration, equation (3) for the FGTA and (4) for the SGTA. Simulating this model, we obtain qualitatively similar results (figure 4); for the FGTA wave expansion is inhibited earlier, whereas for SGTA the wave is able to expand further into the antibiotic gradient. However, the difference in pharmacodynamic function causes a quantitative difference in the spread of the bacterial population. Antibiotics with a Langmuir-like pharmacodynamic function are generally less effective at inhibiting bacterial expansion than antibiotics with a threshold-like pharmacodynamic function, independent of their growth-dependent susceptibility.
It is important to note that the shape of the pharmacodynamic curve may actually be coupled to the growth-dependence drug susceptibility. For example, the model of Greulich et al [33] suggests that, for ribosome-targeting antibiotics, a smoothly decreasing, Langmuir-like inhibition function is intrinsically linked to the FGTA scenario, whereas a threshold-like inhibition function is a property of the SGTA scenario. Therefore, a comparison of the dynamics of a bacterial population in gradients of, for example, streptomycin versus tetracycline, might be complicated (see supplementary material). Understanding more generally the link between growth-dependent susceptibility and the shape of the pharmacodynamic function would be an interesting topic for future research.

For bactericidal antibiotics, more killing is achieved for SGTA than FGTA
We now investigate what happens for bactericidal antibiotics, i.e. those that are able to kill bacteria, rather than just preventing growth. To model a bactericidal antibiotic, we use the pharmacodynamic function φ of equation (6), and implement bacterial death at rate d = −φ when φ < 0. This leads to a qualitatively different outcome in our simulations (figure 5). Focusing initially on the fast-growth targeting antibiotic (FGTA), we observe that bacteria behind the wave front remain alive (figures 5(a)-(c)), although mostly non-growing due to nutrient depletion (figure 5(g)), whereas bacteria at the edge of the wave are dead.
For the slow-growth targeting antibiotic (SGTA), the population wave advances further into the antibiotic gradient than for the FGTA, as is the case for the bacteriostatic antibiotic. However, the population composition within the wave is very different. For the bactericidal slow-growth targeting antibiotic, we observe a 'striped' population structure (figure 5): bacteria at the very edge of the wave are dead, but there is a zone of live bacteria just behind the tip, followed by a zone of dead bacteria and finally, a zone of live, but non-growing bacteria at the very back of the wave. To verify that our key observation of the 'striped' population structure is robust to such discretisation effects, we also simulated a continuum model for the bacteri-cidal antibiotic. Our continuum model consists of the following equations for the population density n(x, t) and nutrient concentration s(x, t): and where D is a diffusion constant, related to the migration rate in our discrete model (see Methods), and is a random variable with mean zero, uniformly distributed in the range ±0.1φ max (c, β). Figure 6 shows the results of our continuum simulations, for the FGTA and SGTA antibiotics. These results have the same qualitative features as those of the discrete model: the population wave advances further into the antibiotic gradient for the SGTA than for the FGTA, but for the SGTA a 'stripe' of dead bacteria trails behind the wave tip. We note that the noise term in equations (10) and (11) is necessary to obtain this behaviour; briefly, this is because the model has a 'tipping point' when the local nutrient concentration is such that the pharmacodynamic function φ = 0. For a bactericidal antibiotic, when the nutrient concentration decreases beyond this point, bacterial death kicks  1 Briefly, in our model, a crossover between growth and death happens at the position i * where φ i (c i , β i ) = 0. Using equations (6) and (8), this implies that i * = α −1 ln [2 + 9(g(s i )/g max )]. However, the nutrient concentration is measured in discrete units of the bacterial yield (i.e. the amount of nutrient needed to produce one bacterium). Starting from the back of the wave, where nutrient is depleted: s i = 0, we encounter a value of i * = α −1 ln [2 + 9(g(0)/g max )] = 141: here bacteria begin to die. Beyond this point, however, the nutrient concentration is higher, and so the bacterial susceptibility is lower: for s i = 1, the value of i * shifts to i * = α −1 ln [2 + 9(g(1)/g max )] = 168. Even further along, the nutrient concentration is even higher, and the susceptibility is even lower: i * (s i = 3) = 191, and so on. Thus, the nutrient concentration in our model undergoes a series of discrete steps close to the wave edge, each of which shifts the threshold between growth and death. Each of these steps is associated with a peak in the number of dead bacteria. The larger size of the first peak, at i = 141, can be explained by migration of live bacteria from the 'live zone' into the 'dead zone'. in, creating the striped population structure. Noise acts to nudge the system over this threshold, while in a purely deterministic version of the model the system simply remains poised at the tipping point (for further discussion of this point, see the supplementary material). While the form of the noise as implemented in equations (10) and (11) is somewhat arbitrary, we expect that stochastic factors that could perform this role would be present in any real scenario.
Returning to the discrete model, our results for the FGTA and SGTA bactericidal antibiotics can be summarised by tracking the total bacterial population size, as well as the percentage of the total population that is alive. Figure 5(i) shows that for the FGTA, the total bacterial population remains quite low (as the population wave does not advance far into the gradient), but almost all the bacteria are alive. In contrast, for the SGTA ( figure 5(j)), the total population size is much higher but a large fraction of these bacteria are dead. This situation suggests a tradeoff: an FGTA antibiotic is optimal where the total bacterial population is the major issue (e.g. to reduce biofilm thickness and hence the potential for inflammatory responses or clogging, or in the case where the biofilm releases toxins), while an SGTA antibiotic may be a better choice for maximising the percentage killing (which is likely to be important for propagation and eradication of infections). , dead bacteria appear at the wave edge once the wave's advance halts. For the SGTA (d)-(f), a 'striped' population structure emerges, in which blocks of live and dead bacteria alternate, with a large 'dead zone' of dead bacteria behind the wave front. Panels (g) and (h) show snapshots of the local growth rate b i for the FGTA and SGTA respectively (in panel (g), the red line lies on top of the orange and green lines). Panels (i) and (j) summarize the situation by tracking the percentage of the population which is alive (green lines) and the total population size (live or dead; purple line), for the FGTA (i) and SGTA (j). Although the total population is larger for the SGTA, the absolute number of live bacteria is similar for the two antibiotics. Therefore there is a trade off between having a small, yet active population (FGTA), or a larger, but more inert one (SGTA).

Discussion
Antibiotics have come to play a key role in modern medicine since their introduction in the first half of the 20th century. Despite this, large gaps remain in our basic understanding of how they work: closing these gaps could lead to better treatment strategies [49][50][51][52][53] and better approaches to counteract the widespread emergence of antimicrobial resistance (AMR) [54,55].
One such gap concerns how bacterial physiology and antibiotic action come together to determine the shape, and the growth-rate dependence, of the pharmacodynamic function [33,34,52,54,56]. Spatial heterogeneity represents another gap: we do not know how the spatial structure of an infection affects its susceptibility to treatment by antibiotics, and its potential for evolving resistant mutant bacteria [19, 22-24, 35, 57, 58]. Here, we bring together these two themes to show that growth-rate dependent antibiotic susceptibility can significantly affect how a bacterial population spreads in a spatial gradient of antibiotic. We show that antibiotics that target fast-growing bacteria are better at preventing bacterial population invasion that antibiotics that target slow-growing bacteria. However, for bactericidal antibiotics, a tradeoff exists: FGTA antibiotics lead to a smaller total population which is mostly alive, while SGTA antibiotics lead to a larger total population, a significant fraction of which are dead. Depending on the situation, either one of these two outcomes might be preferable.
From a fundamental perspective, our work connects with a long-established body of theoretical work on travelling waves [46,59] and more recent work on their application to bacterial systems [23,48,60,61]. Our results can be rationalised in the context of Fisher-KPP wave theory, and the fact that Fisher-KPP waves are 'pulled waves', being governed by what happens at the wave's edge [46,47]. In our model, growth-dependent susceptibility leads to large differences in MIC between the FGTA and SGTA antibiotics at the wave edge, where nutrients are abundant. This in turn leads to differences in wave speed, and in the extent to which the population manages to advance into an antibiotic gradient.
From a practical perspective, our model is simplistic, but its results could have implications for the treatment of biofilm infections [62]. In biofilms, a population of bacteria is embedded in a polymeric matrix and adheres to a surface. Biofilms are notoriously difficult to eradicate with antibiotics, for reasons that remain unclear, but may include lower penetrance of antibiotics, lower growth rate, genetic mutations and phenotypic persistence [57,[62][63][64]. Interestingly, Pamp et al [35] observed that antibiotics such as ciprofloxacin, tetracycline and tobramycin were more effective against the fast growing subpopulation at the surface of a P. aeruginosa biofilm, while colistin, EDTA and SDS were more effective against the non-growing subpopulation within the biofilm. Our model suggests that FGTA antibiotics might be more effective at the biofilm surface, while SGTA antibiotics might be better at eradicating the slow-growing biofilm core.
Our model makes a number of important assumptions. Spatial structure in our model is reduced to a one-dimensional (1D) chain of microhabitats, exposed to an exponentially increasing concentration of antibiotic: features such as the complex structure of a biofilm [65,66], bacterial chemotaxis and fluid flow [67] are neglected [10]. Real infections are subjected to antibiotic gradients in two or three dimensions; it would be interesting to extend our models to more complex geometries. We also make assumptions about bacterial physiology. Perhaps most importantly, we have transferred insights on growth-rate dependent susceptibility obtained from experiments where growth rate is controlled by nutrient richness [68][69][70][71], to the situation where nutrient concentration controls growth rate. In other words, we have assumed that a slow growth rate caused by nutrient depletion has a similar effect on antibiotic susceptibility as a slow growth rate caused by a poor nutrient supplied in abundance. We are not aware of any experimental tests of this hypothesis, although early work by Neidhardt and Magasanik [72] suggests that the RNA/protein ratio (indicative of ribosome abundance) achieved by E. coli in stationary phase on glucose media matches that achieved during steady state growth on very poor nutrients, suggesting that our assumption may be reasonable for ribosome-targeting antibiotics. More data on this topic would be extremely useful. More broadly, we have assumed that the bacterial response to antibiotic is uniquely determined by the nutrient and antibiotic concentrations; in reality, one might expect to see heterogeneity in response among individual bacteria in the population: for example some might be dying while others are still growing, for a given concentration of antibiotic [54,73]. It would be interesting to investigate the consequences of this in future work.
In our model, we have also assumed simple linear relations for the coupling between growth rate and MIC: while this seems to be appropriate for some antibiotics, it is not realistic for others [33]. Our conclusions remain unchanged when we implement non-linear relations between MIC and growth rate; see the supplementary material). We have also assumed that the MIC depends only on the nutrient-dependent component of the growth rate g(s); if it were to depend on the full growth rate g × φ, interesting feedback effects might be obtained. Furthermore, as discussed earlier and in the supplementary material, we have ignored the likely link between growth-rate dependent susceptibility and the shape of the pharmacodynamic function. Because antibiotics with different growth-rate dependent susceptibilities may also have different pharmacodynamic functions, testing exper imentally our model prediction that bacteria exposed to FGTA antibiotics advance less far into an antibiotic gradient than those exposed to SGTA antibiotics may not be a simple task. However, simple experiments might be able to test the prediction of a 'striped' population structure for bactericidal SGTA antibiotics (but not for FGTA bactericidal antibiotics). If confirmed, this live-dead striping phenomenon might be of practical importance in understanding the response of biofilms to antibiotics.
Understanding quantitatively how antibiotics work is important for optimising antibiotic usage and preventing the emergence of AMR. Our work suggests that both growth-dependent susceptibility and spatial structure can play important roles. FGTA antibiotics are predicted to be more effective than SGTA antibiotics at preventing the spatial spread of a bacterial population (e.g. biofilm growth). However, if the goal is to treat an existing population (e.g. a biofilm), bactericidal SGTA antibiotics are predicted to be more effective at killing the slow-growing core.

Model parameters
The number L of microhabitats in our simulations was set to L = 500, following the work of Greulich et al [23], which in turn was loosely based on the experimental protocol of Zhang et al [19]. The migration rate was set to m = 0.1, also following Greulich et al [23], which is small compared to the maximum growth rate b i ≈ 1. As discussed in the main text, this corresponds to a rather low effective bacterial diffusion constant, or to a microfluidic device like that of Zhang et al, in which microhabitats are connected through small channels [19]. The maximal growth rate was set to b max = 1, which matches roughly the growth rate of E. coli on glucose minimal media, if the time units are taken to be hours (see main text). The Monod constant K in the growth equation (2) was set to K = 33 concentration units. For E. coli growing on glucose minimal media, the Monod constant is approximately 1 µm, or equivalently 6 × 10 17 molecules per litre [74,75]. Our units of nutrient concentration are the nutrient needed to make one bacterium, per microhabitat volume. The amount of glucose consumed to make one E. coli bacterium, or the yield, is approximately 1.8 × 10 10 molecules [74,75], while we assume the microhabitat volume to be 1 µl. Therefore K translates into 3.3 × 10 7 bacterial yields per litre, or 33 bacterial yields per microlitre.
In our simulations, the initial nutrient concentration in each microhabitat was set to s max = 500 units, and each bacterial replication event consumed 1 unit of nutrient. Each simulation was initiated with N 0 = 100 bacteria in the first microhabitat.

Simulation algorithm
To simulate our model, we used the Monte Carlo algorithm introduced by Greulich et al [23], which is comparable to the well-known Gillespie algorithm [76,77]. First, a bacterium was chosen at random from the current population (N tot ). Next, the rates for birth (b i ), migration (m) and death (d i , in the case of a bactericidal antibiotic) were calculated for the selected bacterium. Then a random number, r, between 0 and a value R max was generated, where R max > b i + m + d i . In our simulations, we used R max = 1.2 for the bacteriostatic antibiotic and R max = 5.2 for the bactericidal antibiotic. The random number r was used to select a migration event, replication event, a death event or no event, with probability proportional to the rates for these events. Once the chosen event had been executed, the time t was increased by ∆t = 1/N tot R max .

Continuum model for the bactericidal antibiotic
For the bactericidal antibiotic, we simulated the system using the continuous reaction-diffusion population dynamics model defined by equations (10) and (11). The noise was taken to be a random variable with mean zero, uniformly distributed in the range ± 0.1φ max (c, β). The diffusion constant D was chosen to be equivalent to the migration rate in our discrete model. Mapping a discrete 1D diffusion model with lattice size δx and migration rate m onto a continuum model gives D = m(δx) 2 . In our discrete model the microhabitat size is assumed to be δ = 1 mm (see Model section), and m = 0.1 h −1 ; therefore D = 0.01 mm 2 h −1 ≈ 3 µ m 2 s −1 . The model was solved using a simple forwardtime central-space finite difference method, with δt = 0.01 h and δx = 0.1 mm, on a spatial grid of 5000 elements, with zero-flux boundary conditions at both ends.