Analysis of evolved ultrasensitive networks.

<p>(A) Average saturation of enzymes in all of the evolved ultrasensitive networks. The average saturation of enzymes is calculated as the geometric mean of individual Michaelis-Menten constants of the different kinases and phosphatases and their allosteric states normalised by the substrate concentration (kinase, <i>K</i>_<i>1</i>, or phosphatase, <i>K</i>_<i>2</i>). The shape of each data point represents different starting structures to the evolutionary simulations (see <a href="http://www.ploscompbiol.org/article/info:doi/10.1371/journal.pcbi.1004918#pcbi.1004918.s001" target="_blank">S1 Fig</a>). The colours of the data points represent two different evolutionary scenarios; blue: output protein [<i>S</i>_<i>total</i>] = 1000, other signalling proteins (<i>A</i>*) concentrations [<i>A</i>*_<i>total</i>] = 1; red: output protein [<i>S</i>_<i>total</i>] = 10, other signalling proteins concentrations [<i>A</i>*_<i>total</i>] = 10. The blue and red, star-shaped points indicate the average value of the enzyme saturation resulting from these initial concentrations at the start of the evolutionary simulations. Each data point is further labelled with the unique identification number used for each evolutionary simulation. (B) The fraction of different forms of the kinase (<i>y</i>-axis) against the ligand concentration (<i>x</i>-axis) for two different evolved networks (network 20 and 18 in <a href="http://www.ploscompbiol.org/article/info:doi/10.1371/journal.pcbi.1004918#pcbi.1004918.s002" target="_blank">S2 Fig</a>). The fractions of the different forms of the kinase are the substrate-accessible (green), substrate-inaccessible (orange), and substrate-bound (blue) forms. This data is overlaid with the dose-response dynamics; the solid and dashed lines show the steady state concentration of phosphorylated (i.e. response) and unphosphorylated substrate respectively at a given input level. (C) Ratio between <i>K</i>_<i>M</i> values of different conformational states (relaxed “R” state and tensioned “T” states) for kinase (<i>x</i>-axis) and phosphatase (<i>y</i>-axis). The colours, shapes and numbers on the dots are the same as in (A). For enzymes without allosteric regulation the ratio are set to one, so that there are no distinctive conformational differences. (D) The fraction of different forms of the phosphatase (top) and kinase (bottom) (<i>y</i>-axis) against the ligand concentration (<i>x</i>-axis) for two different evolved networks (network 18 and 23 in <a href="http://www.ploscompbiol.org/article/info:doi/10.1371/journal.pcbi.1004918#pcbi.1004918.s002" target="_blank">S2 Fig</a>). The different forms of the enzymes are the different conformational states, relaxed “R” state (green) and tensioned “T” state (orange). These are overlaid with dose-response dynamics; the solid and dashed lines show the steady state concentration of phosphorylated (i.e. response) and unphosphorylated substrate respectively.</p

Plasmids and strains used and the associated literature source.

<p>Plasmids and strains used and the associated literature source.</p

A cartoon diagram of the CheA3-CheA4-CheY6 split kinase system.

<p>The diagram is arranged so to highlight the role of free CheA3 acting as a branching point for the two arms that form competing cycles leading to phosphorylation and dephosphorylation of CheY6. Rate constants are shown on the relevant reactions. In the case of reversible reactions, two rate constants are given (<i>k</i>_forward/<i>k</i>_reverse).</p

Analysis of evolved adaptive networks.

<p>(A) Structure and dynamics of the evolved adaptive network number 1. The upper panel shows a cartoon of the network. The oval shapes represent ligand (top) and the output protein (bottom) (e.g. substrate with a phosphorylation site, <i>S</i>), while rectangles represent all other signalling proteins (e.g. receptors, scaffolds, kinases, or phosphatases). Black lines represent the binding reaction between two sites. Red arrows represent a kinase site (red) phosphorylating a substrate phosphorylation site (purple). Blue arrows represent a phosphatase site (blue) dephosphorylating a substrate site. The green rectangle with a self-pointing green arrow indicates a protein domain whose conformational switching is allosterically regulated. The lower panels show the temporal dynamics of the input ([<i>L</i>], top) and the output protein (black lines in lower panels). The dynamics of the output protein is overlaid with the distribution of the different kinase ([<i>K</i>], middle panels) and phosphatase ([<i>P</i>], bottom panels) complexes: blue for enzyme-substrate complexes, green for free form of the enzymes that are accessible by the substrate, and red for complexes of enzymes with other signalling proteins. (B) Structure and dynamics of the evolved adaptive network number 2. Panels are as in (A).</p

Time-course analyses.

<p>The model is simulated with increasing and decreasing signal levels (<i>k₅</i>) in course of time. <i>k₅</i> is increased from 2 to 6 and decreased in similar fashion at indicated time points (top most, left panel), and changes in each species were measured (as indicated on each panel). The dotted line represents the highest signal level, with equal signal steps on each side of it. The noted asymmetry around this line shows the presence of hysteresis in the system. The x- and y-axis represent time and species concentration respectively, where the latter is normalized by the appropriate total protein levels.</p

Literature source and parameter values used in the analysis of the basic model.

<p>Literature source and parameter values used in the analysis of the basic model.</p

Effects of varying key parameters of the model and addition of different phosphatases.

<p>The x- and y-axis show the signal (<i>k₅</i>) level and the corresponding steady state CheY6-P level respectively. Each panel shows a signal-response analysis for varying model parameters (A–C) or the inclusion of additional phosphatases (D). The results of the basic model are shown in red. Where present, the dark region indicates the region of unstable steady states and hence the presence of bistability. Arrows on panels A, B and C indicate increasing value of the changed parameter. (<b>A</b>) The on rate (<i>k₁</i>) for CheA3:CheA4 complex formation was varied from basic model value [100(µMs)⁻¹] to 10, 1, and 0.208. (<b>B</b>) Concentration of CheA4 was varied from 30 µM, 40 µM (basic model) and 80 µM. (<b>C</b>) The rate of CheA3 mediated dephosphorylation of CheY6-P (<i>k₁₁</i>) was varied from 1 s⁻¹, 2.5 s⁻¹ (basic model) and 5s⁻¹. (<b>D</b>) The basic model has free CheA3 as the sole phosphatase; the effect of having either CheA3-P or CheA3:CheA4 and CheA3:CheA4:ATP as additional phosphatases is shown. See also Figures S1, S2, S3, S4 for additional sensitivity analyses.</p

Designed signalling cycle motif and parameter space for adaptation and ultrasensitivity.

<p>(A) Cartoon showing the designed signalling cycle motif with a sequestering protein. The sequestrating protein (<i>T</i>) binds both the kinase (<i>K</i>) and phosphatase (<i>P</i>), which catalyse the phosphorylation and dephosphorylation of unphosphorylated (<i>S</i>) and phosphorylated (<i>S</i>_<i>p</i>) substrate respectively. (B) The values of key parameters for achieving ultrasensitive (> 0.8) and adaptive (>0.3) responses, when assuming an enzyme-saturated regime ([<i>S</i>_<i>total</i>] = 1, [<i>P</i>_<i>total</i>] = 0.1). Panels on the left show the distribution of Michaelis-Menten constants, for kinase: <i>K</i>_{<i>M</i>,<i>K</i>} = (<i>k</i>₂+<i>k</i>₃)/<i>k</i>₁ (<i>x</i>-axis) and phosphatase <i>K</i>_{<i>M</i>,<i>P</i>} = (<i>k</i>₅+<i>k</i>₆)/<i>k</i>₄ (<i>y</i>-axis). Panels on the right show the distribution of affinities of sequestrating protein <i>T</i> with kinase and phosphatase: <i>K</i>_{<i>D</i>,<i>K</i>} = <i>k</i>₈/<i>k</i>₇ and <i>K</i>_{<i>D</i>,<i>P</i>} = <i>k</i>₁₀/<i>k</i>₉ Note that all four panels are plotted on the same logarithmic range. Each black dot represents a parameter set and the colours shows density of parameters.</p

Nonlinear Dynamics in Gene Regulation Promote Robustness and Evolvability of Gene Expression Levels

<div><p>Cellular phenotypes underpinned by regulatory networks need to respond to evolutionary pressures to allow adaptation, but at the same time be robust to perturbations. This creates a conflict in which mutations affecting regulatory networks must both generate variance but also be tolerated at the phenotype level. Here, we perform mathematical analyses and simulations of regulatory networks to better understand the potential trade-off between robustness and evolvability. Examining the phenotypic effects of mutations, we find an inverse correlation between robustness and evolvability that breaks only with nonlinearity in the network dynamics, through the creation of regions presenting sudden changes in phenotype with small changes in genotype. For genotypes embedding low levels of nonlinearity, robustness and evolvability correlate negatively and almost perfectly. By contrast, genotypes embedding nonlinear dynamics allow expression levels to be robust to small perturbations, while generating high diversity (evolvability) under larger perturbations. Thus, nonlinearity breaks the robustness-evolvability trade-off in gene expression levels by allowing disparate responses to different mutations. Using analytical derivations of robustness and system sensitivity, we show that these findings extend to a large class of gene regulatory network architectures and also hold for experimentally observed parameter regimes. Further, the effect of nonlinearity on the robustness-evolvability trade-off is ensured as long as key parameters of the system display specific relations irrespective of their absolute values. We find that within this parameter regime genotypes display low and noisy expression levels. Examining the phenotypic effects of mutations, we find an inverse correlation between robustness and evolvability that breaks only with nonlinearity in the network dynamics. Our results provide a possible solution to the robustness-evolvability trade-off, suggest an explanation for the ubiquity of nonlinear dynamics in gene expression networks, and generate useful guidelines for the design of synthetic gene circuits.</p></div

<i>In silico</i> simulation under fluctuating environments and adaptation time for genotypes that are robust and evolvable.

<p>(a): Monostable genotypes corresponding to the population mean of parameter values arising from different generations of a single <i>in silico</i> simulation under fluctuating environments, superimposed onto the G-P mapping of circuit I (grey backdrop). Genotypes derived from each generation are color-coded according to generation number, with red indicating the oldest generations. For both the G-P mapping and the <i>in silico</i> simulations, evolvability was computed using parameter perturbations of ±20% (cf. Figure F in <a href="http://www.plosone.org/article/info:doi/10.1371/journal.pone.0153295#pone.0153295.s001" target="_blank">S1 File</a>)), and robustness was calculated using the sensitivity-based measure. See also Figure G in <a href="http://www.plosone.org/article/info:doi/10.1371/journal.pone.0153295#pone.0153295.s001" target="_blank">S1 File</a>, which shows the results obtained for other runs of the evolutionary simulation algorithm. (b): Box plots showing adaptation time for genotypes exhibiting high evolvability (<i>H</i>) and low evolvability (<i>L</i>) under two selection scenarios: selection for a 10% reduction (blue) and for a 10% increase (red) in steady state expression levels. For circuit I, high evolvability (<i>H</i>) was <i>E</i> > 0.7 and low evolvability (<i>L</i>) was <i>E</i> < 0.1 on G-P mappings using robustness and evolvability measures computed from 10% parameter perturbations (cf. <a href="http://www.plosone.org/article/info:doi/10.1371/journal.pone.0153295#pone.0153295.g002" target="_blank">Fig 2c</a>). For circuit II, genotypes from <i>H</i> had <i>E</i> > 1.2 and genotypes from <i>L</i> had <i>E</i> < 0.6 based upon the same measures (cf. <a href="http://www.plosone.org/article/info:doi/10.1371/journal.pone.0153295#pone.0153295.g002" target="_blank">Fig 2f</a>). For both circuits, 11000 genotypes were picked from each respective region in the G-P mapping. For both circuits, adaptation times for the two regions differed significantly (<i>p</i> < 0.01, two-sample Kolmogorov-Smirnov test, n = 11000).</p