Skip to main page content
U.S. flag

An official website of the United States government

Dot gov

The .gov means it’s official.
Federal government websites often end in .gov or .mil. Before sharing sensitive information, make sure you’re on a federal government site.

Https

The site is secure.
The https:// ensures that you are connecting to the official website and that any information you provide is encrypted and transmitted securely.

Access keys NCBI Homepage MyNCBI Homepage Main Content Main Navigation
. 2012;8(9):e1002669.
doi: 10.1371/journal.pcbi.1002669. Epub 2012 Sep 6.

Criticality is an emergent property of genetic networks that exhibit evolvability

Christian Torres-Sosa  1 , Sui Huang, Maximino Aldana
Affiliations

Criticality is an emergent property of genetic networks that exhibit evolvability

Christian Torres-Sosa et al. PLoS Comput Biol. 2012.

Abstract

Accumulating experimental evidence suggests that the gene regulatory networks of living organisms operate in the critical phase, namely, at the transition between ordered and chaotic dynamics. Such critical dynamics of the network permits the coexistence of robustness and flexibility which are necessary to ensure homeostatic stability (of a given phenotype) while allowing for switching between multiple phenotypes (network states) as occurs in development and in response to environmental change. However, the mechanisms through which genetic networks evolve such critical behavior have remained elusive. Here we present an evolutionary model in which criticality naturally emerges from the need to balance between the two essential components of evolvability: phenotype conservation and phenotype innovation under mutations. We simulated the Darwinian evolution of random Boolean networks that mutate gene regulatory interactions and grow by gene duplication. The mutating networks were subjected to selection for networks that both (i) preserve all the already acquired phenotypes (dynamical attractor states) and (ii) generate new ones. Our results show that this interplay between extending the phenotypic landscape (innovation) while conserving the existing phenotypes (conservation) suffices to cause the evolution of all the networks in a population towards criticality. Furthermore, the networks produced by this evolutionary process exhibit structures with hubs (global regulators) similar to the observed topology of real gene regulatory networks. Thus, dynamical criticality and certain elementary topological properties of gene regulatory networks can emerge as a byproduct of the evolvability of the phenotypic landscape.

PubMed Disclaimer

Conflict of interest statement

The authors have declared that no competing interests exist.

Figures

Figure 1
Figure 1. Gene structure.
(A) Schematic representation of the network, showing that each node (circles) “contains” one gene (little bars inside the circles). The arrows represent the regulatory interactions between the genes. (B) Each gene [Image: see text] is composed of a regulatory region and a coding region. The regulatory region contains binding sites which can be added to, or removed from, the regulatory region with the same probability. (C) The binding sites in the regulatory region determine the regulators of [Image: see text] (its input connections). There can be more than one binding site per regulator (as the regulator on the very left), although at the beginning of the evolutionary process the initial networks have only one binding site per regulator. The regulatory region of [Image: see text] determines which other genes it regulates.
Figure 2
Figure 2. Gene expression variability α.
This parameter, which measures whether the genes are frozen in one state, either 0 or 1, or if they more or less switch back and forth between these two states, can be computed in two distinct ways. The first way (horizontal variability [Image: see text]) is to measure [Image: see text] along each attractor state, as shown in A, and then average over all the attractor states and over all the attractors in the attractor landscape. The second way (variability [Image: see text]) shown in B, is to measure the variability for each gene throughout time along the attractor cycle, then average over all the genes in the network and over all the attractors. This is illustrated by the particular example for the same attractor where [Image: see text], whereas [Image: see text].
Figure 3
Figure 3. Evolution towards criticality.
(A) Evolution of the average network sensitivity for four different populations, each initially composed of networks in one of the three dynamical regimes: ordered ([Image: see text], black), critical ([Image: see text], red), and chaotic ([Image: see text], green; and [Image: see text], blue). Under the Darwinian selection given by the ACC and AIC, all the populations quickly become critical ([Image: see text]), regardless of their initial dynamical regime. The inset shows that convergence towards criticality occurs during the first 10000 generation steps. The control curves (in light gray) were obtained by evolving populations without selection, and show that the mutagenic method alone drives the networks into the chaotic regime ([Image: see text]). (B) Distribution of sensitivities at two different generations for the population that started with [Image: see text] chaotic networks. In early generations [Image: see text] is quite broad (dashed line), reflecting a great diversity of networks. However, through evolution all the surviving networks approach criticality and the distribution [Image: see text] narrows down (solid line). The distribution shown here at generation [Image: see text] has [Image: see text].
Figure 4
Figure 4. Genetic variability in the attractors.
(A) Histogram [Image: see text] of the vertical genetic variability [Image: see text] in the attractors of the evolved networks. Note that most of the genetic variability is concentrated around [Image: see text], which indicates that most of the genes in the attractors of the evolved networks switch back and forth between 0 and 1 throughout time. There are almost no frozen genes in these attractors ([Image: see text] is relatively small). (B) Histogram of the vertical variability in the attractors of de novo critical networks that were constructed to be critical by design and did not go through the evolutionary process. Note that in this case most of the genes are frozen in time, as the largest peak at [Image: see text] indicates.
Figure 5
Figure 5. Survival times of the different strains in the population.
(A) Plot of the network labels (strains) that are present in the population at a given generation. Each horizontal line indicates the survival time of a particular strain. The vertical lines indicate the fixation events in which all the networks in the population are relabeled after only one strain was left in the entire population. (B) Distribution [Image: see text] of survival times computed during [Image: see text] generations (black curve). This distribution was computed using logarithmic bins. Only data for [Image: see text] are presented because we checked the existence of strains every 20 generations. The red dashed line is the best fit which corresponds to the power-law [Image: see text]. The inset shows the corresponding cumulative distribution [Image: see text], which better reveals the goodness of the power-law fit. The fact that [Image: see text] has a more or less power-law behavior implies that almost all the strains disappear from the population very quickly, whereas only very few networks are able to survive the Darwinian selection mechanism given by the ACC, the AIC and the α-fitness criterion.
Figure 6
Figure 6. Network structure generated by the evolutionary process.
The top-left network shows the structure of the giant component of the transcription factor interaction network of E. coli according to the RegulonDB . This network has N = 101 nodes and average connectivity K = 2.46. The structure on the top-right corresponds to the typical network that results from our evolutionary algorithm, which in this particular case has N = 100 nodes and average connectivity K = 1.85. Note the existence of global regulators, i.e. nodes with a great number of output connections. The bottom panel presents in a log-log plot the out-degree distribution of these two networks to illustrate their remarkable similarity.
Figure 7
Figure 7. Structural variability within the population.
Although all networks in the final population have the same attractor landscape, they are structurally not completely identical to each other. Here we show three networks randomly chosen from the final population (A, B and C). The image in D is a superposition of all the [Image: see text] networks in that population. The link persistence is defined as [Image: see text], where [Image: see text] is the number of networks in the final population in which the nodes [Image: see text] and [Image: see text] were connected. The links in D have been colored according to their persistence (white for [Image: see text] and red for [Image: see text]). It is apparent that the highly persistent links mostly belong to the global regulators (hubs). This strongly suggests that the global regulators play an important role in determining the phenotypic landscape of the population.
Figure 8
Figure 8. Evolution of the network topology.
(A) The common ancestor network has 10 nodes and one of them (node number 9) is a global regulator that regulates 8 other nodes. (B) Diagram of strain survival times showing the first fixation event at generation g = 6411 (indicated by the red arrow). The common ancestor network is the one that gives rise to the population of the first fixation event. (C) Structure of a randomly chosen network in the final population (generation g = 250000). The initial hub (node number 9) is the one marked with the red circle. Note that at the end this is not a hub anymore, but just another ordinary node of the network. (D) Distribution of the link persistence [Image: see text] for the 10 connections [Image: see text] of the common ancestor. The black and red histograms represent the populations at the first fixation event and at the end of the simulation, respectively. Even after the first fixation event, the links [Image: see text], [Image: see text] and [Image: see text] almost disappear from the population. Furthermore, in the final population none of the links of the initial hub occur at significant frequency. By contrast, link [Image: see text] is present in all the networks of the final population because node 2 became a hub throughout the evolutionary processes. Link [Image: see text] is indicated with the blue bold arrow in (A).
Figure 9
Figure 9. Importance of the gene expression variability as a fitness criterion.
(A) Typical example of the attractors obtained when the evolution of the population is carried out without implementing the [Image: see text]-fitness criterion. In this particular case, instead of the [Image: see text]-fitness we used the S-fitness which assigns a higher replication rate to the networks whose sensitivity is closer to 1. The final attractor landscape consisted of 91 attractors, and only 4 are partially shown here (only the first 45 digits in each attractor state are shown; the remaining 55 digits are all 0's). Note that only the first 10 digits in each attractor state show some activity. There are the 10 genes in the networks of the original population. (B) Plot of the average sensitivity of the population throughout generations, showing that the sensitivity very quickly approaches 1 and remains very close to 1. This is expected since we are explicitly selecting for networks with sensitivities [Image: see text]. (C) Histogram of sensitivities in the final population (generation g = 200000), showing that most networks in the final population have become critical. The inset shows the structure of a typical network in the final population. Note that this network exhibits a more homogeneous random topology with no hubs. This is always the case when [Image: see text]-fitness is not used as selection criterion.
Figure 10
Figure 10. Robustness of the attractor landscape.
Plots of the probability [Image: see text] that a percentage q of the network attractors is conserved after a gene the knockout of one gene for critical networks constructed de novo (A) and the critical networks that result from the evolutionary process (B). The differently colored distributions in B correspond to populations that started in different dynamical regimes (as in Fig. 1). Note the high probability for the de novo networks to lose all their attractors by a gene knockout ([Image: see text]) which does not happen for the evolved networks ([Image: see text] in all cases). Conversely the probability to conserve all attractors is considerably larger for the latter than for the former ([Image: see text] and, [Image: see text] respectively). These data were computed from populations of 1000 networks and 500 attractors per network.
Figure 11
Figure 11. Mutations in the coding region of a gene.
This type of mutation affects how the regulators of a given target gene [Image: see text] (inputs to node [Image: see text]) jointly control its output. (A) Gene [Image: see text] has acquired a new binding site, represented by the leftmost square in green. (B) The new binding site is occupied by one of the already existing regulators of [Image: see text]. (C) The new binding site is occupied by a completely different gene (green circle) that thus becomes a new regulator of [Image: see text]. (D) Another mutation in the regulatory region is the deletion of existing binding sites. Here, the deleted site belongs to a regulator for which there are alternative binding sites left which thus remains a regulator. (E) Conversely, the removal of the binding site can completely break the regulatory interaction between [Image: see text] and the respective regulator. Deletion of binding sites that can leave a given gene with no regulators at all are not allowed.
Figure 12
Figure 12. Mutations in the regulatory region of a gene.
This second type of mutation affects the way in which a given gene [Image: see text] regulates its targets. (A) If [Image: see text] is one of the regulators of a target gene [Image: see text]. Then mutation in the coding region of [Image: see text] may afford gene [Image: see text] the capacity to bind to a new site in the regulatory region of the same target gene [Image: see text] (B), or abrogate the capacity to bind to an existing binding site of that target gene (C). Conversely a change in the coding region of gene [Image: see text] may provide new binding capacity for a site in the regulatory region of a new target, as it is shown in (D) and (E).
See this image and copyright information in PMC

References

    1. Lynch M (2007) The evolution of genetic networks by non-adaptive processes. Nat Rev Genet 8: 803–813. - PubMed
    1. Albert R (2005) Scale-free networks in biology. J Cell Sci 118: 4947–4957. - PubMed
    1. Gama-Castro S, Salgado H, Peralta-Gil M, Santos-Zavaleta A, Muniz-Rascado L, et al. (2010) RegulonDB (version 7.0): transcriptional regulation of Escherichia coli K-12 integrated within genetic sensory response units (Gensor Units). Nucleic Acids Res 39: D98–105. - PMC - PubMed
    1. Serra R, Villani M, Semeria A (2004) Genetic network models and statistical properties of gene expression data in knock-out experiments. J Theor Biol 227: 149–157. - PubMed
    1. Serra R, Villani M, Graudenzi A, Kauffman SA (2007) Why a simple model of genetic regulatory network describes the distribution of avalanches in gene expression data. J Theor Biol 246: 449–460. - PubMed

Publication types

LinkOut - more resources