Regulatory patterns in molecular interaction networks
Abstract
Understanding design principles of molecular interaction networks is an important goal of molecular systems biology. Some insights have been gained into features of their network topology through the discovery of graph theoretic patterns that constrain network dynamics. This paper contributes to the identification of patterns in the mechanisms that govern network dynamics. The control of nodes in gene regulatory, signaling, and metabolic networks is governed by a variety of biochemical mechanisms, with inputs from other network nodes that act additively or synergistically. This paper focuses on a certain type of logical rule that appears frequently as a regulatory pattern. Within the context of the multistate discrete model paradigm, a rule type is introduced that reduces to the concept of nested canalyzing function in the Boolean network case. It is shown that networks that employ this type of multivalued logic exhibit more robust dynamics than random networks, with few attractors and short limit cycles. It is also shown that the majority of regulatory functions in many published models of gene regulatory and signaling networks are nested canalyzing.
Key words: gene regulation; signaling; mathematical model; nested canalyzing function; robustness.
1 Introduction
Elucidating the largescale graph structure of complex molecular interaction networks, from transcriptional networks [2] to proteinprotein interaction networks [31] and metabolic network [7] is an important step toward an understanding of design principles of the cellular architecture. For instance, it has been shown that certain graph theoretic network patterns are much more prevalent in such networks than could be expected. See, e.g., [20]. The next step is to understand cells as complex nonlinear dynamical systems. There now exist many dynamic models of gene regulatory, signaling, and metabolic pathways that provide snapshots of cellular dynamics using a range of modeling platforms. Many of these models represent the interactions of different molecular species as logical rules of some type that describe the combinatorics of how the species combine to regulate others; see, e.g., [21, 18]. The logical rules of Boolean network models are an example of such a description, in which network states are reduced to binary states, with a species either present or absent. It was shown in [15, 16] that a special type of Boolean logical rule which appears frequently in published Boolean network models [11] exhibits robustness properties characteristic of molecular networks. These rules, socalled nested canalyzing functions, capture the spirit of Waddington’s concept of canalyzation in gene regulation [30]. Several other classes of Boolean functions have also been investigated in the search for biologically meaningful rules to describe molecular interactions, including random functions [14], hierarchical canalyzing function [26, 22], chain functions [8], and unate functions [10].
In many cases the regulatory relationships are too complicated to be captured with Boolean logic, and more general models have been developed to represent these. Common other discrete model types, in addition to Boolean networks, are socalled logical models [28], Petri nets [25], and agentbased models [23]. In [29] and [12] it was shown that many of these models can be translated into the rich and general mathematical framework of polynomial dynamical systems over a finite field . (Software to carry out this translation is available at http://dvd.vbi.vt.edu/cgibin/git/adam.pl. Since the algebraic structure of is not relevant for our purposes, we will consider a slightly more general setup. Let be variables, which can take values in finite sets , respectively. Let be the Cartesian product. A dynamical system in the variables is a function
where each coordinate function is a function on a subset of , and takes on values in . Dynamics is generated by iteration of . As an example, if , then each is a Boolean rule and is a Boolean network.
Here, we use this very general framework to give a definition of the notion of nested canalyzing rule, which then applies to all different model types simultaneously. We show through extensive simulations that dynamical systems constructed from such rules as coordinate functions have important dynamic properties characteristic of molecular networks, namely very short limit cycles and very few attractors, compared with the set of all possible functions. Furthermore, we show that many published models use logical interaction rules whose polynomial form is nested canalyzing, thereby providing evidence that general nested canalyzing rules represent a frequently occurring pattern in molecular network regulation.
2 Nested Canalyzing Rules
Here we present the general definition of a nested canalyzing rule in variables with state space .
Definition. Assume that each is totally ordered, that is, its elements can be arranged in linear increasing order. In the Boolean case this could be . Let be subsets that satisfy the property that each is a proper, nonempty subinterval of , that is, every element of that lies between two elements of in the chosen order is also in . Furthermore, we assume that the complement of each is also a subinterval, that is, each can be described by a threshold , with all elements of either larger or smaller than .

The function is a nested canalyzing rule in the variable order with canalyzing input sets and canalyzing output values with if it can be represented in the form:

The function is a nested canalyzing function if it is a nested canalyzing function in some variable order for some permutation on .
It is straightforward to verify that, if for all , then we recover the definition in [15] of a Boolean nested canalyzing rule. As mentioned above, several important classes of multistate discrete models can be represented in the form of a dynamical system , so that the concept of a nested canalyzing rule defined in this way has broad applicability.
3 The dynamics of nested canalyzing networks
Aside from incorporating the biological concept of canalyzation, networks whose nodes are controlled by combinatorial logic expressed by nested canalyzing rules have dynamic properties resembling those of biological networks. In particular, they are robust, due to the fact that they have a small number of attractors, which are therefore large. That is, perturbations are more likely to remain in the same attractor. In addition, limit cycles tend to be very short, compared to random networks, which implies that these networks have very regular behavior. We have performed extensive simulation experiments for this purpose, whose results we report here.
3.1 Computer experiments
We have generated random network topologies with indegree distribution ranging between 2 and 5, i.e., each node depends on at least two inputs and at most on five inputs. This assumption is not unrealistic and is based on the observation that gene regulation networks are sparse [17]. For each network topology we have generated two discrete dynamical systems: one where the update rules are all nested canalyzing and the other where the update rules are randomly chosen (and nondegenerate, in the sense that all inputs indicated in the network topology are realized).
3.1.1 Attractor distribution for nested canalyzing networks versus random networks
We present our results concerning attractor distributions in Figure 1. For each histogram, on the axis we represent the number of attractors, and on the axis we represent the percentage of networks that show the given number of attractors specified on the axis. The parameters , and correspond to the number of nodes, the range for the indegree distribution, and the number of states for each node respectively. For the experiments performed here we have generated the indegree distribution from a uniform distribution, independently for each node and network realization. Figure 1 shows very clearly that the number of attractors in nested canalyzing networks is dramatically smaller than for general networks. Thus, attractor sizes are larger on average than in general networks, leading to more robust behavior under perturbations.
3.1.2 Cycle length distribution for nested canalyzing networks versus random networks
We present the results concerning cycle lengths in Figures 2  5. For each figure, the upper subfigures show the mean number of attractors of length specified on the axis (solid lines) and their standard deviations (dashed lines), the axis of these subfigures are in a logarithmic scale. The bottom left subfigures shows the mean number of attractors of lengths specified on the axis in a loglog plot, and finally the bottom right subfigures shows the percentage of networks that returned a particular cycle of length specified on the axis. The parameters , and correspond to the number of nodes, the range for the indegree distribution, and the number of states for each node, respectively. For all of the experiments, we have generated the indegree distribution from a uniform distribution, independently for each node and network realization. Figures 2  5 show clearly that networks with nested canalyzing rules exhibit significantly smaller cycle lengths, leading to more regular behavior.
4 Nested canalyzing rules are biologically meaningful
We hypothesize that nested canalyzing rules are biologically meaningful. To test this hypothesis we have explored a range of published multistate models as to their frequency of appearance. Table 1 shows that they are indeed very prevalent, providing evidence that the nested canalyzation is indeed a common pattern for the regulatory logic in molecular interaction networks. To illustrate this phenomenon we discuss specific examples. For a complete list of models we have studied see the supporting materials.
4.1 Examples
4.1.1 Lambda Phage Regulation
Thieffry and Thomas [27] built a multistate logical model for the core lambda phage regulatory network. This model encompasses the roles of the regulatory genes CI, CRO, CII, and N. See Figure 6.
The state space for this model is specified by , that is, the first variable has three levels , the second variable has four levels , and the third and fourth variables are still Boolean.
The update rule for , , has inputs and which is nested canalyzing in the variable order , , with canalyzing input sets , and and canalized output values , i.e., (see the supporting materials for complete truth tables)
The update rule for , , is nested canalyzing in the variable order , , with canalyzing input sets , and , and canalized output values , i.e.,
The update rule for , , is nested canalyzing in the variable order , , , with canalyzing input sets , , and , and canalized output values , i.e.,
Finally, the update rule for , , is nested canalyzing in the variable order , , with canalyzing input sets , and , and canalized output values , i.e.,
4.2 Regulation in the p53Mdm2 network
The following model comes from AbouJaude W., Ouattara A., Kauffman M.(2009) [1]. The model represents the interactions of the tumor supressor protein p53 and its negative regulator Mdm2. Here, P, Mn, Mc, and Dam stand for protein p53, nuclear Mdm2, cytoplasmic Mdm2, and DNA damage, respectively.
The state space for this model is specified by , that is, except for the first variable that has three levels , all the other variables are still Boolean.
As shown in Figure 7, acts negatively on . The update rule of , , is nested canalyzing with canalyzing input set and canalized output values , i.e., we can represent as,
Here is the basal value and is the parameter value under the influence of .
Similarly, for , its update rule, , is nested canalyzing with canalyzing input set and canalized output values , i.e., we can represent as,
Here is the basal value and is the parameter value under the influence of .
For , a set of possible parameters for its truth table is given in [1]. We have checked all these cases and found that for each case we either get a nested canalyzing function or a constant function. For example, for the second column of Figure 3 (a) in [1] we get that the update rule for , , is nested canalyzing in the variable order , with canalyzing input sets , and canalized output values , i.e.,
When DNA damage is introduced, it has a negative effect on . From the set of all possible parameters for its truth table given in inequalities (3)(5) at [1], we check that we can always find a nested canalyzing function for its truth table. For example, for the third column of Figure 3 (a) in [1] we get that the update rule for Mn (under DNA damage) is nested canalyzing in the variable order , , , with canalyzing input sets , , and canalized output values , i.e., we can represent as,
Finally, the update rule for DNA damage, , is nested canalyzing in the variable order , , with canalyzing input sets , and canalized output values , i.e.,
5 Discussion
In this paper we have given a definition of a nested canalyzing rule, inspired by the special case of Boolean networks, and we have shown that it appears as a frequent pattern for the regulatory logic of many molecular interaction networks. We have shown that this regulatory pattern leads to networks that have robust and regular dynamics, as a result of having very small numbers of attractors and very short limit cycles, compared to random networks. This behavior is also characteristic of many molecular interaction networks. An important application of this result is to the construction of discrete models, both via a bottomup or a topdown approach. For both approaches, the possibility of restricting the choice of rules to the family of nested canalyzing rules is a significant reduction in the possible model space that is available.
Another interesting aspect of our results is suggested by [13]. There is was shown that in the Boolean case, the class of nested canalyzing Boolean rules is in fact identical to the class of unate cascade functions. These functions have been studied extensively in computer engineering and have been shown [3] to have the important property that they comprise exactly the class of Boolean functions that lead to binary decision diagrams with shortest average path length. Thus, they make good candidates for the representation of efficient information processing. It would be interesting to study this property for general nested canalyzing rules.
The results in [13] were derived by using a special representation of Boolean functions, namely as polynomial functions over the Boolean number field, with arithmetic given by addition and multiplication ”tables,” with the key rule that . Using a parametrization of the family of all nested canalyzing Boolean polynomials, it was shown in [13] that the family of all nested canalyzing polynomials in a given number of variables is in fact identical to the class of unate cascade functions, which, in turn, is equal to the class of Boolean functions that result in binary decision diagrams of shortest average path length [3]. In a paper in preparation we have shown that one can give a similar parameterization of the variety of general nested canalyzing rules, and we use this parameterization to derive a formula for the number of nested canalyzing rules for a given number of variables. It would be interesting to investigate whether this general class of nested canalyzing rules leads to ary decision diagrams that have similar properties to those of Boolean nested canalyzing rules.
Model  References  n ^{a}  % NCF^{b} 

lysislysogeny decision in the lambda phage  [27]  4  100% 
p53Mdm2 regulation  [1]  4  100% 
Signalling pathways controlling Th cell differentiation  [21]  42  92.8% 
Budding yeast exit module  [6, 24]  9  77.7% 
Dorsalventral boundary formation of the Drosophila wing imaginal disc  [9]  24  75% 
Control of Th1/Th2 cell differentiation  [19, 4]  14  71.4% 
Yeast morphogenetic checkpoint  [6, 5]  8  50% 

Number of nodes. Only nodes with indegree are considered, i.e. nonconstant nodes.

Percentage of nodes regulated by nested canalyzing functions.
Figure Legends
References
 [1] Wassim AbouJaoudé, Djomangan A. Ouattara, and Marcelle Kaufman. From structure to dynamics: Frequency tuning in the p53mdm2 network: I. logical approach. Journal of Theoretical Biology, 258(4):561 – 577, 2009.
 [2] Richard Bonneau, Marc T. Facciotti, David J. Reiss, Amy K. Schmid, Min Pan, Amardeep Kaur, Vesteinn Thorsson, Paul Shannon, Michael H. Johnson, J. Christopher Bare, William Longabaugh, Madhavi Vuthoori, Kenia Whitehead, Aviv Madar, Lena Suzuki, Tetsuya Mori, DongEun Chang, Jocelyne DiRuggiero, Carl H. Johnson, Leroy Hood, and Nitin S. Baliga. A predictive model for transcriptional control of physiology in a free living cell. Cell, 131(7):1354 – 1365, 2007.
 [3] Jon T. Butler, Tsutomu Sasao, and Munehiro Matsuura. Average path length of binary decision diagrams. IEEE Transactions on Computers, 54(9):1041–1053, 2005.
 [4] Claudine Chaouiya, Elisabeth Remy, and Denis Thieffry. Qualitative petri net modelling of genetic networks, 2006.
 [5] Andrea Ciliberto, Bela Novak, and John J. Tyson. Mathematical model of the morphogenesis checkpoint in budding yeast. The Journal of Cell Biology, 163(6):1243–54, Dec 2003.
 [6] Adrien Faure, Aurelien Naldi, Fabrice Lopez, Claudine Chaouiya, Andrea Ciliberto, and Denis Thieffry. Modular logical modelling of the budding yeast cell cycle. Mol. BioSyst., 5:1787–1796, 2009.
 [7] Adam M Feist, Christopher S Henry, Jennifer L Reed, Markus Krummenacker, Andrew R Joyce, Peter D Karp, Linda J Broadbelt, Vassily Hatzimanikatis, and Bernhard O Palsson. A genomescale metabolic reconstruction for escherichia coli k12 mg1655 that accounts for 1260 orfs and thermodynamic information. Mol Syst Biol, 3, 06 2007.
 [8] I. GatViks and R. Shamir. Chain functions and scoring functions in genetic networks. Bioinformatics, 19(suppl 1):i108–i117, 2003.
 [9] Aitor Gonzalez, Claudine Chaouiya, and Denis Thieffry. Dynamical analysis of the regulatory network defining the dorsalventral boundary of the drosophila wing imaginal disc. Genetics, 174(3):1625–1634, 2006.
 [10] John Grefenstette, Sohyoung Kim, and Stuart Kauffman. An analysis of the class of gene regulatory functions implied by a biochemical model. Biosystems, 84(2):81 – 90, 2006. Dynamical Modeling of Biological Regulatory Networks.
 [11] Stephen E. Harris, Bruce K. Sawhill, Andrew Wuensche, and Stuart Kauffman. A model of transcriptional regulatory networks based on biases in the observed regulation rules. Complex., 7:23–40, March 2002.
 [12] Franziska Hinkelmann, David Murrugarra, Abdul Jarrah, and Reinhard Laubenbacher. A mathematical framework for agent based models of complex biological networks. Bulletin of Mathematical Biology, pages 1–20, 2010. 10.1007/s1153801095828.
 [13] Abdul Salam Jarrah, Blessilda Raposa, and Reinhard Laubenbacher. Nested canalyzing, unate cascade, and polynomial functions. Physica D: Nonlinear Phenomena, 233(2):167 – 174, 2007.
 [14] S Kauffman. Metabolic stability and epigenesis in randomly constructed genetic nets. Journal of Theoretical Biology, 22(3):437–467, March 1969.
 [15] Stuart Kauffman, Carsten Peterson, Bjrn Samuelsson, and Carl Troein. Random boolean network models and the yeast transcriptional network. Proceedings of the National Academy of Sciences of the United States of America, 100(25):14796–14799, 2003.
 [16] Stuart Kauffman, Carsten Peterson, Bjrn Samuelsson, and Carl Troein. Genetic networks with canalyzing boolean rules are always stable. Proceedings of the National Academy of Sciences of the United States of America, 101(49):17102–17107, 2004.
 [17] Robert D Leclerc. Survival of the sparsest: robust gene networks are parsimonious. Mol Syst Biol, 4, 08 2008.
 [18] Fangting Li, Tao Long, Ying Lu, Qi Ouyang, and Chao Tang. The yeast cellcycle network is robustly designed. Proceedings of the National Academy of Sciences of the United States of America, 101(14):4781–4786, 2004.
 [19] Luis Mendoza. A network model for the control of the differentiation process in th cells. Biosystems, 84(2):101 – 114, 2006. Dynamical Modeling of Biological Regulatory Networks.
 [20] Ron Milo, Shalev Itzkovitz, Nadav Kashtan, Reuven Levitt, Shai ShenOrr, Inbal Ayzenshtat, Michal Sheffer, and Uri Alon. Superfamilies of evolved and designed networks. Science, 303(5663):1538–1542, 2004.
 [21] Aurélien Naldi, Jorge Carneiro, Claudine Chaouiya, and Denis Thieffry. Diversity and plasticity of th cell types predicted from regulatory network modelling. PLoS Comput Biol, 6(9):e1000912, 09 2010.
 [22] S. Nikolajewa, M. Friedel, and T. Wilhelm. Boolean networks with biologically relevant rules show ordered behavior. Biosystems, 90(1):40 – 47, 2007.
 [23] Mark Pogson, Rod Smallwood, Eva Qwarnstrom, and Mike Holcombe. Formal agentbased modelling of intracellular chemical interactions. Biosystems, 85(1):37 – 45, 2006.
 [24] Ethel Queralt, Chris Lehane, Bela Novak, and Frank Uhlmann. Downregulation of pp2acdc55 phosphatase by separase initiates mitotic exit in budding yeast. Cell, 125(4):719 – 732, 2006.
 [25] L. Jason Steggles, Richard Banks, Oliver Shaw, and Anil Wipat. Qualitatively modelling and analysing genetic regulatory networks: a Petri net approach. Bioinformatics, 23(3):336–343, 2007.
 [26] Z. Szallasi and S. Liang. Modelling the normal and neoplastic cell cycle with ”realistic genetic networks”: their application for understanding carcinogenesis and assessing therapeutic strategies. Pac. Symp. Biocomput., pages 66–76., 1998.
 [27] Denis Thieffry and René Thomas. Dynamical behaviour of biological regulatory networksII. immunity control in bacteriophage lambda. Bulletin of Mathematical Biology, 57:277–297, 1995. 10.1007/BF02460619.
 [28] René Thomas and Richard D’Ari. Biological Feedback. CRC Press, Boca Raton, Florida., first edition, 1989.
 [29] Alan VelizCuba, Abdul Salam Jarrah, and Reinhard Laubenbacher. Polynomial algebra of discrete models in systems biology. Bioinformatics, 26(13):1637–1643, 2010.
 [30] C.H. Waddington. Canalization of development and the inheritance of acquired characters. Nature, 150:563–565, November 1942.
 [31] Haiyuan Yu, Pascal Braun, Muhammed A. Yildirim, Irma Lemmens, Kavitha Venkatesan, Julie Sahalie, Tomoko HirozaneKishikawa, Fana Gebreab, Na Li, Nicolas Simonis, Tong Hao, JeanFranÃ§ois Rual, AmÃ©lie Dricot, Alexei Vazquez, Ryan R. Murray, Christophe Simon, Leah Tardivo, Stanley Tam, Nenad Svrzikapa, Changyu Fan, AnneSophie de Smet, Adriana Motyl, Michael E. Hudson, Juyong Park, Xiaofeng Xin, Michael E. Cusick, Troy Moore, Charlie Boone, Michael Snyder, Frederick P. Roth, AlbertLaszlo Barabasi, Jan Tavernier, David E. Hill, and Marc Vidal. Highquality binary protein interaction map of the yeast interactome network. Science, 322(5898):104–110, 2008.