Fluxbased classification of reactions reveals a functional bowtie organization of complex metabolic networks
Abstract
Unraveling the structure of complex biological networks and relating it to their functional role is an important task in systems biology. Here we attempt to characterize the functional organization of the largescale metabolic networks of three microorganisms. We apply flux balance analysis to study the optimal growth states of these organisms in different environments. By investigating the differential usage of reactions across flux patterns for different environments, we observe a striking bimodal distribution in the activity of reactions. Motivated by this, we propose a simple algorithm to decompose the metabolic network into three subnetworks. It turns out that our reaction classifier which is blind to the biochemical role of pathways leads to three functionally relevant subnetworks that correspond to input, output and intermediate parts of the metabolic network with distinct structural characteristics. Our decomposition method unveils a functional bowtie organization of metabolic networks that is different from the bowtie structure determined by graphtheoretic methods that do not incorporate functionality.
pacs:
82.39.Rt 87.18.Vf 87.18.hI Introduction
Biological systems provide many examples of the intricate relationship between the structure and functionality of complex networks Hartwell et al. (1999); Bornholdt et al. (2003); Barabási and Oltvai (2004); Wagner (2005); Sneppen and Zocchi (2005); Alon (2006); Kaneko (2006). Cellular metabolism is a complex biochemical network of several hundred metabolites that are processed and interconverted by enzymecatalyzed reactions Heinrich and Schuster (1996); Jeong et al. (2000); Wagner and Fell (2001); Ma and Zeng (2003); Csete and Doyle (2004); Palsson (2006). Metabolic networks have a dynamic flexibility that enables organisms to survive under diverse environmental conditions. A key goal of systems biology is to unveil the functional organization of metabolic networks explaining their systemlevel response to different environments. To this end, we have attempted to decompose metabolic networks into functionally relevant subnetworks. Flux balance analysis (FBA) has been widely used to harness the knowledge of largescale metabolic networks and investigate genotypephenotype relationships Price et al. (2004); Feist and Palsson (2008); Oberhardt et al. (2009). FBA has been successful in predicting the growth and deletion phenotypes of organisms Edwards et al. (2001); Ibarra et al. (2002); Segre et al. (2002). Reaction fluxes carry information about the flows on metabolic networks and, as such, describe the functional use of the network by the organism. In this paper, we have used this information to decompose the network into functionally relevant subnetworks.
The paper is organized as follows: In section II we describe the modelling framework in which we study metabolic networks. In section III we discuss the classification of active reactions in metabolic networks into three categories by an algorithm that is blind to their biochemical roles. Section IV shows that the three categories are functionally relevant for the organism. In section V we compare the bowtie architecture obtained by our functional classification of reactions with that obtained by graphtheoretic methods that do not employ functional information. In the last section we conclude with a summary.
Property  E. coli  S. cerevisiae  S. aureus 

Number of metabolites  761  1061  648 
Number of reactions in the model  931  1149  641 
Number of onesided reactions in the equivalent network  1167  1576  863 
Number of external metabolites  143  116  84 
Number of organic external metabolites (carbon sources)  131  107  68 
Number of biomass metabolites  49  42  56 
Number of feasible minimal environments  89  43  27 
Number of active reactions  585  482  418 
Number of reactions in category I  185  89  84 
Number of reactions in category IIa  147  117  194 
Number of reactions in category IIb  42  46  28 
Number of reactions in category III  211  230  112 
Ii The modelling framework
ii.1 Flux balance analysis (FBA)
Flux balance analysis (FBA) is a computational approach widely used to analyze the capabilities of genomescale metabolic networks Price et al. (2004); Feist and Palsson (2008); Oberhardt et al. (2009). The stoichiometric matrix encapsulates the stoichiometric coefficients of different metabolites involved in various reactions of the metabolic network. FBA primarily uses structural information of the metabolic network contained in this matrix to predict the possible steady state flux distribution of all reactions and the maximum growth rate of an organism. In any metabolic steady state, the metabolites achieve a dynamic mass balance wherein the vector of fluxes through the reactions satisfies the following equation representing the stoichiometric and mass balance constraints:
(1) 
Equation 1 is an underdetermined linear system of equations relating various reaction fluxes in genomescale metabolic networks leading to a large solution space of allowable fluxes. The space of allowable solutions can be reduced by incorporating thermodynamic and enzyme capacity constraints. To obtain a particular solution, linear programming is used to find a set of flux values  a particular flux vector  that maximizes a biologically relevant linear objective function . The linear programming formulation of the FBA problem can be written as:
(2) 
where vectors and contain the lower and upper bounds of different fluxes in and the vector corresponds to the coefficients of the objective function . In FBA, the objective function is usually taken to be the growth rate of the organism. The environment, or medium, is defined in this approach by the components of and corresponding to the transport reactions, which determine, in particular, the set of metabolites whose uptake is allowed.
ii.2 Largescale metabolic networks
In this work, we have analyzed the largescale metabolic networks of three microorganisms: Escherichia coli (version iJR904 Reed et al. (2003)), Saccharomyces cerevisiae (version iND750 Duarte et al. (2004)) and Staphylococcus aureus (version iSB619 Becker and Palsson (2005)). Table 1 gives the number of metabolites and reactions in the metabolic networks of these three organisms. The metabolic networks contain internal and transport reactions. Internal reactions occur within the cell boundary. Transport reactions represent processes involving import or export of metabolites across the cell boundary. Each model also contains a pseudo biomass reaction that simulates the drain of various biomass precursor metabolites for growth in the specific organism. Starting from the published metabolic network, we obtain an equivalent reaction network as follows: Every reversible reaction in the network is converted into two onesided (irreversible) reactions so that all reaction fluxes in the equivalent system are positive. A few reactions appear in duplicate in these networks, and only a single copy of each reaction is kept in the equivalent network. The equivalent metabolic network is a reaction set consisting of unique onesided reactions where is 1167, 1576 and 863 for E. coli, S. cerevisiae and S. aureus, respectively (cf. Table 1).
ii.3 Feasible minimal environments and associated flux vectors
In this work, we have considered minimal environments each characterized by the presence of a limited amount of one organic carbon source (external nutrient metabolite) along with unlimited amounts of the following inorganic metabolites: ammonia, iron, potassium, protons, pyrophosphate, sodium, sulfate, water and oxygen. The number of environments we consider for each organism thus coincides with the number of organic external metabolites (carbon sources) in its metabolic network (cf. Table 1). We used FBA to determine the set of minimal environments supporting growth in the metabolic networks of E. coli, S. cerevisiae and S. aureus. A minimal environment is termed as feasible if the growth rate predicted by FBA is nonzero for that carbon source. The number of feasible minimal environments in E. coli, S. cerevisiae and S. aureus was obtained to be 89, 43 and 27, respectively (cf. Table 1) Samal et al. (2006). For each organism, and for each feasible minimal environment for that organism, we obtained an dimensional optimal flux vector maximizing growth in the metabolic network of the organism, and whose component gives the flux of reaction . For every organism, this led to a set of flux vectors corresponding to feasible minimal environments which were stored in the form of a matrix () of dimensions where the rows (=1,2,,) correspond to different reactions in network and columns (=1,2,,) to different feasible minimal environments. is defined as the flux of reaction in the optimal flux vector obtained for environment .
ii.4 Active reactions
A given reaction is termed as active in an environment if . The activity of a reaction denotes the number of minimal environments in which the reaction is active. The activity for a reaction ranges from 0 to with equal to 89, 43 and 27 for E. coli, S. cerevisiae and S. aureus, respectively. A reaction is termed as active in an organism if 1 (i.e., if it is active in at least one feasible minimal environment for that organism). The number of active reactions in E. coli, S. cerevisiae and S. aureus was obtained to be 585, 482 and 418, respectively (cf. Table 1). This paper primarily focuses on decomposing this set of active reactions into functionally relevant subnetworks.
Iii Classification of active reactions
We ask the question: How does the activity of a reaction vary across different environments? To address this question, we determine the frequency distribution of the activity of reactions in an organism. Fig. 1 shows the histogram of the activity of reactions in the E. coli metabolic network. The distribution is bimodal. Most reactions in E. coli are either onceactive (=1) or always active (=89); the number of reactions for any given intermediate activity in the range 189 is small. Thus, the largest number of active reactions in the metabolic network are used in either one environment or in all environments. The histograms of activity of reactions in S. cerevisiae and S. aureus also have a pattern similar to that in E. coli (cf. Fig. 1). The frequency distribution of activity of reactions in the three organisms suggests a natural classification of active reactions into three categories:

Category I reactions or onceactive reactions (=1)

Category II reactions or always active reactions (=)

Category III reactions with intermediate activity (1)
iii.1 Subclassification based on correlation of reaction fluxes
Clustering of gene expression data using the correlation coefficient has been successful in predicting regulatory modules associated with a biological function across diverse conditions Eisen et al. (1998). We used the correlation coefficient to identify sets of reactions whose fluxes are correlated across different environments. We used the set of flux vectors corresponding to feasible minimal environments contained in matrix to obtain the matrix where is the correlation coefficient between two active reactions and and is given by:
(3)  
If then reactions and are perfectly correlated with each other in the given set of environments. Perfect clusters in metabolic networks are maximal sets of reactions that are perfectly correlated to each other pairwise. Perfect clusters are similar to enzyme subsets Pfeiffer et al. (1999); Stelling et al. (2002), correlated reaction sets Papin et al. (2002); Reed and Palsson (2004) or fully coupled sets Burgard et al. (2004) which have been used to detect modules in metabolic networks.
We use Eq. 3 to identify perfect clusters in metabolic networks of E. coli, S. cerevisiae and S. aureus. In particular, a large perfect cluster of 147 reactions was found in E. coli that is a subset of category II reactions. We refer to this subset of perfectly correlated reactions within category II as category IIa reactions. The remaining 42 category II reactions that are always active but not perfectly clustered with category IIa reactions are part of category IIb. Similarly, large prefect clusters of sizes 117 and 194 were found in category II reactions of S. cerevisiae and S. aureus, respectively. In Fig. 1, category IIa and IIb reactions are shown in pink and blue colours, respectively. We have shown elsewhere that perfect clusters are metabolic modules that can be explained by studying the connectivity of their constituent metabolites Samal et al. (2006).
Note that we have used a single optimal flux vector obtained using FBA for each of the feasible minimal environments to determine the activity of a reaction and the set of active reactions in the metabolic network of an organism. However, it is well known that there exist multiple flux vectors or alternate optimal solutions in most largescale metabolic networks that maximize growth in a given environment Lee et al. (2000); Mahadevan et al. (2003); Reed and Palsson (2004); Samal (2008). In principle, due to the presence of alternate optima, the set of active reactions can change depending on the choice of the flux vectors. In Appendix A, we show the robustness of our reaction categories to the presence of alternate optima.
iii.2 Scatter plot of standard deviation versus mean flux of reactions across environments discriminates the three categories
For each active reaction, following Almaas et al Almaas et al. (2004), we have calculated the mean flux and the standard deviation around this mean by averaging the flux of the reaction over feasible minimal environments. Fig. 2 shows the scatter plot of versus for active reactions in E. coli. From Fig. 2, we can distinguish between categories I, II and III, respectively, as they show up quite distinctly (upper line, category I; lower line, category IIa; with category IIb and category III largely in between the two lines). The upper line in Fig. 2 is the expected curve for category I reactions and the lower line is the curve , where is obtained via best fit of data for category IIa reactions. Appendix B gives the derivation of the relation between and for category I and IIa reactions. Thus, we find that the three categories of reactions are distinct from each other by virtue of their statistical properties.
Iv Functional relevance of the three categories of reactions
Until now our classification of active reactions into the three categories was solely motivated by the activity of reactions in E. coli, S. cerevisiae and S. aureus with two very prominent peaks for onceactive and always active reactions (cf. Fig. 1). However, we now show that our three categories I, II, and III obtained using a computational algorithm blind to the biochemical nature of pathways corresponds to the input, output and intermediate subnetworks, respectively. Thus, each category of reactions is a subnetwork with a distinct functional role in metabolism.
iv.1 Category I: Fanin of input pathways
Fig. 3 shows the subnetwork of all 185 category I reactions in E. coli. The figure shows a number of essentially linear paths of one to about five reactions starting from an external nutrient metabolite, often converging to some other metabolite. These are the input pathways of those metabolites, typically starting from their transport reaction that brings them into the cell, and subsequent catabolic reactions that break them down into a smaller set of metabolites. Input pathways of 86 out of the 89 external nutrient metabolites (carbon sources) characterizing different feasible minimal environments are contained in category I, thereby implying that category I essentially covers all the input pathways of metabolism. Similarly, we find that category I reactions in S. cerevisiae and S. aureus contain input pathways for most external nutrient metabolites characterizing different feasible minimal environments. Thus, category I essentially corresponds to input part of the metabolic network.
Fig. 4 shows a portion of category I reactions belonging to sugar input pathways in E. coli where several external sugar metabolites converge downstream into a few intermediate metabolites. Thus, the input pathways in category I exhibit the fanin property whereby diverse external nutrient metabolites are first catabolized into a smaller set of intermediate metabolites before being drawn into the interior of the metabolic network. Usually the external nutrients whose input pathways converge to a common metabolite belong to the same biochemical class (cf. Figures 3 and 4). Fig. 3 contains a number of disconnected subgraphs each describing the input pathways of one or more biochemically similar metabolites; these disconnected paths get connected to the larger metabolic network via further downstream reactions that belong to other categories and are not shown in Fig. 3.
iv.2 Category II: Output biosynthetic pathways
A key biological function of the metabolic network is to convert nutrient metabolites in the environment into biomass metabolites required for growth and maintenance of the cell. The biomass metabolites, which include all the amino acids, nucleotides, lipids and certain cofactors, may be considered to be the output of the metabolic network. Category II reactions are alwaysactive and have a nonzero flux for any feasible minimal environment. We found that the category II subnetwork has biosynthetic pathways for 30 out of the 49 biomass metabolites in E. coli. These pathways are typically the sole production pathways of those biomass metabolites in E. coli Samal et al. (2006). Thus, this subnetwork is at the output end of the metabolism.
Of the 189 category II reactions in E. coli, 147 reactions belong to category IIa, whose fluxes are perfectly correlated across the different minimal environments. Fig. 5 shows the graph of the category IIa subnetwork in E. coli, which is the single largest perfect cluster of reactions. The remaining 42 reactions in category II constitute the category IIb, which are always active but not perfectly correlated with category IIa reactions and with each other. Thus, the fluxes of category IIb reactions vary in a more complicated manner across minimal environments. Categories IIa and IIb exist with similar properties in the metabolic networks of the other two organisms (cf. Table 1). In our previous work, we have shown that most of the category II reactions are essential for growth irrespective of the environment Samal et al. (2006). The set of category II reactions is a superset of reactions in the activity core found earlier by Almaas et al Almaas et al. (2005) which are reactions always used across minimal as well as rich environments.
iv.3 Category III: Intermediate pathways between input and output
Fig. 6 shows the subnetwork of category III reactions in E. coli, which are neither onceactive nor always active; the activity of these reactions depends on the availability of nutrients in a more complicated manner. Category III reactions may be considered to constitute the intermediate part of the network. By comparing the structures of the three categories, it is evident that category III has a highly reticulate and complex architecture compared to category I and II. There is a functional reason for the observed complexity in category III subnetwork. The biomass metabolites collectively contain several different types of chemical structures (moieties), and the E. coli metabolic network is capable of producing these biomass metabolites from different minimal environments, each containing a different (and single) carbon source. A typical external carbon source has one or a few moieties with different nutrients containing different subsets of moieties. Category I reactions transport the carbon sources into the cell and break it down into a small set of moieties. The function of category III reactions is to start with a small set of moieties and produce all the moieties required for biomass production. This requires a complex set of internal transformations and the exact set of transformations required depends on the nature of the input moieties. Thus, the activity of category III transforming reactions depends upon the biochemical nature of available nutrients in different minimal environments. We find that category III contains most of the reactions in central metabolism such as the citric acid cycle. A similar architecture of category III subnetwork was found in the metabolic networks of the other two organisms as well. Some of the biomass metabolites are produced in category III itself. For the other biomass metabolites category III produces precursors which are then taken up in the biosynthetic pathways of category II to produce the biomass metabolites.
V Comparison of functional bowtie decomposition with graphtheoretic bowtie decomposition
Ma and Zeng Ma and Zeng (2003); Ma et al. (2004) have used graphtheoretic measures to reveal a bowtie architecture of metabolic networks similar to that seen in World Wide Web (WWW) Broder et al. (2000), wherein the network can be decomposed into an incomponent, outcomponent and a giant strong component. Given a directed graph, a strong component is a maximal subgraph such that for any pair of nodes and in the set there exists a directed path from to and from to within the subgraph. In general, a directed graph can have many strong components, and the strong component with the largest number of nodes is designated as the giant strong component (GSC). The associated incomponent consists of nodes which have access to GSC nodes via some directed path, but cannot be reached from any GSC node via a directed path. The outcomponent consists of nodes which can be reached from the GSC nodes via some directed path, but lack access to any GSC node via a directed path.
In this work, we have decomposed the metabolic network into three categories using a simple algorithm based on activity patterns of reactions across different minimal environments. Our categorization reveals a functional bowtie architecture wherein the input pathways (category I reactions) fan into intermediate metabolism (category III reactions) which forms the knot of a bowtie and from where the output pathways (category II reactions) for various biomass components fan out.
In our functional bowtie decomposition, the three categories I, II and III of reactions discussed above broadly correspond to the incomponent, outcomponent and GSC, respectively, of the graphtheoretic bowtie decomposition by Ma and Zeng Ma and Zeng (2003); Ma et al. (2004). However, the corresponding sets of reactions in the two decompositions differ in detail. For example, we find that the end products of several (and often long) chains of reactions in the category II subnetwork are recycled resulting in feedback loops. Such feedback loops in category II subnetwork presumably minimize wastage and could be instrumental in producing the biomass metabolites in the desired ratios. An example of such a feedback loop in category II subnetwork is the one involving metabolite 5mdr1p (which can be seen in the electronic version of Fig. 5 upon zooming). The biosynthetic pathways involved in such feedback loops appropriately belong to the output part of metabolism because they connect the precursor metabolites to the outputs. However, the graphtheoretic bowtie decomposition would classify such category II reactions in feedback loops into the GSC. Thus, our functional bowtie decomposition based on fluxes of reactions across different environments gives better a insight and is biochemically more realistic. The picture of the metabolic network our decomposition reveals is similar in spirit to the one envisioned by Csete and Doyle Csete and Doyle (2004). Further, it is important to note that our fluxbased categorization does not involve the a priori exclusion of high degree currency metabolites as was needed in the graphtheoretic bowtie decomposition of the metabolic network Ma and Zeng (2003); Ma et al. (2004).
Vi Conclusions
In this paper, we have performed flux balance analysis (FBA) for the metabolic networks of three microorganisms: E. coli, S. cerevisiae and S. aureus to obtain fluxes of reactions in the network under diverse environmental conditions. We have followed a purely algorithmic approach leveraging on the predicted fluxes of reactions across different minimal environments to decompose the metabolic network into functionally relevant subnetworks. We find that the activity of a reaction given by the number of minimal environments for which it has a nonzero flux is an important indicator of the functional role of a reaction. We have classified the reactions into three functional categories based on their activity. Category I contains onceactive reactions which are used in only one minimal environment. Most reactions belonging to the category I subnetwork are uptake pathways for external nutrients in feasible minimal environments, and the primary function of these reactions is to catabolize external nutrients into simpler metabolites which can be further processed by intermediary metabolism. Category II contains always active reactions which are used in all minimal environments. The category II subnetwork is critical for the survival of the organism and accounts for the majority of the biosynthetic pathways for the production of the biomass metabolites at the output end of metabolic network. Category III contains reactions which are used in an intermediate number of minimal environments, and is responsible for generating the ‘precursor’ molecules that are eventually converted into biomass metabolites by Category II reactions. We find that while category I and II subnetworks are dominated by simple linear pathways, the structure of the category III subnetwork is highly reticulate. In summary, our decomposition method for largescale metabolic networks based on activity of reactions captures the proposed functional bowtie organization by Csete and Doyle: the input pathways (category I reactions) for nutrients in the environment fan into intermediate metabolism (category III reactions) which forms the knot of bowtie from where the output biosynthetic pathways (category II reactions) for biomass components fan out. Our results are valid for metabolic networks of three phylogenetically different organisms (two distinct prokaryotes and a eukaryote), which suggests that the observed functional bowtie organization could be quite common in living systems.
Appendix A Robustness of categorization of reactions to alternate optimal solutions
In this work, flux balance analysis (FBA) was used to obtain a particular flux vector or optimal solution that maximizes the objective function taken as the growth rate in a given minimal environment. However, for largescale metabolic networks, there exist multiple flux vectors or alternate optimal solutions that maximize growth in a given minimal environment, i.e., there are many flux vectors with exactly the same value of the objective function but use different alternate pathways in the network Lee et al. (2000); Mahadevan et al. (2003); Reed and Palsson (2004); Samal (2008). FBA finds one of many possible alternate optima for a given minimal environment that maximizes growth. In the main text, we have used a single optimal flux vector for each of the feasible minimal environments to determine the activity of a reaction and the set of active reactions in the metabolic network of an organism. Since, in principle, the activity of a reaction can change depending on the particular flux vector considered, we study the robustness of our categorization of reactions to the presence of alternate optima.
Flux variability analysis (FVA) Mahadevan et al. (2003) can be used to determine the set of reactions whose fluxes vary across alternate optima for a given minimal environment. Specifically, FVA determines the maximum and minimum flux value that each reaction can take across alternate optima for a given minimal environment. FVA involves the following steps:

Determine using FBA the maximum value of the objective function or growth rate in a given minimal environment .

Fix the flux of the biomass reaction equal to .

Change the objective function to be the flux of a reaction .

Using linear programming determine the maximum flux value of reaction in the minimal environment , constraining the biomass reaction to have a flux equal to .

Using linear programming determine the minimum flux value of reaction in the minimal environment , constraining the biomass reaction to have a flux equal to .

The range to gives the variability of flux of reaction across different alternate optima.

The above steps c, d, e and f can be repeated for every reaction in the metabolic network to determine the flux variability of each reaction across alternate optima for a given minimal environment .
We have used FVA to determine and for each reaction and for each feasible minimal environment in the E. coli metabolic network. A reaction is designated as blocked if =0 for all feasible minimal environments Schuster and Schuster (1991); Burgard et al. (2004). We found 329 blocked reactions in the E. coli metabolic network. The remaining 838 reactions, for which 0 for at least some environment are designated as potentially active reactions. This set includes the 585 active reactions considered in the main text. We define a reaction as essential for a given minimal environment if 0. 484 reactions were found to be essential for some in the E. coli metabolic network which are a subset of the 585 active reactions considered in the main text. We now classify these 484 reactions into the following three categories:

Essential category I: Reactions which satisfy 0 for exactly one minimal environment. We found 162 reactions in the E. coli metabolic network to be in Essential category I. Of these, 153 reactions belong to category I of the main text.

Essential category II: Reactions which satisfy 0 for all minimal environments. We found 171 reactions in the E. coli metabolic network to be in Essential category II. All of these belong to category II of the main text.

Essential category III: Reactions which satisfy 0 for minimal environments where 1. We found 151 reactions in the E. coli metabolic network to be in Essential category III. Of these, 145 belong to category III of the main text.
Thus we find that the classification discussed in the main text which uses a particular flux vector correctly predicts the essential category I, II or III of 469 out of the 484 essential reactions.
Appendix B Relation between standard deviation and mean flux for category I and category IIa reactions
In Fig. 2, we plot the standard deviation versus mean flux for active reactions in a metabolic network across its feasible minimal environments. Here, we derive the relation between mean flux and standard deviation for reactions in category I and category IIa shown as upper and lower lines, respectively, in Fig. 2.
b.1 Category I reactions
In a given organism any reaction belonging to category I has activity =1, and is active for a single environment (say ). The mean flux of a category I reaction across feasible environments is given by:
(4)  
where is the flux of reaction in the environment . is the flux of reaction in the only feasible minimal environment where the reaction has nonzero value and in all other feasible minimal environments the flux of reaction is 0.
Thus, the standard deviation for a category I reaction is given by:
(5)  
where we have used the result in Eq. 4.
b.2 Category IIa reactions
For two reactions and in category IIa, their flux and in a given environment are perfectly correlated to each other. The fluxes of category IIa reactions are proportional to each other having the same proportionality constant for all minimal environments. For a minimal environment , we can write the flux of category IIa reaction as:
(6) 
where is a constant for the minimal environment and is some number. For any two reactions and in category IIa with fluxes correlated across minimal environments, we have:
(7)  
where and are two different feasible minimal environments for the organism.
The mean flux of reaction is:
(8)  
where is the mean of across the set of feasible minimal environments.
The standard deviation for category IIa reaction is given by:
(9)  
where we have used the result in Eq. 8.
Acknowledgements.
SS and VG acknowledge support from University Grants Commission (UGC), AS from Council for Scientific and Industrial Research (CSIR), and SJ from Department of Biotechnology (DBT), India.References
 L. Hartwell, J. Hopfield, S. Leibler, and A. Murray, Nature 402, C47 (1999).
 S. Bornholdt, H. Schuster, and J. Wiley, Handbook of graphs and networks, Vol. 2 (Wiley Online Library, 2003).
 A. Barabási and Z. Oltvai, Nature Reviews Genetics 5, 101 (2004).
 A. Wagner, Robustness and evolvability in living systems (Princeton University Press Princeton, NJ:, 2005).
 K. Sneppen and G. Zocchi, Physics in molecular biology (Cambridge University Press, 2005).
 U. Alon, An introduction to systems biology: design principles of biological circuits, Vol. 10 (Chapman & Hall/CRC, 2006).
 K. Kaneko, Life: An introduction to complex systems biology, Vol. 171 (Springer Heidelberg, Germany:, 2006).
 R. Heinrich and S. Schuster, The regulation of cellular systems, Vol. 416 (Chapman & Hall New York, 1996).
 H. Jeong, B. Tombor, R. Albert, Z. Oltvai, and A. Barabási, Nature 407, 651 (2000).
 A. Wagner and D. Fell, Proceedings of the Royal Society of London. Series B: Biological Sciences 268, 1803 (2001).
 H. Ma and A. Zeng, Bioinformatics 19, 1423 (2003).
 M. Csete and J. Doyle, Trends in Biotechnology 22, 446 (2004).
 B. Palsson, Systems biology: properties of reconstructed networks (Cambridge University Press, 2006).
 N. Price, J. Reed, and B. Palsson, Nature Reviews Microbiology 2, 886 (2004).
 A. Feist and B. Palsson, Nature biotechnology 26, 659 (2008).
 M. Oberhardt, B. Palsson, and J. Papin, Molecular Systems Biology 5 (2009).
 J. Edwards, R. Ibarra, B. Palsson, et al., Nature Biotechnology 19, 125 (2001).
 R. Ibarra, J. Edwards, and B. Palsson, Nature 420, 186 (2002).
 D. Segre, D. Vitkup, and G. Church, Proceedings of the National Academy of Sciences 99, 15112 (2002).
 J. Reed, T. Vo, C. Schilling, B. Palsson, et al., Genome Biol 4, R54 (2003).
 N. Duarte, M. Herrgård, and B. Palsson, Genome Research 14, 1298 (2004).
 S. Becker and B. Palsson, BMC Microbiology 5, 8 (2005).
 A. Samal, S. Singh, V. Giri, S. Krishna, N. Raghuram, and S. Jain, BMC bioinformatics 7, 118 (2006).
 M. Eisen, P. Spellman, P. Brown, and D. Botstein, Proceedings of the National Academy of Sciences 95, 14863 (1998).
 T. Pfeiffer, F. Montero, S. Schuster, et al., Bioinformatics 15, 251 (1999).
 J. Stelling, S. Klamt, K. Bettenbrock, S. Schuster, E. Gilles, et al., Nature 420, 190 (2002).
 J. Papin, N. Price, and B. Palsson, Genome Research 12, 1889 (2002).
 J. Reed and B. Palsson, Genome Research 14, 1797 (2004).
 A. Burgard, E. Nikolaev, C. Schilling, and C. Maranas, Genome Research 14, 301 (2004).
 S. Lee, C. Phalakornkule, M. Domach, and I. Grossmann, Computers & Chemical Engineering 24, 711 (2000).
 R. Mahadevan, C. Schilling, et al., Metabolic engineering 5, 264 (2003).
 A. Samal, Systems and synthetic biology 2, 83 (2008).
 E. Almaas, B. Kovacs, T. Vicsek, Z. Oltvai, and A. Barabási, Nature 427, 839 (2004).
 J. Ellson, E. Gansner, L. Koutsofios, S. North, and G. Woodhull, in Graph Drawing (Springer, 2002) pp. 594–597.
 E. Almaas, Z. Oltvai, and A. Barabási, PLoS Computational Biology 1, e68 (2005).
 H. Ma, X. Zhao, Y. Yuan, and A. Zeng, Bioinformatics 20, 1870 (2004).
 A. Broder, R. Kumar, F. Maghoul, P. Raghavan, S. Rajagopalan, R. Stata, A. Tomkins, and J. Wiener, Computer networks 33, 309 (2000).
 S. Schuster and R. Schuster, Journal of Mathematical Chemistry 6, 17 (1991).