Network reconstruction is a long and arduous process that involves building a large network in a step-by-step fashion by identifying one reaction at a time. It is foundational to the bottom-up approach to biology. A reconstruction collects all the available biochemical, genetic, and genomic (BiGG) information that is available on a cellular process of interest, and then organizes it in a formal, mathematical fashion that is consistent with the corresponding fundamental chemical and genetic properties. In this chapter we illustrate this process by looking at the familiar glycolytic pathway. Then we introduce the module-by-module nature of the network reconstruction process. We then show how detailed information about the enzymes, genetic information, and structural properties are incorporated in a reconstruction. Next, we will detail the reconstruction of the central metabolic pathways in Escherichia coli and how a knowledge base is formed from the reconstruction process. We close the chapter by discussing the features of genome-scale reconstructions and the computational models formed from them.

Many Reactions and Their Stoichiometry

Networks are composed of compounds (nodes) and reactions (links). When many reactions are known that share reactants and products, they can be graphically linked together. More and more reactions can be added to a graphical representation as one grows the scope of the network under consideration, and Figure 2.1A shows the first few steps in glycolysis as an example. Such reaction maps were formed historically as more metabolic reactions were discovered and they were eventually joined together.

Such a map can be represented mathematically (Figure 2.1B). Such mathematical representation is based on the stoichiometric coefficients that count the molecules that are consumed and produced by a biochemical reaction. All such coefficients for all the reactions in a network can be organized in a matrix format; a mathematical matrix akin to a table. This information is exact, quantitative, and forms the basis for mathematical characterization and assessment of integrative biochemical properties of a network as a whole.

Chemical reactions link cellular components together to form a network. Although we can specify the chemical properties of links in biological networks, it is the way in which a multitude of such links form networks that determines phenotypic functions. These integrated network functions are also called functional states, and they correspond to the observed biological functions or phenotypic states that networks create. A functional state may be viewed as the outcome of the execution of the genetic program written in the DNA. In this chapter we detail the concept of a functional state of a genome-scale network and how it represents a physiologically observable state. The following chapters then describe the framework for computing functional states using the constraint-based approach.

Components vs. Systems

Components come and go Biological components all have a finite turnover time. Most metabolites turn over within a minute in a cell, mRNA molecules typically have two-hour half-lives in human cells [463], 3% of the extracellular matrix in cardiac muscle is turned over daily, and so forth. So a cell that you observe today, compared with the same cell yesterday, may only contain a small fraction of the same molecules.

Similarly, cells have finite lifetimes. The cellularity of the human bone marrow turns over every two to three days. The renewal rate of skin is on the order of five days to a couple of weeks. The lining of the gut epithelium has a turnover time of about five to seven days. Slower tissues, like the liver, turnover their cellularity approximately once a year. So a mammal that you observe today may only contain a small fraction of the same cells as the same mammal observed a year ago. Thus, the components of a biological system come and go, and their turnover takes place on multiple time scales.

However, the system remains Most of the cells that are contained in an individual today were not there just a few years ago.

We now turn our attention to examining how the genome-scale reconstruction process for metabolic networks has been applied to particular organisms. The reconstruction process is enabled by genome sequencing. The first full genomic sequence appeared for Haemophilus influenzae in 1995. Four years later, the genome-scale metabolic reconstruction for H. influenzae appeared, the first of its kind [105]. The second genome-scale reconstruction to appear was for E. coli K-12 MG1655 in 2000 [106]. Due to a wealth of bibliomic data available for E. coli, there have been several subsequent updates and expansions of this reconstruction. This chapter discusses the history of the reconstruction of the genome-scale metabolic network in E. coli. We then discuss how this network reconstruction can be converted to a computational model and illustrate the range of basic scientific questions that can be addressed and describe its practical uses.

Some Basic Facts about E. coli

E. coli is a Gram-negative bacterium, a prokaryote belonging to the Enterobacteria family. It is perhaps the best-characterized organism as evidenced by its high Species Knowledge Index (SKI) [184, 295]. Non-pathogenic strains have proven to be the workhorse of the biotechnology industry. E. coli is normally harmless, but needs only to acquire a combination of a few mobile genetic elements to become a highly adapted pathogen [92]. Pathogenic strains cause a range of diseases, from gastroenteritis to extra-intestinal infections of the urinary tract, bloodstream, and central nervous system (see Figure 4.1). The worldwide burden of these diseases is high. With the low cost of sequencing, there are now many strains of E. coli that have been sequenced [251].

E. coli is a remarkable and diverse organism. It can grow on a wide variety of different carbon and nitrogen sources both aerobically and anaerobically. Fundamentally, its life cycle consists of five key shifts in growth conditions: a heat shock (i.e., from ambient temperature to 37°C), a pH shock (from neutral to gastric pH of about 2), oxygen deprivation, nutritional richness and community competition, and finally a cold shock back to ambient temperature.

It is becoming increasingly clear that making non-trivial predictions of phenotypic functions from genomic information requires a thorough understanding of the molecular processes involved in the expression of genes into functional proteins, the organization of these proteins with lipids and nucleic acids in spatiotemporal structures, and the assessment of the role of such structures in carrying out metabolic processes, creating ion currents involved in membrane excitability, maintaining electro-neutrality and osmotic homeostasis, and performing other cellular functions. This statement represents a grand challenge for biological and medical research for the twenty-first century. The COBRA methods covered in this text can address some of these issues. They represent a first step in addressing this grand challenge.

The Brief History of COBRA

Constraint-based reconstruction and analysis has been applied to biochemical reaction networks for over two decades. Research articles that utilize COBRA methods for interpreting and predicting biological phenotypes appearing from 1986 to mid-2013 have been collected and the cumulative rate of progress in the field can be graphed (Figure 29.1). The analysis of this published literature shows that the history of deployment of COBRA methods can be divided into four phases [52].

Initial studies (1986–1998) Early interest in constraint-based models was directed towards the determination of theoretical pathway yields and metabolite overflows [121, 258]. Then, experimental metabolic fluxes and growth rates were found to be consistent with computation based on optimization of physiologically meaningful objective functions, including minimal production of reactive oxygen species for hybridoma cells [373] and maximal growth rate for laboratory strains of E. coli [438]. The quantitative match between model predictions and measured cellular behavior opened up the possibility of predicting cellular phenotypes from a biochemically reconstructed network.

A reconstructed network can have many functional states even in the same environment. Constraint-based optimization problems can have multiple equivalent solutions that yield the same numerical value of the objective that is being optimized. This issue is an important one in characterizing networks and is a common occurrence in large-scale networks. In this chapter, we will cover the methods that have been developed to study what are called alternative optimal solutions (AOS). We demonstrate them by applying them to the core E. coli model, and discuss illustrative studies and issues that arise with finding AOS at the genome-scale.

Equivalent Ways to Reach a Network Objective

An historical note Even for explorations of network properties, such as those illustrated in the last chapter, alternative optimal solutions with the same optimal value of the objective function are found. This issue arose in the early days of the development of FBA approaches when two solutions for the formation of E4P were discovered (see Figure 20.1). The optimization that was being performed was to compute the maximal yield of the biosynthetic precursor E4P from glucose as a substrate. Two different solutions were computed that gave a yield of 1.33 molecules E4P per molecule glucose, or an 88.7% carbon conversion. The two solutions differ in the utilization of the TCA cycle or the glyoxalate shunt to achieve the same value of the objective.

The existence of AOS is a fundamental issue in COBRA methods and understanding their biological implications and relevance is important. Thus, a fair amount of effort has been devoted towards computing and studying AOS.

Characterizing alternative optimal states As outlined in the last chapter, LP can be used to find single optimal solutions. The solution to this optimization problem may not be unique, as illustrated in Figure 20.1, and AOS may exist. Various methods have been developed to characterize AOS.

The stoichiometric matrix is a connectivity matrix. Elementary topological properties of the network it represents can be computed directly from the individual elements of S. Direct topological studies are interesting from a variety of standpoints. They focus on relatively easy to understand and intuitive properties of the structure of the network. Elementary topological properties relate to how connected a network is, and how its components participate in forming the connectivity properties of the network. There may be many functional states for a given network structure (see Chapter 16). Topological properties are thus global and less specific than functional states of networks. Some of the differences between functional states and network topology are covered in Part III.

The Binary Form of S

The elementary topological properties are determined based on the non-zero elements in the stoichiometric matrix. Thus, we define the elements of a new matrix ŝ as

ŝij = 0 if sij = 0

ŝij = 1 if sij ≠ 0

which is the binary form of S. This matrix is composed of only zeros and ones. If ŝij is unity, it means that compound i participates in reaction j. Note that in the rare case where a homodimer is formed, i.e., in a reaction of the type 2A → A2, the stoichiometric coefficient of two becomes unity in the binary form of S.

S is a sparse matrix A number of genome-scale stoichiometric matrices have been reconstructed (see [2]). As there are typically only a handful of compounds that participate in a reaction out of hundreds of compounds participating in a network, the stoichiometric matrix is sparse. A sparse matrix is mostly composed of zero elements. For instance, if there are on average 3 compounds that participate in a reaction, but there are m compounds in the network, then the fraction of non-zero elements in the matrix is 3/m. If m is 300, then only 1% of the elements are non-zero and the matrix is sparse.

Pathways as basis vectors for the null space are useful for studying the capabilities of a network and to determine network properties. An alternative approach to characterizing the contents of solution spaces is uniform random sampling. This approach involves obtaining a statistically meaningful number of solutions uniformly distributed throughout the entire solution space and then studying their properties. Randomized sampling of candidate states throughout an entire solution space gives an unbiased assessment of its properties and can result in significant biological insights and understanding.

The Basics

A simple flux split Uniform random sampling can be illustrated by looking at a simple flux split (see Figure 14.1). The null space can be shown in two dimensions using its two conical basis vectors:

b1 = (1, 0, 1) and b2 = (0, 1, 1)

Note that the maximum weight that can be placed on these vectors is α1,max = 6 and α2,max = 8, whereas the minimum values in both cases is zero. There is an additional constraint association with reaction 3 that states that

α1,max + α2,max≤ α3,max = 14

thus both α1,max and α2,max cannot be attained simultaneously.

The null space can be enclosed with a parallelepiped by ignoring the constraint in Equation (14.2), which is a parallelogram in two-dimensions.

Ask not what you can do for your reconstruction, but what your reconstruction can do for you – with apologies to JFK

Part I described the network reconstruction process, Part II detailed its conversion into a mathematical format, and Part III discussed the characterization of network properties using constraint-based methods. What is all this effort good for? A well-curated genome-scale model, in principle, can be used to study all the phenotypic functions that can be produced from the genome of the target organism. Thus, the possible range of applications of GEMs is broad.

In June 2013, 645 papers had appeared that used COBRA tools to explain existing data or predict biological functions [52]. Analysis of these publications showed that the history of GEMs and their uses can be divided into three phases. Shortly after the appearance of the first GEMs in 1999 and 2000, there was a period of creativity around algorithms and analysis method development. In the mid-2000s, experimental validation studies began to accumulate resulting in the availability of well-curated and validated GEMs. Around 2010, a series of studies begun to appear that demonstrated that a variety of predictions of biological functions could be made using GEMs.

In this part of the book we describe the use of genome-scale networks. Chapters 22 and 23 describe how environmental and genetic parameters are represented and studied with genome-scale models. GEMs have proven to be remarkably useful for studying these two types of parameters. Media composition and growth requirements can be studied productively using GEMs. The prediction of the outcome of phenotypic screens of KO strain collections using GEMs with tens of thousands of outcomes represent perhaps the largest-scale effort for predicting biological functions.

Then we move onto specific application areas in a series of four chapters. First, we discuss how GEMs give a useful context for the analysis of omics data sets. Mapping such data against a known background proves to increase the resolution and use of the data.

Since the late 1990s, there has been an explosion in the development of technologies that measure cellular content on a genome scale. These methods generate large amounts of data, generally referred to as omic data; sometimes referred to as content. The quality and coverage of omic data sets has improved steadily with time. Individual data points in omic data sets are treated as independent variables. They can then be correlated statistically to find patterns in the presence of cellular components. Such statistical correlations do not confer causality. GEMs, given their mechanism-based construction, can be used as a context for the analysis of polyomic data sets [358]. Such contextualization can lead to the establishment of causation. Furthermore, GEMs can integrate mechanistically multiple omic data sets, thus aiding in the determination of how various components come together to produce cellular functions and phenotypic states. For some, this big data to knowledge conversion represents a grand challenge in biology. Here we will cover the basic principles of omic data analysis using GEMs. More detailed reviews have appeared describing these data-mapping methods and their use [41, 175, 209, 232, 253].

Context for Content

Over the past 10–15 years, many ingenious methods have been developed to profile the molecular content of cells and to determine the interactions between components. There are vast repositories of the resulting omic data available on various websites, and there are many sources that describe the omics data types, how they are generated, and their availability. We will not repeat these here, but instead turn our attention to the use of COBRA methods to analyze omics data sets.

Networks as a backdrop or context for omics data-mapping Two main network approaches are used to extract biological insights from omic data sets [348]: inference based and knowledge-based. Both approaches use an interconnected network of biological components to interpret omic data sets.

Some 60 years ago, the promise of molecular biology held that if we knew and understood the function of the molecules that comprise cells, then we could understand cells and their functions. Although this was true in principle (and in practice in a few cases), the sheer number of molecules made it very difficult to comprehend so many simultaneous functions. The simultaneous measurement of the majority of these molecules became possible over the last 10–15 years through the development of many ingenious technologies. As a result, we now have a growing number of data sets that give us the composition of particular cells and organisms under certain conditions. The chemical interactions between many of these components are now known and this knowledge gives rise to reconstructed biochemical reaction networks on a genome-scale that underlie various cellular functions. Thus, enter (molecular) systems biology.

Systems biology is not necessarily focused on the components themselves, but on the nature of the links that connect them and on the functional states of the biochemical networks that result from the collection of all such links. These functional states of networks correspond to observable physiological or homeostatic states. Completing the relationship between all the chemical components of a cell, with their genetic bases, and its physiological functions is the promise of (molecular) systems biology. This undertaking represents the de facto construction of a mechanistic genotype– phenotype relationship.

The Genotype–Phenotype Relationship

The concept Through breeding experiments, Gregor Mendel discovered that there are discrete quanta of information passed from one generation to the next that determine the form and function of an organism. These quanta, or packets, of information are now generally referred to as genes. The collection of all the genes and the particular version of them found in a genome of an individual organism is referred to as its genotype. The form and function of an organism is referred to as its phenotype. How the phenotype is related to the genotype represents the fundamental relationship of biology.

The E. coli metabolic reconstruction and its GEM are highly developed. There are now a number of other microorganisms for which GEMs have been built and put to use. A significant feature of microbial life and biology revolves around metabolism and growth, and thus GEMs of microbial metabolism have been used to address a number of important issues in industrial, environmental, and medical microbiology. We will describe some of these accomplishments in this chapter. It is likely that the number of such reconstructions and their applications will grow significantly over the coming years given the rapidly growing amount of sequenced genomes and semi-automated reconstruction tools.

State of The Field

Phylogenetic coverage Bacteria display an astonishing spectrum of metabolic capabilities that reflect the wide spectrum of microenvironmental niches in which they grow. Metabolism in some of the branches of the phylogenetic tree, such as enterobacteria, is generally well characterized, whereas little molecular detail is known about members of other branches. In fact, a large number of microorganisms cannot even be cultured in the laboratory, but their genomes can be sequenced and a metabolic reconstruction can be performed. Interestingly, metabolic network reconstructions can then be used to formulate growth media [441].

The phylogenetic coverage of existing genome-scale metabolic reconstructions is shown in Figure 5.1A. It shows that the coverage of the tree is biased by well-understood model species. There is clearly a need to get better uniform coverage of the phylogenetic tree so that we can have a broad view of the metabolic capabilities of microbes that inhabit the various microenvironments on this planet.

The number of manually curated reconstructions has grown steadily since the publication of the first GEM in 1999 (Figure 5.1B). The iterative nature of the process also shows up in this summary, as metabolic reconstructions for some organisms have gone through multiple iterations (see Figure 3.11 for more details).

Genetic parameters can be represented explicitly in genome-scale models enabling the assessment of their effects on phenotypic functions. There are three issues associated with genetic parameters discussed in this chapter. First, the analysis of the consequences of the loss of function (LOF) or gain of gene function (GOF) will be explored. For essential genes, it is important to understand how the loss of a gene product may lead to cell death. Further, sophisticated methods have been developed to assess the minimal perturbation of phenotypic states resulting from non-lethal gene knock-outs. Second, the simultaneous removal of two genes allows one to address issues related to synthetic lethality and epistasis in general. Such understanding has implications for many studies, such as finding drug targets. Third, sequence variations in genes or their gene dosage, i.e., gene copy number, can determine their activity level. Gene dosage may be of interest in many situations, ranging from imprinted genes in humans (where one copy is active), to aneuploidy in cancer, to cell line engineering (where gene dosage can be manipulated). The coordinated functions of genes can be assessed through reaction co-sets, and genetic variations in a co-set may have similar consequences on network functions, and thus phenotypic states. The following contains some examples of how to describe genetic parameters within GEMs; there will undoubtedly be many more applications developed in the future.

Single Gene Knock-outs

23.1.1 Concept

The network function of a gene is traced through its gene–protein–reaction (GPR) association to the reactions with which the gene product is involved. Thus, if a gene is deleted, then flux constraints for the reaction(s) affected by the loss of a gene product are put to zero and thus no flux can take place through the corresponding reaction.

With the publication of the first full genome sequence in the mid 1990s, it became possible, in principle, to identify all the gene products involved in complex biological processes in a single organism. In practice, almost 17 years later, this has proven difficult to accomplish using just sequence information. A toolbox of molecular biology methods that are implemented on a genome-scale are now available and it allows us to measure many properties of a genome with unprecedented resolution. A workflow that systematically integrates these data types can lead to the definition of the transcription unit (TU) architecture of a genome. This information, in turn, enables the reconstruction of transcriptional regulatory networks (TRNs) and other features. The full application and integrative analysis of genome-wide measurements led to the concept of a metastructure of a prokaryotic genome as there are many more attributes to a genome than its sequence and 3D arrangement.

The Concept of a Metastructure

With the plethora of genome-scale measurement methods, genomes can be characterized at multiple different organizational levels. Genomes used to be thought of in terms of their base sequence and three-dimensional structure. It is now clear that there are many more dimensions to genome organization, use, and information content (Figure 8.1). This realization has led to the concept of a metastructure of a bacterial genome [74].

Detailed examination of genome-wide data sets results in the definition of structural, operational, and functional annotations [74, 349]. Structural genome annotation provides the foundation for further operational and functional annotation and consists of coding (open reading frames, or ORFs) and non-coding genes, as well as intergenic regions. Elucidating the precise structural genome annotation subsequently allows the decoding of the operational genome annotation, which consists of operons and TUs. As a higher level of genome organization, the operon structure is a key to deciphering the flow of information encoded in the genome. A functional genome annotation assigns a function to an ORF and describes the biochemical properties of the gene products.

Organisms used for bioprocessing are engineered to achieve the production of the desired chemical compound. Such compounds can be natural or non-natural. In the former case, the existing gene portfolio of the organism is altered in its function to change the production phenotype of the cell. In the latter case, a series of heterologous genes are introduced to confer the capability on the host to carry out a new biochemical function. Once such a function is established, the host functions are modified to achieve an optimal function of the pathway incorporated. Metabolic engineering and microbial factory design represent the most ambitious efforts at directed and systematic construction of a new phenotypic function. Over the past decade, GEMs have played an increasing role in this process [207, 223, 224, 328].

Historical Background

Many compounds can be made biologically The chemical and pharmaceutical industries produce products with sales in excess of $3.5T per year, representing one of the world's largest industries. A broad variety of chemical compounds are manufactured in large quantities. Most of these products are derived currently from petroleum as a feedstock. With technological advances in the life sciences, it has become clear in recent years that a fair fraction of these compounds can be made through biological means (Figure 27.1). Given the significant economic impetus for producing chemicals from sustainable feedstocks, major effort is going into designing cells to produce industrially, commercially, and pharmaceutically valuable compounds.

Cell factories The notion of a cell factory is illustrated at the center of Figure 27.1. Substates are fed to an engineered organism that can produce and secrete chemical compounds that it does not produce naturally. Substrates are typically sugars (both five- and six-carbon sugars) derived from renewable biomass, although there is an increasing interest in the use of one-carbon substrates (methane, carbon dioxide, and synthetic gas).