Bacterial persistence is essential for susceptible cell survival in indirect resistance, mainly for lower cell densities

Antibiotic-susceptible bacteria may survive bactericidal antibiotics if other co-inhabiting bacteria detoxify the medium through antibiotic degradation or modification, a phenomenon denominated as indirect resistance. However, it is unclear how susceptible cells survive while the medium is still toxic. One explanation relies on the speed of detoxification, and another, non-exclusive explanation, relies on persistence, a state of bacterial dormancy where cells with low metabolic activity and growth rates are phenotypically tolerant to antibiotics and other cytotoxic substances. Here we simulated the fate of susceptible cells in laboratory experiments in the context of indirect resistance to understand whether persistence is necessary to explain the survival of susceptible cells. Depending on the strain and experimental conditions, the decay of persister populations may follow an exponential or a power-law distribution. Therefore, we studied the impact of both distributions in the simulations. Moreover, we studied the impact of considering that persister cells have a mechanism to sense the presence of a toxic substance–a mechanism that would enable cells to leave the dormant state when the medium becomes nontoxic. The simulations show that surviving susceptible cells under indirect resistance may originate both from persister and non-persister populations if the density of detoxifying cells is high. However, persistence was necessary when the initial density of detoxifying cells was low, although persister cells remained in that dormancy state for just a few hours. Finally, the results of our simulations are consistent both with exponential and power-law decay of the persistence population. Whether indirect resistance involves persistence should impact antibiotic treatments.


Introduction
Susceptible bacterial populations do not perish instantaneously in the presence of bactericidal antibiotics. Instead, they decay exponentially, typically for a few hours for wild-type bacterial strains. After this first period, a second population has a significantly lower death rate [1,2]. These cells are in the persistent state and usually account for less than 1% of the original bacterial community. The persistent cells do not harbor resistance genes, but they thrive in the presence of a drug or other harsh environments by lowering their metabolic activity and growth rate. Importantly, they are responsible for recurrent and chronic infections because of their ability to resume growth following antibiotic therapy [2][3][4]. The ubiquitous distribution of persistence among bacteria, fungi, and cancer cells [2,5,6], together with its impact on antibiotic resistance development among bacteria [7,8], highlights the need for a better understanding of the role of bacterial persistence in the survival of pathogenic bacteria. Persistence is expected to be involved in indirect resistance, mainly for longer distance between detoxifying and susceptible cells [9]. In fact, during indirect resistance, susceptible cells are protected against a bactericidal antibiotic because other co-inhabiting bacterial cells detoxify the medium through antibiotic degradation or modification [10,11]. Once the environment becomes nontoxic, cells that leave the persistent state survive and thrive [9]. Indirect pathogenicity is an alternative name for indirect resistance because, in many cases, it involves antibiotic-susceptible pathogenic bacteria and cells from a non-pathogenic bacterial species that detoxify the environment enabling the growth of the pathogens (see, for example, refs. [10,12,13].
In the context of indirect resistance, the chances of survival of susceptible cells depend on several factors, not just entering into the dormant state of persistence. For example, medium detoxification certainly takes some time to be completed. Moreover, the survival of susceptible cells depends on detoxifying cells' density and the total cell density. High population density implies that susceptible and resistant cells tend to be close neighbors, increasing the odds of susceptible cells [14][15][16]. Furthermore, the survival of susceptible cells should depend on the death rate of non-persister cells in the presence of antibiotics and, importantly, on the persistent cells' behavior and death rate.
While the bactericidal antibiotic is still present, most cells returning to growth dies. As mentioned above, in the presence of a bactericidal antibiotic, the non-persistent population of wild-type strains seems to decline exponentially (that is, according to exp(-k.t), where t is time, and k is a constant). The exponential decay is a direct consequence of the homogeneity of bacterial clonal populations and involve many independent entities (cells), each having the same constant probability per unit of time of starting growth. After some time of decaying exponentially and fast, a second phase begins. In this phase, only persister cells are alive. This population also decays, although at a lower rate [1], because persister cells also tend to die when returning to growth if the medium is still toxic.
Several studies suggest that the decay of the persister population follows an exponential (see, for example, the case of the hipA7 mutant, where the surviving population follows a descending straight line in a log-linear scale [17]). However, other works have suggested that the persistence state results from different faults and errors in cell division rather than an evolved genetic program [18][19][20][21]. In that case, the persistent population is physiologically heterogeneous, comprising several sub-populations, each with its proper exponential decay. A recent paper mathematically shows that the sum of all these negative exponentials results in a power-law curve (i.e., proportional to t β , where t is time and β is a negative exponent), which is not an exponential curve [22]. Importantly, this reasoning also tells us that the exponent of the power-law decay, β, should be close to -2. In other words, the population of persister cells decayed proportionally to close to t -2 , where t is time. However, for the cases where persister populations decay exponentially, there is no way to calculate the exponential decay rate (that is, there is no way to compute the value of the constant k in exp(-k.t)). This prediction was experimentally corroborated in the same paper by showing that the persister population of E. coli K12 NCM3722 ΔmotA strain declined close to t -2 [22]. Fig 1 shows what happens to a clonal drug-sensitive bacterial population exposed to a cytotoxic drug, as well as the relative importance of exponential and power-law decays. In both Fig 1A and 1B, the first decaying phase consists of exponential decay but, after some time (instant t = τ o ), the decline is much slower. However, even if the exponential decay is prolonged, it crosses the powerlaw decay curve sooner or later (Fig 1B). A power-law decline may allow the persistent population's survival for much longer. The long tail of power-law distributions is relevant for indirect resistance because resistant (detoxifying) cells may take a long time to detoxify the medium. This paper has two main objectives. First, we aimed to understand whether persistence is necessary for the survival of sensitive cells in the context of indirect resistance (Fig 2). Second, we aimed at understanding whether the decay of the persistent population is better explained by a negative exponential or by power-law with a negative exponent. To achieve both objectives, we took advantage of previous experiments performed in our laboratory, where we measured the degree of protection of susceptible cells when co-cultured together with β-lactamase-producing cells and in the presence of the β-lactam antibiotic ampicillin [14]. In the present paper, we performed computer simulations to understand how many non-persister and persister cells contributed to the survival of susceptible cells, identified which parameters explain such The horizontal axis, representing time t (minutes), is linear, while the vertical axis is on a logarithmic scale and represents the proportion of the population that is still alive. Descending full lines represent the exponential decay of the nonpersistent population, according to exp(-0.04 t). When t = τ 0 , only persisters are alive. The broken lines represent decay according to the power-law 1/t 2 . The dotted lines represent exponential decay, i.e., according to exp(-k t). A: Persisters decay according to constant k = 0.005 min -1 . B: Persisters decay according to constant k = 0.0025 min -1 .
https://doi.org/10.1371/journal.pone.0246500.g001 contributions, and directly compared results with those obtained experimentally by Domingues et al (2017) [14]. In the end, we hoped to understand better indirect resistance and the role of persistence, which is critical for developing medical strategies against pathogenic bacteria.

Methods
In the context of indirect resistance, we intend to understand how susceptible bacteria survive while the medium is still toxic. For that, we compared experimental results from our previous experiments [14] with simulations performed in the present work. In this work, we use the mean value of the triplicated experiments performed in ref. [14] with six different initial conditions and six final values, a total of twelve time-points. In the simulations, sensitive cells decay exponentially with a decay rate of k 1 (that is, according to e À k 1 :t ) until time t = τ 0 , after which decays slower, either according to another exponential with decay rate k 2 (that is, e À k 2 :t ) or according to a power-law with the exponent β (that is, according to t −β ). We ran several simulations by varying these parameters to find the set of parameters that better explain the twelve experimental data from ref. [14].
The method is detailed in the following sections.

Previous experimental data used in this study
In the experiments, we used bacterial cells of Escherichia coli K12 MG1655 to measure the degree of protection of susceptible cells in the presence of ampicillin when co-cultured with cells encoding a β-lactamase (resistant cells). The basic experimental setup was to initiate the co-culture with a specific total initial density and frequency of resistant:sensitive cells in plates with rich medium supplemented with ampicillin at the concentration of 100 μg/mL (about ten-fold its minimal inhibitory concentration for this strain). After incubating for 24h, we quantified the density of both susceptible and resistant cells. These cells were resistant to ampicillin because they harbored the natural isolated R1 plasmid, which encodes a β-lactamase that detoxifies the medium by breaking the β-lactam ring through hydrolyzation. This plasmid is conjugative, so we also quantified transconjugants (here defined as cells that received the plasmid plus their descendants). However, the frequency of transconjugants remained very low [14], which is a consequence of the low conjugation rate of the R1 plasmid in the E. coli strain used in the experiments [23][24][25][26].
To develop our study, we used the experimental data for (i) two initial total cell densitiesapproximately 10 6 to 10 7 cfu/mL and 10 4 to 10 5 cfu/mL, henceforth denominated as high and low density respectively; (ii) three proportions between resistant (R) and susceptible (S) cells-1R:99S, 50R:50S and 99R:1S (where, e.g., 1R:99S means a frequency of a resistant cell for 99 susceptible cells); and (iii) the six final number of susceptible cells in each conditions, a total of twelve data points. The relevant information about the initial experimental conditions and final results from the Domingues et al. study (ref. [14]) are in Table 1, where we can see the average of three replicates.
As a control, Domingues et al. (ref. [14]) noted that, if alone, no sensitive cells survived after 24 h in the presence of ampicillin, which agrees with a previous study where Wiuff et al. observed between 0 and 50 cells after 24 h with ampicillin with an E. coli strain [27].

Computational model-flow of the simulation
Here we describe the algorithm of the simulation process. Table 2 and Fig 3 show the respective pseudocode and flowchart. All code is available on GitHub (https://github.com/jrebelo27/ Simulation-code-of-persistence).
We simulated the spread of resistant and sensitive cells with a given total density and at specific proportions (1R:99S, 50R:50S, or 99R:1S) in a medium plate, and we saved in a file the distances between each susceptible cell and the nearest resistant cell.
As in the experiments by Domingues et al. (ref. [14]), we assumed that the plate contains nutrients and the antibiotic ampicillin. Susceptible cells are, by definition, sensitive to this antibiotic, and resistant cells can detoxify their surroundings, clearing up the cytotoxic antibiotic. A decreasing ampicillin concentration gradient is generated from the inside out by the diffusion of the β-lactamase enzyme that degraded the antibiotic.
The simulation is composed of several cycles, as much as the number of generations completed by the resistant cells in the experiments by Domingues et al. (ref. [14]). The following i. Non-persister cells: the computer program randomly takes a percentage of susceptible bacteria, defined according to the exponential decrease (Fig 1), and tests whether its site is already nontoxic (that is, if the distance to the nearest resistant cell is less than the total detoxified radius around the resistant cell). If yes, the susceptible bacterium survives. If not, that susceptible cell dies; in practice, the program removes that cell from the simulation's next steps. Decay of susceptible cells in the cycles before t = τ 0 For each cycle between t = 0 and t = τ 0 , do (obtain the proportion of bacteria that start dividing (function A 1 . exp(-k 1 .t)); retrieve as many distances as the number of bacteria that resume growth from the list of distances between all susceptible to the closest producer; obtain the number of susceptible cells that are at a shorter distance than the radius of the detoxified area; save this number as surviving non-persistent bacteria).
Decay of susceptible cells in the cycle containing t = τ 0 For the cycle that includes τ 0 , do (obtain the proportion of bacteria that start dividing (function A 1 .exp(-k 1 .t)+A 2 . t β or A 1 .exp(-k 1 .t)+ A 2 .exp(-k 2 .t)); retrieve as many distances as the number of bacteria that resume growth from the list of distances between all susceptible to the closest producer; obtain the number of susceptible cells that are at a shorter distance than the radius of the detoxified area; calculate the proportion of nonpersistent and persistent bacteria (function A 1 .exp(-k 1 .t) +A 2 .t β or A 1 .exp(-k 1 .t)+A 2 .exp(-k 2 .t)); save the number of surviving non-persistent bacteria; save the number of surviving persistent bacteria).
Decay of susceptible cells in the cycles after t = τ 0, when bacteria leave the dormant state when the medium is detoxified For each cycle after t = τ 0 , do (obtain the proportion of bacteria that start dividing (function A 2 .t β or A 2 .exp(-k 2 . t))); retrieve as many distances as the number of bacteria that resume growth from the list of distances between all susceptible to the closest producer; obtain the number of susceptible cells that are at a shorter distance than the radius of the detoxified area; all persistent cells in a detoxified area resume growth; save the sum of these numbers as surviving persistent bacteria).
Decay of susceptible cells in the cycles after t = τ 0 when bacteria do not leave the dormant state when the medium is detoxified For each cycle between after t = τ 0 , do (obtain the proportion of bacteria that start dividing ((function A 2 .t β or A 2 .exp(-k 2 .t))); retrieve as many distances as the number of bacteria that resume growth from the list of distances between all susceptible to the closest producer; obtain the number of susceptible cells that are at a shorter distance than the radius of the detoxified area; save this number as survival persistent bacteria).
� The program code was implemented in R programming language. https://doi.org/10.1371/journal.pone.0246500.t002 ii. For persister cells, the computer program may follow two different approaches, depending on the biological assumptions. Either (a) persister cells leave the dormant state stochastically only according to the function considered (exponential-law or power-law), or (b) persister cells resume growth (leaving the dormant state) stochastically according to an exponential distribution or whenever their site becomes detoxified. If the biological assumption is that persister cells leave the dormant state only stochastically (independently of the antibiotic's presence/absence), the computer program follows the approach (a). In this case, a percentage of susceptible bacteria is randomly taken (according to the powerlaw or the exponential-law) (Fig 1). The program tests whether each susceptible cell's site is already nontoxic (that is, if the distance of each susceptible (persister) cell to the nearest resistant cell is less than the total detoxified radius around the resistant cell). If yes, the susceptible cell survives. If not, that susceptible cell dies, which means that the program removes this cell in the simulation's next steps. If the biological assumption is that persister cells leave the dormant state as soon as their site becomes nontoxic, the computer program follows the approach (b). First, we calculate the difference between the number of susceptible bacteria in the dormant state in the simulations and the one predicted by the exponential function (in the next paragraph, we explain why, in this case, we do not consider the case of a power-law function). If that difference is positive, we randomly chose that number of susceptible cells to leave the dormant state. The program tests if the nearest resistant cell's distance is less than the total detoxified radius around it. If yes, that cell survives. Otherwise, the program removes this cell from the simulation (the cell dies because the antibiotic is still present). Then we look for all persister cells present in the detoxified area; these cells leave the dormant state, resuming growth. However, if the difference is negative, i.e., if the number of susceptible bacteria in the dormant state in the simulations is lower than the one predicted by the power-law function or the exponential function, we look for all persister cells present in the detoxified area. These cells leave the dormant state, resuming growth.
Note that, returning to growth immediately after detoxification implies a mechanism, suggesting that persistence is an evolved genetic system, not the result of inadvertent metabolic and cell replication problems. In that way, it would be contradictory to assume simultaneously that persistent cells leave the dormant state as soon as the medium is nontoxic and that the persistent population decays according to power-law.
The population of genetically susceptible cells that resume growth while t < = τ 0 are, by definition, in the non-persister state, while those susceptible cells that resume growth when t > τ 0 are persister cells. Each cycle represents a generation time. Here we assume that one generation time is 30 minutes. Time does not flow continuously in the simulations, but rather in intervals of 30 to 30 minutes. Given the division of time into these intervals of 30 minutes, the interval containing τ 0 has both persister and non-persister bacteria. Using the decay curve of the population of susceptible bacteria, we calculate the percentage of persister and non-persister bacteria in this period. For example, if τ 0 = 70 mins, the simulation performs 60 mins (two cycles, each representing 30 mins) plus 10 mins decaying as non-persisters, and the remaining 20 mins decaying as persisters. Therefore, 1/3 (= 10/30) of the remaining genetically susceptible cells resume growth as non-persisters (i.e., according to the exponential) and 2/3 (= 20/30) resume growth as persisters.
At the end of a simulation, we have gathered information on how many persister and nonpersister cells generated the population of susceptible cells observed after 24h. Moreover, we can also know how many persister and non-persister cells have survived in each generation.
By performing simulations with several combinations of parameters, we can find those that better explain experimental results.

Details of the computational model
Simulating the medium plate and bacterial cells in the plate. The main procedure was to simulate the experiments performed in Domingues et al. (ref. [14]), where susceptible and resistant bacteria were mixed and cultured in agar plates. Escherichia coli cells are rod-shaped cells about 2 μm long and 0.5 μm diameter, hence occupying a 2-dimensional area of about 1 μm 2 and a plate has a diameter of 9 cm = 90000 μm. Therefore, we considered that each point in the agar plate, computationally defined by two integers (coordinates x and y), is the center of a square with an area of 1 μm 2 . The computer program's first step was to simulate the random distribution of cells in the plate, assigning random coordinates to all cells. Then, we calculated the distances between each susceptible cell and the nearest resistant cell, saving the values in a file.
Calculating the number of generations. In each experimental setup, it is possible to calculate how many generations were completed by the resistant population (resistant cells are not affected by the antibiotic). The appropriate mathematical expression is: Number of generations = Log 2 [Final number of resistant cells/Initial number of resistant cells].
In the simulations, both susceptible and resistant cells are otherwise isogenic and consume the same nutrients. We further assumed that there was no resistance cost, i.e., in the absence of antibiotics, resistant and susceptible cells replicate at the same speed. Therefore, if there were no antibiotics and given that both strains are of the same species, they would complete the same number of generations.
Simulating the radial spread of β-lactamase around resistant cells. We simulated the spread of β-lactamase as an expanding circle centered in each resistant bacterium. At any time, these circles represent an antibiotic-free area. According to the Einstein equation for the Brownian motion, the mean displacement of a small particle diffusing in a medium is proportional to the root square of the time elapsed. Therefore, the circle radius grows proportionally to the square root of time, Sqrt(time). Counting the time in bacterial generations, we may express this as R = C.Sqrt(number of generations), where R is the circle's radius, and C is a constant that depends on the diffusion constant, which may depend on the medium conditions (e.g., the agar concentration). Henceforth, we name this constant C as the "spreading parameter". Note that in the initial moment (generation 0), the value of R is 0. Therefore, all susceptible bacteria that start dividing at that moment dies.
Comparison of results between experiments and simulations. The parameters to adjust were k 1 , τ 0 , β or k 2 , and C. As explained above, the parameters A 1 and A 2 depend on k 1 , τ 0 , and β or k 2 . The program ran as many generations as those completed by resistant cells in the experiments performed in Domingues et al. [14]. Therefore, the final number of resistant cells should be the same in both experiments and simulations.
We ran several simulations by varying the parameters k 1 , τ 0 , β or k 2 , and C, to find the set of parameters that better explain the experimental results found in ref. [14] ( Table 1). In these comparisons between experiments and computer simulations, we considered that experiments had an associated experimental error. For instance, agar thickness and other physical conditions of the agar plates that may influence the spreading parameter may constitute a variance source. Furthermore, experiences are also subject to unknown errors. For these reasons, we accept our results to deviate from experimental results. We calculated the lower and upper limits of the intervals according to the following: Lower limit = Final number of susceptible bacteria obtained experimentally / Margin of error Upper limit = Final number of susceptible bacteria obtained experimentally � Margin of error We tested different margins of error and we only obtained results considering margins of error equal to or greater than 4 (a margin of error of 3.9 was not enough).
As explained above, we studied two initial cell densities and three initial frequencies of susceptible to resistant cells. In the simulations, we combined all experimental cases with our parameters. For each combination, we performed three repetitions. If one repetition result is contained in an interval, we consider that the simulated parameters explain the set experimental results for that margin of error.

Results
In this work, we took advantage of the experimental results previously obtained by our research group [14]. The authors spread resistant cells (producers of the detoxifying enzyme βlactamase) and susceptible cells in a nutrient-rich medium plate with ampicillin (a β-lactam antibiotic), followed by the quantification of susceptible (and resistant) cells after one day. This was done in one of the three frequencies (for a specific initial total density), namely, 99% of susceptible cells and 1% of resistance cells (denominated as 1R:99S), the reverse (99R:1S), and also 50% of each (50R:50S). Resistant cells can produce β-lactamase because they harbor the R1 plasmid encoding the enzyme. This naturally isolated plasmid is transferable by conjugation, so later, we check the impact of conjugation on the survival of susceptible cells.

The encounter probability of resistant and susceptible cells does not explain the survival of susceptible cells
We started by addressing the hypothesis that surviving susceptible cells are those that were very close to β-lactamase-producing cells. According to this hypothesis, the probability of encounter between resistant cells and susceptible cells would be the main factor for the survival of susceptible cells. This hypothesis is unlikely valid because the probability that susceptible cells find themselves in detoxified microenvironments is higher for higher numbers of resistant cells, that is, higher cell density or when the proportion of detoxifying cells is high (99R:1S frequency). Nevertheless, it is instructive to analyze the importance of the encounter probability between resistant and sensitive cells when spread in the agar plate and contrast with experimental observations. If the hypothesis was valid, the number of surviving susceptible cells (and their descendants after 24h) should be the same for the 99R:1S and 1R:99S frequencies. The encounter probability of a resistant and a susceptible cell is proportional to 99/100x1/100 (for the case 99R:1S) = 1/100x99/100 (for the case 1R:99S). If the encounter probabilities were the same, the number of surviving susceptible cells would be similar. However, they differed considerably (Table 1). In the high-density case, the final number of susceptible cells for the frequency 1R:99S was 3.27x10 4 , whereas for 99R:1S was 1.23x10 7 , hence differing by more than three-hundred-fold ( Table 1).
The encounter probability for the 50R:50S frequency is proportional to 50/100x50/100. This probability is approximately 25-fold higher than the encounter probability for the 99R:1S and 1R:99S frequencies. Therefore, the above hypothesis predicts that the number of surviving cells in the 50R:50S frequency should be 25 fold higher than in the 99R:1S and 1R:99S frequencies. This prediction is also far from experimental observations ( Table 1). For example, for the high-density case, the final number of susceptible cells for the 50R:50S frequency was 1.08x10 6 , which is about 10-fold less, not 25-fold higher than the 1.23x10 7 cells observed for the 99R:1S frequency ( Table 1).
As expected, these results suggest that the encounter probability is not an essential factor for the survival of susceptible cells in the indirect resistance phenomenon.

Persistence is required for susceptible cells survival
After the inoculation of susceptible and resistant cells, the latter replicate for several generations until resources present in the plate are over. The number of generations completed by the resistant cells can be calculated (see the Methods section). Assuming that the resistance cost is negligible and that all susceptible cells start replicating at the same time as resistant cells, we can also estimate how many susceptible cells should have survived when inoculated to explain their final number. Table 3 shows these estimations for the six conditions.
In two cases shown in Table 3 (high density, frequencies 50R:50S and 99R:1S), the estimated number of surviving susceptible cells is higher than one cell, but it was lower than one cell in the other four cases (high density, frequency 1R:99S, and the three frequencies when density was low). These four cases of less than one cell are important in this work because, at least in these cases, one or more bacteria have entered the persistence state and started duplicating only some hours after inoculation. In this state, susceptible bacteria can survive in the presence of ampicillin because they are not replicating, and resistant bacteria continue to produce and release β-lactamase into the culturing medium.

Persistence in indirect resistance
The time that each bacterium remains in the persistence state varies from one bacterium to another. When should persistent cells leave the dormant state? We have analyzed three possibilities. The subpopulation of persister cells resumes growth, either according to a power-law or to an exponential-law distribution. For the latter case, dormant cells may or may not resume growth as soon as the medium is nontoxic.
It is impossible to determine the number of persister cells needed to give rise to the final number of susceptible cells observed experimentally. Suppose we observed exactly four surviving susceptible cells at the end of an experiment. We wouldn't know whether: (i) the four bacteria were in a dormant state all the time; (ii) two cells were in the persistence state most of the time but resumed growth (replicating once) about 30 minutes before the end of the experiment; or (iii) one cell was dormant most of the time but resumed growth about 60 minutes before the end of the experiment. Any of these three scenarios would explain four susceptible cells at the end of the experiment. Increasing the number of final susceptible cells would sharply increase the number of possible scenarios. Therefore, we performed simulations, varying several parameters (more details in the next section), to estimate the number of persister and non-persister cells necessary to explain the experimental number of surviving susceptible cells.

Simulations to estimate the growth of susceptible cells
We had to consider the spreading of β-lactamase produced by the resistant bacteria and the decline in the susceptible population while exposed to the β-lactam antibiotic. We have seen that this decaying period has two main phases (Fig 1). The non-persistent population decays exponentially until t = τ 0 . At t = τ 0, only persistent cells survived. They resume growth and die if the antibiotic is still present. In that case, we tested two alternative possibilities for the decay of the persistent population: according to a power-law distribution or according to another exponential distribution.
We used different parameters to describe the population decay: (i) k 1 , the rate constant in the first exponential decay, is the decay rate of the non-persistent population; (ii) τ 0 is the time from which only persister cells are alive, which is when the probability distribution changes from the exponential decay to the power-law or the second exponential decay; (iii) β, the power-law exponent or k 2 , the rate constant in the second exponential decay. Therefore, in our simulations, we considered these three parameters together (k 1 , τ 0 , and β or k 1 , τ 0 , and k 2 ). We used a fourth parameter (spreading constant) representing the rate increase of the detoxified area-this area increases around each resistant bacterium due to the detoxifying enzyme's diffusion.
To find the parameters that best fit the experimental results [14], we combined the following parameters:  3 = 4913 combinations (assuming that the persistent population decays according to power-law) or 2635 combinations (assuming exponential decay), since k 2 has to be lower than k 1 . For each combination of these parameters, we tested spreading rates C 2 {0.2, 0.4, 0.6, 0.8, 1, 1.2, 1.4, 1.6, 1.8, 2, 4, 6, 8, 10, 12, 14, 16, 18, 20}. Then, by comparing the computational results with the experimental ones, we obtained a set of parameters that explained the experimental results of the three frequencies (1R:99S, 50R:50S, 99R:1S) in low and high density. We included the possibility of experimental error; that is, we allowed our results to differ from experimental results according to a certain error-margin (see the Methods section).
Assuming an error-margin between 1 (no error) and 3.9, we could not find any combination of parameters (k 1 , τ 0 , β, nor k 1 , τ 0, and k 2 ) explaining all the experimental results. With an error margin of four, we found three solutions sets of parameters explaining the experimental data from ref. [14]. (Table 4). The detoxified area's spreading parameter varied considerably in these combinations, probably due to different experimental conditions (see Discussion).
In the simulations presented until now, persister cells were resuming growth as soon as the medium became detoxified. These are the "yes" cases shown in Table 4. As explained before, one should expect a power-law decay of the persister population if there is no genetic system for persistence. Therefore, as explained in the methods section, we should not assume that persister cells decay like a power-law resume growth as soon as the medium became detoxified because that would imply a genetic system to leave the dormant state.
We performed further simulations, this time assuming that persister cells stay in the dormant state even when the medium becomes nontoxic. Therefore, in this case, persister cells resume growth stochastically, independently of the antibiotic's presence in the medium, and, as in the previous case, they do that according to the power-law distribution or the second exponential distribution. We found one combination with an error margin of 4 if the persistent population decays according to a power-law distribution (Table 4) and five combinations if the persistent population decays exponentially (Table 4).
During simulations, we also quantified the number of persistent and non-persistent bacteria throughout generations. Therefore, we can analyze how many susceptible cells in the final population originated from non-persister and how many from persister cells. We have done this analysis for each combination of parameters presented in Table 4. The S1 to S9 Tables show the results of the analyses. S1 to S9 Figs show the corresponding number of resistant, non-persister cells, persister cells, and their respective descendants through time using the parameter found in Table 4.
All susceptible bacteria observed at the end of the experiments in the low-density cases descend from persistent bacteria (S1 to S9 Tables). When the density was high, some non-persistent bacteria also survived in the early generations. As non-persister bacteria duplicate in each generation (contrary to dormant persister cells), they may become strongly represented at the end of the experiment (this is observed for the cases of high density, frequencies 1R:99S and 50R:50S), even if they were a minority of the surviving cells (S1 to S9 Tables).
When the initial number of β-lactamase-producing cells was too low to detoxify the agarplate fully, persister cells maintained their state until the end of the experiment (i.e., until 24h later when cells were finally plated in a medium without antibiotic for quantification). Such permanence of bacteria in the persister state occurred when the initial number of β-lactamaseproducer cells was low, when the initial total cell density was low (for the three frequencies), or when the initial cell density was high but the initial frequency of β-lactamase-producing cells was low (1R:99S) or intermediate (50R:50S) (S1 to S9 Tables). In these cases, resistant cells spent all resources before the detoxification of the agar-plate.

The impact of plasmid transfer in susceptible cells survival is negligible
As explained before, we used experimental data obtained with resistant cells that were encoding the detoxifying enzyme in a transferable (conjugative) plasmid, the R1 plasmid. Therefore, plasmids may move (by replication) into susceptible cells and form transconjugants, here broadly defined as cells that received the R1 plasmid and their descendants. Transconjugants become producers of β-lactamase, hence able to detoxify the environment.
Transconjugants represent a small percentage (between 0% and 1%) of the susceptible cells in the experiments performed with the R1 plasmid by Domingues et al. (ref. [14]) ( Table 1).
Because we estimated the number of generations completed in the agar-plate, it is possible to assess transconjugants' impact on indirect resistance. If there are T transconjugants at the end of the experiment, and assuming that the contribution to the detoxification is the highest if all transconjugants formed at the end of the first generation, there would be T/2ng transfer events, where ng is the total number of generations (we are assuming that, at the end of the first generation, resistant cells have already replicated once). In interestack Table 5, we show the number of generations, the final number of transconjugants observed in the agar-plate, the estimated number of plasmid transfers, and the proportion of resistant cells that are transconjugants. In all cases, the proportion of transconjugants among all cells capable of detoxifying the medium was extremely low (0.0011% or less). Therefore, the impact of transconjugants on the detoxification of the medium must have been shallow.

Mathematical description of the persistent sub-population and biological implications
Although close to -2, the exponent found by Simsek and Kim (ref. [22]) was -2.1, and the one found with our simulations was 2.2. In both cases, the power-law decays slightly faster than according to 1/t 2 . Therefore, it is relevant to understand how heterogeneous populations should decay. We argue in this section that, if the persistent population is heterogeneous and persistence is not genetically controlled, it should decay according to a distribution close to a power-law but not precisely according to this distribution.
Following the argument by Simsek and Kim (2019), consider a homogeneous population of antibiotic-susceptible cells in the presence of a bactericidal antibiotic. If a non-growing cell rejuvenates (here defined as resuming growth, see below), it dies due to the antibiotic. Therefore, the number of cells still alive at a given time t decreases according to dnðtÞ dt ¼ À k 1 :n t ð Þ, where k 1 is the rejuvenation probability constant. The solution of this differential equation is n

(t) = n(0).exp(−k 1 .t). The rejuvenation probability refers to the number of cells resuming growth in the time interval, which is proportional to the number of cells still alive:
NðtÞ:Dt / À DnðtÞ ¼ À k 1 :expðÀ k 1 :tÞ:Dt: On the other hand, a subpopulation of the cells with various problems in the metabolism, in the cell replication cycle, or even the cell's response to these problems, stop dividing for some time [21,22]. Each bacterium may have a different issue from a big group of possible problems. Therefore, these bacteria should present a wide range of rejuvenation constants [22]. This bacterial population is heterogeneous, with many different k constants. The number of cells resuming growth at a particular time t, N(t), is proportional to: If the population decays between time τ 0 and t Max , then the integral's limits are a > 1/τ 0 and b < 1/t Max . This lower limit b is close to zero because t Max is high-the persistent population endures a long time [22].
Note that, until now, we only know that the upper limit of the integral, a, has to be higher than 1/τ 0 . We now argue that this upper limit has to be lower than k 1 . This limit arises from the fact that cells in this heterogeneous sub-population rejuvenate later than the non-persister cells-their rejuvenation constant should be lower than that of the non-persister cells. Therefore, the integral becomes: where k 1 is the rejuvenation constant of the non-persister population. In general, 1À expðÀ k 1 :tÞð1þk 1 :tÞ t 2 ≲t À 2 . Experimental results from ref. [22] have shown that the persistent population starts decaying after about 93 min and k 1 is about 0.063 min -1 . Therefore, k 1 .t ⋍ 5.859 or higher and increases in time, so the numerator in Eq 1 is 0.98 (that is 1−exp(−k 1 .t)(1+k 1 .t) = 0.98). Therefore, in general 1À expðÀ k 1 :tÞð1þk 1 :tÞ t 2 ≲t À 2 . This result may explain why Simsek and Kim (ref. [22]) derived an exponent from their experiments of -2.1, which is slightly lower than their theoretical prediction of -2. However, when t increases, the numerator of Eq 1 converges to 1, which means that the power-law t β should converge to t -2 when t increases.
Likewise, our results for the power-law decay (Table 4) suggest that the persistent population starts decaying after τ o = 60 min, and k 1 is about 0.07 min -1 . Therefore, when t = τ o , k 1 .t ⋍ 4.2 and increases in time (because t increases), so the numerator in Eq 1 is close to 0.92. Again, our results suggest that 1À expðÀ k 1 :tÞð1þk 1 :tÞ t 2 ≲t À 2 and also explains why we obtained exponents slightly lower than -2 (if the persister population indeed decays according to a power-law).

Fitting a power-law in the decay of persister populations from previous works
Given the intriguing similarity of the power law exponent found in this paper and the one predicted and experimentally obtained by Simsek and Kim [22], we looked at a few previous studies of bacterial persistence and estimated a putative power-law exponent. We note, however, that the following is neither a systematic review nor a meta-analysis. In a log-log plot, a powerlaw α.t β curve (where time t is the independent variable and α and β are two constants) follow a straight line with a slope of β. Following this procedure, experimental results involving persisters' decay of E. coli under different concentrations of ampicillin (all above MIC) from refs. [17,[27][28][29][30] present slopes with a median of -2.25 and IQR = -3.11 to -1.85 (S10 Table).

Discussion
This work aimed at addressing if persistence explains sensitive cells' survival in the context of indirect resistance with different initial cell densities and resistant/susceptible cell ratios. Our results strongly suggest that persistence supports the growth of sensitive cells in the presence of antibiotics, mainly when the cell density is low. These results help to understand our prevous work, where we observed that susceptible cells grew even when resistant cells weren't able to grow and detoxify the medium [14]. Interestingly, some of those persister cells were in a dormant stage just for a few hours. However, we could not achieve our second objective, namely deciding whether persister populations decay exponentially or like a power-law.
To understand persisters' behavior, we started by asking whether they were responsible for the survival of susceptible cells in the context of indirect resistance. For that, we carried out simulations to mimic the experiments that we have performed in a previous work where we spread a mixture of susceptible and β-lactamase-producing cells in agar-plates supplemented with a β-lactam antibiotic [14]. We simulated the behavior of persister cells in three different ways: (i) in the presence of a bactericidal antibiotic, the persistent population decays according to an exponential-law versus according to a power-law; (ii) persister cells leave the dormant state as soon as the medium becomes detoxified (just for the exponential decay) versus independently of the medium detoxification, hence merely according to the probability mentioned above. Our simulations suggest that persister cells and their descendants were a part, or even all, of the surviving susceptible population, irrespectively of the three alternative behavior models of the persister cells implemented in the simulations. Persisters were the only survivors in the indirect resistance phenomenon when the initial cell density was low.
Given persistent cells' involvement, we used the results to go more in-depth and understand their nature. Based on the several molecular mechanisms underlying prokaryotic or eukaryotic persistence states, the prevalent view is that persistence is an evolved characteristic [2,5,6]. If genetically encoded, the expectation would be that the persistent population is homogeneous and decays exponentially [17]. Instead, a few recent works have proposed that persistence may also be an accidental consequence of inadvertent cell problems and errors [20,21]. In this case, the persistent populations should be heterogeneous because cells would have different reasons for showing low metabolism, and the consequent theoretical prediction is that the persistent population should decay, not exponentially, but according to a power-law with the exponent of -2 [22] or slightly lower than -2 (this paper).
The exponential decays are direct consequences of first-order kinetics. The exponential declines occur in various situations, from radioactive decay to the drop of atmospheric pressure with increasing height above sea level. And, of course, the non-persister bacterial population also decays exponentially in time because the bacterial population is large, homogenous, and the law of large numbers holds. To our knowledge, there is no theoretical prediction for the decay rate if the persistent population declines exponentially. Our simulations show that an exponential decline of persisters is possible only for shallow values of the decay constantthis allows the survival of persister cells for several hours in the experiments. However, Simsek and Kim [22] were able to mathematically predict the exponent in the power-law case, namely that it should be -2. Likewise, we found an exponent close to -2 in the simulations where we assume that persisters' decay follows a power-law (Table 4).
It is relevant to emphasize that, despite the similarity of the exponent values found here (based on the experiments from ref. [14]) and in the Simsek and Kim' study [22], the experimental methods of these two studies were significantly different. While Simsek and Kim [22] studied the decay of the susceptible population in a liquid and well-mixed medium without resistant cells, the experiments simulated here (based in ref. [14]) were performed in agarplates where some susceptible cells die due to the antibiotic and others survive thanks to persistence or the effect of indirect resistance. The similarity of the exponents found with two different experimental methods, with the one predicted theoretically [22], is, therefore, intriguing.
On the other hand, we also found parameter sets for exponential decay of the persister population. With the exponential decay rate found with our simulations, the persistent population decays faster than with a power-law. With the exponential decays, there should be no persister cells after 24 h, a prediction that is consistent with a control experiment performed in ref. [14] with pure cultures of susceptible cells in ampicillin where no persister cells were found after 24 h. The absence of persisters after 24 h in these conditions is also consistent with other works (e.g. ref. [27]). The absence of persister cells after 24 h may result from the fact that the culture of susceptible cells was in stationary phase for just a few hours, while Simsek and Kim work [22], cultures were in the stationary phase for three days. In this respect, our experimental data are more consistent with an exponential decay of the persistent sub-population than with a power-law decay.
We evaluated the impact of persister cells resuming growth as soon as the medium is nontoxic versus resuming growth stochastically, independently of the antibiotic's presence in the medium. We found reasonable sets of parameters using both behavioral models. Therefore, strictly speaking, we could not conclude whether persister cells leave the dormant state and resume growth when the medium is nontoxic. It is interesting to find that the error margin is the same in both models. Possibly, the reason is the following. When the initial cell density is low, detoxifying cells are insufficient to detoxify the agar plate thoroughly. When the initial density is high, detoxification is complete. However, in these cases, detoxification is so fast that the origin of most surviving susceptible cells is the non-persister subpopulation. Such fast detoxification allows the duplication of these surviving cells for several generations; such exponential growth allows their numbers to reach very high values. So, in most cases with high density, persisters leaving the dormant stage are not relevant among the final number of surviving susceptible cells. Therefore, the mechanistic way of leaving the dormant stage (simply stochastically or according to the medium detoxification) becomes irrelevant.
The similarity of the exponent values in this work and Simsek and Kim [22] work is impressive, but there should be an explanation for the discrepancy from the theoretical prediction of -2. We have shown that heterogeneous populations should decay according to 1À expðÀ k 1 :tÞð1þk 1 :tÞ t 2 , which is close to but slightly lower than 1 t 2 . Such discrepancy may explain why our simulations and Simsek and Kim's experiments point to exponents slightly lower than -2. Both theoretical predictions assumed that several sub-populations of cells constitute the persistent population. The difference between the two theoretical predictions is that our derivation assumes that no hypothetic subpopulations are decaying faster than non-persister cells. This assumption implies fewer cells alive in the persistent state than predicted before [22]. As time passes, the two mathematical predictions converge because even if we were including the subpopulations decaying faster than the non-persistent population, those cells would already be dead.
Given that we simulated bacteria in the agar-plate, we had to consider the radial spreading of β-lactamase around their producers (resistant cells) and the subsequent decrease in antibiotic concentration. The system has some complexity because, in some simulations, the initial number of resistant cells can be high and because there is undoubtedly diffusion of β-lactamase from each resistant bacterium and of the antibiotic towards each resistant bacterium. It is even possible that the detoxifying area increases as a diffusion wave. Moreover, resistant cells duplicate every half an hour, probably increasing the β-lactamase enzyme production outwards the resistant colony. Therefore, we assumed that the detoxified area increases monotonically in time. Future studies should scrutinize the relevance of this assumption. We had to consider a wide range of values for the detoxified area's speed of increase to fit the experimental results. This range may have several causes. For example, although the agar concentration was the same in all experiments, some plates could be more dried than others, eventually facilitating or hampering the detoxifying enzyme molecules and antibiotic molecules' movement. Table 5 shows that transconjugants' participation in the detoxification of the agar-plate must have been low. This result agrees with previous works showing that the transfer rate of the R1 plasmid is low [14,23,24,31,32].
Several resistance determinants, including genes and chromosomal mutations, are responsible for the burden of antibiotic resistance. This burden is tremendous. Just in the European Economic Area, antibiotic resistance is responsible for 33000 deaths/year and 874000 disability-adjusted life-years [33]. Unfortunately, to survive bactericidal antibiotics, bacteria do not even need to harbor resistance determinants. Susceptible bacteria may rely on indirect resistance and bacterial persistence, as we have seen. Therefore, this work's conclusion that persistence is involved in indirect resistance is worrying.
The power-law distribution has a long tail (that is, longer than that of an exponential), which means that, at least theoretically, some susceptible bacteria could survive for several weeks, eventually after the end of antibiotic usage by the patient. Long-lived persisters may dictate treatments' failure because some of these cells may leave the dormant state and reinitiate their pathogenic effects. This risk goes in line with the reports on persistence being a significant cause for recurrent and chronic infections, dictating the patients' disease progression and outcome [29,34]. However, in the experimental system studied here [14], simulations suggest that only short-term persister cells were involved in indirect resistance.
This study has limitations. In the computer model, we hypothesized that the decay rate of both sub-populations (non-persisters and persisters) was constant in time. Several previous works would predict that it is a reasonable assumption when sensitive cells are alone and in a liquid media, but, to our knowledge, the decay of these sub-populations in structured environments with companion resistant cells has not been studied. Likewise, we assumed that the persistent population decays either like a power-law or exponentially. If there is a genetic system coding for persistence, the resulting persister population may result from a mixture of cells, some of them expressing the genes involved in the persistence (decaying exponentially), while other cells would be in the persistence state due to errors and faults (following a power-law close to t -2 ). In structured habitat like agar-plates used in ref. [14], biochemical signals may vary locally, resulting in different responses to the antibiotic, conceivably leading to a mixture of both decays. Another limitation concerns the spread of detoxified areas around each resistant colony. We assumed that those circular areas increase linearly in time, but the system is more complex because, on the one hand, the number of detoxifying cells in the colony increases in time; on the other hand, one should expect some inward diffusion of ampicillin. A third limitation is our assumption that all cells replicate at the same speed. As argued above, that may be considered a strong assumption in structured habitats given the expected variance across each agar plate. These assumptions of the model and the exponential growth of bacteria may lead to large margins of error. Such margins implied our inability to determine if the persister subpopulation decays exponentially or like a power-law.
In conclusion, this work shows the involvement of persister cells in indirect resistance and that the way the persister population decays is consistent with both hypotheses studied: the persistent population decays according to a power-law with an exponent slightly lower than -2, or that the persistent population decays exponentially. As Simsek and Kim (ref. [22]) argued, a power-law decay with an exponent close to -2 could mean that persistence in the E. coli strain in those specific experimental conditions is the consequence of accidental problems involving replication and metabolism, instead of being an evolved character (see also [20,21]). If confirmed, the implication is that, in the case of the experimental conditions in ref. [14], persistence is maladaptive, despite its frequent dramatic medical consequences. A strategy to find anti-persistent drugs should perhaps be different if persisters are moribund cells versus the result of an evolved genetic program.  Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 30, k 1 = 0.1, k 2 = 0.01. Results of simulations when we assumed that the persister population decays according to an exponential and that persister cells leave the dormant state as soon as the medium becomes detoxified. (DOCX) S2 Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 50, k 1 = 0.055, k 2 = 0.01. Results of simulations when we assumed that the persister population decays according to an exponential and that persister cells leave the dormant state as soon as the medium becomes detoxified. (DOCX) S3 Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 60, k 1 = 0.045, k 2 = 0.01. Results of simulations when we assumed that the persister population decays according to an exponential and that persister cells leave the dormant state as soon as the medium becomes detoxified. (DOCX) S4 Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 60, k 1 = 0.07, β = -2.2. Results of simulations when we assumed that the persister population decays according to a power-law and that persister cells do not leave the dormant state as soon as the medium becomes detoxified. (DOCX) S5 Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 30, k 1 = 0.09, k 2 = 0.01. Results of simulations when we assumed that the persister population decays according to an exponential and that persister cells do not leave the dormant state as soon as the medium becomes detoxified. (DOCX) S6 Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 30, k 1 = 0.095, k 2 = 0.01. Results of simulations when we assumed that the persister population decays according to an exponential and that persister cells do not leave the dormant state as soon as the medium becomes detoxified. (DOCX) S7 Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 30, k 1 = 0.1, k 2 = 0.01. Results of simulations when we assumed that the persister population decays according to an exponential and that persister cells do not leave the dormant state as soon as the medium becomes detoxified. (DOCX) S8 Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 50, k 1 = 0.055, k 2 = 0.01. Results of simulations when we assumed that the persister population decays according to an exponential and that persister cells do not leave the dormant state as soon as the medium becomes detoxified. (DOCX) S9 Table. Persister and non-persister cells that originated the final susceptible population considering τ 0 = 60, k 1 = 0.045, k 2 = 0.01. Results of simulations when we assumed that the persister population decays according to an exponential and that persister cells do not leave the dormant state as soon as the medium becomes detoxified. (DOCX) S10 Table. Slope of the decay of the persister population assuming a power-law. The exponent of a putative power-law describing the decay of the persister population. (DOCX)