When is gene expression noise advantageous?

The variability of gene expression levels, also known as gene expression noise, is an evolvable trait subject to selection. While gene expression noise is detrimental in constant environments where the expression level is under stabilizing selection, it may be beneficial in changing environments when the phenotype is far from the optimum. However, expression noise propagates along the gene network, making the evolution of connected genes interdependent. Here, we explore how their position in the gene network constrains the evolution of genes under selection using an in silico evolution experiment. We simulate the evolution of populations of model gene regulatory networks under directional and fluctuating selection on the gene expression level while allowing the basal expression level and expression noise level to mutate. We find that expression noise is only transiently favoured under directional selection, but high levels of noise can be maintained in a fluctuating selection regime. Furthermore, target genes, regulated by other genes, were more likely to increase their gene-specific expression noise than regulator genes. These findings suggest that both the mean and variance of gene expression levels respond to selection due to changing environments – and do so in a network-dependent manner. They further point at gene expression noise as a putative mechanism for populations to escape extinction when facing environmental changes.


Introduction
Gene expression levels directly a ect the viability and tness of the organism.An important component of the gene expression pro le, in addition to the mean gene expression level, is the expression variability around the mean.While gene expression has been extensively measured using bulk transcriptomics, averaging expression levels over multiple cells, the advent of singlecell biology permits quantifying the intrinsic variability of expression levels between isogenic cells, rmly establishing gene expression as a stochastic process (Kaern et al., ).The cell-tocell variance in expression level, or in other words, the (un)predictability of expression level, is termed gene expression noise and stems from the typically low number of molecular players, such as transcription factors, within each cell.It has been shown to be an evolvable trait independent of expression mean.For example, stabilizing selection on gene expression level reduces gene expression noise (Lehner, ).However, the response of expression noise to other selection scenarios, such as directional and uctuating selection, has not yet been su ciently explored.
Expression noise has been suggested to be bene cial during the adaptation of the mean expression level to a new expression level optimum in genes under directional selection (Duveau et al., ).Increasing expression noise as a bet-hedging strategy for adapting to uctuating environments was experimentally demonstrated in previous studies.Bet-hedging is the phenomenon of increasing phenotypic variation, which decreases the short-term average tness of the population, but increases the long-term tness in changing environments, thereby allowing it to survive environmental uctuations (Grimbergen et al., ).In the same way, expression noise can create phenotypic heterogeneity in a clonal population and improve its capacity to survive in changing environments.For instance, the increased cell-to-cell variability of a signal transduction system in Escherichia coli permits growth in case of rapid oxygen availability uctuations (Carey et al., ).
It was shown that expression noise of the transcription factor comK drives cell fate determination in Bacillus subtilis and enables competency in a proportion of the cell population (Maamar et al., ), and increasing the noise increases the response range of the competency circuit (Mugler et al., ).A recent study has shown that the increase of gene expression noise at low growth rates confers a tness advantage to microbes in unpredictable environments and allows them to explore new phenotypes to adapt to unfavourable conditions (de Groot et al., ).Also, several observed properties of gene regulatory networks have been attributed to uctuating selection in gene networks (Tsuda and Kawata, ).
The response of expression noise to di erent selection scenarios has been studied in single genes, but not in the context of genetic networks.Since expression noise has been demonstrated to propagate in gene networks (Pedraza, ), it is important to take into account the network background of a gene under selection when investigating the evolution of expression noise.Changing the expression noise of a gene changes the expression noise of downstream elements because expression noise propagates between genes in the network.It has been argued that most genes are under stabilizing selection on expression level and, consequently, it can be assumed that low expression noise will be maintained.Increasing expression noise in some parts of the network as a response to uctuating or directional selection would increase noise in neighbouring parts of the network, which might con ict with stabilizing selection acting on other genes and decrease the e ciency of selection.Therefore, the adaptation to directional or uctuating selection on the gene expression level of individual genes might be constrained by their position in the gene network.
In this study, we used a computational gene regulatory network model to simulate the evolution of populations of model gene regulatory networks in two di erent selection scenarios, directional and uctuating selection acting on gene expression levels, to investigate the adaptability constraints imposed on individual genes by their gene network background.We found that higher gene expression noise is transiently bene cial under directional selection, while the mean expression level is evolving towards the new expression optimum.Higher gene expression noise is also bene cial under uctuating selection.Also, the evolvability of gene expression mean and noise is not independent of the network background.Namely, regulator genes are less likely to adapt to directional selection than non-regulator genes because the change to their mean expression level impacts the mean expression level of downstream genes.Regulators are also less likely to adapt to uctuating selection by increasing their expression noise as part of a bet-hedging strategy.

Material and Methods
To study the evolution of gene expression noise and mean expression level in gene regulatory networks in changing environments, we simulated the evolution of populations of model gene regulatory networks under directional and uctuating selection. .

Gene regulatory network model with evolvable gene expression mean and noise level
We used the formalism introduced by Wagner (Wagner, ) to model networks of genes using a regulatory network matrix = 1≤ ≤ , 1≤ ≤ , in which entries determines the presence, strength and sign of the e ect of gene on the expression level of gene .The sign of the elements of the regulatory network matrix determines whether the interaction represents an activation or repression of the downstream expression level, and the value determines the strength of the interaction (Figure A).The equilibrium expression level of each gene, if any, can be obtained and is considered to be the phenotype of the network.A tness value is then calculated as a distance to a given optimal expression level, the tness of the full network being a function of the tness of each individual gene.The model was extended to account for the intrinsic noise of each gene via gene-speci c parameters { int } 1≤ ≤ (Puzović et al., ) and basal expression level of each gene The genotype is realized into the phenotype through a set of di erence equations, updating the expression level ( ) of each gene in every time step using the following rule: where { basal } 1≤ ≤ represents the basal expression level of each gene, that is, the constitutive expression level of each gene that is present regardless of the input from regulatory genes.In this stochastic version, the expression level of each gene in each time step is drawn from a normal distribution, with a mean equal to the sum of the activation rate in the previous time step ( ) and the basal expression level, and a standard deviation equal to the intrinsic noise parameter.
The activation rate is de ned as the sum of the e ects of all regulators: If the drawn expression level value is below 0 or above = 100, its value is set to 0 or , respectively, thereby constraining the expression level range to an interval ( , ).The expression levels of each gene are synchronously updated in every time step for time steps ( = 50 in this study).The expression level vector { } 1≤ ≤ at the nal time step 50 is taken as the phenotype of the individual network.An individual genotype consists of the regulatory network matrix , the vector of intrinsic noise parameters { int } 1≤ ≤ , and the vector of basal expression levels .

Forward-in-time simulation of expression mean and noise evolution
Populations of model gene regulatory networks were evolved for = 10, 000 generations using an evolutionary algorithm consisting of repeated cycles of phenotype realization, reproduction, mutation and recombination (Figure B).First, all individuals in the population had their genotype realized into a phenotype.Next, the tness of each individual was calculated as a function of the distance of gene expression levels from an optimal gene expression level vector, weighted by the tness contribution of each gene, noted 1≤ ≤ : ( ) The tness contribution de nes the magnitude of the tness cost resulting from the deviance from the optimum expression level.In this study, we set ∀ , = 1, de ning an equally strong selective pressure on all genes.Individuals in the next generation were randomly sampled from the population at the current generation, with genotype probabilities proportional to their respective tness values.A constant population size was maintained.Mutations were introduced in the intrinsic noise vector and the basal expression level vectors with a probability of = 0.005 and = 0.005 per gene, respectively.The mutation values for intrinsic noise mutations were drawn from a uniform distribution  (0, 200), and the mutation values for basal expression level mutations were drawn from a uniform distribution  (0, 100 .

Simulation pipeline
We generated , unique -gene random (Erdős-Rényi) network topology samples with a network density of = 0.05 using the igraph package (Csardi and Nepusz, ) in R (version . ., Team ( )).Each unique network topology sample was used to simulate the evolution of one population of networks, which share the same network con guration but may di er in their basal expression level or expression noise.
Gene regulatory networks generated by randomly assigning the strength of the interactions usually had a network-wide gene expression pro le; many genes were expressed at the maximum value or not expressed at all.For a more realistic expression level pro le with intermediate gene expression levels, a network establishment step was performed through evolutionary simulations similar to the main noise evolution simulations.In this rst step, the strength of the interactions was allowed to mutate, but the overall network topology was xed; that is, non-null regulatory strengths cannot become zero, and regulatory strengths equal to zero cannot change value.The intrinsic noise of all genes was set to zero, and the values for basal expression levels were drawn from a uniform distribution  (0, ) to form the starting population.Both intrin-sic noise and basal expression levels were immutable.The target (optimal) expression level for all genes in this step was { } 1≤ ≤ = { 2 , ..., 2 }.For each of the , network topologies, a population of = 1, 000 individuals was evolved for = 5, 000 generations.The values for the strength of regulatory interactions were drawn from a uniform distribution  (−3, 3) with a mutation probability of = 0.001 per link.Networks which had oscillatory gene expression trajectories were removed.The ttest network at the end of the network establishment process had more genes with intermediate expression levels than the network in the starting population, and its equilibrium expression level values were used to set the optimal expression level in the rest of the simulation.The genotype of the ttest network was used as a founding network to create a uniform population for the next step, the burn-in.
Stabilizing selection was applied to all genes in each network population for = 5, 000 generations until the mutation-selection-drift balance was reached (burn-in step).In this step, the intrinsic noise was set to zero for all genes in the population of = 1, 000 individuals at the beginning of the simulation.The intrinsic noise and basal expression levels were allowed to mutate.Intrinsic noise mutations were drawn from a uniform distribution  (0, 200) with a mutation rate of = 0.01 per gene, and basal expression level mutations were drawn from a uniform distribution  (0, 100) with a mutation rate of = 0.01 per gene.The burn-in phase ensured that all genes in each network had already been evolving under stabilizing selection for long enough that mutation-selection-drift balance had been reached, and any selection response observed after an environmental shift would be due to directional or uctuating instead of stabilizing selection.The heterogeneous population of genotypes at the end of the burn-in phase was the starting point for the main simulation, where an environmental shift was applied.
One gene in each network was randomly chosen to undergo an environmental shift, in which its optimum expression level was increased or decreased by % of relative to its previous value.The remaining genes remained under stabilizing selection, with their respective optimum expression level values unchanged.Two selective scenarios were simulated: directional selection, in which the optimum expression level was either increased or decreased by % of , and uctuating selection, in which the optimum expression level was oscillating between an increase and a decrease of % of every other generation.The intrinsic noise and basal expression levels were mutable, and the network topology was immutable.The populations were evolved for = 10, 000 generations.The evolutionary simulation was replicated times for each selective scenario, and the resulting data (changes in phenotypes and genotypes) was averaged over the replicates for each gene in each scenario.
To understand the network background's speci c e ect on the evolvability of gene expression in changing environments, we also simulated the evolution of single, isolated genes under directional and uctuating selection.A dataset of , genes with a random basal expression level value drawn from  (20, 80) was used to generate populations of = 1, 000 individuals.
Their evolution was simulated using the previously described pipeline: each population evolved for = 5, 000 generations under stabilizing selection, after which directional or uctuating selection was applied for an additional = 30, 000 generations.As a control, we also let the populations evolve under stabilizing selection.To distinguish the response of the mean expression level and intrinsic noise, and their co-evolutionary dynamics, three scenarios were simulated: i) with immutable basal expression level and mutable intrinsic noise; ii) with both basal expression level and intrinsic noise mutable; iii) with mutable basal expression level and immutable intrinsic noise.

Results
To study the evolution of gene expression noise in regulatory networks when the optimal expression level under selection is not constant, we simulated model networks evolving under scenarios where the gene expression level is under directional or uctuating selection.In each network, one focal gene was randomly chosen to have its optimal expression level shifted, while the other genes in the network remained under stabilizing selection.The populations of model networks were evolved with mutable gene intrinsic noise and basal expression level, and a xed network topology.
To better understand the e ect of the network structure, we rst analyse the evolution of single, isolated genes and then add interacting genes.In order to distinguish the response of the mean expression level, intrinsic noise, and their co-evolutionary dynamics, we simulated three scenarios: (i) with mutable intrinsic noise, but xed basal expression level, (ii) with mutable intrinsic noise and basal expression level, and (iii) with mutable basal expression level, but xed intrinsic noise.We then studied the phenotypic response (mean expression level, and variance of expression level) and genotypic response (intrinsic noise parameters, and basal expression level parameters of the genotypes) of the populations evolved under directional or uctuating selection on the gene expression level. .

Gene expression noise transiently increases as a response to directional selection
To understand how gene expression level distributions respond to an environmental shift that imposes directional selection on gene expression level, we simulated the evolution of , isolated genes.The evolutionary trajectory of one example gene under directional selection in the rst scenario (mutable intrinsic noise, but xed basal expression level) is shown in Fig. A-D.
The mean expression level did not change as the basal expression level was immutable (Fig. A,   D).However, the expression variance and intrinsic noise increased (Fig. B, C), indicating that noise was bene cial in genes whose expression level was far from the optimum.In the dataset of , simulated genes, the average increase of the intrinsic noise was signi cantly higher in genes under directional selection than genes which remained under stabilizing selection (p-value < 2.2 × 10 −16 , Wilcoxon's test, Fig. E).
In the second scenario (mutable intrinsic noise and mutable basal expression level) we observed two distinct evolutionary phases.The adaptive phase, during which the basal expression level was evolving towards the new optimum expression level, and the postadaptive phase, af- These results showed that gene expression noise has a tness bene t if the mean expression level is not near the optimal level, i.e. the gene is under directional selection.We simulated the evolution of genes under directional selection, by either increasing or decreasing the optimal ex-pression level.We observed consistent results in the two cases and for concision only reported the results from directional selection with an increased optimal expression level in the main text, while the results from directional selection with a decreased optimal expression level can be found in Section of the Supplementary Information.We observed that expression noise remains elevated under directional selection if the mean expression level is constrained and cannot change.To assess the e ect of the interaction strength on the adaptation propensity, we tted linear models with two node-level network metrics (absolute node instrength and absolute node outstrength), the results of which are summarized in Table .Absolute instrength and absolute outstrength had a signi cant e ect on the adaptation probability under directional selection.Absolute instrength, a metric of how strongly a gene is being regulated by other genes, had a signi cant positive e ect on the adaptation probability (Table ).Contrary, absolute outstrength, a metric of how strongly a gene is regulating other genes, had a signi cant negative e ect on the adaptation probability.There is no signi cant correlation between the two network metrics (Spearman's = −0.02,p-value = 0.28), but there is a relationship in which genes with high values of instrength have an outstrength of zero, and vice versa, i.e. the distributions of both network metrics are zero-in ated (Supplementary Information, Section ).Therefore, we tted the same linear models in subsets of each gene regulatory category.In target genes, which have a nonzero value of absolute instrength, and an absolute outstrength of zero, the absolute instrength did not have a signi cant e ect on the probability of adapting to selection (Table ).In intermediate transcription factors, which have non-zero values for both absolute instrength and outstrength, absolute instrength had a signi cant positive e ect on the probability of adaptation, while absolute outstrength had a signi cant negative e ect on the probability of adaptation, echoing the results from the joint dataset.In master transcription factors, genes which have non-zero values for outstrength, but a zero value for instrength, the absolute outstrength had a signi cant negative e ect on the probability of adaptation.Notably, the negative e ect size of the absolute outstrength on the adaptation probability was consistently higher than the positive e ect size of the absolute instrength.
A similar pattern was observed for the speed of adaptation.Again, the absolute instrength had a signi cant positive e ect on the number of generations until the new optimum was reached (Table ), while absolute outstrength had a signi cant negative e ect.Absolute outstrength also had a signi cant negative e ect on adaptation speed among transcription factors.In target genes and intermediate transcription factors the two network metrics did not have a signi cant e ect on adaptation speed.
Lastly, we investigated the response of the expression variance (phenotypic noise), intrinsic expression noise and basal expression level during the adaptation and postadaptation phases in the genes that managed to adapt to the new expression level optimum.As in the case of isolated genes, genes in networks under directional selection had signi cantly increased their expression variance (p-value < 2.2 × 10 −16 , Wilcoxon's test) and intrinsic noise (p-value < 2.2 × 10 −16 , Wilcoxon's test) during the adaptation phase, relative to their pre-optimum shift and relative to genes that remained under stabilizing selection.A signi cant increase in expression variance and intrinsic noise was observed in all three gene regulatory categories (Fig G , H).After adaptation, the expression variance and intrinsic noise decreased to similar values than before the optimum shift.Genes that remained under stabilizing selection after the directional selection regime was applied to a random gene in their network also showed a statistically signi cant, but extremely small, increase in average expression variance and intrinsic noise, but not in their basal expression level.This small average increase in expression variance and intrinsic can be explained by their connection to genes under directional selection, whose response to directional selection would a ect local genes through noise propagation and in ate the expression variance of the network. .

Gene expression noise is increased as a response to uctuating selection
To investigate whether elevated gene expression noise could be maintained in constantly changing environments, we simulated the evolution of isolated genes under uctuating selection on the gene expression level.This uctuating selection was implemented by imposing switches of the optimum expression level every other generation so that the population did not have time to adapt to the new optimum before a new one occurred.As in the directional selection case, we rst report results from the simulations of single, isolated genes and then study the case where the gene is part of a regulatory network.For the isolated gene evolution simulations, we also simulated three Since uctuating selection imposes a new optimal expression level every other generation, it is not bene cial for the population to evolve towards any single optimum.Because the evolution of the mean expression level is constrained, the increased expression noise as a response to uctuating selection is maintained, as in the case of single genes under directional selection with immutable mean expression level (Fig. E).In the dataset of , genes, the average change of intrinsic noise relative to pre-selection values was signi cantly higher in genes under uctuating selection than in genes that remained under stabilizing selection (p-value < 2.2 × 10 −16 , Wilcoxon's test) (Fig. J).
The average change of basal expression level was signi cantly di erent between genes under uctuating and stabilizing selection (Fig. K), though the di erence is negligibly small.
In the third scenario, in which the basal expression level was mutable, but the intrinsic noise was not, the mean expression level and basal expression level showed large uctuations between the two expression level optima (Fig. L, O).However, the expression variance Fig. M) was not increased as in the second scenario, in which the mean and noise could evolve.This suggests that the increase of expression variance in the second scenario (mutable intrinsic noise and basal expression level) cannot be attributed to the population heterogeneity of basal expression level genotypes while the population was evolving towards a new peak.Instead, the increase of expression variance in the second scenario re ects the increase of intrinsic noise as a tness bene t in response to uctuating selection on the gene expression level.
The increase of intrinsic noise after the new selection regime was applied could also be explained by a relaxation of purifying selection.Namely, the drastic environmental shift or alternating environmental shifts throw the population far from its adaptive peak and drastically reduce the average tness of the population.If the tness landscape is signi cantly rugged, the population may be trapped in regions of very low tness and selection may be impeded.More deleterious mutations may accumulate, i.e. higher intrinsic noise variants would be tolerated and not removed by selection, and the population-wide intrinsic noise would increase.However, if that were the case, these mutations would be expected to segregate at low frequency, resulting in a population-wide heterogeneity with a mixture of low-and high-noise variants, which is not observed in our study.After the environmental shift happens, expression variance and intrinsic noise are consistently increased (Fig. E, J) and the populations after uctuating selection are quite homogeneous in high intrinsic noise genotypes.Therefore, the noise increase is a signal of selective advantage, not a spurious signal from a relaxation of selection and more pronounced neutral evolution.
These results demonstrate that gene expression noise is selectively favoured in genes under uctuating selection on gene expression level, which can be understood as the evolutionary emergence of a bet-hedging strategy.Lastly, we looked into how the network background of the gene under uctuating selection a ects its evolutionary propensity.
. Target genes in gene regulatory networks respond more strongly to uctuating selection than non-target genes We investigated whether the gene regulatory network background has an e ect on the evolvabil- The adaptation of gene expression in the case of uctuating selection is not an adaptation of mean expression level, but of expression noise, since the mean expression level is prevented from evolving to any optimum by its frequent switches.Therefore, we classi ed genes as "adapted" to uctuating selection if their relative change of intrinsic noise parameters was higher than , which represents a twofold increase of intrinsic noise relative to pre-selection switch values.The proportion of genes di ered signi cantly between the three gene regulatory categories (p-value < 2.2 × 10 −16 , proportion test).Out of all genes put under uctuating selection, % of target genes (TGs) increased their intrinsic noise more than twofold, as opposed to % of intermediate transcription factors (interTF) and % of master transcription factors (masterTF To study the e ect of the interaction strength on the noise adaptation probability, we tted linear models with absolute node instrength and absolute node outstrength as node-level network metrics.The results are summarized in Table .In the joint dataset, the absolute instrength had a signi cant positive e ect on the adaptation probability, while absolute outstrength had a signi cant negative e ect on adaptation probability.The absolute outstrength also had a signi -cant negative e ect on adaptation probability in intermediate transcription factors and in master transcription factors (Table ).
Finally, we looked into the evolutionary trajectory of the expression variance (phenotypic noise), intrinsic noise and basal expression level of genes in gene regulatory networks undergoing uctuating selection.In the entire dataset, consisting of , genes from , -gene network topologies, genes under uctuating selection had a signi cantly higher increase of expression Notably, many genes that remained under stabilizing selection also signi cantly increased their expression variance and intrinsic noise after another gene in the same network was put under uctuating selection.In fact, the average change of expression variance and intrinsic noise for genes that remained under stabilizing selection was signi cantly higher than zero, although still signi cantly lower than for genes under uctuating selection.This e ect can be explained by noise propagation within the gene regulatory networks; that is, genes propagate their expression noise downstream, and increasing the intrinsic noise of one gene would increase noise propa-gated to downstream elements in the network.If a gene is a transcription factor and responds to selection by increasing its intrinsic noise, its increased expression variance would propagate to the downstream genes, increasing their expression variance, as well as the overall expression variance in the network.Similarly, if a target gene is under uctuating selection, it can increase its expression variance by increasing its intrinsic noise, or by increasing the intrinsic noise of any upstream element, which would then propagate the increased noise to the focal gene.These cases showcase that the evolution of gene expression noise, and gene expression level, is not independent among genes in a gene network.Furthermore, an adaptive property of gene expression, such as expression noise, for a gene under selection can come not from the intrinsic gene properties (cis regulatory elements), but from a gene it is connected to in the network (trans regulatory elements).On the other hand, the same mechanism can prevent a transcription factor from increasing noise to respond to uctuating selection, as it would propagate noise to other genes in the network, which would be deleterious other parts of the network are under stabilizing selection.Target genes do not have the same constraint -as they do not have downstream elements -and are free to adapt to uctuating selection by increasing their intrinsic expression noise.

Discussion
We investigated how natural selection in changing environments a ects the gene expression mean and noise levels in gene regulatory networks.We hypothesized (i) that high levels of expression noise may be favoured after an environmental shift because of the increased proportion of t phenotypes in the population and (ii) that genes di er in adaptability because of the di erential constraints resulting from the structure of the regulatory network.To test these hypotheses, we developed a gene regulatory network evolution model with evolvable basal expression level and intrinsic expression noise.We simulated the evolution of gene expression mean and noise in populations under directional and uctuating selection.We found that expression noise was transiently increased during adaptation but became counter-selected as the population reached the new optimum and the selection regime switched back to stabilizing selection.Maintaining higher levels of expression noise requires a regime of uctuating selection where optima are frequently changing, before the system enters a stabilizing selection phase.Importantly, regulator genes had a lower probability of responding to selection and responded less strongly than non-regulator genes, indicating a constraining e ect of the gene regulatory network on the adaptability of its constituent genes.In the following, we further discuss the implications of these results.

Expression noise and adaptive walks
Our model implements two types of mutations, a ecting the basal expression level and each gene's expression noise, respectively.In this framework, the two types of mutations are independent.
We showed that both mutation types are favoured when the expression level is away from the optimum: mutations a ecting the basal expression level bring the population closer to the tness peak, while mutations increasing the expression noise may increase the proportion of t phenotypes generated by chance.However, the putative advantage of the second mutation type decreases as the genotypes in the population get closer to the optimal expression levels.High noise levels, therefore, can only persist if the mean expression level is kept away from the optimum level, as we showed with our toy model where the basal expression level was prevented from mutating.Such a situation also arises in our uctuating selection regime, where optimal expression levels alternate every generation, preventing the mean expression level from adapting.These results explain why genes involved in the immune system typically display high levels of expression noise (Shalek et al., ), which may result from the optimal expression level being determined by constantly changing pathogens.More generally, we posit that high levels of expression noise should be generally expected in any gene involved in an evolutionary arms race.
While higher expression noise levels do not bring the populations closer to the tness peak and may only be transiently advantageous, mutations increasing the noise level may still prove crucial to the adaptation process.High levels of expression noise increase the phenotypic heterogeneity in the population, ensuring that some proportion of the population has a t phenotype, whichever the environmental conditions it nds itself in; an evolutionary strategy known as bethedging (Grimbergen et al., ). Bet-hedging may permit the population's survival by allowing the production of individuals with adaptive phenotypes by chance, triggering an evolutionary rescue mechanism where a declining population manages to survive an environmental shift and recover.In this study's simulation framework, the population size is kept constant; the individuals are reproduced into the next generation by sampling until the xed population size is reached.
Consequently, the population never goes extinct, even if the individuals have extremely low tness, and the population will have time to accumulate potentially bene cial mutations and adapt.
Conversely, natural populations may go extinct because of mutation load due to the population being far from its tness optimum.The evolution of bet-hedging might be more apparent if the populations that did not evolve higher levels of expression noise went extinct instead of being propagated to a constant population size.The simulation framework used here can be modi ed to include the possibility of extinction events by adding a variable population size and a viability threshold.Extending the model in such a way would shed light on the putative role of expression noise as an evolutionary rescue mechanism.

Network-driven epistasis
Adaptive processes have been analyzed using the framework of adaptive landscapes, rst introduced by Sewall Wright (Wright, ), where populations of genotypes "walk" through a landscape consisting of tness peaks and valleys, corresponding to genotype con gurations that are more or less t to the given environment.The adaptability of a population, that is, its capacity to reach a certain tness peak and the rate at which it does so, depends on many factors, such as population genetics parameters (e.g. e ective population size), the initial frequency of adaptive alleles, and the genetic architecture of the selected trait (Olson-Manning et al., ).A known factor that may slow adaptation is the non-additivity of mutations' e ects, or epistasis, which generates so-called rugged landscapes (Bank, ).One epistasis component results from the gene network, where the expression levels of (directly or indirectly) connected genes are not in-dependent.Here, we show that network-driven epistasis leads to a di erential adaptability of connected genes.Genes with a large regulatory output are less likely to get closer to a new optimum value as their expression may pull multiple other genes away from their respective optima.
The potential bene t of increasing their expression noise may also not compensate for the cost of ter the basal expression level has reached the optimum expression level and the selection regime switched back to stabilizing selection.The evolutionary trajectory of one example gene is shown in Fig. F-I.A co-evolutionary pattern was observed -expression level variance, determined by the intrinsic noise, was elevated while the basal expression level was evolving towards the new optimum (Fig. F, G).Expression noise and intrinsic noise were reduced once the mean expression level reached the optimum expression level (Fig. G, H).The average increase of the intrinsic noise during the adaptive phase was signi cantly higher in genes under directional selection than genes which remained under stabilizing selection (p-value = 1.25×10 −10 , Wilcoxon's test) or genes under directional selection in the postadaptive phase (p-value = 8.19 × 10 −07 , Wilcoxon's test) (Fig. J).As expected for genes under directional selection, the average change of the basal expression level was higher in genes under directional selection than genes which remained under stabilizing selection in both adaptive (p-value < 2.2 × 10 −16 , Wilcoxon's test) and postadaptive phase (p-value < 2.2 × 10 −16 , Wilcoxon's test) (Fig. K).In the third scenario (mutable basal expression level, but xed intrinsic noise) we again observed the adaptive and postadaptive phases of expression mean evolution.The evolutionary trajectory of an example gene is shown in Fig.L-O.There is no increase of expression level variance during the adaptive phase (Fig.M), contrary to the increase of expression level variance due to increased intrinsic noise in the second scenario (Fig.G, H).Importantly, the lack of a signi cant increase in expression level variance during the adaptive phase when the intrinsic noise cannot evolve indicates that the basal expression level variants segregating in the population did not cause the elevated expression variance signal in the scenario in which both mean and noise could evolve (Fig.G).Therefore, the elevated expression variance signal stems from the selectively favoured increase of intrinsic noise during the period in which the mean expression level is evolving towards a new expression optimum.
However, if the mean expression level can evolve towards a new optimum, a switch in the selection regime occurs as the population gets closer to the new optimum -the expression level is under stabilizing selection again, and expression noise becomes deleterious.Therefore, expression noise is bene cial transiently while the expression mean is evolving towards a new optimum, but becomes deleterious once it has reached the new peak.If the expression mean cannot evolve, expression noise has a constant tness bene t.This may happen if the evolution of the mean expression level is constrained by the position of the gene in the regulatory network, which is what we explore in the next section. .Regulator genes in gene regulatory networks are less adaptable to directional selection than non-regulator genes Next, we investigated whether the network background has an e ect on the evolvability of genes connected in a regulatory network.The evolutionary trajectory of a gene under directional selection from an example -gene network is shown in (Fig A-D).Before a directional selection regime was applied, the population was at mutation-selection-drift equilibrium and the expression level of the focal gene was under stabilizing selection.Consequently, the population-wide mean expression level showed little variation (Fig A, MSD balance phase).After the optimal expression level changed, the population mean expression increased toward the new optimum.We categorized genes as adapted if their mean expression level reached ±10% .As in the single gene case, we distinguish two evolutionary phases in genes whose expression level responded to directional selection: i) adaptation, during which the expression mean is evolving towards the new optimum, and ii) postadaptation, after the expression mean has converged toward a new op-timum.During these two phases, we track the changes in the expression variance or phenotypic noise (Fig B), the intrinsic noise genotypes (Fig C) and the basal expression level genotypes (Fig D).In the depicted example gene, the gene expression level converged toward the new optimum.During the adaptation phase, the phenotypic and intrinsic noise increased.Conversely, in the postadaptation phase, the phenotypic and intrinsic noise decreased to their average values before the optimum shift.Out of genes under directional selection, ( %) adapted after , generations.The regulatory function of the selected gene had a signi cant e ect on the probability of adaptation and the time to reach the new optimum.We categorized genes into three regulatory categories based on their regulatory interactions: i) master transcription factors (masterTFs), which regulate one or more downstream genes but are not regulated by any upstream gene; ii) intermediate transcription factors (interTFs), which are both regulated and regulate one or more genes; and iii) target genes (TGs), which are regulated by one or more upstream genes, but do not regulate any downstream gene themselves.A higher proportion of target genes adapted to directional selection than intermediate transcription factors and master transcription factors (p-value < 2.2 × 10 −16 , proportion test).Out of all genes put under directional selection, % of target genes (TG) managed to adapt to the new expression level optimum, contrary to % of intermediate transcription factors (interTF) and % of master transcription factors (masterTF), indicating a constraint on the evolvability of transcription factor expression levels (Fig E).The time to adaptation, de ned as the number of generations until the new optimum expression level is reached, also di ered between the three regulatory categories.Adaptation was signi cantly slower in master transcription factors (p-value < 2.2 × 10 −16 , Wilcoxon's test) and intermediate transcription factors (p-value < 2.2 × 10 −16 , Wilcoxon's test) than in target genes (Fig F).
As expected, genes under directional selection increased their basal expression level (Fig I), as opposed to genes under stabilizing selection.These results show that the position of the gene in the gene regulatory network imposes constraints on the evolvability of its expression level in the face of one-time environmental shifts that would impose directional selection.Namely, transcription factors are less likely to adapt to a new optimal expression level, and do so more slowly.Next, we explore what happens in genes experiencing periodic environmental shifts, which impose a uctuating selection regime on the gene expression level.
scenarios: i) with mutable intrinsic noise only, ii) with mutable both intrinsic noise and basal expression level, and iii) with mutable basal expression level only.The evolutionary trajectory of one example gene under uctuating selection in the rst scenario is shown in Fig. A-D.The optimal expression level alternates every second generation, and the mean expression level does not adapt to any single optimum (Fig. A), as the basal expression level is xed in this scenario (Fig. D).However, the expression variance (Fig. B) and intrinsic noise (Fig. C) are signi cantly increased after uctuating selection is applied.In the entire dataset of , simulated genes, the average change of intrinsic noise relative to pre-selection values was signi cantly higher in genes under uctuating selection than in genes that remained under stabilizing selection (p-value < 2.2 × 10 −16 , Wilcoxon's test) (Fig. E).Allowing the basal expression level to evolve (second scenario, mutable intrinsic noise and mutable basal expression level) did not alter the observed pattern of increased intrinsic noise after the uctuating selection regime was applied.The evolutionary trajectory of one example gene is shown in Fig. F-I.The mean and basal expression levels did not signi cantly change under uctuating selection (Fig. F, I), even though the basal expression level was mutable in this scenario.
ity of genes under a uctuating selection regime.The evolutionary trajectory of a gene under uctuating selection from an example -gene network is shown in Fig A-D.In the mutationselection-drift balance phase, before uctuating selection was applied, the expression level of the focal gene was under stabilizing selection.Consequently, the population-wide mean expression level shows little variation (Fig A, MSD balance).The mean expression level persisted after uctuating selection was applied.However, the phenotypic noise (Fig B) and the intrinsic noise (Fig C) increased as a response to uctuating selection, as in the case of an isolated gene under uctuating selection.
the increased propagated noise to target genes under stabilising selection.Reducing the number and strength of interactions between genes, i.e. reducing the connectivity of the gene regulatory network, would reduce epistasis and thereby increase adaptability to potential selective pressure in the future.Real biological gene regulatory networks are, indeed, sparse (Leclerc,), and network sparsity was shown to be an emergent property resulting from optimizing the explorability of new phenotypes (Busiello et al., ).In a recent study, Burban et al. ( ) used a Wagner gene network evolution model and showed that domestication, here implemented as a strong directional selection combined with a bottleneck, led to strong rewiring of the gene regulatory network.More speci cally, they observed that domestication resulted in an increase in the number of connections.Several factors may explain the apparent discrepancy between our conclusions and the results of Burban et al. ( ).First, the model of Burban et al. ( ) does not implement intrinsic expression noise while the network structure is xed in the model we introduced.In particular, the Burban et al. ( ) model did not implement any cost to increased complexity, while we showed that noise propagation imposes a strong constraint on gene evolution, in a network centrality-dependent manner.Adding expression noise to the Burban et al. ( ) o ers an interesting perspective for understanding the impact of domestication (and more generally, adaptation) on the evolution of regulatory interactions.Ultimately, while being computationally challenging, an extended model where basal expression levels, intrinsic expression noise, and gene interaction evolve would provide insights into the evolution of network complexity in relation to adaptation.
variance than genes in the rest of the gene network that remained under stabilizing selection (pvalue < 2.2 × 10 −16 , Wilcoxon's test).Genes under uctuating selection also had a higher increase of intrinsic noise than genes under stabilizing selection (p-value < 2.2 × 10 −16 , Wilcoxon's test).The increased expression variance and intrinsic noise trend held across all three gene regulatory categories (Fig F-H).However, the degree of increase of expression variance and intrinsic noise under uctuating selection was dependent on the position of the gene in the gene regulatory network.Speci cally, target genes under uctuating selection had a signi cantly higher increase