Study of the Structure and Dynamics of Complex Biological Networks

The metabolic network represents the set of biochemical reactions that are responsible for the uptake of food molecules or nutrients from the external environment and converting them into other molecules that are building blocks required for the growth and maintenance of the cell. The latter include ATP, the energy currency of the cell, other nucleotides, amino acids, lipid molecules and other molecules. These are sometimes called ‘biomass’ metabolites and constitute the output of the metabolic network. The inputs are the food molecules such as sugars as well as other organic molecules varying from organism to organism and inorganic molecules such as water, hydrogen ions, oxygen and sources of nitrogen, phosphorus, sulphur, iron, sodium, potassium, etc. These typically enter the cell through its membrane. The bulk of the metabolic network are the chemical reactions (ranging from several hundred to more than a thousand reactions in different organisms) which transform the input molecules into the output molecules. The various reactions in the metabolic network are catalyzed by enzymes. Enzymes are proteins which are coded by genes. Metabolites are the reactants or products of various reactions.

The metabolic network can be represented as a bipartite graph consisting of two types of nodes: metabolites and reactions. An example of a directed bipartite graph is shown in Fig. 1.1. In the directed bipartite metabolic graph, there is a link from a metabolite pointing to a reaction node if the metabolite is a reactant of the reaction, and a link from a reaction pointing to a metabolite node if the metabolite is a product of the reaction. In the bipartite metabolic graph, there are no direct links between two metabolites or two reactions. In chapter 2 of this thesis, we have studied the metabolic networks inside three organisms as a directed bipartite graph.

The metabolic networks are also sometimes represented as unipartite graphs (which could be directed or undirected) in which there is only one type of node (metabolite nodes or reaction nodes). In the undirected unipartite metabolite graph, for example, the nodes of the network represent metabolites, with edges between two metabolites if they participate in a single reaction. Similarly, in the undirected unipartite reaction graph, the nodes of the network represent reactions, with edges between two reactions if they share the same metabolite. In the directed unipartite reaction graph, there is a link pointing from one reaction node to another if the former produces a metabolite that is consumed by the latter.

The bipartite representation of the metabolic network has more information compared to the above unipartite representations. A more detailed description of the metabolic network would require three types of nodes: metabolites, reactions and enzymes. Such a tripartite network could additionally account for the correspondence between enzymes and reactions and regulatory interactions between metabolites and enzymes inside the cell.

Genes are segments of DNA that code for proteins inside the cell. Transcription is the process by which an enzyme, RNA polymerase, reads the sequence of bases on a gene and constructs an mRNA molecule from that sequence. Translation is the process in which a ribosome, a macromolecular assembly, reads the information contained in the mRNA molecule and synthesizes a protein molecule from the sequence on the mRNA molecule. Thus, each protein molecule is a product of the gene that codes for it. In turn, proteins are responsible for carrying out various functions inside the cell, including catalyzing reactions of the metabolic network as enzymes, building various cellular structures, etc.

The transcription of a gene to an mRNA molecule is also regulated by proteins. The proteins that regulate the expression of genes inside cells are referred to as transcription factors. A transcription factor may activate or inhibit the expression of a gene inside the cell by binding to regions upstream or downstream of the gene on the DNA molecule. This process may in turn facilitate or prevent RNA polymerase and rest of the transcription machinery from binding and initiating the transcription of the gene. Thus, the genes inside cells interact amongst each other via intermediate transcription factors to influence each other’s expression. This network of interacting genes inside the cell is referred to as the transcriptional regulatory network. This network may be represented by a graph in which every gene is represented by a node, and where an arrow from one gene node to another means that the former codes for a transcription factor that regulates the transcription of the latter gene.

A more detailed description of the transcriptional regulatory network would require us to characterize the regulatory links in the network as activating or repressing. Further, an input function needs to be specified for each gene node when more than one transcription factor regulates the expression of that gene in the network. In chapters 3 and 4 of this thesis, we have studied the part of the transcriptional regulatory network in E.coli that controls metabolism.

Inside the cell, the proteins interact with each other to influence each other’s activity. Extracellular signals are mediated to the inside of a cell by protein-protein interactions of signaling molecules. Proteins that interact for a long time with each other structurally can form part of a protein complex. A protein may also act as a carrier for another protein. A protein may interact briefly with another protein leading to a transfer of a phosphate group. These protein-protein interactions are crucial for living systems. We can represent various protein-protein interactions inside the cell by a network with nodes as proteins and links between two nodes representing a physical interaction between two proteins. Such networks are referred to as protein-protein interaction networks. The links between two nodes can be of various different types corresponding to different types of interactions between proteins.

It is important to emphasize that the metabolic, transcriptional regulatory and protein-protein interaction networks mentioned above are not independent of each other inside the cell. The state of the genes in the transcriptional regulatory network determines the activity of the metabolic network. The concentration of metabolites in the metabolic network determines the activity of transcription factors or proteins which regulate the expression of genes in the regulatory network. The protein-protein interactions determine the activity of various proteins inside the cell. Thus, the above mentioned and other biochemical networks together with the interactions at the interface of these networks form a ‘network of networks’ inside the cell that determines the overall behaviour of the organism.

The knowledge of enzymes along with their functional assignments have led to a reconstruction of nearly complete lists of organism specific metabolic reactions [22, 23, 24]. Jeong et al [25] studied the structure of the metabolic networks inside 43 different organisms. They represented the metabolic network as a bipartite graph with two types of nodes: metabolites and reactions. The degree of a node is defined as the number of links attached to that node in the graph [26, 27]. Since the metabolic network is a directed graph, each metabolite in the network has an in-degree and an out-degree. The in-degree of a node denotes the number of links that the node has to other nodes in the directed graph. The out-degree of a node denotes the number of links that start from the node to other nodes in the directed graph. The degree distribution of a graph, $P(k)$, gives the probability that a randomly selected node has exactly $k$ links in the graph [26, 27]. Jeong et al found both the in-degree and out-degree distribution of metabolites to approximate a power law form $P(k)\sim k^{-\gamma}$ for the metabolic networks inside 43 organisms. $\gamma$ is the degree exponent which they found to be universal and close to 2.2 for all the organisms studied, both for the in-degree and out-degree distribution.

Independently, Wagner and Fell [28] studied the large scale metabolic network of E. coli. They represented the E. coli metabolic network as two different unipartite graphs: metabolite graph and reaction graph [28]. Wagner and Fell also found the connectivity of the metabolites in the network to follow a power law distribution. The two independent studies by Jeong et al [25] and Wagner and Fell [28] show that most metabolites participate in a few reactions while there are a few metabolites which participate in many reactions in metabolite networks. The ubiquitous metabolites such as ATP that participate in several reactions and have a high degree are also referred to as hubs of the network.

The distance or shortest path between nodes $i$ and $j$ in a graph is defined as the minimum number of links that have to be traversed to reach from node $i$ to $j$. The average path length of a graph is defined as the average over the shortest paths between all pairs of nodes in the network [26, 27]. The diameter of a graph is defined as the supremum of the shortest paths between all pairs of nodes in the network [26, 27]. The studies by Jeong et al [25] and Wagner and Fell [28] found the metabolic networks inside organisms to have the small-world [1] property, i.e., any two nodes in the system can be connected by relatively short paths along existing links. Jeong et al showed that the sequential removal of the high degree nodes or hubs from metabolic networks results in a sharp rise of network diameter as the network disintegrates into small isolated clusters. On the other hand, when they removed a set of randomly chosen metabolite nodes from the network, the average path length between the remaining nodes was not affected [25]. This observation led them to conclude that the hubs of the metabolic network are crucial for maintaining functionality of metabolic networks.

Ma and Zeng [29] have further explored the global connectivity structure of metabolic networks by classifying nodes in the metabolite graph into four subsets based on their mutual connectivity properties and location in the network: a giant strong component, in-component, out-component and an isolated subset. Two nodes $i$ and $j$ are said to strongly connected in a directed graph, if there is a path from $i$ to $j$ and from $j$ to $i$ in the graph. A strong component is a subset of nodes in the directed graph such that for any pair of nodes $i$ and $j$ in the subset there is a path from $i$ to $j$ and $j$ to $i$ in the directed graph. The largest strong component of a directed network is referred to as the giant strong component. The set of nodes that are not in the giant strong component but from which the nodes in the giant component can be reached forms the in-component. The set of nodes that are not in the giant strong component but that can be reached from the nodes in the giant component forms the out-component. The set of nodes that have no path to the nodes in the giant strong component forms the isolated subset. The decomposition of the metabolite nodes into the above mentioned four connected components revealed a ‘bow-tie’ macroscopic structure of the metabolic network [29]. The bow-tie structure of the metabolic network was similar to that observed by Broder et al for the World Wide Web [30]. In uncovering the bow-tie structure, Ma and Zeng removed the connections through the ubiquitous currency metabolites in the metabolic network. Csete and Doyle have argued that the bow tie architecture of the metabolic network with a conserved core and plug-and-play modularity around core can contribute toward robustness and evolvability of the system [31].

Ma and Zeng [29] also tried to account for the preferred directionality of reactions in the graph for the metabolic network, and found the average path length between metabolites to be almost double of that observed by Jeong et al. The average path length between nodes in the giant strong component was found to determine the average path length of the whole network. An alternative study by Arita [32] tried to account for the actual structural changes in connecting metabolites in the graph of the E. coli metabolic network, and found again the average path length to be almost double of that observed by Jeong et al. Thus, accounting for the directionality of reactions, activity of reactions and functional transfer of biochemical groups in a more biologically meaningful graph representation of the metabolic network gives an average path length larger than that obtained by Jeong et al.

The above mentioned studies of the structure of the metabolic networks have revealed a large variation in the metabolite connectivity inside these networks. It has been suggested that one of the important consequences of power law degree distribution is the vulnerability of the network to selective attack on hubs while being robust to random deletion of nodes from the network as most nodes are of low degree and their deletion does not change the average path length between remaining nodes in the network [33]. For the protein-protein interaction network of S. cerevisiae, it was shown that the essentiality of a protein is correlated with its degree in the network [34]. This observation has been suggested as evidence for the importance of the hubs in maintaining the overall structure and function of cellular networks. Although the role of high degree metabolites or hubs in maintaining the overall structure of the metabolic networks has been well emphasized in the literature, the role of low degree metabolites has attracted little or no attention.

In chapter 2 of this thesis, we show that certain low degree metabolites introduce fragility for flows in metabolic network. There we have used a computational method to determine essential reactions for growth in the metabolic networks inside three organisms. A reaction is designated as ‘essential’ if its knockout from the metabolic network renders the organism unviable. We show that the low degree metabolites as opposed to high degree metabolites explain essential reactions in metabolic networks [35]. It is the low degree metabolites that are critical from the point of view of functional robustness of the system.

In chapter 2, we also show that certain low degree metabolites lead to clusters of reactions with highly correlated reaction fluxes in the metabolic network. We then show that genes corresponding to reactions of such clusters predict regulatory modules in E. coli. Thus, the modularity observed by us at the metabolic level is also reflected at the genetic level. Our work therefore shows that low degree metabolites play a role in two hitherto unconnected properties of biological networks: on the one hand they cause certain reactions to become essential for the viability of the organism, and on the other hand they contribute to the modularity of biological networks.

The available information regarding the target genes of transcription factors has led to a reconstruction of transcriptional regulatory networks inside model organisms like E. coli and S. cerevisiae [36, 37, 38, 39, 40]. The presently available transcriptional regulatory maps are highly incomplete due to ongoing annotation of the fully sequenced genomes. It has been estimated [41] that the coverage of the known transcriptional regulatory network of E. coli is only about $25\%$ of the actual network inside the organism.

The out-degree distribution for the known transcriptional regulatory network of E. coli and S. cerevisiae was shown to approximate a power law. However, the in-degree distribution for the two networks followed a restricted exponential function [37, 38]. This example shows that not all biological networks are characterized by a power law degree distribution of nodes. The out-degree of a node in the transcriptional regulatory network represents the number of target genes regulated by a transcription factor. The in-degree of a node in the transcriptional regulatory network represents the number of transcription factors regulating a target gene. The exponential in-degree distribution for the transcriptional regulatory network suggests that a very large promoter region required for the combinatorial regulation of a target gene by many transcription factors is highly unlikely inside the cell.

‘Network motifs’ have been defined as patterns of interconnections or subgraphs that are over-represented in a real network compared to randomized versions of the same network with similar local connectivity [37, 42]. Network motifs can be detected by algorithms that compare the patterns found in the real network to those found in suitably randomized networks. The method is analogous to detection of sequence motifs in genomes as recurring sequences that are very rare in random sequences. By studying the E. coli transcriptional regulatory network, Alon and colleagues found that the ‘feed-forward loop’ (FFL) is a motif in the regulatory network [37]. The structure of a feed-forward loop (FFL) motif is defined by a transcription factor X that regulates a second transcription factor Y, such that both X and Y jointly regulate a gene or operon Z. Other motifs found in the E. coli transcriptional regulatory network include the ‘single-input module’ (SIM) and ‘dense overlapping regulons’ (DOR) [37]. The structure of single-input module (SIM) is defined by a set of genes or operons that are controlled by a single transcription factor. The structure of dense overlapping regulons (DOR) is defined by a layer of overlapping interactions between genes or operons and a group of input transcription factors. Later, the motifs found in the E. coli transcriptional regulatory network were also found in the S. cerevisiae transcriptional regulatory network [38, 42].

By studying the dynamics of various motifs found in the regulatory networks, it has been shown that these motifs may perform important information processing tasks. The coherent feed-forward loop motif has been shown to filter out noise or spurious signals in the network [37, 16]. The single-input module motif was shown to help generate temporal programs of gene expression [37, 43]. Recent studies have shown that the appearance of same motifs in regulatory networks of E. coli and S. cerevisiae does not necessarily imply that these motifs are evolutionarily conserved [44, 45]. Evolution seems to have converged on the same patterns of interconnections in different organisms after a lot of tinkering perhaps due to the specific information processing tasks these motifs perform inside a cell [9, 45]. An objective of the exercise of detecting motifs in different biological networks is to create a library of motifs and their possible functions.

Ma et al [46] tried to decompose the E. coli transcriptional regulatory network into various connected components. They found that there were no strong components in the E. coli transcriptional regulatory network. This was consistent with earlier observations by Shen-Orr et al [37] that the network had no cycles of length $\geq$ 2. There were only autoregulatory loops in the presently known transcriptional regulatory network of E. coli [37]. Ma et al [46] found the structure of the E. coli transcriptional regulatory network to be hierarchical. A similar lack of strong components or cycles was also observed for the transcriptional regulatory network of S. cerevisiae [39]. The transcriptional regulatory network of E. coli had a five layered hierarchical architecture with genes at the top layer having no incoming links. Balaszi et al found that the genes at the top layer are regulated by distinct environmental signals in the transcriptional regulatory network of E. coli [47].

In chapters 3 and 4 of this thesis, we have studied the large scale structure and system level dynamics of the transcriptional regulatory network controlling metabolism in E. coli. Our study reinforces the hierarchical, essentially acyclic structure with environmental control of the genes belonging to the top layer of the regulatory network of E. coli, also found by previous studies mentioned above. Further, our dynamical study of regulatory network of E. coli metabolism elucidates the functional consequences of this observed architecture of the network [48]. We show that the regulatory network of E. coli metabolism exhibits two types of robustness. One, the regulatory network of E. coli metabolism exhibits an insensitivity to perturbations of gene configurations for a fixed environment, i.e., the system returns to the same attractor when gene configurations are perturbed for a fixed environment. Two, the regulatory network of E. coli metabolism exhibits a flexible response to changed environments, i.e., the system moves to a new attractor that enables it to maintain its key functionality when it encounters a changed environment in a sustained manner. The hierarchical acyclic architecture of the regulatory network of E. coli metabolism with control variables as external metabolites explains the observed robust dynamics of the system. Further, we observe a highly disconnected and modular architecture at the intermediate level of the hierarchical graph of the regulatory network. We find that the modules at the intermediate level of this hierarchical graph are regulated by different sets of environmental signals, and the modules interact only at the lowest level of the graph contributing to the robust response of the system to changed environments. This modular architecture of the regulatory network may also contribute towards the evolvability of the system. Thus, our study sheds new light on how structural design features of the regulatory network of E. coli contribute towards robustness and modularity of the system.

The list of reactions in a reconstructed metabolic network for an organism includes both reversible and irreversible reactions. Starting from the reconstructed metabolic network of an organism, we prepare a list of metabolic reactions that has each reversible reaction in the original network replaced by two unidirectional reactions (one reaction each for the forward and the backward direction). From this list of unidirectional metabolic reactions, construct a matrix ${\bf A}$ = $(A_{ij})$ of dimensions $m\times n$, where $m$ is the number of internal metabolites and $n$ is the number of reactions in the network. The rows of matrix ${\bf A}$ correspond to metabolites and the columns correspond to reactions in the metabolic network. The matrix element $A_{ij}$ is set equal to -1, if metabolite $i$ is consumed in reaction $j$, +1 if metabolite $i$ is produced in reaction $j$, and 0 if metabolite $i$ does not participate in reaction $j$. The matrix ${\bf A}$ is a compact representation of the bipartite metabolic graph.

A UP(UC) metabolite has in-degree (out-degree) equal to unity in the bipartite metabolic graph. A metabolite $i$ in the network is UP(UC), if the $i^{th}$ row of matrix ${\bf A}$ has exactly one entry that equals -1 (+1). A UP-UC metabolite has in-degree and out-degree equal to unity in the bipartite metabolic graph. A metabolite $i$ in the network is UP-UC, if the $i^{th}$ row of matrix ${\bf A}$ has exactly one entry that equals -1, one entry that equals +1 and has all other entries 0. A reaction $j$ is UP(UC), if the $j^{th}$ column of matrix ${\bf A}$ has at least one entry +1 (-1) such that the row corresponding to the entry +1 (-1) in column $j$ has no other entry equal to +1 (-1).

As mentioned earlier, we do not consider the external metabolites while determining the set of UP(UC) metabolites. This amounts to excluding the rows corresponding to external metabolites from the matrix ${\bf A}$ for the computation of UP(UC) metabolites. Thus, the matrix ${\bf A}$ does not contain any rows corresponding to external metabolites and the rows of matrix ${\bf A}$ correspond only to internal metabolites in the network. The biomass reaction which defines the ratios of various metabolic precursors that are required for the unit biomass production of the organism is also included in the list of reactions while determining the UP(UC) metabolites. Usually, the last column of the bipartite matrix ${\bf A}$ corresponds to the biomass reaction.

The databases were downloaded from the website [22]. The E. coli metabolic network iJR904 accounts for 761 metabolites participating in 931 reactions (686 irreversible and 245 reversible reactions). The S. cerevisiae metabolic network iND750 accounts for 1061 metabolites participating in 1149 reactions (719 irreversible and 430 reversible reactions). The S. aureus metabolic network iSB619 accounts for 645 metabolites participating in 644 reactions (423 irreversible and 221 reversible reactions).. Following the steps outlined in section 2.1.1, the bipartite matrix ${\bf A}$ was constructed for the three metabolic networks. We found the dimensions ($m$,$n$) of matrix ${\bf A}$ to be (618,1177), (945,1580) and (561,866) for the metabolic networks of E. coli, S. cerevisiae and S. aureus, respectively. In obtaining the bipartite matrix ${\bf A}$ for the three metabolic networks, each reversible reaction was converted into two one sided reactions. Further, the biomass reaction was added to the list of metabolic reactions. The number of external metabolites in the metabolic networks of E. coli, S. cerevisiae and S. aureus was 143, 116 and 84, respectively. The rows corresponding to these external metabolites were not included in the matrix ${\bf A}$ for determining UP or UC metabolites.

Using the bipartite matrix ${\bf A}$ for each one of the three organisms, we determined UP(UC) metabolites and reactions in the metabolic networks of E. coli, S. cerevisiae and S. aureus. We found the number of UP(UC) metabolites in the metabolic networks of E. coli, S. cerevisiae and S. aureus to be 291 (285), 395 (376) and 282 (237), respectively. We found the number of UP-UC metabolites in the metabolic networks of E. coli, S. cerevisiae and S. aureus to be 185, 178 and 145, respectively. Examples of UP-UC metabolites in E. coli and S. aureus networks are shown in Fig. 2.1. We found the number of UP(UC) reactions in the metabolic networks of E. coli, S. cerevisiae and S. aureus to be 289 (272), 391 (370) and 277 (218), respectively. The number of reactions that were either UP or UC or both (the ‘UP/UC reactions’) in the metabolic networks of E. coli, S. cerevisiae and S. aureus were found to be 418, 583 and 376, respectively.

We now describe in detail the algorithm to determine UP-UC clusters in any metabolic network.

1. 1.

   Starting from the reconstructed metabolic network of an organism, construct the bipartite matrix ${\bf A}$ of dimensions $m\times n$, where $m$ is the number of internal metabolites and $n$ is the number of reactions in the network as described in section 2.1.1.

2. 2.

   In the matrix ${\bf A}$, determine the rows corresponding to UP-UC metabolites as described in section 2.1.1.

3. 3.

   Obtain a matrix ${\bf B}$ from matrix ${\bf A}$ by setting every entry of each row in ${\bf A}$ that corresponds to a non UP-UC metabolite equal to zero, i.e., delete all links in the graph except those going into or out of UP-UC metabolites.

4. 4.

   From the matrix ${\bf B}$, construct a reaction-reaction graph in which each node corresponds to a reaction. The $n\times n$ adjacency matrix ${\bf C}$ = $C_{jk}$ of this graph is defined as $C_{jk}$ = 1 if $B_{ij}$ = 1 and $B_{ik}$ = -1, else $C_{jk}$ = 0. The matrix ${\bf C}$ represents a directed graph.

5. 5.

   The weak components of size $\geq$ 2 of the graph ${\bf C}$ are the various UP-UC clusters. These are obtained as follows: First convert the directed graph ${\bf C}$ into the associated undirected graph ${\bf\tilde{C}}$ by dropping all the directions of the arrows, i.e., $\tilde{C}_{jk}$ = 1 if $C_{jk}$ = 1 or $C_{kj}$ = 1 or both, else $\tilde{C}_{jk}$ = 0. Two nodes $j$ and $k$ in ${\bf C}$ are said to be weakly connected if there exists a path between them in the associated undirected graph ${\bf\tilde{C}}$. A weak component is a maximal set of nodes that are weakly connected to each other.

   By construction, if two reaction nodes $j$ and $k$ are adjacent (i.e., connected by a link) in ${\bf\tilde{C}}$, there exists a UP-UC metabolite that is produced in one of those reactions and consumed in the other. Thus, the fluxes of those two reactions will have a constant ratio in all steady states. This logic extends to entire connected cluster in ${\bf\tilde{C}}$ to which those reactions belong.

In practice, we cannot directly knockout a reaction in the metabolic network via molecular biology techniques in order to experimentally determine essential reactions for the survival of an organism. However, it is possible to knockout genes via molecular biology techniques and observe the effect of the knockout on the viability of the organism for different growth media in the wet lab. The knockout of a gene in effect results in the coded protein or enzyme being eliminated from the network. This would in effect result in the elimination of one or more reactions from the metabolic network that are catalyzed by the enzyme. If the knockout of a gene renders the organism unviable for certain growth medium then the gene is deemed ‘essential’ for that growth medium. The overall process of determining essential genes using in vivo experiments is very tedious and time consuming. Also, there is not always a one to one correspondence between genes and reactions in the network, and it may be difficult to determine essentiality of certain reactions using wet lab experiments. This is because the knockout of a single gene may in turn eliminate an enzyme that may catalyze multiple reactions in the metabolic network which are all removed from the network as a result of the knockout of the gene.

We have used a here computational method to determine essential reactions in metabolic networks. This method relies on the technique of flux balance analysis (FBA) [79, 80, 81, 82, 83, 84, 53, 85, 86, 87, 54, 55, 56, 88, 89, 66, 90, 91, 67, 92, 11]. Flux balance analysis (FBA) is a computational modelling technique which can be used to obtain the maximal growth rate of the organism supported by the metabolic network for any given nutrient medium. It also gives the steady state fluxes of all reactions in the metabolic network for any medium (for a detailed description of FBA see Appendix B). We have used FBA to study the metabolic networks of E. coli (version iJR904), S. cerevisiae (version iND750) and S. aureus (version iSB619).

FBA technique was used to compute fluxes of all metabolic reactions and optimal growth rate of E. coli for all possible aerobic minimal media. Any aerobic minimal media is characterized by availability of a single carbon source and key inorganic sources (ammonium, Fe²⁺, oxygen, phosphate, potassium, proton, sodium, sulfate and water) for the uptake of E. coli. If for a particular medium, the optimal growth rate obtained using the technique of FBA is zero, then under that condition the metabolic network does not permit the organism to grow. Using FBA, it was found that the E. coli metabolic network iJR904 supports nonzero growth under 89 aerobic minimal media [93]. Further, a reaction in the E. coli metabolic network was designated as ‘active’ if it has a nonzero flux value for at least one of the 89 minimal media, and ‘inactive’ otherwise. Similarly, using FBA, it was found that the metabolic networks for S. cerevisiae and S. aureus supported nonzero growth under 43 and 27 aerobic minimal media respectively. These 89, 43 and 27 aerobic minimal media that supported growth in E. coli, S. cerevisiae and S. aureus, respectively, were designated as ‘feasible minimal media’.

We used FBA to computationally determine essential reactions in the metabolic networks of E. coli, S. cerevisiae and S. aureus. We checked the effect of ‘switching off’ or removal of a reaction one by one from the metabolic network on the optimal growth rate of the organism obtained using FBA for different feasible minimal media. In FBA, a reaction can be switched off or removed from the network by setting the maximum flux (a parameter input in FBA) through the reaction equal to zero. We designated a reaction as ‘essential’ for a particular minimal medium, if switching the reaction off resulted in a zero optimal growth rate for that medium. We designated a reaction as ‘globally essential’ for an organism, if it was essential for all its feasible minimal media under aerobic conditions. A program to determine essential reactions in the E. coli metabolic network iJR904 is contained in Appendix C.

The number of essential reactions for each of the 89 minimal media varied between 200 and 240 and the number of globally essential reactions was 164 for the E. coli metabolic network. The number of essential reactions for each of the 43 minimal media varied between 165 and 187 reactions and the number of globally essential reactions was 127 for the S. cerevisiae metabolic network. The number of essential reactions for each of the 27 minimal media varied between 222 and 256 reactions and the number of globally essential reactions was 196 for the S. aureus metabolic network.

We found a set of 164 metabolic reactions to be globally essential in the E. coli metabolic network. Similarly, the number of globally essential reactions in S. cerevisiae (S. aureus) metabolic network was found to be 127 (196). We then tried to understand why certain reactions happen to be globally essential in terms of the underlying structure of the metabolic network. Notice that if a UP or UC metabolite is an essential intermediate in the production of a metabolite that are part of the biomass reaction then that reaction responsible for the production or consumption of that UP or UC metabolite becomes essential for the growth of the organism. Of the 164 globally essential reactions in the E. coli metabolic network, 133 were found to be either UP or UC. Similarly, a high fraction of globally essential reactions in the metabolic networks of S. cerevisiae and S. aureus were found to be UP or UC (see Table 2.2). This explains why the subset of 133, 86 and 157 reactions are globally essential in E. coli, S. cerevisiae and S. aureus, respectively, namely, there is simply no other path around these reactions in the entire network to produce or consume some metabolite that is presumably required for the eventual production of biomass.

The probability of such a high overlap between the set of globally essential reactions and set of UP/UC reactions occurring by pure chance is very small. To quantify this, the result was compared to a null model in which the two sets corresponding to globally essential reactions and UP/UC reactions were considered to be independent of each other. The total number of reactions in the E. coli metabolic network is 1176. The number of globally essential reactions is 164 and the number of UP/UC reactions is 417 for the E. coli network. The probability that out of a set of 1176 reactions in E. coli, two independently chosen subsets of size 417 and 164 will have an intersection of 133 or greater is $p<10^{-37}$ (any one or both of the subsets is chosen randomly). Similarly, we obtained very small $p$ values for the metabolic networks of S. cerevisiae and S. aureus (see Table 2.2).

In an earlier paper [94]  Mahadevan and Palsson had determined the ‘lethality fraction’ for each metabolite in the network. The lethality fraction of a metabolite was defined as the fraction of reactions in which the metabolite is involved that are essential. Mahadevan and Palsson had observed that this lethality fraction of the low degree metabolites is on average comparable to high degree metabolites. In particular, they found that some metabolites with in-degree and out-degree unity (that we have designated as UP-UC metabolites) have lethality fraction unity. We have presented here a stronger result regarding the role of low degree metabolites: most globally essential reactions involve at least one UP or UC metabolite. The essential reactions may involve other metabolites of higher degree, but their essentiality is due to their uniqueness in producing or consuming a UP or UC metabolite.

It was shown above that 133 out of 164 globally essential reactions were associated with a UP or UC metabolite in the E. coli metabolic network. To understand the remaining globally essential reactions, a reduced or pruned version of the E. coli network was considered.

The databases of reconstructed metabolic networks that have been used here contain certain reactions that can only have a zero flux value under any steady state due to stoichiometric reasons. Such reactions have been referred to as ‘strictly detailed balanced’ reactions [95, 96] or ‘blocked’ reactions [75] (see section B.2 in Appendix B for more details). The blocked reactions can be removed from the reconstructed metabolic networks for any steady state analysis. We used a previously described algorithm by Burgard et al [75] to determine the blocked reactions in the metabolic networks of E. coli, S. cerevisiae and S. aureus. This algorithm to determine blocked reactions has been described in detail in section B.2.1 in Appendix B. A program to determine blocked reactions in the E. coli metabolic network is contained in Appendix C. The number of blocked reactions in the metabolic networks of E. coli, S. cerevisiae and S. aureus were found to be 290, 800 and 294, respectively (for details see Table 2.3). We removed the 290 blocked reactions in the E. coli metabolic network from the list of 1176 reactions to obtain the ‘reduced network’ of 886 reactions in E. coli. We emphasize that the removal of blocked reactions from the original network does not affect our results obtained using FBA for any of the 89 minimal media considered for E. coli. The 290 blocked reactions removed from the original E. coli metabolic network of 1176 reactions had a zero flux value for all 89 minimal media studied here. Further, since the blocked reactions are guaranteed to have a zero flux value for any minimal media, they can be never essential. So, the set of essential reactions obtained by implementing FBA on the reduced network is exactly the same as that obtained from the original network for any minimal media. Thus, the set of 164 globally essential reactions is the same for the original network of 1176 reactions and the reduced network of 886 reactions in E. coli. Similarly, the reduced networks for S. cerevisiae and S. aureus were obtained by removing the corresponding blocked reactions from the original network (see Table 2.3).

The set of UP(UC) metabolites was then determined for the reduced metabolic network obtained after removing blocked reactions from the original network. A metabolite that was not UP(UC) in the original metabolic network can become UP(UC) in the reduced network. This may happen due to the removal of blocked reactions which may be contributing to the degree of a metabolite in the original network. The set of UP(UC) metabolites and reactions were obtained for the reduced networks of E. coli, S. cerevisiae and S. aureus. We found 352, 306 and 276 reactions to be UP/UC in the reduced networks of E. coli, S. cerevisiae and S. aureus, respectively (see Table 2.3). The set of UP/UC reactions in the reduced network turns out to smaller than that for original network. This is so because several reactions that were UP/UC in the original network happen to be blocked and have been removed from the reduced network. Also, some metabolite that was earlier not UP(UC) can now become UP(UC) in the reduced network due to removal of blocked reactions which adds new reactions to the UP/UC set. However, the number of reactions added turns out to be smaller than the number of reactions removed (see Table 2.3). The new UP(UC) metabolites have, by definition, their in-degree (out-degree) unity in the reduced network. Even in the original network the set of UP(UC) metabolites have a low degree. In E. coli, the average in-degree (out-degree) of the set of UP(UC) metabolites in the reduced network was found to be 1.31(1.33) in the original network. We remark that UP(UC) reactions in the reduced network are uniquely determined starting from the original network.

We found 156 out of the 164 globally essential reactions (95$\%$) to be UP or UC in the E. coli reduced network ($p<10^{-62}$). Similarly, it was found that almost all globally essential reactions in S. cerevisiae and S. aureus were either UP or UC in the reduced network (92$\%$ and 93$\%$, respectively). Thus, we have shown here that nodes with a low degree of connectivity (i.e., UP(UC) metabolites) play an ‘essential’ role in metabolism (see Table 2.2).

The results obtained here provide some insight into the structural or topological origin of essential reactions in metabolic networks. It is, of course, obvious that if a metabolite is an essential intermediate for the production of some biomass metabolite, and if this metabolite is uniquely produced or uniquely consumed, then the corresponding production or consumption reaction will be essential for the growth of the cell. However, the converse of this statement that all globally essential reactions in the network should be UP/UC is far from obvious. The finding that about 5-8 $\%$ of globally essential reactions do not satisfy the UP/UC property proves that the converse statement is indeed false. Thus, the fact that the overwhelming majority (92-95 $\%$) of globally essential reactions satisfy the UP/UC topological property in the network is a characterization of the nature of metabolic networks found in organisms. We remark that we do not as yet understand why the remaining globally essential reactions happen to be essential.

We found that there were 352 UP/UC metabolic reactions in the E. coli reduced network. 156 of these 352 reactions were essential for all 89 minimal media or globally essential in E. coli. It was found that there were 400 reactions in the E. coli network which were essential for at least one of the 89 minimal media. These 400 reactions that were found to be essential for some of the 89 minimal media may be designated as ‘conditionally essential’ for the E. coli metabolic network. We found 288 of the 352 UP/UC reactions (82 %) in the E. coli reduced network to be conditionally essential. Such a large overlap is very unlikely ($p<10^{-74}$) between the two sets of UP/UC reactions and conditionally essential reactions. Some of these conditionally essential UP/UC reactions were part of the input pathways of only one carbon source and hence, they were essential only for that minimal media. In S. cerevisiae and S. aureus, we found the number of conditionally essential reactions to be 269 and 331, respectively. In S. cerevisiae, we found 170 out of 306 UP/UC reactions (56 %) in the reduced network to be conditionally essential, while in S. aureus, 257 out of 276 (93 %) were found to be conditionally essential. We found the $p$ values for such large overlaps between the set of UP/UC reactions and conditionally essential reactions in S. cerevisiae and S. aureus to be $p<10^{-22}$ and $p<10^{-67}$, respectively. In E. coli and S. aureus more than 80% of the UP/UC reactions were conditionally essential while in S. cerevisiae only about 56% of the UP/UC reactions were conditionally essential. The substantial difference in the fraction of UP/UC reactions which are conditionally essential in yeast as opposed to that for the two bacteria may reflect a more evolved metabolic structure in the eukaryote. Yeast may have many more alternative pathways and may be more robust to random failures as compared to the two bacteria. This observation needs to be further investigated.

We had mentioned earlier that it is not possible to always experimentally determine essential reactions for an organism due to lack of one to one correspondence between genes and reactions in the metabolic network. Experimentally one can only determine essential genes for a given medium. Gerdes et al [97] have determined experimentally the list of essential genes in E. coli for a rich medium. To compare the computationally predicted results with experimental data of Gerdes et al, we implemented FBA under rich medium for the E. coli metabolic network iJR904.

A set of 95 reactions were found to be essential under rich medium for the E. coli metabolic network. 89 of these 95 reactions were found to be either UP or UC in the E. coli reduced network. Of the 95 essential reactions in rich medium, information about the corresponding genes (coding for the enzymes) was available for only 85 reactions in the database iJR904 for the E. coli metabolic network. Of these 85 reactions, 14 reactions had known isozymes, i.e, multiple enzymes catalyzing a single reaction in the network. Hence, the genes corresponding to isozymes are not expected to be essential for these 14 reactions as the effect of the knockout of a single gene can be compensated by another gene. Of the remaining 71 reactions, 5 had associated genes whose essentiality was undetermined in the database by Gerdes et al [97]. Of the remaining 66 reactions, a fairly high fraction of 38 reactions had associated genes that were found to be essential in the database.

Conversely, of the 618 genes determined to be essential for E. coli by Gerdes et al in rich medium, 158 genes were also part of the E. coli metabolic network iJR904. Note that there are 904 genes coding for various enzymes catalyzing reactions in the E. coli metabolic network. We found that 103 of the 158 essential genes had their enzymes catalyzing only a single reaction in the E. coli metabolic network. Of these 103 essential genes, 62 were associated with a UP or UC reaction in the original E. coli metabolic network. Further, it was found that 73 of the 103 essential genes were associated with a UP or UC reaction in the E. coli reduced network. A possible reason for the difference between theoretical prediction and experimental data could be reconciled on the basis of incompleteness of the reconstructed metabolic network iJR904 used for our study. There may be alternative pathways in the real organism that are missing in the reconstructed network iJR904. Further, for certain reactions in the metabolic network there may be isozymes that are presently not included in the network.

In this study, we have considered external environment corresponding to buffered minimal media characterized by $m_{i}$ that are constant in time. For buffered media, we have

independent of time $t$ which represents a constant external environment for the cell. For each buffered medium considered, the components of m corresponding to the metabolites present in the external environment were set to unity and the remaining components were set to zero. The E. coli metabolic network iJR904 is capable of transporting 143 metabolites into and out of the cell. The 143 external metabolites in the metabolic network iJR904 include 131 organic molecules and 12 inorganic molecules. Further, 96 of the 143 external metabolites in the metabolic network iJR904 are also included in the regulatory network part of the database iMC1010^v1. The 96 external metabolites accounted in the regulatory network can be divided into 86 organic molecules and 10 inorganic molecules. In this study, we have considered the following classes of minimal media:

• (a)

  List of 93 minimal media: Using the 96 external metabolites (86 organic and 10 inorganic external metabolites) contained in the regulatory network, we can construct in principle 86 aerobic and 86 anaerobic minimal media (i.e, in total 172 minimal media). A minimal medium is characterized by the availability of single organic source of carbon and the ions of ammonium, sulphate, phosphate, hydrogen, iron, potassium and sodium. The components of m corresponding to these metabolites were set to unity and others were set to zero in a given minimal medium. Oxygen was set to unity in the aerobic media and to zero in anaerobic media. FBA was performed for each of the 172 possible minimal media for the E. coli metabolic network iJR904 to obtain the optimal growth rate for each of the 172 minimal media (see Appendix B for details on FBA) . Using FBA, it was found that only 62 aerobic and 31 anaerobic minimal media supported nonzero growth. The list of 93 minimal media (62 aerobic and 31 anaerobic) is provided in Table LABEL:minimalmedialist. Note that in obtaining this list of 93 minimal media, we have used the FBA model without incorporating any regulatory constraints. We have used this list of 93 minimal media for obtaining most our results.

• (b)

  Library of 109732 media: Following Barrett et al [111], we obtained a list of 109732 media. The 143 external metabolites in the E. coli metabolic network were divided into the four different sets: carbon sources, nitrogen sources, sulphur sources and phosphorus sources. Note that the four sets are not mutually exclusive. For example, an amino acid is a carbon source as well as a nitrogen source. A large library of 109732 different media was obtained by picking a single metabolite from each of the four sets corresponding to carbon sources, nitrogen sources, sulphur sources and phosphorus sources. In addition, each media may or may not contain oxygen. Some of the work reported here uses this larger library of 109732 media.

The 19 Boolean variables $c_{i}(t)$ corresponding to conditions can be divided as follows:

• •

  15 of the 19 Boolean variables $c_{i}(t)$ depend upon the configuration of a subset of transcription factors and external metabolites at time $t$, i.e.,

  $$c_{i}(t)=C_{i}({\bf t}(t),{\bf m(t)}),\quad\quad i=1,2,\ldots,15,$$

  (3.9)

  where the Boolean functions $C_{i}$ are given by the database iMC1010^v1. The functions given by Eq. 3.9 can be substituted in Eq. 3.3 which eliminates these 15 $c_{i}$ variables from the dynamical system at the expense of more complicated effective dependence of $g_{i}(t)$ on ${\bf t}(t)$ and ${\bf m}$.

• •

  Another condition variable $c_{i}(t)$ corresponds to the growth of the cell. This variable representing the growth of the cell is fixed to unity as most media considered in this study give a nonzero growth rate for the cell.

• •

  Another condition variable $c_{i}(t)$ corresponds to the pH of the external environment. In this study, pH is taken to be between 5.5 and 7, i.e., weakly acidic. For example, one of the habitats of E. coli is the human gut where the pH is weakly acidic. The pH condition affects only 3 genes in the network. For two of the genes that are affected by pH condition, the operative regulatory clause is ‘pH $<$ 4’. The Boolean variable $c_{i}$ corresponding to pH is fixed to zero (false) for these two genes. For the third gene the clause is ‘pH $<$ 7’. For this gene, the Boolean variable $c_{i}$ corresponding to pH is fixed to unity (true).

• •

  The remaining two condition variables $c_{i}(t)$ correspond to ‘surplus FDP’ and ‘surplus PYR’ in the database iMC1010^v1. They represent whether surplus amounts of fructose 1,6-bisphosphate and pyruvate are being produced in the cell. These two conditions depend upon the values of some of the internal fluxes $v_{i}$ and the presence of an external metabolite, fructose, through specified Boolean functions in the database. The variable corresponding to external fructose is treated as unity, if the minimal medium includes fructose and zero otherwise. The treatment of the internal fluxes $v_{i}$ is discussed next.

The state of some of the genes in the regulatory network is influenced by the state of 21 components of vector ${\bf v}$ representing fluxes of internal metabolic reactions. As described in Covert et al [112], the fluxes of internal reactions in the database are surrogate for internal metabolite concentrations inside the cell that determine the state of genes in the regulatory network. The internal metabolite concentrations have been approximated by internal fluxes in the database iMC1010^v1 as the database is intended to model the integrated metabolic and regulatory system f E. coli using the framework of FBA. Since FBA cannot determine the internal metabolite concentrations, the internal metabolite concentrations have been approximated by internal fluxes of reactions in the network. In this study, we have treated the variables corresponding to the internal fluxes in two distinct ways leading to two slightly different dynamical systems as discussed below.

• (a)

  System A: In the first approach, for a given external environment ${\bf m}$, we determined whether each of the 21 internal reactions in the E. coli metabolic network whose fluxes correspond to the 21 components of vector ${\bf v}$ was ‘blocked’ or not [95, 75, 35]. A reaction is said to be blocked in a particular environmental condition, if under that medium no steady state flux is possible through it [95, 96, 75]. The algorithm to determine blocked reactions for a given environment ${\bf m}$ is described in detail in section B.2 in Appendix B. An internal flux variable $v_{i}(t)$ was set equal to zero for a medium (specified by vector ${\bf m}$), if the reaction was found to be blocked for that medium, and unity otherwise. Thus, in this approach, the internal fluxes $v_{i}$ were not dynamical variables, but rather fixed parameters (albeit fixed with an eye on self-consistency).

• (b)

  System B: In the second approach, the Boolean variables corresponding to internal fluxes $v_{i}$ were allowed to be dynamical by making a simplifying assumption about their dynamics. In the metabolic network, the flux of a reaction is determined by the concentrations of the participating metabolites and the catalyzing enzymes. The enzyme concentrations are determined by the activity or state of their respective genes. In a discrete-time approximation, an enzyme is present at time $t$ if the genes coding for it were active at the previous time $t-1$. Thus, the internal flux $v_{i}(t)$ is set equal to unity if the genes coding for the enzyme catalyzing that metabolic reaction were active at the previous time $t-1$, and zero otherwise. This could be done for a subset of 10 out of 21 internal reaction fluxes, since the genes of their catalyzing enzymes were part of the 583 genes in the database. Genes coding for the enzymes of the remaining 11 internal reaction fluxes were not part of the regulatory network database and hence the corresponding $v_{i}$ could not be made dynamical variables in this fashion. These latter $v_{i}$ were fixed as in part (a) for each medium by checking their blocked status for the medium under consideration. The approach (b) introduces feedbacks in the genetic regulatory network.

Note that the choices (a) and (b) for the ${\bf v}$ variables yield different dynamical systems for Eq. 3.1 which we denote as system A and system B, respectively. In system B, 6 out of 583 genes have additional links from other genes in the set compared to system A. Programs implementing the two dynamical systems A and B are discussed in Appendix C.

In chapter 2, we have introduced the notion of ‘uniquely produced’ (‘UP’) and ‘uniquely consumed’ (‘UC’) metabolites in the metabolic network. We have designated a metabolite that is both UP and UC as a ‘UP-UC’ metabolite in the network. A UP-UC metabolite has only one reaction producing it and one reaction consuming it in the network. We have designated a set of reactions connected by UP-UC metabolites as a ‘UP-UC cluster’ in the metabolic network. Under any steady state, the flux of the reaction producing a UP-UC metabolite is proportional to the flux of the reaction consuming it. Hence, in any steady state, the fluxes of all reactions that are part of a single UP-UC cluster in the metabolic network are always proportional to each other. We have determined the list of UP(UC) metabolites and UP-UC clusters of reactions in the metabolic networks of E. coli, S. cerevisiae and S. aureus. Further, we have used the computational technique of flux balance analysis (FBA) to determine essential reactions for growth in the metabolic networks of E. coli, S. cerevisiae and S. aureus.

In chapter 2, we have shown that UP-UC metabolites in the metabolic network lead to correlated UP-UC clusters of reactions in the network. We found that large UP-UC clusters are over-represented in the real metabolic network as opposed to randomized networks with same local connectivity as the real network. We then showed that UP-UC clusters at the metabolic level correspond, with a high probability, to sets of genes forming expression modules at the regulatory level in E. coli. We have used here only the operon information in determining the correspondence between UP-UC clusters in the metabolic network and sets of coexpressed genes in regulatory network. It is possible for genes that do not belong to the same operon to be coexpressed inside the cell. Thus, an extended analysis of E. coli gene expression data may find a better correspondence between UP-UC clusters at the metabolic level and sets of coexpressed genes at the regulatory level. Further, the analysis of gene expression data to check coregulation of genes corresponding to UP-UC clusters of reactions should be extended to other organisms where operon information is not available.

Metabolic networks inside different organisms have been shown to follow a power law degree distribution [25, 28]. It has been suggested that the power law degree distribution of real networks contributes towards its robustness [33]. In particular, it has been emphasized that networks with power law degree distribution are vulnerable to selective attack on its high degree nodes or hubs, while removal of many low degree nodes from these networks have negligible effect on the functionality of the system [33]. For example, in the case of the internet, the removal of high degree nodes corresponding to routers with many connections can turn out to be fatal for the communication system. Further, for the S. cerevisiae protein-protein interaction network, the essentiality of a protein was also found to be correlated with the degree of the protein in the network [34]. However, in case of the metabolic network, the metabolites participate in reaction processes where they are produced or consumed. We can only control the reaction process through its catalyzing enzyme that is a product of some gene. It is unclear how a biological process can lead to removal of metabolites from the network. A removal of high degree metabolite from the metabolic network would require elimination of all reactions in which the metabolite participates or knockout of all genes associated with reactions in which the metabolite is involved. Further, genetic mutations effectively correspond to removal of biochemical reactions from the metabolic network. Thus, in case of the metabolic network, one is interested in determining the effect of removing a reaction rather than the effect of removing a metabolite. In chapter 2, we have shown that most globally essential reactions in the metabolic networks of E. coli, S. cerevisiae and S. aureus either produce or consume a low degree UP or UC metabolite. The essential reactions may involve other metabolites of higher degree, but their essentiality is due to their uniqueness in producing or consuming an intermediate low degree UP or UC metabolite that is needed for the eventual production of biomass. Thus, we have shown here that, in a consideration of robustness or fragility of metabolic networks to naturally occurring perturbations, it is the role of low degree metabolites that needs to be considered rather than high degree metabolites. This is an example of how the study of dynamical properties (in this case, flows) and functionality at the system level can alter one’s perspective on what is significant for robustness or fragility of real world complex systems. We remark that the importance of low degree nodes in determining the robustness or fragility of complex networks has also been observed elsewhere [115, 116]. There the crashes of an evolving catalytic network were found to be due to ‘core-shifts’ caused by changes in low degree nodes of the network.

The existence of essential reactions in the metabolic network is an indicator of the fragility in the system. Even though the metabolic networks we have studied have many reaction nodes, a small perturbation such as the removal of a single essential reaction node from the network, destroys the functionality of the complete metabolic network. One way of dealing with this fragility is by introducing redundancy at the level of genetic network by associating genes coding for isozymes with essential reactions. Isozymes are multiple enzymes catalyzing a single reaction in the metabolic network. We have to knockout multiple genes at the same time in order to remove an essential reaction with associated isozymes from the metabolic network. Hence, we may expect essential reactions in the metabolic network to have a greater probability of being associated with isozymes than non-essential reactions. However, in an earlier paper [117], Papp et al showed that both essential and non-essential reactions in the metabolic network of S. cerevisiae are equally likely to be associated with isozymes. In that paper, isozymes were found to be associated with greater probability with reactions carrying high flux rather than with essential reactions in the metabolic network [117]. This finding suggests that organisms in the course of evolution have developed redundancies in the genetic network predominantly to associate isozymes with metabolic reactions carrying high flux rather than with essential reactions in the network.

This raises the question: why has evolution tolerated the fragility associated with essential reactions? The finding that essential reactions can be tagged by UP or UC metabolites may provide some insight into this. UP or UC metabolites that participate in very few reactions perhaps do so in part because some feature of their chemical structure prohibits ready association with other molecules in nature. If so, the low degree of UP or UC metabolites is a consequence of constraints coming from chemistry. Then evolution seems to tolerate essential reactions that produce or consume such metabolites because chemistry has left it with no other choice.

Alternatively, the fragility associated with low degree metabolites may be a byproduct of some other desirable property that contributes to evolvability or robustness of the system, such as modularity. In chapter 2, we have shown that UP-UC metabolites which have low degree in the network lead to correlated clusters of reactions in the metabolic network. This raises the question: if UP-UC metabolites contribute to modularity, could it be that the evolutionary advantages of that have outweighed the disadvantage of the above mentioned fragility caused by the same low degree metabolites? Is it the case that evolution has preferred ‘chemically constrained’ low degree metabolites in spite of the fragility they cause because they contribute to modularity?

A goal in biology is to understand highly evolved biological organization in terms of simpler and more inevitable structures [118]. In chapter 2, we have presented evidence that certain regulatory modules, in particular certain operons, mirror the UP-UC structure of the metabolites whose production and consumption they regulate in E. coli. This could be an example of how the origin of certain regulatory structure can be traced to simple chemical constraints. In future, we hope to shed light on this by investigating metabolic and regulatory networks inside many other organisms.

The fragility caused by low degree metabolites in metabolic networks can have potential applications in medicine. We have observed that essential reactions are explained by UP/UC structure in three organisms. Thus, it is likely that UP/UC structure explains essential reactions in other organisms that are pathogens for humans. This generates candidate targets for therapeutic intervention. It is conceivable that drugs could be found that incapacitate the enzymes of one or the other essential reactions in pathogens.

We now discuss some limitations and caveats associated with our results reported in chapter 2. The results have been obtained by using the metabolic network databases iJR904 [66], iND750 [67] and iSB619 [68] for E. coli, S. cerevisiae and S. aureus, respectively. These databases were reconstructed using the annotated genome sequences for the three organisms and the available biochemical information in the literature. Since the annotation of the fully sequenced genomes of E. coli, S. cerevisiae and S. aureus is not yet complete, the list of reactions contained in the three metabolic network databases iJR904, iND750 and iSB619 is still incomplete.

We may expect additions to the list of biochemical reactions in the reconstructed metabolic network databases in future based on more complete annotation information. The addition of reactions to the reconstructed metabolic network databases can introduce alternate pathways for certain crucial metabolites required for the eventual production of biomass metabolites in the databases for the three organisms studied here. Thus, the introduction of alternate pathways may affect the set of essential metabolic reactions obtained using FBA for the three organisms. Similarly, the inclusion of additional reactions in the metabolic network databases would render some of the presently UP(UC) metabolites non-UP(non-UC).

In this context, it is worth noting that for the metabolic network databases as they stand, the present definition of UP(UC) metabolites allows us to establish a connection between distinct properties of the metabolic network. In chapter 2, we have shown that UP-UC metabolites lead to UP-UC clusters which are correlated sets of reactions in the metabolic network. The genes corresponding to UP-UC clusters were shown to predict regulatory modules in E. coli. We have also shown that most essential metabolic reactions are explained by their association with a UP or UC metabolite. Also, most of our results have been shown to hold for metabolic networks inside three distinct organisms. This suggests that our definition of UP(UC) metabolites and reactions does capture a certain useful pattern in the metabolic network. We remark that a reduction in the number of UP(UC) metabolites and essential reactions would only improve the correspondence of UP-UC clusters and regulatory modules described in section 2.3 and between the theoretically predicted and experimentally observed essential genes discussed in section 2.6.4.

In chapter 3, we have studied the genetic network controlling E. coli metabolism as represented in the database iMC1010^v1 [41] as a Boolean dynamical system. We have constructed an effective Boolean dynamical system of genes and external metabolites using the information contained in the database iMC1010^v1. We have studied the dependence of the attractors of the Boolean dynamical system representing the genetic network controlling E. coli metabolism to changes in initial conditions of genes and configuration of the external environment. Robustness is an inherent property of living systems which enables them to maintain their functionalities in the face of diverse perturbations [119, 120, 121]. Biological systems encounter externally induced perturbations in the form of changes in the environment and internal perturbations such as fluctuations in gene configurations. In chapter 3, we have shown that the genetic network controlling E. coli metabolism for a fixed external environment has essentially a unique attractor for the configuration of the genes regardless of the initial conditions or perturbations of gene configurations. So, we found the attractors of the dynamical system to be fixed points or low period cycles for any given fixed environmental condition and the system exhibits the property of homeostasis. However, robustness is a more general property than homeostasis which concerns itself with maintaining system functionality rather than system state. In chapter 3, we have then shown that when the dynamical system encounters a changed external environment, the system again flows and stabilizes to another unique attractor state regardless of the initial configuration of genes; however, the attractor for the changed environment may be widely separated from the attractor for the previous environment. We found the genetic network controlling E. coli metabolism as a dynamical system to be flexible to different environmental conditions as their attractors show a wide variation. In chapter 3, using the technique of FBA, we have further shown that the attractors of the genetic network for most environmental conditions allow optimal functioning of the underlying metabolic network. Thus, we have shown that the robustness of the genetic network controlling E. coli metabolism manifests itself in two different ways. When the system encounters internal perturbations in form of changes in genes’ configuration for a fixed external environment, the system returns to the present attractor. When the system encounters external perturbations in form of changed environment, the system moves to a new attractor that maintains the systems functionality in the sense of maintaining growth rate of the organism. Further, we have shown that the robust behaviour of the genetic network in moving to a different stable attractor state in response to a changed external environment also ensures metabolic efficiency for the organism.

In chapter 4, we have studied in detail the design features of the genetic network controlling E. coli metabolism in order to understand the origin of the observed system level dynamical properties of homeostasis and response flexibility. We have found the genetic network controlling E. coli metabolism as represented in the database iMC1010^v1 to be hierarchical and an largely acyclic graph. Cycles where they did exist were found to be of small length and localized. The nodes with no incoming links in the genetic network are referred to as the root nodes of the hierarchy which correspond to different external metabolites. The nodes with no outgoing links in the genetic network are referred to as the leaves of the hierarchy which correspond to genes coding for metabolic enzymes. The nodes at the intermediate level of the hierarchical acyclic graph have both incoming and outgoing links in the genetic network and correspond to genes coding for transcription factors. The observed hierarchical architecture of the genetic network controlling E. coli metabolism with root control in the hands of external metabolites endows the system with the twin dynamical properties of homeostasis and response flexibility. Thus, the robust behaviour of the dynamical system has its origin in the underlying architecture or topology of the genetic network.

As mentioned above, the genetic network controlling E.coli metabolism as represented in the database iMC1010^v1 lacks feedback from genes to other genes via transcription factors. The models studied originally by Kauffman [57, 59] were random Boolean networks. Those networks had substantial feedbacks between genes and hence had more complicated attractors and dynamics. One of our main results is that the genetic network controlling E.coli metabolism as a real biological system is structured (and hence departs from random networks) in such a way that it has simple attractors and dynamics. Thus, while at the abstract level the Boolean dynamical system studied by us is not very different from Kauffman’s (apart from the inclusion of the external environment), our dynamical results are quite different because the underlying network has a very different structure from the one considered by Kauffman.

In chapter 4, we have shown that removing the leaf nodes at the bottom of the hierarchy corresponding to the genes coding for metabolic enzymes along with their links from the full graph shown in Fig. 3.1 leads to a subgraph shown in Fig. 4.2 with many disconnected components. The various disconnected components shown in Fig. 4.2 are dynamically independent of each other and may be regarded as modules of the genetic network. The different disconnected components shown in Fig. 4.2 are regulated by different external metabolites. Thus, for any given environmental condition, only a subset of the disconnected components shown in Fig. 4.2 get switched on. The highly disconnected structure at the level of genes coding for transcription factors compared to the highly connected structure of the full network implies that the leaf nodes corresponding to the genes coding for metabolic enzymes get incoming links from many different disconnected components shown in Fig. 4.2. Thus, although there is a lack of crosstalk between disconnected components shown in Fig. 4.2, there is a lot of crosstalk between the modules at the leaf level, that of genes coding for metabolic enzymes in the genetic network.

In chapter 4, we have argued that the above mentioned disconnected structure at the level of genes coding for transcription factors can contribute towards the evolvability of the genetic network. Since the different disconnected components are regulated by different subsets of external metabolites in the environment, if a population of such cells encounters an environment where a particular set of external metabolites is consistently absent, then the species can discard the disconnected components regulated by these external metabolites from the genetic network without a loss of functionality with respect to other external metabolites. Conversely, if the population encounters a new set of external metabolites in a sustained manner, it can develop a new module with pathways controlled by this new set without affecting the existing functionality, thereby adapting itself to the new niche.

The results reported in chapters 3 and 4 have been obtained using the integrated regulatory and metabolic network database iMC1010^v1 [41] for E. coli and there are limitations that stem from the reconstructed database itself. The database iMC1010^v1 for E. coli has been reconstructed using various literature sources. The database iMC1010^v1 describes the regulation of 583 genes in E. coli. Of the 583 regulated genes in the database iMC1010^v1, 479 genes code for metabolic enzymes in E. coli. These 479 enzyme coding genes are a subset of the 904 enzyme coding genes associated with various reactions in the E. coli metabolic network database iJR904. Thus, the database iMC1010^v1 does not describe the regulation of a large fraction of genes coding for metabolic enzymes in E. coli. In the database iMC1010^v1, the set of 583 genes are regulated by a set of 96 external metabolites. In the E. coli metabolic network iJR904, there are 143 external metabolites which can be transported across the cell boundary. The set of 96 external metabolites accounted in the database iMC1010^v1 are a subset of the 143 external metabolites contained in the metabolic network database iJR904. Thus, the database iMC1010^v1 does not completely account for the effect of all possible external metabolites that the E. coli can uptake or excrete on the state of genes in the regulatory network. Further, the present incomplete genetic network database iMC1010^v1 of 583 genes and 96 external metabolites could have many connections as false positives or false negatives, especially the latter. It is also possible that some of the Boolean rules in the database iMC1010^v1 may be incorrect or incomplete due to lack of detailed experimental data for different possible environmental conditions. Hence, we expect modifications and expansion of the present database iMC1010^v1 describing the genetic network controlling E. coli metabolism in the near future. This could modify our results on the dynamics of the genetic network controlling E. coli metabolism reported in this thesis.

We next ask the question: How would the future introduction of additional gene nodes and connections into the existing network database iMC1010^v1 affect the results reported in this thesis? If new nodes and connections corresponding to genes coding for metabolic enzymes are added to the present database iMC1010^v1, it is unlikely to affect our qualitative conclusions about the nature of attractors of the genetic network controlling E. coli metabolism significantly. The reason is that most of the genes coding for metabolic enzymes are likely to be leaves of the genetic network like the nodes at the bottom of Fig. 4.1, in which case they would have no outgoing links and would not affect the dynamics of other gene nodes in the network. However, the inclusion of genes coding for metabolic enzymes as well as additional connections to existing genes in the genetic network would add to the constraints on FBA. It would be then interesting to see the extent to which regulatory dynamics enhances metabolic efficiency in different environmental conditions for an enlarged genetic network accounting for the regulation of most genes coding for enzymes in the E. coli metabolic network.

If new nodes and connections corresponding to genes coding for transcription factors are added to the present database iMC1010^v1, it can affect our qualitative conclusions about the dynamics and the nature of the attractors of the genetic network controlling E. coli metabolism. In particular, the introduction of feedback loops or cycles between genes coding for transcription factors could lead to longer cycles as attractors of the dynamical system. This is a question to be ultimately settled by more detailed knowledge of the empirical facts. Several regulatory systems with feedbacks are known, e.g., the Yeast cell-cycle network and the Drosophila segment polarity network [61, 62]. Several genes in E. coli are known to have autoregulatory self-loops [36] that are not included in the present database iMC1010^v1. The introduction of autoregulatory loops in the present network could produce two-cycles at the individual nodes even at constant input. Earlier works have observed that apart from self-loops, the transcriptional regulatory network of E. coli is largely acyclic [37, 46, 47] and has a small depth of about 5. We remark that the lack of feedback from genes to other genes via transcription factors is not an assumption on our part, rather it reflects the way this biological system actually is as captured in the present database iMC1010^v1 and also in other studies [37, 46, 47].

In addition to the feedback from genes to other genes via transcription factors, discussed above, there can be another kind of feedback in the genetic network controlling E.coli metabolism. This is the feedback of the concentrations of internal metabolites in the metabolic network into the genetic network. The feedback of internal metabolite concentrations could affect both the genes coding for enzymes and genes coding for transcription factors. In the database iMC1010^v1, the effect of internal metabolite concentrations on the state of the genes has been approximated via fluxes of internal reactions in the metabolic network. In the present database iMC1010^v1, there are 21 internal fluxes affecting the state of genes in the network. The feedback of these 21 internal fluxes from the metabolic network affect the state of only 5 genes coding for transcription factors and 16 genes downstream of these coding for metabolic enzymes in the present network. In chapter 3, we found that these 21 genes in the present network which get affected by metabolic feedback via internal fluxes may have their states undetermined across the different attractors for a given fixed environment. The state of these 21 genes can oscillate back and forth between zero and one in the two-cycle attractors obtained for any given fixed environment. This set of 21 genes in the genetic network controlling E.coli metabolism forms the twinkling island.

We may expect introduction of additional feedback from internal metabolite concentrations into the genetic network during the future network expansion. If the additional feedback from internal metabolite concentrations affect enzyme coding genes in the genetic network, then it is unlikely to change our qualitative conclusions about the observed dynamics of the system. Since the enzyme coding genes happen to be the leaf nodes of the hierarchical acyclic graph, they do not affect other nodes in the genetic network controlling E.coli metabolism. Thus, the effect of feedback from internal metabolite concentrations on to an enzyme coding gene is unlikely to affect the dynamics of other nodes in the genetic network. If the additional feedback from internal metabolite concentrations affect genes coding for transcription factors in the genetic network, then it can change our qualitative conclusions about the observed dynamics of the system. It is possible that there is a larger set of genes forming the twinkling island, and possibly longer cycles as attractors of the genetic network controlling E. coli metabolism. In chapter 4, we have shown that the genetic network controlling E. coli metabolism has a highly disconnected structure at the level of genes coding for transcription factors. This disconnected architecture of the network at the level of genes coding for transcription factors as well as the overall small depth ($\sim$4) of that graph suggests that the introduction of addition feedback loops at the level of genes coding for transcription factors may only lead to cyclic attractors which are of low period and localized. We further showed that most Boolean functions in the present database iMC1010^v1 satisfy the canalyzing property. These canalyzing functions may provide additional stability to the network when additional feedback loops get introduced.

The metabolic network needs to function at all times as it is responsible for utilizing the nutrients available in the external environment to produce the key molecules required for growth and maintenance of the cell. The lack of feedback between genes coding for transcription factors in the genetic network controlling E. coli metabolism observed here may be expected. However, it is known that there are many feedback loops at the level of the metabolic network where the concentration of an internal metabolite that is the end product of a metabolic pathway regulates the activity of the enzyme or the protein catalyzing the reaction at the start of the pathway. Such processes where the activity of the enzymes or proteins that are products of genes at the leaf level in the acyclic graph representing the genetic network controlling E.coli metabolism are controlled by internal metabolite concentrations are not included in the database iMC1010^v1. Further, the regulation of enzyme activity by internal metabolite concentrations occurs at a much faster time scale than transcription processes inside the cell. E. coli is known to double its population in about 20 minutes in rich medium while transcriptional process inside cells occur in the order of minutes to hours. Thus, it is possible that metabolism being a functionality that needs to be active whenever food is available is largely regulated without cycles at the genetic level, with feedbacks typically entering at the level of metabolites regulating enzymes to ensure efficient functioning on a faster time scale. Feedback loops or cycles are expected in systems exhibiting explicit temporal phenomena such as cell cycle or circadian rhythms. Nevertheless it would be important to explore these questions with an enlarged database representing the genetic network controlling E.coli metabolism. We further need to check the universality of the observed architecture for the genetic network controlling E. coli metabolism in other organisms as and when an integrated metabolic and regulatory network database like iMC1010^v1 becomes available for them.

In chapter 3, we have approximated the nonzero flux of an internal reaction by the activity of the gene coding for the enzyme catalyzing the reaction in the metabolic network. Hence, we have converted the metabolic feedback via internal fluxes into effective feedbacks from enzyme coding genes to other genes. This was done in order to have a simplified dynamics and a closed system of the genes alone (along with external metabolites). We remark that there exist in the literature alternative ways of treating metabolic feedback on regulation and in particular the flux variables. These include the regulatory FBA [112, 41, 65] and dynamic FBA [122, 123] in which the fluxes and the genes are dynamically coupled to each other. However, in these methods one makes an arbitrary choice of the flux vector out of many alternative flux vectors satisfying the constraints. Another method, SR-FBA [124], has been proposed that systematically accounts for multiple optimal metabolic steady states given a regulatory state. However, SR-FBA cannot be used for dynamical simulations since it only yields the various steady states for the metabolic-regulatory system. In chapter 3, our treatment of the internal fluxes is simpler compared to the above mentioned methods in that it eliminates the flux variables in favour of an effective feedback of the genes upon other genes. In the context of the present database, we believe that our broad conclusions would not change significantly because of our simplified approach to the treatment of the internal fluxes since the feedbacks from the fluxes affect only 5 genes coding for transcription factors and 16 genes downstream of these coding for enzymes, thus affecting only 21 genes out of 583. A better treatment of the feedbacks from internal metabolites than is achieved by our approach and the other approaches mentioned above requires metabolite concentrations which are difficult to compute at the present time due to the paucity of kinetic information for such large scale networks.

We end this chapter with a comment relating our results presented in chapters 3 and 4 to earlier works by Kauffman and a speculation. Kauffman [57, 59] has found biologically motivated random Boolean networks to possess multiple attractors that he has interpreted as different cell types of a multicellular organism. In chapters 3 and 4, we have studied the nature of the attractors of the genetic network regulating E. coli metabolism. Since E. coli is a prokaryote, perhaps it is not surprising that we get a much simpler picture of its genetic network exhibiting a much higher degree of order in the dynamics than the systems investigated by Kauffman. While we have shown in chapter 3 that the configuration of the genetic network regulating E. coli metabolism can settle down into different attractors, yet, unlike Kauffman, for whom different attractors were realized via different initial conditions of the genes, in our case the different attractors are realized when the control variables corresponding to external metabolites have different configurations. In chapter 3, we have shown that when the control variables are held fixed, there is none (or very little) multiplicity of attractors irrespective of the initial conditions of genes. The architecture and dynamics we have found in chapters 3 and 4 is probably quite suitable for prokaryotic lifestyle and evolution. The question remains open whether for eukaryotes and especially multicellular ones, the hypothesis that associates different cell types with different attractors of the regulatory network is valid. While that hypothesis remains an enticing possibility, it is worth noting that the simple architecture of the genetic network regulating E. coli metabolism described in detail in chapter 4 could have its uses in the eukaryotic case as well. Environmental control through the root nodes corresponding to external metabolites of the cellular attractors can itself cause a cell to differentiate into another type, the environment being determined in the multicellular case by the state of other cells in the organism. The modular structure of the the genetic network regulating E. coli metabolism at the intermediate level of the essentially acyclic graph described in chapter 4 would even permit a cell to hop through several attractors in the course of development of the organism as the environmental cues to this cell change. Such an architecture could thus contribute to developmental flexibility and, potentially, evolvability of eukaryotes as well. The multiplicity of internal attractor basins pointed out by Kauffman would be an asset in keeping the cell in its new attractor after a transient environmental cue has caused it to shift from one basin to another. It would be interesting to investigate these questions when a database similar to iMC1010^v1 becomes available for a multicellular organism.