Gene expression noise can promote the fixation of beneficial mutations in fluctuating environments

Nongenetic phenotypic variation can either speed up or slow down adaptive evolution. We show that it can speed up evolution in environments where available carbon and energy sources change over time. To this end, we use an experimentally validated model of Escherichia coli growth on two alternative carbon sources, glucose and acetate. On the superior carbon source (glucose), all cells achieve high growth rates, while on the inferior carbon source (acetate) only a small fraction of the population manages to initiate growth. Consequently, populations experience a bottleneck when the environment changes from the superior to the inferior carbon source. Growth on the inferior carbon source depends on a circuit under the control of a transcription factor that is repressed in the presence of the superior carbon source. We show that noise in the expression of this transcription factor can increase the probability that cells start growing on the inferior carbon source. In doing so, it can decrease the severity of the bottleneck and increase mean population fitness whenever this fitness is low. A modest amount of noise can also enhance the fitness effects of a beneficial allele that increases the fraction of a population initiating growth on acetate. Additionally, noise can protect this allele from extinction, accelerate its spread, and increase its likelihood of going to fixation. Central to the adaptation-enhancing principle we identify is the ability of noise to mitigate population bottlenecks, particularly in environments that fluctuate periodically. Because such bottlenecks are frequent in fluctuating environments, and because periodically fluctuating environments themselves are common, this principle may apply to a broad range of environments and organisms.


Code implementation
We performed all modelling in python 2.7. During simulations, we observed that excessively small concentrations could lead to negative concentration values of some metabolites. To avoid negative concentrations, we set any metabolite concentration to zero if there is less than 1 molecule either in the growth environment or within a given cell. To calculate the number of molecules for any one molecular species in a cell, the concentration of these molecules is multiplied with the cell mass for intracellular metabolites or with the volume of the environment for extracellular metabolites, and divided by Avogadro's number. We performed all statistical analyses in R 3.2.3 [1]. All figures are plotted in matplotlib 2.1.0 [2].

Additional model details
In this part we describe details of the system of ordinary differential equations that model the change in both extracellular and intracellular metabolite concentrations, given the rates of the metabolic reactions occurring in every cell. We also discuss the process of cell efflux, where cells are flushed out of the growth environment at random.

Extracellular metabolite concentrations
We simulate growth in an environment with constant flow of medium entering and leaving the system at a "dilution" rate D h −1 . The concentrations of glucose and acetate in the incoming medium are denoted by glc in and ac in mM respectively. A population of cells takes up these carbon source at rates glc uptake and ac uptake mM h −1 , respectively. These considerations give rise to two ODEs for the concentrations of glucose glc and acetate ac in the environment (in mM).
Because the uptake rates of glucose and acetate, J Gi and J Ai , vary between cells, we obtain the overall rate of a population's carbon source consumption by summing over all N individuals in the population.
Because intracellular metabolite concentrations are given per unit biomass, to calculate the change in concentration in the growth environment we need to convert the the concentrations of the intracellular metabolites in mmol g −1 h −1 to those of the extracellular metabolites in mM h −1 using the biomass of each cell B i and the volume of the growth environment V . In this way, we arrive at the following expressions for both carbon source uptake rates

Reaction rates
The rate of each chemical reaction is calculated for each cell using the intra-and extracellular concentrations of the relevant metabolites, the number of enzymes catalysing the reaction in the cell, as well as the cell mass at any given point in time. The turnover number of every enzyme kcat is converted from s −1 to mmol h −1 using equation 5 [3].
kcat = kcat s · 3600 6.023 × 10 20 (5) because there are 3600 s h −1 and 6.023 × 10 20 mmol −1 . All enzymatic steps follow Michaelis-Menten kinetics except for anabolism An, which converts fbp into biomass and follows Monod-Wyman-Changeux kinetics [4]. All reaction rates are divided by the current cell mass B to account for an enzyme's dilution during cell growth. In sum, reaction fluxes are given by Because both protein and metabolite concentrations are given per unit mass (g −1 ), metabolite fluxes are given in mmol g −1 h −1 .

Cell efflux
A constant flow of nutrients through the simulated growth environment also leads to cells being flushed out at random, a process that we model stochastically as follows. Given an initial populations size of N 0 cells, the number of survivors after a given time interval ∆t follows a binomial distribution where p is the probability of survival and decays exponentially over time at a constant rate set by the

A deterministic version of the model is bistable
In this section, we use a combination of analytical and numerical approaches to confirm that our model generates two stable states on acetate that correspond to a growing and a nongrowing subpopulation.
We use the term stability here in the sense that once metabolism settles in a given growth state, it remains there unless perturbed from outside and will return to that state if perturbed. The following analysis also allows us to show precisely how the components of our model interact to form this bistable behaviour. In order to do so, we use a simpler, deterministic variant of our cell growth model.
where c is a conversion factor between the anabolism reaction rate J An and the growth rate µ (see Materials and Methods in main text). The ODE system contains only three variables, namely the intracellular concentrations pep and f bp, and the concentration of the acetate incorporation protein Ai.
As we show below, the external acetate concentration ac plays an important role in determining the equilibrium behaviour of this ODE system.
We start the analysis by finding out where the ODE system is at equilibrium, that is, for which values of the three variables (and given the extracellular concentration of acetate) there is no change in the system. To find these points, one can determine where in the variable space each differential equation is equal to zero (Fig 1). For each ODE, this is a line through variable space that is called the nullcline.
Points where all nullclines of the ODE system intersect are equilibrium points. In order to identify these nullclines, it is easiest to start with the differential equation governing the concentration of the acetate incorporation enzyme (dAi/dt), whose nullcline we will call the Ai nullcline for convenience. This differential equation can be factored into two parts; the growth rate c · J An (pep, f bp) and the difference between the current protein concentration and the protein production rate P Ai (f bp) − Ai. Consequently, the enzyme concentration Ai remains constant (dz/dt = 0) if (i) cells are not growing or (ii) the protein number per cell is equal to the production rate (Fig 1a). Because An important point is that when f bp > 0 or ac > 0, dAi/dt is zero only if Ai = P Ai (f bp), and thus the Ai nullcline is entirely determined by the fbp concentration f bp (Fig 1a). Given that any equilibrium point must lie on the Ai nullcline, this also means that the equilibrium valueÂi is determined byf bp, and thus the deterministic system can be simplified at equilibrium by eliminating Ai. In other words, we can analyse the reduced system and proceed simply by identifying the nullclines of these two equations.
The nullcline of dpep/dt (the pep nullcline) is determined by the rate of acetate incorporation J Ai (P Ai (f bp), ac), which is determined by both the fbp concentration f bp and the external acetate concentration ac. Keeping f bp constant, the nullcline with respect to ac lies on a line determined by the Michaelis-Menten kinetics of the acetate uptake rate J Ai (P Ai (f bp), ac) (Fig 1b). Conversely, when ac is constant, the pep nullcline is determined by the Hill kinetics of P Ai (f bp) (Fig 1c). As the external acetate concentration ac increases, the pep nullcline in pep-f bp space moves increasingly towards the right on the pep axis, i.e. to higher pep concentrations, until the cell reaches a maximum rate of Ai production, and the rate of acetate uptake saturates.
We next consider the fbp nullcline (df bp/dt = 0). For each value of pep there is a unique value of f bp where f bp does not change. The converse is not true, however. That is, for each given value of f bp, there may be up to three concentrations of pep for which f bp does not change (Fig 1d). This asymmetry is a consequence of the MWC kinetics of J An (pep, f bp); at low pep, the rate of fbp production from pep J Lg (pep) increases faster with with greater pep concentrations than the rate of anabolism In this parameter regime there are thus three intersections, the two 'outer' intersections moving along the fbp nullcline towards higher values of pep, and one 'inner' intersection between them that is moving back towards lower pep and higher f bp values (Fig 2c). This latter intersection eventually joins with the original intersection still travelling down the fbp nullcline. The pep nullcline leaves the fbp nullcline at this point, leaving the intersection further down the fbp nullcline to continue onwards to higher pep concentrations (Fig 2d).
We can get some intuition about the stability of these intersections (which are equilibrium points) by visually inspecting how the ODE system behaves around these points (Fig 3a-b). If the ODE system moves towards an equilibrium point, that point is stable, while if it moves away from the point, the  The two outer equilibrium points are stable, and the point in between is unstable. The vector field shows the behaviour of the two-equation ODE system close to the three equilibria. Equilibria to which all arrows in the immediate neighbourhood point are equilibria to which the system converges. These equilibria are stable, because the system will return to these equilibria after perturbation. In contrast, if some arrows point away from an equilibrium (as is the case for the equilibrium in the middle), then the equilibrium is unstable and the system will diverge from it if perturbed.
equilibrium is unstable. If we examine the ODE system at an acetate concentration where there are three equilibrium points, we observe that the two outer points are stable and the inner point is unstable.
Thus, at this acetate concentration the ODE system is bistable.
To determine the number and stability of equilibrium points for a range of acetate concentrations, we returned to the full ODE system with three differential equations. Using the python package sympy (version 1.1.1), we created expressions for the three nullclines and the Jacobian matrix of the full ODE system. Sympy was not able to find equilibrium points by solving the ODE system at a given acetate concentration. We thus used the expressions of the nullclines to find the equilibrium points numerically  (B) The system shows hysteresis in growth rates to changes in acetate concentration. For some acetate concentrations (grey rectangle), two stable equilibrium points exist. In this region the ODE system shows hysteresis. That is, when acetate concentrations are slightly changed, the new equilibrium point depends on the system's previous acetate concentration: As acetate concentrations increase from a low value, the stable growth rate will jump to a high value at a relatively high acetate concentration (6.2 mM). Conversely, when the acetate concentration is lowered from a high value, the system will jump to a lower growth rate at a lower acetate concentration of 5.3 mM.
of pep production reaches a maximum, and consequently the pep nullcline does not advance further along the pep axis. The results confirmed the qualitative analysis above. In the pep-f bp plane (Fig 4a), two regions with stable equilibria are interrupted by a region containing unstable equilibria. This latter region occurs where the fbp nullcline decreases along the f bp axis while increasing along the pep axis. One of the points in this unstable region is the unstable equilibrium point described in the qualitative analysis.
Taken together, these observations confirm our intuition that allosteric activation of the anabolic enzyme An is necessary to create bistability on acetate: If this allosteric activation of An by pep were not present, the fbp concentration would not decrease, and Cra activity would remain low, thus excluding growth on acetate.
When the equilibria are plotted in the plane described by the growth rate (c · J An (pep, f bp)) and the acetate concentration (Fig 4b) the system displays dynamics reminiscent of classic hysteresis. Nonlinear systems with hysteresis have memory effects: That is, the current system state depends on its previous state. This behaviour is apparent in the response of the growth rate to changing acetate concentrations.
Because cells will avoid unstable growth rates (red circles in Fig 4b) the acetate concentration at which a cell jumps between a high and a low growth rate depends on whether the cell is currently experiencing slow or fast growth; the jump from fast to slow occurs at a lower acetate concentration than the jump from slow to fast growth. We note that when acetate concentrations change, the concentration of the acetate incorporation enzyme Ai in the full deterministic model with three ODEs does not travel along the Ai nullcline determined by its production rate P Ai (f bp). Consequently, the trajectory of the threeequation ODE system in response to a change in extracellular acetate concentration will follow a more direct path than that given by the equilibrium points in the acetate-growth rate plane.
7 Searching parameter space: How likely is bistability, and what influences it?
In the previous section, we used a numerical analysis to confirm that a deterministic version of our model is bistable within a given range of acetate concentrations. However, on its own this observation might still mean that bistability of this kind is a rare and simply a peculiarity of our choice of parameter values.
In other words, the bistability we simulate would not be robust to change in these parameters. If this were the case, then the mechanism we assume to cause bistability in carbon metabolism might also not persist over evolutionary time, because mutational changes might alter the biochemical parameters that allow bistability. We therefore searched the parameter space for combinations of parameters that lead to bistability in our model. The fraction of parameter samples that result in bistability is an estimate of the volume of parameter space that allows bistable behaviour, and is closely associated with its robustness to perturbations [5].
We used a simple, 'brute-force' approach to search the parameter space of our model. We  ) we sampled values as 10 x , where x was drawn from a uniform distribution. We did so because preliminary results indicated that low values preferentially led to bistability. The MWC allosteric constant L An , is a ratio between the relaxed and tense forms of an enzyme [4] that determines the enzyme activity and the enzyme's response to a ligand that changes enzyme activity. Enzymes in their relaxed form are active, whereas enzymes in the tense form are inactive. Allosteric activators or inhibitors act by shifting the ratio between these two forms either towards the relaxed or the tense form respectively. The MWC allosteric constant is the ratio between these two states in the absence of the ligand, where a value of one means that both forms are equally likely [4]. We therefore sampled x from a uniform distribution in order to explore cases where either one or the other form predominates in absence of the An allosteric activator pep. For the MWC cooperativity constant n f bp we sampled integers values, because the MWC cooperativity constant relates directly to the number of binding sites [4].
To determine if a parameter set leads to bistable behaviour, we used the same numerical approach as in the previous section, and again applied it to the deterministic model of three ODEs (ODE system 12).
For each sampled parameter combination, we investigated 500 pep concentrations sampled at regular intervals on a logarithmic scale from 10 −10 mmol g −1 to 10 2 mmol g −1 . For each pep concentration, we then solved the system of ODEs for the corresponding fbp and Ai concentrations where the system lies at an equilibrium point (i.e. that all three ODEs are equal to zero). We then evaluated whether the equilibrium point thus identified was stable. Specifically, we used sympy to numerically determine whether the real parts of the eigenvalues of the Jacobian matrix at the equilibrium are negative. In this way, we obtained 500 equilibria along the pep concentration axis together with the stability of each equilibrium. If we found two regions with stable equilibria seperated by a region of unstable equilibria along the pep concentration axis (as in the previous section, Fig 4a), we concluded that the parameter set in question leads to bistable behaviour.
Overall, our sampling analysis showed that 1.5% of all parameter sets led to bistable behaviour. This percentage is probably an underestimate, as we examined only a modest range of pep concentrations.
We observed that bistability can be found in a wide range of the parameter space we examined (Fig 5).
For L An , bistability becomes more likely as L An increases because this means that allosteric activation by pep makes a greater impact on the activity of An (Fig 5g). A larger L An means that in the absence of pep, the enzyme is less active. This result is consistent with a less formal analysis performed by Kotte et al. [8]. We also note that increased cooperativity of the anabolic enzyme An increases the probability of observing bistability (Fig 5l). A greater degree of cooperativity n f bp makes the response of An to changes in the concentration of its activator pep more sigmoidal. Larger values of n f bp make the response almost binary ('on-off'). This switch-like behaviour is conductive to bistability because it reinforces the difference between the growing (high rate of fbp removal, 'on') and nongrowing (low rate of fbp removal, 'off') state. When we examined the presence or absence of skew in the probability of bistability, we observed that bistability is most likely in those regions of parameter space where both the growing and nongrowing state are stable and as distinct from each other as possible. These results are consistent with our analysis in the previous section and the intuitions laid out in the main text. To illustrate how parameter values can influence the stability of either one or both growth states, we will consider parameters affecting the activity of the three metabolic enzymes required for growth on acetate (Ai, Lg and An). For each of these enzymes, activity can be affected by enzyme concentration, by how tight the enzyme binds to its substrate, and by the binding of allosteric ligands [9].
For instance, bistability is more likely if the dissociation constant of the acetate incorporation enzyme Ai for acetate is high, i.e. if Ai-acetate binding is loose (Fig 5b). A higher dissociation constant means that, everything else being equal, a cell starts growing only at higher acetate concentrations. This insensitivity to low acetate concentrations stabilizes the nongrowing state. At the same time, bistability is more likely if the turnover number (Fig 5h) or the maximum expression level (Fig 5n) of Ai is higher.
Both of these lead to a faster rate of acetate incorporation when acetate concentrations are high, and helps to stabilize the growing state.
Within a cell at steady-state, lowering the activity of an enzyme can increase the intracellular concentration of its substrate, and vice versa [9]. For example, in the enzyme Lg (Fig 5c), which converts pep into fbp, loose Lg-pep binding facilitates the accumulation of pep and thus the activation of the anabolic enzyme An. This process is required for the growing state to be stable once acetate flows into a cell. A low Lg activity can be also be achieved in other ways. For example, bistability also becomes more likely when the amount (Fig 5o) or turnover number (Fig 5i) of Lg is low.
Another example is the anabolic enzyme An (Fig 5d). Bistability is more likely when it has a lower affinity for fbp; lower affinity increases fbp concentrations and favours the nongrowing state. This low affinity is complemented, however, by a strong binding affinity of An for pep (Fig 5e), so that if pep is present, the rate of the reaction draining fbp from the cell increases, and cells enter the growing state.
We note that neither the amount (Fig 5p) nor the turnover number (Fig 5j) of An has a strong effect on the probability of observing bistability. This is because both affect the maximum rate of anabolism, while bistability is primarily caused by the difference in the reaction rate of anabolism between the growing and nongrowing state.

The number of Cra proteins affects the probability of growth on acetate
To make sure that the number of Cra proteins per cell affects the probability that a cell resides in the growing and non-growing state, we investigated the consequences of changing the number of Cra  (Fig 6a), and it increases the steady-state growth rate at all acetate concentrations. A similar effect is visible in the stochastic model (Fig 6b).
Specifically, as the number of Cra molecules per unit biomass increases, the distribution of growth rates shifts from unimodal low growth rates (less than 0.1 h −1 with 50 Cra molecules per cell) to unimodal high growth rates (more than 0.2 h −1 with 100 Cra molecules per cell), with two modes clearly visible for Cra molecule numbers between 70 and 80 molecules.

Protein degradation does not remove the fitness benefit of noise when Cra-fbp binding is tight
Here we investigate the consequences of relaxing the assumption that proteins are not degraded. Proteins with a high turnover respond faster to changes in their production rate, and thus their concentrations are more sensitive to changes in gene expression [10]. Consequently, active degradation of the acetate incorporation protein Ai can make this protein more sensitive to variation in Cra activity, so that (i) cells may respond faster after the environment shifts to acetate and (ii) cell growth rate may become more vulnerable to random fluctuations in Cra numbers due to noise. Active protein degradation of Ai therefore has the potential to modify the fitness effects of Cra expression noise.
Active degradation of Ai may affect the circuit regulating a cell's response to an environment with only acetate, changing the probability at which cells occupy the growing and the nongrowing state. We simulated four degradation rates, with each 1000 cell lineages growing in an environment with a fixed concentration of acetate (Fig 7). We kept the average protein number per cell constant by simultaneously increasing the protein production rate k 1 (see Materials and Methods in main text). The distribution of the acetate incorporation protein Ai per cell shows two modes, with the high and low number of proteins corresponding to the growing and nongrowing state respectively. Increasing the degradation rate does not affect the positions of the two modes in the number of Ai proteins per cell. However, it does shift the population towards the nongrowing state. This decrease in the probability of growth occurs because fluctuations in Cra activity occur faster at high than at low growth rates, so that the growing state is more vulnerable to a random decrease in Cra activity than the nongrowing state is to a random increase.
We next considered how protein degradation may influence the bottleneck that populations experience after a switch to acetate. We chose a degradation rate at the upper end of the rates reported in E. coli (0.7 h −1 , BNID 109921, [11] and considered two Cra-fbp binding strengths (dissociation constants of 0.05 and 0.1 mmol g −1 , tighter and looser binding respectively). As in the main text, we modelled four levels of Cra expression noise. For each combination of Cra expression noise and the Cra-fbp dissociation constant, we simulated a population of cells growing for four days in an environment with constant nutrient flow. The environment switched from glucose to acetate after two days. We found that when Cra-fbp binding is tight (Fig 8a), noise in the expression of Cra makes a considerable difference in the trajectory of the population size. Instead of falling to a minimum and then recovering in population size, as happens in the two quieter populations (Cra expression noise η 2 = 0.01 and 0.1, smallest population size for η 2 = 0.1 is shown in Fig 8b), the populations with noisier Cra expression (η 2 = 1.0 and 10.0, the smallest population size for η 2 = 10.0 is shown in Fig 8b) decelerate their decline in population size much earlier, but continue falling until they reach the same cell number as the populations with quieter Cra expression (only shown for η 2 = 10 in Fig 8a). The result is that the two populations with higher Cra expression noise do not experience a bottleneck. The smallest population sizes in populations with low Cra expression noise are about 500 cells smaller than in populations with higher noise (Fig 8b).
Differences in the smallest population size are statistically significant between populations with different levels of Cra expression noise (pairwise Wilcoxon rank sum test with Holm's correction, p < 0.001) except between the two noisiest populations (p = 0.15). In populations with looser Cra-fbp binding, population sizes do not fall beneath the carrying capacity in acetate at any level of noise (Fig 8c).
Population decline temporarily stopped, and in the populations with the lowest Cra expression noise the population size even increases momentarily, before falling slowly towards the carrying capacity. As a result, the smallest population sizes are more indicative of the carrying capacity than a population bottleneck (Fig 8d). The mean smallest population size lies between 1820-1830 cells for all populations except for the noisiest (pairwise Wilcoxon rank sum test with Holm's correction, p = 0.8 for comparisons between all populations except with the noisiest population). In contrast, the noisiest population has a mean smallest population size of 1760 cells (± 30 cells), which is significantly lower (p < 0.001 for all pairwise comparisons with the noisiest population). This decrease in carrying capacity (also visible in Fig 8c) is a consequence of the greater sensitivity to random fluctuations in Cra activity, which reduces the probability of maintaining growth on acetate (Fig 7). This increased sensitivity of Ai amounts to Cra activity also drives a reduction in the carrying capacity of populations with tight Cra-fbp binding (Fig 8a) compared to simulations without protein degradation (Fig 2a in  We also investigated how protein degradation affects the relation between Cra expression noise and fitness at these two Cra-fbp binding strengths (Cra-fbp dissociation constants of 0.05 and 0.1 mmol g −1 ).
We again estimated fitness through competition with a reference population, which also underwent degra- noise occurs despite the fact that these populations resume growth the fastest (Fig 8a). However, this initial advantage is offset by the fact that these populations also have a lower carrying capacity (Fig 8d).
Fitness differences between populations with different Cra expression noise are an order of magnitude lower than in populations with tighter Cra-fbp binding and protein degradation (Fig 8e) or in the absence of protein degradation (Fig 3a in the main text).
We therefore conclude that when Cra-fbp binding is sufficiently tight and Cra expression noise is low to intermediate, populations will undergo a bottleneck after a switch to acetate regardless of whether the acetate incorporation enzyme is actively degraded or not. In addition, increasing Cra expression noise results in a fitness benefit. Because nongenetic variation is predicted to be beneficial only in populations with low mean fitness, and tight Cra-fbp binding leads to low fitness (Fig 3b in the main text), we can explore the effects of Cra expression noise in the absence of active protein degradation without loss of generality.

Invasion probabilities of a neutral allele
To find out whether the impact of noise on the spread of a beneficial allele is caused by selection, we simulated the introduction of a neutral allele into a population. We also used these simulations to investigate how Cra expression noise influences the process of drift, because random nongenetic variation can affect the strength of drift [12,13]. We simulated a population of 2000 cells, of which one carried a different allele at a neutral locus. We simulated growth in an environment that switches between glucose and acetate every two days, and continued the simulations until one of the two alleles had gone to fixation. The fixation of a neutral allele in an asexual haploid population occurs when the cell that carried that allele at the beginning of the simulation has become the ancestor of all surviving cells. Assuming that noise has no effect on the probability of fixation, such fixation occurs with a probability that is the inverse of the population size (N=2000) [14], and thus equals 0.0005. This probability of fixation is indeed close to the fraction of simulation runs in which the allele in question went to fixation (Fig 9a). The number of runs in which the neutral allele reaches fixation did not differ significantly between populations with different Cra expression noise (two-sided Fisher's exact test, p = 0.598, 50,000 replicates).
Greater nongenetic variation strengthens the influence of genetic drift on a population [12,13], and therefore neutral alleles that do reach fixation should do so faster in populations with higher Cra expression noise. We observed that increasing Cra expression noise does indeed shorten fixation times. The average time to fixation of a neutral allele (Fig 9b)  in the frequency of the neutral allele (Fig 9c). The standard deviations in allele frequency after four days are, from low to high noise, 0.004, 0.006, 0.010 and 0.015, an almost fourfold increase between the highest and the lowest level of noise. In contrast, the average allele frequency is close to 0.0005 (for the highest Cra expression noise, it is about 0.0006). The distribution of frequencies is significantly different between the simulations with the highest level of Cra expression noise and the other levels of noise, but not among the latter (Wilcoxon rank sum test with Holm's correction, p = 0.009 for pairwise comparison between η 2 = 10 and 0.01, p = 0.005 for comparisons between η 2 = 10 and 0.1 or 1.0, and p = 1 otherwise).
We therefore concluded that the differences in fixation probability of the beneficial allele that we observed in the main text are indeed due to its fitness difference to the wild-type allele. Also, we observed that while Cra expression noise does not affect the outcome of drift in terms of the fixation probability of a neutral allele, the fixation times are primarily driven by gene expression noise, and not by demographic stochasticity arising from small population sizes after the shift to acetate.
11 The benefits of Cra expression noise are less pronounced in environments that fluctuate randomly Here we investigate whether Cra expression noise also results in faster fixation of a beneficial allele in environments that are fluctuating randomly rather than periodically. We simulated an environment that switches between glucose and acetate on average every 48 hours and a standard deviation of also 48 hours, with the duration of each period drawn at random from a gamma distribution. We used the same four levels of Cra expression noise (0.01 to 10) and the same alleles as in the main text, that is, a beneficial allele that loosens Cra-fbp binding (Cra-fbp dissociation constant of 0.07 mmol g −1 ) and thus decreases lag times, and its competing wild-type allele (0.05mmol g −1 ). Each simulation run started with an initial population size of 2000 cells and half of the population carrying the beneficial allele. The simulation continued until one of the competing alleles reached fixation.
Due to the increased environmental stochasticity, we observed a greater number of simulations where either the beneficial allele was lost and/or the population went extinct (Fig. 10a). In fact, because we allowed the simulation to continue until the population had finished the current glucose-acetate cycle, we observed several replicates in which one of the alleles reached fixation and the population subsequently went extinct. For example, in populations with the lowest Cra expression noise, four replicates lead to population extinction after the fixation of the beneficial allele. Consistent with our results in periodically fluctuating environments, the probability that the beneficial allele fixes in a population that survives until the end of the simulation increases with noise. Specifically, it increases from 86% to 100% from the least noisy to the noisiest populations.
We also observed that the effect of Cra expression noise on the fixation times of the beneficial allele (if the allele reached fixation) is less pronounced than in environments that fluctuate periodically (Fig. 10b).
The beneficial allele did reach fixation the fastest in populations that had the second highest level of Cra expression noise (10 days ± 5 days standard deviation), three days faster than in the populations with the highest (13 days ± 5 days) and lowest noise (13 days ± 9 days), and a day faster than the populations with the second-lowest level of noise (11 days ± 7 days). The differences in fixation time were still significant between the highest level of noise and the second and third highest levels of noise (pairwise Wilcoxon's rank sum test with Holm's correction, p = 0.03 for the pairwise comparison of η 2 = 0.1 with η 2 = 10 and p = 0.007 for the comparison between η 2 = 1 and η 2 = 10). In sum, expression noise is less able to accelerate the spread of a beneficial mutation in highly unpredictable environments.