## Abstract

To better understand the dynamic regulation of optimality in metabolic networks under perturbed conditions, we reconstruct the energetic‐metabolic network in mammalian myocardia using dynamic flux balance analysis (DFBA). Additionally, we modified the optimal objective from the maximization of ATP production to the minimal fluctuation of the profile of metabolite concentration under ischemic conditions, extending the hypothesis of original minimization of metabolic adjustment to create a composite modeling approach called M‐DFBA. The simulation results are more consistent with experimental data than are those of the DFBA model, particularly the retentive predominant contribution of fatty acid to oxidative ATP synthesis, the exact mechanism of which has not been elucidated and seems to be unpredictable by the DFBA model. These results suggest that the systemic states of metabolic networks do not always remain optimal, but may become suboptimal when a transient perturbation occurs. This finding supports the relevance of our hypothesis and could contribute to the further exploration of the underlying mechanism of dynamic regulation in metabolic networks.

## Introduction

Producing a modeling method that permits the formulation of hypotheses on missing components and unidentified interactions could lead to an understanding of systemic behaviors with limited input of dynamics information (Palsson, 2000; Ideker *et al*, 2001; Kitano, 2002; Collins *et al*, 2003; Nicholson *et al*, 2004). Metabolic flux analysis is a significant effort towards this goal (Fell, 1996; Schilling *et al*, 1999). Flux balance analysis (FBA), a method that assumes that optimal biological functions are acquired through an evolutionary process (Ibarra *et al*, 2002), has been especially promising (Varma and Palsson, 1994). Able to model an optimal systemic state based on stoichiometry by integrating a restricted set of parameters in a system, it has been widely used to rebuild and analyze behaviors of semicharacterized metabolic networks (Edwards and Palsson, 2000; Covert *et al*, 2001; Ibarra *et al*, 2002). Dynamic flux balance analysis (DFBA), a further development of FBA, has been used to represent dynamic process of metabolic networks based on the theory of dynamic optimal control. With the ability to predict the metabolite concentration, it has been successfully applied to the analysis of diauxic growth in *Escherichia coli* (Mahadevan *et al*, 2002).

However, as FBA and DFBA ignore the possibility that, under abnormal conditions, metabolic networks may not immediately regulate towards the optimal objective, the method of minimization of metabolic adjustment (MOMA) has been designed based on the hypothesis that perturbed metabolic fluxes undergo a minimal redistribution (Segrè *et al*, 2002). MOMA is considered to be an important improvement over static FBA for the simulation of perturbed metabolic networks; since its development, it has been used frequently in modeling mutated metabolic networks (Kauffman *et al*, 2003; Stephanopoulos *et al*, 2004).

Myocardial ischemia can be defined as insufficient myocardial blood flow for a given myocardial energetic demand. Investigation of the energetic‐metabolic network under ischemic conditions is important for understanding the pathology of coronary artery disease, which frequently results in myocardial ischemia. The primary effect of ischemia is impaired mitochondrial metabolism, resulting in less ATP formation, and accumulation of lactate in the cell. The oxidation of fatty acids provides 60–90% of the energy for ATP synthesis in a normal, healthy heart, with the balance coming from carbohydrates (glucose and lactate). Studies in mammalians show that under mild and moderate ischemic conditions, the dominant contribution of fatty acid to oxidative ATP synthesis does not alter, even though fatty acid has a lower efficiency in ATP synthesis given limited oxygen supply (Stanley *et al*, 1997). The mechanisms responsible for aggravated insufficient energy production are still unclear. Previous mechanistic models simulating the energy metabolism in ischemic myocardia (Salem *et al*, 2002) represented this phenomenon under the condition of 40% blood flow, but could not give any relevant explanation.

In this study, we modify the DFBA method with the MOMA hypothesis, naming the composite modeling approach M‐DFBA. We reconstruct the M‐DFBA model of energy metabolism in ischemic mammalian myocardia to explore the dynamic process of systemic adjustment under abnormal conditions. The original MOMA method based on static FBA does not integrate the information on metabolite concentrations, which represents the characteristics of a systemic state (Raamsdonk *et al*, 2001). The transition of the profile of metabolite concentrations may cause a series of dynamic processes in the cell, resulting in regulation at a genetic or metabolic level. Consequently, for the M‐DFBA model, we extend the MOMA hypothesis from a minimal redistribution of metabolic fluxes to a minimal fluctuation of the profile of metabolite concentrations in perturbed metabolic networks. The resulting model is able to represent the dynamic regulation of utilization of metabolic substrates for energy production. The simulations of the M‐DFBA model show a closer correlation to experimental data than do those of the original DFBA model, in which maximal ATP production is set as the optimal objective under ischemic conditions. Furthermore, the M‐DFBA model describes the phenomenon of predominant contribution of fatty acids to oxidative ATP production under mild and moderate ischemic conditions, which the DFBA model is not able to do. The consistency between M‐DFBA and experimental reports indicates that the metabolic system does not keep the optimal objective of maximization of ATP production under transient perturbation, but rather adopts a minimal regulation of its systemic state. This finding may reflect a conservative survival response of life forms under sudden perturbation prior to their making any adaptations that may be required for long‐term evolutionary survival.

Integrating the MOMA method with DFBA for the first time, this work promotes the significance of the MOMA hypothesis for predicting the reactions of transiently perturbed systems, which may help deepen understanding of the dynamic regulations of optimality in myocardial metabolic networks undergoing environmental perturbations.

## Results and discussion

### Rebuilding the metabolic network by M‐DFBA method

The following pathways are considered in the model of myocardial energy metabolism, including glycolysis, fatty acid oxidation, glycogen oxidation, phosphocreatine synthesis and breakdown, and the TCA cycle. Substrates of energetic metabolism, such as glucose, fatty acid, lactate and oxygen, are transported from blood flow into myocardia; glycogen and phosphocreatine act as endogenic supplements.

Based on the premise that under normal conditions the energetic‐metabolic network in myocardia maximizes ATP production (Cairns *et al*, 1998), this complex network was simplified according to these rules: priority is placed on the pathways that contribute most to ATP production; anaerobic metabolism, as a rescue mechanism, is also considered, even though it produces less ATP in normal myocardia than does aerobic metabolism; finally, metabolic intermediates, such as pyruvate and acetyl‐CoA, are overlooked in accordance with the modeling rules of DFBA. Thus, we obtained a simplified network including eight crucial pathways involving seven important metabolites (Figure 1 ).

We reconstruct the myocardial energetic‐metabolic network using the dynamic optimization approach (DOA), one kind of formulations within DFBA. As mentioned above, the maximal ATP production is supposed to be the systemic optimal objective of normal myocardia (Cairns *et al*, 1998). While under ischemic conditions, we modify the optimal objective by minimizing the fluctuation of the profile of metabolite concentrations, according to the hypothesis of M‐DFBA. To describe the energy demand in myocardia, a lower limit constraint of ATP synthesis velocity is required under ischemic conditions. The solution to the dynamic optimization problem is described in the Materials and methods section, and a simple example of how this algorithm works is provided in Supplementary information.

### The M‐DFBA modeling results

This M‐DFBA model describes the dynamic behaviors of the metabolic network under normal (*F*=1) and ischemic (*F*<1) conditions over 10 min, where *F* represents the blood flow. According to the extent of ischemia, the modeling conditions are defined as mild (*F*=0.8), moderate (*F*=0.6 and 0.4), and severe (*F*=0.2) ischemia, respectively. The inputs of fuels to the myocardia, such as glucose, free fatty acids, lactate, and oxygen, are assumed to be infinite, and the concentrations of them in the blood are assumed to be constant throughout the simulation (Salem *et al*, 2002).

Figure 2 shows the principal simulation results. In general, the metabolite concentration remains almost constant in normal myocardia; in contrast, during ischemic myocardia the consumption of both endogenic and exogenic fuels significantly changes, and less ATP is produced (Supplementary Figure 2). It is noticeable that there are slight local fluctuations in those curves, which might be attributed to the presence of a large null space in the process of searching for the optimal objective. We attend to the transition trend of the curves in the present study.

The major modeling results that were consistent with previous publications are summarized as follows:

Under ischemic conditions, endogenic glycogen and phosphocreatine are consumed (Figure 2C and D) (Stanley

*et al*, 1997). The utilization of fatty acid decreases (Figure 2F) and it accumulates in myocardia (Figure 2B) (Stanley, 2000; Lee*et al*, 2004).Lactate, an exogenic fuel, increases in intracellular concentration (Figure 2A), reversing from an input flux to an output flux during the aggravation of ischemia (Figure 2G) (Pantely

*et al*, 1990; Stanley*et al*, 1997; Salem*et al*, 2002); this occurrence could alleviate the toxicity of lactate accumulation caused by raised anaerobic metabolism under ischemic conditions.Under moderate ischemic conditions, glucose uptake increases (Figure 2E) and oxygen uptake decreases (Figure 2H) (Salem

*et al*, 2002). However, when ischemia becomes severe, glucose uptake begins to decrease, while still remaining at a higher level than normal (Figure 2E) (Stanley*et al*, 1997). It seems that although glucose combustion produces more ATP, it does not continue to increase to prevent lactate accumulation.

It is noticeable that the M‐DFBA model predicts several phenomena that are more consistent with experimental data than are the predictions of the DFBA model (see Supplementary information for details):

(1) Under the condition of 40% blood flow, or moderate ischemia, the value of lactate concentration in the M‐DFBA model increases up to around 7 μmol g^{−1} in 10 min, which is much closer to the average experimental data of about 5 μmol g^{−1} (Arai *et al*, 1991) than is the approximate value of 20 μmol g^{−1} in the DFBA model (see Supplementary Figure 4C and D). Moreover, the velocity of lactate uptake in the M‐DFBA model decreases to −0.6∼−0.7 μmol g^{−1} min^{−1}, which is fairly consistent with the experimental data of −0.2∼−0.7 μmol g^{−1} min^{−1} (Stanley *et al*, 1994) compared to the value of −4 μmol g^{−1} min^{−1} in the DFBA model (see Supplementary Figure 4A and B). Considering the reliability of M‐DFBA modeling results, we propose that the metabolic network regulates in a mild way with respect to the flux configuration and metabolite concentration distribution rather than in an extreme way to meet the optimization.

(2) In theory, glucose is regarded as preferable to fatty acid for ATP production under ischemic conditions, given that when 1 mol O_{2} is utilized, glucose combustion produces 6.0 mol ATP, while palmitic acid combustion produces 5.2 mol ATP. Consequently, it is puzzling that mild or moderate ischemia does not actually alter the relative contributions of fatty acid and carbohydrates to oxidative ATP synthesis (Stanley, 2000). The DFBA model cannot simulate this phenomenon; the consumption of fatty acid in the DFBA model decreases sharply under ischemic conditions and loses its predominant contribution to energy production (see Supplementary Figure 1D–F). However, the M‐DFBA model is able to predict the actual dynamics under ischemic conditions (Figure 2D–F) since the increase in consumption of both glucose and glycogen and the decrease in consumption of fatty acid are all much less than that found in the DFBA model. According to the hypothesis of M‐DFBA, we propose that under ischemic conditions the energetic‐metabolic system does not function according to the optimal objective of ATP production maximization; rather, without immediately regulating mechanisms, the system is transferred to a suboptimal state. It is interesting that an alternative treatment of ischemia is to inhibit myocardial fatty acid oxidation and promote carbohydrate oxidation by which the insufficiency of energy production is mitigated (Stanley, 2001; Lee *et al*, 2004).

### Conclusion

In the present work, we reconstructed the energetic‐metabolic network in mammalian myocardia under normal and ischemic conditions using the M‐DFBA method. The simulation results are basically supported by the previously published experimental data. We propose that the modification of the optimal objective, from the maximization of ATP production in the DFBA model to the minimization of metabolite concentration fluctuation under ischemic conditions in the M‐DFBA model, is the key point of this work. Actually, the biological relevance of the selection of systemic optimality has always been questioned, especially for a system under perturbations. Our results show that the metabolic network in myocardia does not follow the optimal objective of ATP production maximization under ischemic conditions, but instead reaches a suboptimal level of energy production with minimal adjustment with respect to metabolite concentration distribution, which extends the MOMA hypothesis. Also, these results may help us to understand the underlying mechanisms involved in the dynamic regulation of optimality of myocardial metabolic networks under ischemic conditions.

## Materials and methods

### Modeling approach

#### FBA and DFBA methods.

As was described previously (Varma and Palsson, 1994; Kauffman *et al*, 2003), the FBA modeling method constrains the capability of a metabolic network by balancing the metabolic fluxes. When the network is in a steady state, the mass balances can be described as follows:

where *S* is the stoichiometric matrix and *V* is the flux vector. By imposing the flux constraints, systemic behaviors can be restricted to an enclosed solution space. Finally, a solution is obtained by optimizing an objective function, such as the maximum of biomass, or the minimum production of toxin, using linear programming.

The classical FBA method mentioned above does not take into account the factors of time and molecular concentration, making it difficult to represent the dynamic process of the biological system. Hence, DFBA based on optimal control theory was developed to remedy these shortcomings (Mahadevan *et al*, 2002). Currently, DFBA includes two different formulations: DOA and static optimization‐based DFBA approach (SOA). DOA performs with greater precision than SOA, while SOA is scalable to larger metabolic networks, in which the number of variables does not increase exponentially with the scale of a network (Mahadevan *et al*, 2002). In need of precision, we chose the DOA to reconstruct the myocardial energetic‐metabolic model. A general dynamic optimization problem is represented as follows:

where Objective is the objective function, *L(V,X)* is the linear constraints, *C(V,X)* is the nonlinear constraints, and *X*_{0} is the original concentration of a particular metabolite. The solution to the dynamic optimization problem with the constraints of a differential‐algebraic equation (DEA) is to parameterize the dynamic equations by an orthogonal collocation of finite elements (Cuthrell and Biegler, 1987). Then, the dynamic optimization problem can be transformed into a nonlinear programming problem, and solved by the *fmincon* function in MATLAB^{®}6.5 (The MathWorks Inc., Natwick, MA). Mahadevan *et al* (2002) introduced the details of this solution process in their published paper.

#### The formulation of M‐DFBA model.

As myocardia continuously exchange metabolic substances with blood, such as glucose, free fatty acids, and lactate, we used the following equation to describe the rate of change of metabolite concentration:

where *X* is a vector of metabolite concentration, *t* is the time, *S* is the stoichiometric matrix, *V* is a vector composed of the values of all reactions and transport fluxes, *F* is blood flow into myocardia, Out_{x} is the concentration of a particular fuel in arterial blood, &_{x} is the blood–tissue partition coefficient, In_{x} is the fuel concentration in cardiac myocytes, and *F*(Out_{x}−&_{x}In_{x}) is the flux of a particular fuel from blood flow to cardiac myocytes.

According to the DOA, the M‐DFBA formulation of this model is described in equation (4).

As the concentration of various interior metabolites represents the characteristic of a systemic state (Raamsdonk *et al*, 2001), we define the minimization of the fluctuation of concentration vector *X* to be the optimal objective under ischemic conditions, extending the MOMA hypothesis. Where *N* is the number of metabolites in the network, *x*_{i,j} represents the value of the concentration of metabolite *i* on the time point of orthogonal root *j*, represents the Euclidean distance between two adjacent orthogonal roots (see Supplementary Figure 3), δ(*t*) is the Dirac‐delta function, and *M* is the number of orthogonal roots. The goal is to find the vector *X* such that the integral of Euclidean distances is minimized.

The stoichiometic matrix of the M‐DFBA model is received from the simplified network shown in Figure 1. The adjusted stoichiometric coefficients are calculated based on the data set of flux values (Salem *et al*, 2002; Supplementary Table 5). The constraints of the DEA describe the rate of change of metabolite concentration, where *S*_{Gluc}, *S*_{Gly}, *S*_{FA}, *S*_{Lac}, ^{S}o_{2}, *S*_{PC}, and *S*_{ATP} indicate rows of the stoichiometric matrix associated with glucose, glycogen, fatty acid, lactate, oxygen, phosphocreatine, and ATP, respectively. *X*_{0} is the vector of original concentration. In addition, three linear constraints of particular fluxes and molecular concentrations are integrated. The first one is qualitatively estimated to represent the capability of myocardia to act as an energetic buffer during blood flow reduction, namely, the concentration of a cellular metabolite will not fall directly to zero under ischemic conditions. The other two constraints are the maximum flux values of glucose glycolysis and oxygen utilization, respectively. Moreover, considering the feedback restraint of pyruvate dehydrogenase by acetyl‐CoA, we constrained the proportion of pyruvate‐to‐acetyl‐CoA flux to acetyl‐CoA‐to‐CO_{2} flux. This constraint indicates that the increase in combustion of fatty acids, which causes the accumulation of acetyl‐CoA will inhibit the oxidation of carbohydrate (glucose and lactate) in myocardia. Finally, as the normal optimal objective of ATP production is replaced under ischemic conditions by minimization of metabolite concentration fluctuation, a new constraint is required for the purpose of describing the unaltered requirement of energy in ischemic myocardia. We made a qualitative estimate of the lower limit constraint of ATP synthesis velocity to express this energetic demand, where α represents the normal velocity of ATP consumption (the values of the parameters used in equation (4) are listed in Supplementary Table 1–5):

## Acknowledgements

We are grateful to Yike Guo and Qiang Lu for their valuable comments on this work. The National High‐Tech Research and Development Program of China (2002AA234051, 2003AA231011) provided support for this work, and the National ‘973’ Basic Research Program of China (2004CB518606).

## Supplementary Information

Supplementary Figure 1

The modeling results of the DFBA model. (A) lactate concentration increased to 10∼35 µmol g‐1 during moderate ischemia and continued to rise during severe ischemia Lactate concentration, it is severe higher than the experimental data. (B) Fatty acid accumulated. (C, D) The endogenic phosphocreatine and glycogen were consumed to produce ATP under ischemic conditions. (E) During moderate ischemia, the uptake of glucose increased, while during severe ischemia it decreased, but still showed a higher level than normal. (F, H) The uptake of fatty acid and oxygen decreased. But under mild and moderate ischemia, considering the ATP production decreases, the predominant contribution of fatty acid to oxidative ATP synthesis compared to carbohydrates is actually altered. This prediction is not conformable with previous findings. (G)The input flux of lactate converted to an output flux during moderate and severe ischemia, however it is severe higher than the experimental data. [msb4100071-sup-0001.tiff]

Supplementary Figure 2

The slope of the curve at each point represents the velocity of ATP production at particular time. (A) In DFBA model, under normal condition (F=1), the velocity of ATP production approximately equals to 34.9 µ mol g wet wt‐1 min‐1, the value of the normal velocity of ATP consumption; while under ischemic conditions (F<1), it decreases. (B) In M‐DFBA model, the velocity of ATP production approximately equals to that in DFBA model under normal condition; while under ischemic conditions, it is always lower than those in DFBA model. [msb4100071-sup-0002.tiff]

Supplementary Figure 3

The explanation of how calculating the Euclidean distance described in objective function. The upper figure shows the process of dividing time period into a finite number of intervals (finite elements). In this example, the number of intervals is 4. Then each variable was parameterized at the roots of an orthogonal polynomial within each finite element. In this example, the number of orthogonal points in each finite element is 2. The lower figure shows the Euclidean distance of the vector of metabolite concentrations between two adjacent orthogonal roots (j and j+1). [msb4100071-sup-0003.tiff]

Supplementary Figure 4

Modeling results of lactate uptake and lactate concentration in DFBA and M‐DFBA model respectively. The experimental data at the original time point are measured under normal condition (F=1); while the data marked on the time point of 10 minutes in the figures are experimentally measured under the condition of F=0.4 in 10 minutes or longer time. (A,C) The modeling results in DFBA model show that the velocity of lactate uptake decreases approximately to −4 µmol g‐1 min‐1 under the conditions of F=0.4, which is much lower than the experimental data of −0.2 ∼ −0.7µmol g‐1 min‐1 (Stanley *et al*., 1994); while the value of myocardial lactate concentration increases approximately to 20 µmol g‐1 under the conditions of F=0.4, which is much higher than the average experimental data of 5 µmol g‐1 (Arai *et al*., 1991). (B,D) The M‐DFBA modeling curves of lactate uptake and myocradial lactate concentration under the condition of F=0.4 are more consistent with experimental data compared with those in DFBA model. [msb4100071-sup-0004.tiff]

Supplementary Table 1 [msb4100071-sup-0005.pdf]

Supplementary Table 2 [msb4100071-sup-0006.pdf]

Supplementary Table 3 [msb4100071-sup-0007.pdf]

Supplementary Table 4 [msb4100071-sup-0008.pdf]

Supplementary Table 5 [msb4100071-sup-0009.pdf]

Supplementary Information 1 [msb4100071-sup-0010.pdf]

Supplementary Information 2 [msb4100071-sup-0011.pdf]

Supplementary Information 3 [msb4100071-sup-0012.pdf]

Supplementaryfigure and table legends [msb4100071-sup-0013.pdf]

## References

- Copyright © 2006 EMBO and Nature Publishing Group