Metabolic capabilities of Escherichia coli: I. synthesis of biosynthetic precursors and cofactors

Metabolism of living cells converts substrates into metabolic energy, redox potential and metabolic end products that are essential to maintain cellular function. The flux distribution among the various biochemical pathways is determined by the kinetic properties of enzymes which are subject to strict regulatory control. Simulation of metabolic behavior therefore requires the complete knowledge of biochemical pathways, enzyme kinetics as well as their regulation. Unfortunately, complete kinetic and regulatory information is not available for microbial cells, thus preventing accurate dynamic simulation of their metabolic behavior. However, it is possible to define wider limits on metabolic behavior based solely on flux balances of biochemical pathways. We present here comprehensive information about the catabolic pathways of the bacterium Escherichia coli. Using this biochemical database, we formulate a stoichiometric model of the bacterial network of fueling reactions. After logical structural reduction, the network consists of 53 metabolic fluxes and 30 metabolites. The solution space of this under-determined system of equations presents the bounds of metabolic flux distribution that the bacterial cell can achieve. We use specific objective functions and linear optimization to investigate the capability of E. coli catabolism to maximally produce the 12 biosynthetic precursors and three key cofactors within this solution space. For the three cofactors, the maximum yields are calculated to be 18.67 ATP, 11.6 NADH and 11 NADPH per glucose molecule, respectively. The yields of NADH and NADPH are less than 12 owing to the energy costs of importing glucose. These constraints are made explicit by the interpretation of shadow prices. The optimal yields of the 12 biosynthetic precursors are computed. Four of the 12 precursors (3-phosphoglycerate, phosphoenolpyruvate, pyruvate and oxaloacetate) can be made by E. coli with complete carbon conversion. Conversely, none of the sugar monophosphates can be made with 100% carbon conversion and analysis of the shadow prices reveals that this conversion is constrained by the energy cost of importing glucose. Three of the 12 precursors (acetyl-coA, α-ketoglutarate, and succinyl-coA) cannot be made with full carbon conversion owing to stoichiometric constraints; there is no route to these compounds without carrying out a decarboxylation reaction. Metabolite flux balances and linear optimization have thus been used to determine the catabolic capabilities of E. coli .     

High-throughput metabolic flux analysis based on gas chromatography-mass spectrometry derived 13C constraints

13C-constrained flux balancing analysis based on gas chromatography-mass spectrometry data is presented here as a simple and robust method for the estimation of intracellular carbon fluxes. In this approach, the underdetermined system of metabolite balances deduced from stoichiometric relations and measured extracellular rates is complemented with 13C constraints from metabolic flux ratio analysis. Fluxes in central carbon metabolism of exponentially growing Escherichia coli were estimated by 13C-constrained flux balancing from three different 13C-labeled glucose experiments. The best resolution of the network was achieved using 13C constraints derived from [U-13C]glucose and [1-13C]glucose experiments. The corresponding flux estimate was in excellent agreement with a solution that was independently obtained with a comprehensive isotopomer model. This new methodology was also demonstrated to faithfully capture the intracellular flux distribution in E. coli shake flasks and 1-ml deep-well microtiter plates. Due to its simplicity, speed, and robustness, 13C-constrained metabolic flux balancing is promising for routine and high-throughput analysis on a miniaturized scale.     

Modeling with a view to target identification in metabolic engineering: A critical evaluation of the available tools

The state of the art tools for modeling metabolism, typically used in the domain of metabolic engineering, were reviewed. The tools considered are stoichiometric network analysis (elementary modes and extreme pathways), stoichiometric modeling (metabolic flux analysis, flux balance analysis, and carbon modeling), mechanistic and approximative modeling, cybernetic modeling, and multivariate statistics. In the context of metabolic engineering, one should be aware that the usefulness of these tools to optimize microbial metabolism for overproducing a target compound depends predominantly on the characteristic properties of that compound. Because of their shortcomings not all tools are suitable for every kind of optimization; issues like the dependence of the target compound's synthesis on severe (redox) constraints, the characteristics of its formation pathway, and the achievable/desired flux towards the target compound should play a role when choosing the optimization strategy. © 2010 American Institute of Chemical Engineers Biotechnol. Prog., 2010     

Optimal stoichiometric designs of ATP-producing systems as determined by an evolutionary algorithm.

The design of metabolic pathways is thought to be the result of an optimization process such that the structure of contemporary metabolic routes maximizes a particular objective function. Recently, it has been shown that some essential stoichiometric properties of glycolysis can be explained on the basis of the requirement for a high ATP production rate. Because the number of stoichiometrically feasible designs increases strongly with the number of reactions involved, a systematic analysis of all the possibilities turns out to be inaccessible beyond a certain system size. We present, therefore, an alternative approach to compute in a more efficient way the optimal design of glycolysis interacting with an external ATP-consuming reaction. The algorithm is based on the laws of evolution by natural selection, and may be viewed as a particular version of evolutionary algorithms. The following conclusions are derived: (a) evolutionary algorithms are very useful search strategies in determining optimal stoichiometries of metabolic pathways. (b) Essential topological features of the glycolytic network may be explained on the basis of flux optimization. (c) There is a strong interrelation between the optimal stoichiometries and the thermodynamic and kinetic properties of the participating reactions. (d) Some subsequences of reactions in optimal pathways are strongly conserved at variation of system parameters, which may be understood by applying principles of metabolic control analysis.     

Profiling of dynamic changes in hypermetabolic livers

The liver plays an important role in the overall negative nitrogen balance leading to muscle wasting commonly observed in patients following many conditions, including severe injury, cancer, and diabetes. In order to study changes in liver metabolism during the establishment of such catabolic states, we used a rat skin burn injury model that induces hypermetabolism and muscle wasting. At various times during the first week following the injury, livers were isolated and perfused in a recirculating system under well-defined conditions. We applied a steady-state metabolic flux analysis model of liver metabolism and then used k-means clustering to objectively group together reaction flux time profiles. We identified six distinct groups of reactions that were differentially responsive: (1) pentose phosphate pathway (PPP); (2) amino acid oxidation reactions leading to the formation of tricarboxylic acid (TCA) cycle intermediates; (3) gluconeogenesis; (4) TCA-cycle and mitochondrial oxidation; (5) lipolysis, beta-oxidation, and ketone body formation; and (6) urea-cycle. Burn injury sequentially upregulated the urea-cycle, the PPP, and the TCA-cycle, in order, while beta-oxidation and gluconeogenesis remained unchanged. The upregulation of the PPP was transient, whereas the rise in urea- and TCA-cycle fluxes was sustained. An ATP balance predicted an increased production of ATP and energy expenditure starting on day 3 post-burn, which correlated with the induction of the oxidative phosphorylation uncoupler uncoupling protein-2. We conclude that metabolic profiling using flux analysis and clustering analysis is a useful methodology to characterize the differential activation of metabolic pathways in perfused organs and to identify specific key pathways that are sensitive to a stimulus or insult without making a priori assumptions.     

Interplay between Constraints,Objectives, and Optimality for Genome-Scale Stoichiometric Models

High-throughput data generation and genome-scale stoichiometric models have greatly facilitated the comprehensive study of metabolic networks. The computation of all feasible metabolic routes with these models, given stoichiometric, thermodynamic, and steady-state constraints, provides important insights into the metabolic capacities of a cell. How the feasible metabolic routes emerge from the interplay between flux constraints, optimality objectives, and the entire metabolic network of a cell is, however, only partially understood. We show how optimal metabolic routes, resulting from flux balance analysis computations, arise out of elementary flux modes, constraints, and optimization objectives. We illustrate our findings with a genome-scale stoichiometric model of Escherichia coli metabolism. In the case of one flux constraint, all feasible optimal flux routes can be derived from elementary flux modes alone. We found up to 120 million of such optimal elementary flux modes. We introduce a new computational method to compute the corner points of the optimal solution space fast and efficiently. Optimal flux routes no longer depend exclusively on elementary flux modes when we impose additional constraints; new optimal metabolic routes arise out of combinations of elementary flux modes. The solution space of feasible metabolic routes shrinks enormously when additional objectives---e.g. those related to pathway expression costs or pathway length---are introduced. In many cases, only a single metabolic route remains that is both feasible and optimal. This paper contributes to reaching a complete topological understanding of the metabolic capacity of organisms in terms of metabolic flux routes, one that is most natural to biochemists and biotechnologists studying and engineering metabolism.     

Metabolic flux estimation--a self-adaptive evolutionary algorithm with singular value decomposition

Metabolic flux analysis is important for metabolic system regulation and intracellular pathway identification. A popular approach for intracellular flux estimation involves using ¹³C tracer experiments to label states that can be measured by nuclear magnetic resonance spectrometry or gas chromatography mass spectrometry. However, the bilinear balance equations derived from ¹³C tracer experiments and the noisy measurements require a nonlinear optimization approach to obtain the optimal solution. In this paper, the flux quantification problem is formulated as an error-minimization problem with equality and inequality constraints through the ¹³C balance and stoichiometric equations. The stoichiometric constraints are transformed to a null space by singular value decomposition. Self-adaptive evolutionary algorithms are then introduced for flux quantification. The performance of the evolutionary algorithm is compared with ordinary least squares estimation by the simulation of the central pentose phosphate pathway. The proposed algorithm is also applied to the central metabolism of Corynebacterium glutamicum under lysine-producing conditions. A comparison between the results from the proposed algorithm and data from the literature is given. The complexity of a metabolic system with bidirectional reactions is also investigated by analyzing the fluctuations in the flux estimates when available measurements are varied.     

Network analysis of intermediary metabolism using linear optimization. I. Development of mathematical formalism.

Analysis of metabolic networks using linear optimization theory allows one to quantify and understand the limitations imposed on the cell by its metabolic stoichiometry, and to understand how the flux through each pathway influences the overall behavior of metabolism. A stoichiometric matrix accounting for the major pathways involved in energy and mass transformations in the cell was used in our analysis. The auxiliary parameters of linear optimization, the so-called shadow prices, identify the intermediates and cofactors that cause the growth to be limited on each nutrient. This formalism was used to examine how well the cell balances its needs for carbon, nitrogen, and energy during growth on different substrates. The relative values of glucose and glutamine as nutrients were compared by varying the ratio of rates of glucose to glutamine uptakes, and calculating the maximum growth rate. The optimum value of this ratio is between 2-7, similar to experimentally observed ratios. The theoretical maximum growth rate was calculated for growth on each amino acid, and the amino acids catabolized directly to glutamate were found to be the optimal nutrients. The importance of each reaction in the network can be examined both by selectively limiting the flux through the reaction, and by the value of the reduced cost for that reaction. Some reactions, such as malic enzyme and glutamate dehydrogenase, may be inhibited or deleted with little or no adverse effect on the calculated cell growth rate.     

13C tracer experiments and metabolite balancing for metabolic flux analysis: comparing two approaches

Conventional metabolic flux analysis uses the information gained from determination of measurable fluxes and a steady-state assumption for intracellular metabolites to calculate the metabolic fluxes in a given metabolic network. The determination of intracellular fluxes depends heavily on the correctness of the assumed stoichiometry including the presence of all reactions with a noticeable impact on the model metabolite balances. Determination of fluxes in complex metabolic networks often requires the inclusion of NADH and NADPH balances, which are subject to controversial debate. Transhydrogenation reactions that transfer reduction equivalents from NADH to NADPH or vice versa can usually not be included in the stoichiometric model, because they result in singularities in the stoichiometric matrix. However, it is the NADPH balance that, to a large extent, determines the calculated flux through the pentose phosphate pathway. Hence, wrong assumptions on the presence or activity of transhydrogenation reactions will result in wrong estimations of the intracellular flux distribution. Using 13C tracer experiments and NMR analysis, flux analysis can be performed on the basis of only well established stoichiometric equations and measurements of the labeling state of intracellular metabolites. Neither NADH/NADPH balancing nor assumptions on energy yields need to be included to determine the intracellular fluxes. Because metabolite balancing methods and the use of 13C labeling measurements are two different approaches to the determination of intracellular fluxes, both methods can be used to verify each other or to discuss the origin and significance of deviations in the results. Flux analysis based entirely on metabolite balancing and flux analysis, including labeling information, have been performed independently for a wild-type strain of Aspergillus oryzae producing alpha-amylase. Two different nitrogen sources, NH4+ and NO3-, have been used to investigate the influence of the NADPH requirements on the intracellular flux distribution. The two different approaches to the calculation of fluxes are compared and deviations in the results are discussed. Copyright 1998 John Wiley & Sons, Inc.     

Investigation of metabolic objectives in cultured hepatocytes

Using optimization based methods to predict fluxes in metabolic flux balance models has been a successful approach for some microorganisms, enabling construction of in silico models and even inference of some regulatory motifs. However, this success has not been translated to mammalian cells. The lack of knowledge about metabolic objectives in mammalian cells is a major obstacle that prevents utilization of various metabolic engineering tools and methods for tissue engineering and biomedical purposes. In this work, we investigate and identify possible metabolic objectives for hepatocytes cultured in vitro. To achieve this goal, we present a special data-mining procedure for identifying metabolic objective functions in mammalian cells. This multi-level optimization based algorithm enables identifying the major fluxes in the metabolic objective from MFA data in the absence of information about critical active constraints of the system. Further, once the objective is determined, active flux constraints can also be identified and analyzed. This information can be potentially used in a predictive manner to improve cell culture results or clinical metabolic outcomes. As a result of the application of this method, it was found that in vitro cultured hepatocytes maximize oxygen uptake, coupling of urea and TCA cycles, and synthesis of serine and urea. Selection of these fluxes as the metabolic objective enables accurate prediction of the flux distribution in the system given a limited amount of flux data; thus presenting a workable in silico model for cultured hepatocytes. It is observed that an overall homeostasis picture is also emergent in the findings.     

Identification of flux control in metabolic networks using non-equilibrium thermodynamics

A method is presented to identify flux controlling reactions in metabolic networks using experimentally determined flux distributions. The method is based on the application of Ziegler's principle for the maximization of entropy production. According to this principle a metabolic network tends to maximize the entropy production rate while satisfying mass balances and maximal rate constraints. Experimental flux data corresponding to four different metabolic states of Saccharomyces cerevisiae were used to identify the corresponding flux controlling reactions. The bottleneck nature of several of the identified reactions was confirmed by earlier studies on over-expression of the identified target genes. The method also explains the failure of all the previous trials of increasing the glycolysis rate by direct over-expression of several glycolytic enzymes. These findings point to a wider use of the method for identification of novel targets for metabolic engineering of microorganisms used for sustainable production of fuels and chemicals.     

Metabolic network reconstruction and flux variability analysis of storage synthesis in developing oilseed rape (Brassica napus L.) embryos

Computational simulation of large‐scale biochemical networks can be used to analyze and predict the metabolic behavior of an organism, such as a developing seed. Based on the biochemical literature, pathways databases and decision rules defining reaction directionality we reconstructed bna572, a stoichiometric metabolic network model representing Brassica napus seed storage metabolism. In the highly compartmentalized network about 25% of the 572 reactions are transport reactions interconnecting nine subcellular compartments and the environment. According to known physiological capabilities of developing B. napus embryos, four nutritional conditions were defined to simulate heterotrophy or photoheterotrophy, each in combination with the availability of inorganic nitrogen (ammonia, nitrate) or amino acids as nitrogen sources. Based on mathematical linear optimization the optimal solution space was comprehensively explored by flux variability analysis, thereby identifying for each reaction the range of flux values allowable under optimality. The range and variability of flux values was then categorized into flux variability types. Across the four nutritional conditions, approximately 13% of the reactions have variable flux values and 10–11% are substitutable (can be inactive), both indicating metabolic redundancy given, for example, by isoenzymes, subcellular compartmentalization or the presence of alternative pathways. About one‐third of the reactions are never used and are associated with pathways that are suboptimal for storage synthesis. Fifty‐seven reactions change flux variability type among the different nutritional conditions, indicating their function in metabolic adjustments. This predictive modeling framework allows analysis and quantitative exploration of storage metabolism of a developing B. napus oilseed.     

Pathway analysis of liver metabolism under stressed condition

Pathway analysis is a useful tool which reveals important metabolic network properties. However, the big challenge is to propose an objective function for estimating active pathways, which represent the actual state of network. In order to provide weight values for all possible pathways within the metabolic network, this study presents different approaches, considering the structural and physiological properties of the metabolic network, aiming at a unique decomposition of the flux vector into pathways. These methods were used to analyze the hepatic metabolism considering available data sets obtained from the perfused livers of fasted rats receiving burn injury. Utilizing unique decomposition techniques and different fluxes revealed that higher weights were always attributed to short pathways. Specific pathways, including pyruvate, glutamate and oxaloacetate pools, and urea production from arginine, were found to be important or essential in all methods and experimental conditions. Moreover the pathways, including serine production from glycine and conversion between acetoacetate and B—OH-butyrate, were assigned higher weights. Pathway analysis was also used to identify the main sources for the production of certain products in the hepatic metabolic network to gain a better understanding of the effects of burn injury on liver metabolism.     

A metabolic network analysis & NMR experiment design tool with user interface-driven model construction for depth-first search analysis

A Windows program for metabolic engineering analysis and experimental design has been developed. A graphical user interface enables the pictorial, "on-screen" construction of a metabolic network. Once a model is composed, balance equations are automatically generated. Model construction, modification and information exchange between different users is thus considerably simplified. For a given model, the program can then be used to predict all the extreme point flux distributions that optimize an objective function while satisfying balances and constraints by using a depth-first search strategy. One can also find the minimum reaction set that satisfies different conditions. Based on the identified flux distributions or linear combinations, the user can simulate the NMR and GC/MS spectra of selected signal molecules. Alternately, spectra vectorization allows for the automated optimization of labeling experiments that are intended to distinguish between different, yet plausible flux extreme point distributions. The example provided entails predicting the flux distributions associated with deleting pyruvate kinase and designing 13C NMR experiments that can maximally discriminate between the flux distributions.     

Reconstruction and flux analysis of coupling between metabolic pathways of astrocytes and neurons: application to cerebral hypoxia

Background

It is a daunting task to identify all the metabolic pathways of brain energy metabolism and develop a dynamic simulation environment that will cover a time scale ranging from seconds to hours. To simplify this task and make it more practicable, we undertook stoichiometric modeling of brain energy metabolism with the major aim of including the main interacting pathways in and between astrocytes and neurons.

Model

The constructed model includes central metabolism (glycolysis, pentose phosphate pathway, TCA cycle), lipid metabolism, reactive oxygen species (ROS) detoxification, amino acid metabolism (synthesis and catabolism), the well-known glutamate-glutamine cycle, other coupling reactions between astrocytes and neurons, and neurotransmitter metabolism. This is, to our knowledge, the most comprehensive attempt at stoichiometric modeling of brain metabolism to date in terms of its coverage of a wide range of metabolic pathways. We then attempted to model the basal physiological behaviour and hypoxic behaviour of the brain cells where astrocytes and neurons are tightly coupled.

Results

The reconstructed stoichiometric reaction model included 217 reactions (184 internal, 33 exchange) and 216 metabolites (183 internal, 33 external) distributed in and between astrocytes and neurons. Flux balance analysis (FBA) techniques were applied to the reconstructed model to elucidate the underlying cellular principles of neuron-astrocyte coupling. Simulation of resting conditions under the constraints of maximization of glutamate/glutamine/GABA cycle fluxes between the two cell types with subsequent minimization of Euclidean norm of fluxes resulted in a flux distribution in accordance with literature-based findings. As a further validation of our model, the effect of oxygen deprivation (hypoxia) on fluxes was simulated using an FBA-derivative approach, known as minimization of metabolic adjustment (MOMA). The results show the power of the constructed model to simulate disease behaviour on the flux level, and its potential to analyze cellular metabolic behaviour in silico.

Conclusion

The predictive power of the constructed model for the key flux distributions, especially central carbon metabolism and glutamate-glutamine cycle fluxes, and its application to hypoxia is promising. The resultant acceptable predictions strengthen the power of such stoichiometric models in the analysis of mammalian cell metabolism.     

Metabolite-balancing techniques vs. 13C tracer experiments to determine metabolic fluxes in hybridoma cells

The estimation of intracellular fluxes of mammalian cells using only mass balances of the relevant metabolites is not possible because the set of linear equations defined by these mass balances is underdetermined. In order to quantify fluxes in cyclic pathways the mass balance equations can be complemented with several constraints: (1) the mass balances of co-metabolites, such as ATP or NAD(P)H, (2) linear objective functions, (3) flux data obtained by isotopic-tracer experiments. Here, these three methods are compared for the analysis of fluxes in the primary metabolism of continuously cultured hybridoma cells. The significance of different theoretical constraints and different objective functions is discussed after comparing their resulting flux distributions to the fluxes determined using 13CO2 and 13C-lactate measurements of 1 - 13C-glucose-fed hybridoma cells. Metabolic fluxes estimated using the objective functions "maximize ATP" and "maximize NADH" are relatively similar to the experimentally determined fluxes. This is consistent with the observation that cancer cells, such as hybridomas, are metabolically hyperactive, and produce ATP and NADH regardless of the need for these cofactors.     

Thermodynamics-based metabolic flux analysis

A new form of metabolic flux analysis (MFA) called thermodynamics-based metabolic flux analysis (TMFA) is introduced with the capability of generating thermodynamically feasible flux and metabolite activity profiles on a genome scale. TMFA involves the use of a set of linear thermodynamic constraints in addition to the mass balance constraints typically used in MFA. TMFA produces flux distributions that do not contain any thermodynamically infeasible reactions or pathways, and it provides information about the free energy change of reactions and the range of metabolite activities in addition to reaction fluxes. TMFA is applied to study the thermodynamically feasible ranges for the fluxes and the Gibbs free energy change, Delta(r)G', of the reactions and the activities of the metabolites in the genome-scale metabolic model of Escherichia coli developed by Palsson and co-workers. In the TMFA of the genome scale model, the metabolite activities and reaction Delta(r)G' are able to achieve a wide range of values at optimal growth. The reaction dihydroorotase is identified as a possible thermodynamic bottleneck in E. coli metabolism with a Delta(r)G' constrained close to zero while numerous reactions are identified throughout metabolism for which Delta(r)G' is always highly negative regardless of metabolite concentrations. As it has been proposed previously, these reactions with exclusively negative Delta(r)G' might be candidates for cell regulation, and we find that a significant number of these reactions appear to be the first steps in the linear portion of numerous biosynthesis pathways. The thermodynamically feasible ranges for the concentration ratios ATP/ADP, NAD(P)/NAD(P)H, and H(extracellular)(+)/H(intracellular)(+) are also determined and found to encompass the values observed experimentally in every case. Further, we find that the NAD/NADH and NADP/NADPH ratios maintained in the cell are close to the minimum feasible ratio and maximum feasible ratio, respectively.     

Combining pathway analysis with flux balance analysis for the comprehensive study of metabolic systems

The elucidation of organism-scale metabolic networks necessitates the development of integrative methods to analyze and interpret the systemic properties of cellular metabolism. A shift in emphasis from single metabolic reactions to systemically defined pathways is one consequence of such an integrative analysis of metabolic systems. The constraints of systemic stoichiometry, and limited thermodynamics have led to the definition of the flux space within the context of convex analysis. The flux space of the metabolic system, containing all allowable flux distributions, is constrained to a convex polyhedral cone in a high-dimensional space. From metabolic pathway analysis, the edges of the high-dimensional flux cone are vectors that correspond to systemically defined "extreme pathways" spanning the capabilities of the system. The addition of maximum flux capacities of individual metabolic reactions serves to further constrain the flux space and has led to the development of flux balance analysis using linear optimization to calculate optimal flux distributions. Here we provide the precise theoretical connections between pathway analysis and flux balance analysis allowing for their combined application to study integrated metabolic function. Shifts in metabolic behavior are calculated using linear optimization and are then interpreted using the extreme pathways to demonstrate the concept of pathway utilization. Changes to the reaction network, such as the removal of a reaction, can lead to the generation of suboptimal phenotypes that can be directly attributed to the loss of pathway function and capabilities. Optimal growth phenotypes are calculated as a function of environmental variables, such as the availability of substrate and oxygen, leading to the definition of phenotypic phase planes. It is illustrated how optimality properties of the computed flux distributions can be interpreted in terms of the extreme pathways. Together these developments are applied to an example network and to core metabolism of Escherichia coli demonstrating the connections between the extreme pathways, optimal flux distributions, and phenotypic phase planes. The consequences of changing environmental and internal conditions of the network are examined for growth on glucose and succinate in the face of a variety of gene deletions. The convergence of the calculation of optimal phenotypes through linear programming and the definition of extreme pathways establishes a different perspective for the understanding of how a defined metabolic network is best used under different environmental and internal conditions or, in other words, a pathway basis for the interpretation of the metabolic reaction norm.     

Combining flux and energy balance analysis to model large-scale biochemical networks

Stoichiometric Network Theory is a constraints-based, optimization approach for quantitative analysis of the phenotypes of large-scale biochemical networks that avoids the use of detailed kinetics. This approach uses the reaction stoichiometric matrix in conjunction with constraints provided by flux balance and energy balance to guarantee mass conserved and thermodynamically allowable predictions. However, the flux and energy balance constraints have not been effectively applied simultaneously on the genome scale because optimization under the combined constraints is non-linear. In this paper, a sequential quadratic programming algorithm that solves the non-linear optimization problem is introduced. A simple example and the system of fermentation in Saccharomyces cerevisiae are used to illustrate the new method. The algorithm allows the use of non-linear objective functions. As a result, we suggest a novel optimization with respect to the heat dissipation rate of a system. We also emphasize the importance of incorporating interactions between a model network and its surroundings.