Data from 31P-nuclear magnetic resonance spectroscopy of human forearm flexor muscle were analyzed based on a previously developed model of mitochondrial oxidative phosphorylation (PLoS Comp Bio 1: e36, 2005) to test the hypothesis that substrate level (concentrations of ADP and inorganic phosphate) represents the primary signal governing the rate of mitochondrial ATP synthesis and maintaining the cellular ATP hydrolysis potential in skeletal muscle. Model-based predictions of cytoplasmic concentrations of phosphate metabolites (ATP, ADP, and Pi) matched data obtained from 20 healthy volunteers and indicated that as work rate is varied from rest to submaximal exercise commensurate increases in the rate of mitochondrial ATP synthesis are effected by changes in concentrations of available ADP and Pi. Additional data from patients with a defect of complex I of the respiratory chain and a patient with a deficiency in the mitochondrial adenine nucleotide translocase were also predicted the by the model by making the appropriate adjustments to the activities of the affected proteins associates with the defects, providing both further validation of the biophysical model of the control of oxidative phosphorylation and insight into the impact of these diseases on the ability of the cell to maintain its energetic state.
- computational model
- cellular energetics
- oxidative phosphorylation
- 31P-NMR spectroscopy
mitochondrial oxidative ADP phosphorylation is the primary source of ATP in skeletal muscle during aerobic exercise. Thus, to maintain the free energy state of the cytoplasmic phosphoenergetic compounds ATP, ADP, and Pi, oxidative phosphorylation is modulated to match the rate of ATP utilization during exercise. It has recently been shown through computational model-based analysis of data obtained from 31P-NMR spectroscopy of working in vivo dog hearts that the primary control mechanism operating in cardiomyocytes is feedback of substrate concentrations for ATP synthesis (5). In other words, changes in the concentrations of the products generated by the utilization of ATP in the cell, ADP and Pi, effect changes in the rate at which mitochondria utilize those products to resynthesize ATP (5).
Here the question of whether this same mechanism can explain the observed data on the control of oxidative metabolism in skeletal muscle is investigated. Previous analyses of 31P-NMR spectroscopy (31P-MRS) data on energy balance in exercising skeletal muscle have mainly focused on testing ADP feedback control of mitochondrial ATP synthesis using black box descriptions of the mitochondrial ATP synthetic pathway (8, 14–16, 28), Pi acceptor control (7), and thermodynamic control involving quasi-linear relations between cytoplasmic Gibbs free energy of ATP hydrolysis and mitochondrial ATP synthesis flux (13, 18, 31). Yet, to date, these 31P-MRS data have not been adequately explained based on a detailed mechanistic model of oxidative phosphorylation and cellular energetics.
To analyze and interpret data from skeletal muscle, our previously published model of oxidative ATP synthesis and metabolism in cardiomyocytes (5) is adapted to skeletal muscle by setting intracellular concentration pools of creatine and phosphate to appropriate values (based on measured data) and appropriately adjusting the cellular mitochondrial content to match the available morphometric data. Data on cytoplasmic ADP and Pi concentrations as a function of work rate in human forearm flexor muscle from 20 untrained healthy subjects (13), 6 subjects with complex I deficiency in skeletal muscle (11, 22), and a single patient with deficiency in adenine nucleotide translocase (ANT) (2, 3) were analyzed based on the computational model.
Results of this analysis indicate that the existing computational model of the kinetics of mitochondrial oxidative phosphorylation accurately captures the in vivo kinetics of oxidative ATP synthesis and transport of phosphate metabolites between mitochondria and cytoplasm in skeletal muscle. The mechanism of control of oxidative phosphorylation is through feedback of substrates for ATP synthesis. No additional regulatory mechanisms, such as feed-forward control of certain enzymes via cytosolic calcium levels (9) or functional coupling between mitochondrial creatine kinase and ANT (23, 26, 27), are necessary to explain the majority of the observed data. In addition, the impact of specific protein deficiencies on the relationship between oxidative phosphorylation and cytoplasmic ADP and Pi is successfully explained by making the appropriate modifications to the mitochondrial enzymes altered in the diseases.
Overview of Computational Model
ATP utilization, cytoplasmic phosphoenergetic buffers, and oxidative ATP synthesis are simulated in a model of skeletal muscle energetics illustrated in Fig. 1. The cell is divided into cytoplasmic and mitochondrial compartments; the variables simulated within the compartments are listed in Table 1, with a brief description of the variables and units associated with each variable. The computational model for cellular energetics and oxidative phosphorylation is derived from recently published computational models developed for cardiac mitochondria (4) and cardiomyocytes (5). A complete description of the computational model is provided in the appendix.
Model Parameter Values
Model parameter descriptions and assigned values are listed in Table 2. With the exception of the total pool of exchangeable phosphate (TPP) in the cell, all parameter values in the model are fixed at values justified by previous studies. The total exchangeable phosphate pool is computed from the equation (1) where Vcyto and Vmito are the volume densities of cytoplasm and mitochondria in the myocyte model (in units of volume cytoplasm or mitochondria per cell volume); Wx, and Wi are the matrix and intermembrane space water volumes (in units of volume of water per volume of mitochondria), respectively; and Wc is the water fraction of the cytoplasmic space. By comparing simulation predictions with experimental data (see results), values of TPP of 36.8, 36.3, and 30.3 mM are chosen as the values most consistent with the experimental observations for healthy subjects, complex I-deficient subjects, and ANT-deficient subjects. To adjust the TPP value, the value of [Pi]c used as an initial condition in model simulations is adjusted to obtain optimal model fits to the observed data.
The parameter values listed in Table 3 are organized structure/volume parameters, physicochemical parameters, mitochondrial model parameters, fixed concentrations and concentration pools, and binding constants.
All values for concentrations of pooled metabolites are set according to values reported in previous studies, with the exception of TPP, which is estimated below. Binding constants are obtained from the literature; enzyme activities for reactions maintained near equilibrium are set to arbitrarily high values.
31P-NMR spectroscopic measurements of phosphate metabolites in the ulnar finger flexor muscle in the human forearm were acquired at 1.5 T at rest and during voluntary ramp exercise using 1H imaging-guided localization in all subjects and analyzed according to the methods described in detail elsewhere (12).
Analysis of Data from Healthy Subjects
Model predictions of steady-state concentrations of cytoplasmic ADP and Pi as functions of the work rate (rate of cellular ATP consumption) are compared with experimental measures obtained from healthy subjects in Fig. 2. The model predicts that ADP concentration increases from ∼17 μM at rest to ∼110 μM at the maximum observed work rate and that Pi concentration increases from 0.3 mM at rest to ∼18 mM at the maximum observed work rate.
The experimental data plotted in Fig. 2 are obtained from the previous study of Jeneson et al. (13). Data were collected during steady-state exercise, during which the contribution to ATP synthesis from anaerobic glycolysis was negligible (13). Cellular pH stayed within 0.1 pH units of 7.0 in all subjects. Jeneson et al. (13) reported measured phosphate metabolite levels as a function of the work rate measured as a fraction of the maximal exercise rate. These fractional work rates were translated to absolute ATP synthesis fluxes based on the observations that the average rates of ATP hydrolysis in human forearm flexor muscle at rest and 65% of maximal exercise are 0.008 mmol (liter cell water)−1·s−1 (6) and 0.22 mmol (liter cell water)−1·s−1, respectively (14), and the assumption that power output and rate of ATP hydrolysis are linearly proportional (21).
The agreement between the model simulations and the observed data is striking considering that a single adjustable parameter (TPP) was varied to match model simulations to the data. In fact, the model-predicted [ADP]c values are not sensitive to the value of TPP. A 10% increase in the value of this parameter results no significant change in the model predicted [ADP]c and a >100% increase in the mean-squared difference between model predictions of [Pi]c and the observed data. Thus the nature of the relationship between work rate and [ADP]c does not depend on the value of TPP. The predicted [ADP]c and [Pi]c as a functions of workload at the value of TPP = 40.5 mM (10% greater than the optimal value) is plotted as a dashed lines in Fig. 2. It is apparent that the higher value of TPP results in an improvement in the model fit to the [Pi]c data at low work rates, but an overall agreement between the model predictions and experimental data that is worse than for the optimal value.
To investigate the factors controlling this relationship, we analyzed the sensitivity of the model predictions of Fig. 2A to the parameter values used in the kinetic model for the ANT flux of Eq. A10. This expression invokes three parameters, XANT, θ, and Km,ADP (assumed values for these parameters are listed in Table 2). To quantify the impact of variation in these parameters on the work-ADP relationship, we fit the predicted data to the function V = Vmax([ADP]c − xo)/([ADP]c − xo + Km), where V represents the ATPase flux and Vmax, xo, and Km are fitting parameters. The predictions plotted in Fig. 2A are well represented by this function for Vmax = 0.44 mmols−1 (l cell)−1, xo = 16.5 μM, and Km = 0.11 mM. Sensitivity analysis reveals that the apparent Km is highly sensitive to the value of θ (sensitivity coefficient |∂ ln Km/∂ln θ| ≈5), and less sensitive to XANT (sensitivity coefficient | ∂ ln Km /∂ ln XANT | ≈ 0.5). The Vmax is most sensitive to the value of XANT (sensitivity coefficient |∂ ln Vmax /∂ ln XANT|≈ 0.6) and less sensitive to θ (sensitivity coefficient |∂ ln Vmax/∂ ln θ| ≈ 0.3). Neither Km nor Vmax is sensitive to the value of Km,ADP. Thus, in terms of the ANT transporter model, the apparent Km for the relationship between ADP and work is primarily controlled by the parameter θ.
Analysis of Data from Patients with Complex I Deficiency
Model predictions and experimental data from six patients (three of which were first- and second-degree relatives) who were diagnosed with mitochondrial complex I deficiency at the isolated mitochondria level are plotted in Fig. 3. These data were collected under the same protocol as for the healthy subjects; data were first published in Ref. 11. The experimental data show steeper relationships between [ADP]c (or [Pi]c) and ATP hydrolysis rate than are observed in healthy subjects, resulting from an impaired capacity of the mitochondria to synthesize ATP and transport it to the cytoplasm as levels of ADP and Pi increase.
To match the observed data on complex I-deficient patients, the mitochondrial model was modified in two ways. First, based on observations that whole-body resting oxygen consumption is increased in the patients with complex I deficiency compared with healthy subjects (22), it was assumed that no protons are pumped by complex I, altering the stoichiometry of the reaction model. Changes to the governing equations to account for this phenomenon are described in the appendix; see Eq. A19 for details. This change in the proton stoichiometry results in a 47% increase in the predicted rate of resting oxygen consumption in the muscle compared with the normal (healthy) case. Measured whole-body resting oxygen consumption in three patients with this deficiency was 28 ± 14% greater compared with healthy subjects (22). Second, the activity of complex I was reduced to ∼1/1,000 of the normal value to match the observed data. Thus, to obtain the model predictions illustrated in Fig. 3, the activity of complex I was reduced by approximate 1,000-fold compared with the normal case.
Here the model shows excellent agreement to the observed data on both [ADP]c and [Pi]c. As is the case for the data from healthy subjects, the predictions of [Pi]c are sensitive to the value of TPP, whereas the predictions of [ADP]c are not. Thus the relationship between work rate and substrates for ATP synthesis is explained by a drastic reduction in the activity and a loss in proton pumping of mitochondrial complex I.
Analysis of Data from Patient with ANT Deficiency
Figure 4 shows data on [ADP]c and [Pi]c measured in a single patient characterized as having an ANT deficiency in muscle (2, 3). The ANT transporter exchanges mitochondrial ATP for cytoplasmic ADP. Thus impairment in the activity of ANT results in a reduction in the ability of the mitochondrion to deliver ATP to the cytoplasm. Data were collected under the protocol described for healthy and complex I-deficient subjects; a subset of these data was first published in Ref. 11.
Western blot analysis revealed that in the patient the ANT protein was present at 25% of the concentration found in healthy subjects (3). Thus to account for the deficiency in the model, the activity of the ANT transporter was reduced to 25% of the normal value for healthy subjects.
Figure 4 shows that ADP was much higher at rest and increased more rapidly with exercise in this patient than in healthy subjects and in the complex I-deficient patients. The ADP concentration of 150 μM that was measured at the modest work rate of 0.081 mmol ATP consumed per second per liter cell was greater than any value measured in the healthy or complex I deficient subjects. Model predictions are similar to the measured data; the model predicts [ADP]c is 42 μM at rest (compared with the measured value of 64 μM) and increases sharply with work rate.
Maintenance of free energy of ATP hydrolysis.
Figure 5 shows model-predicted free energy of ATP hydrolysis vs. ATP hydrolysis rate, computed using the parameters obtained for the healthy subjects, complex I-deficient patients, and the ANT deficient patient. Since cellular pH values are 7.0 ± 0.1 for all the subjects, the free energy of ATP hydrolysis is computed by the relationship (2) where ΔGATP is the standard free energy of ATP hydrolysis at pH 7.0, R is the universal gas constant, and T is temperature in degrees Kelvin (for parameter values, see Table 2). [fADP]c and [fATP]c denote magnesium unbound ADP and ATP concentration in the cytoplasm, respectively. For the normal subjects, Jeneson et al. (9) observed quasi-linear relationship between the free energy of ATP hydrolysis and power output; the model predictions verify this observation except at low work rates where the predicted values of [Pi]c tend to be lower than the observed values.
This work introduces an integrated computational model for skeletal muscle oxidative phosphorylation and fluxes of ATP, ADP, AMP, CrP, and Pi in the cytoplasm. Although the central component of the model—mitochondrial oxidative phosphorylation—is based on a mitochondrial model previously developed to match data on isolated mitochondria from rat heart, the integrated model matches a rich set of data on in vivo phosphate compounds from human skeletal muscle in healthy and complex I deficient individuals. The model also produces reasonable predictions for the ANT deficient subject, although the data available for comparison are sparse.
The analysis predicts that the rate of oxidative phosphorylation is primarily regulated through concentrations of the substrates for ATP synthesis (ADP and Pi), since no additional control mechanisms, such as feed-forward control of certain enzymes via cytosolic calcium levels (9) and functional coupling between mitochondrial creatine kinase and ANT (23, 26, 27) that have been proposed to operate in the heart, were incorporated into the model. The current analysis does not rule out the possibility that ancillary control mechanisms are active in skeletal muscle (16, 28); however, it shows that major contributions of such mechanisms to the overall regulation of the mitochondrial ATP synthetic pathway are not necessary to explain the thrust of the observed data.
Although the model predicts cytoplasmic ADP, Pi, and ATP (not shown) concentrations that agree well with observed data, the present model systematically underpredicts Pi concentration in the resting state. As illustrated in Fig. 2B, 31P MRS measurements indicate that inorganic phosphate concentrations are ∼3 mM at rest, whereas the model predicts resting [Pi]c to be only 0.3 mM. This may be due to in part to the fact that the mitochondrial model was constructed to match data obtained from mitochondria isolated from cardiomyocytes. Inorganic phosphate concentrations in heart are significantly lower than in mixed fiber type skeletal muscle such as human skeletal muscle (19) if not undetectable at low work rates (32).
We used the empirical fitting function V = Vmax([ADP]c − xo)/([ADP]c − xo + Km) to capture key trends in the model prediction for the control case and to characterize the sensitivity of the model predictions to key model parameters in terms of sensitivity coefficients for the fitting parameters. Sensitivity analysis revealed that the parameter θ for ANT flux in the model has a significant impact on the apparent affinity for ADP in the relationship between cytoplasmic ADP concentration and rate of oxidative metabolism. Future studies including a full-scale sensitivity and metabolic control analysis of the current model as well as next-generation models incorporating ancillary biochemical detail will be necessary to further improve agreement of predicted data with empirical knowledge.
The control of oxidative phosphorylation by substrate concentrations allows the mitochondria to maintain a free energy of ATP hydrolysis of less than −55 kJ/mol over the observed range of work rates in human forearm muscle of healthy subjects. However, as shown in Fig. 5, the magnitude of ΔGATP drops more quickly with increasing work in the ANT-deficient patients than in normal subjects, and the predicted magnitude of ΔGATP in the complex I-deficient patients at rest is significantly lower than that of the other two sets of subjects. These abnormalities in the complex I- and ANT-deficient subjects result in reduced capacity to do work.
Data on complex I-deficient patients are fit by reducing the activity of complex I compared with the normal case, and assuming that complex I pumps no protons in the complex I deficient patients. However, while the closest fit to the observed data on ADP and Pi is obtained by reducing the complex I activity by a factor of ∼1/1,000 compared with the normal case, it is necessary to be cautious in interpreting the scaling factor of ∼1/1,000. Since the complex I activity for the normal case was not identified with significant sensitivity (4), the ratio between activities in normal and complex I deficient patients also cannot be determined sensitively. Yet, regardless of the sensitivity of the estimate, it was clear from oxygen polarographic studies on isolated mitochondria from patient leg muscle biopsies (H. R. Scholte, unpublished observations) that the activity of complex I is significantly diminished in the complex I deficient patients.
To analyze data from the patient with a deficiency of ANT, it was not necessary to introduce an arbitrary scaling factor to fit the measured data. Since the level of ANT expressed in mitochondria of the patient was directly assayed and found to be 25% of that in healthy subjects, it was possible to incorporate this measurement directly in the model by scaling the healthy ANT activity by 0.25. The fact that differences observed in cytoplasmic phosphoenergetic compounds between healthy subjects and this patient are explained based on this single adjustment to the model provides independent validation of the mitochondrial model and the conclusions drawn from its behavior.
The computational model for skeletal muscle energetics and oxidative phosphorylation is derived from a previously published model applied to cardiac tissue (5). For the current application, the cardiac model has been modified in two ways to adapt the model to analyze data from human skeletal muscle. First, the oxygen transport component of the previous model has been removed. It is assumed that the skeletal muscle remains normoxic during the experiments and the cellular oxygen concentration, [O2], which is a variable in the model of Beard (5), appears as a fixed parameter in the current model. Second, the mitochondrial volume of the muscle cell, which is ∼30% of cell volume in cardiomyocytes (29), is set to Vmito = 0.056 (l mito) (l cell) −1, a value measured from biopsies of human vastus lateralis muscle (30).
The model is expressed in terms of the following set of differential equations: (A1)
In the above set of equations, the subscripts “x”, “i”, and “c”, denote mitochondrial matrix, intermembrane space, and cytoplasm, respectively. All of the variables in this set of equations are defined in Table 1.
In addition to the state variables treated in Eq. (A1), the concentrations of several species are computed (A2) where NADtot, Qtot, cytCtot, and Atot, are the total concentrations of NAD(H), ubiquinol, cytochrome c, and adenine nucleotide in the matrix, respectively, and CRtot is the total creatine plus creatine phosphate concentration in the cytoplasm.
Parameters that appear in the above equations are described in detail below. The fluxes that appear on the right-hand side of the governing equations are tabulated in Table 2. For mitochondrial species, the governing equations follow from Ref. 4. For cytoplasmic species, the reactions modeled are ATP consumption, creatine kinase reaction, adenylate kinase reaction, and transport between the cytoplasm and the mitochondrial intermembrane space.
Mathematical Expressions for Mitochondrial Fluxes
The expressions for the mitochondrial fluxes in the model are described in detail in Ref. 4, and are listed here without detailed explanations. Definitions of the variables and parameters that appear in the following expressions are listed in Tables 2 and 3.
Dehyhdrogenase flux: (A3) Complex I flux (A4) Where ΔGH = FΔΨ + RT ln([H+]c/[H+]x).
Complex III flux (A5) Complex IV flux (A6) where [O2] is the O2 concentration in the cell, which is set at the fixed constant of 3.48 × 10−5 M.
F1F0-ATPase flux (A7) Magnesium binding fluxes (A8) where [fATP]x, [fADP]x, [fATP]i, and [fADP]i denote magnesium unbound ATP in the matrix, ADP in the matrix, ATP in the intermembrane space, and ADP in the intermembrane space, respectively.
Substrate transport fluxes (A9) Adenine nucleotide translocase (ANT) flux (A10) where θ is an empirical parameter with value set to 0.35.
The phosphate-hydrogen cotransporter flux (A11) where (A12)
Mitochondrial adenylate kinase flux (A13) Proton leak flux (A14) Potassium-hydrogen ion exchange (A15)
Mathematical expressions for cytoplasmic reaction fluxes
Four biochemical processes are modeled in the cytoplasm-the adenylate kinase reaction, the creatine kinase reaction, ATP hydrolysis, and binding of magnesium ions to ADP and ATP.
The binding of magnesium to ATP and ADP in the cytoplasm takes the same form as the binding fluxes in the mitochondria (A16) where [fATP]c and [fADP]c denote magnesium unbound ATP and ADP in the cytoplasm. Similarly, the cytoplasmic adenylate kinase is analogous to the mitochondrial reaction (A17)
In Eq. A18, KAK is the equilibrium constant for the reaction 2ADP ⇋ATP + AMP, and XAK is the enzyme activity, which is set to a large enough value so that the reaction is effectively maintained in equilibrium.
The creatine kinase flux is modeled using the expression (A18) where the activity XCK is set to a large enough value so that the equilibrium KCK = ([ATP]c[Cr]c/[ADP]c[CrP]c[H+]c)eq is maintained during simulations. The value of the apparent equilibrium constant assumed here (see Table 2) is calculated to account for the intracellular ionic strength and magnesium ion concentration (25).
The flux JAtC is defined as the flux through the reaction ATP → ADP + Pi. Mathematical models for the ATP consumption flux are considered in the Results section.
Modification of the Model for Patients with Complex I Deficiency
Equation A1 assumes that 4 protons are pumped from the matrix to the intermembrane space for each pair of electrons transferred via the reaction at complex I. For patients with complex I deficiency, it is assumed that no protons are pumped at complex I, and the equations for [H+]x and ΔΨ are modified as follows (A19)
This work was supported by the National Institutes of Health Grants HL-072011 and EB-005825.
We thank Dr. H. R. Scholte for sharing oxygen polarographic data of isolated mitochondria of the complex I-deficient patients.
The costs of publication of this article were defrayed in part by the payment of page charges. The article must therefore be hereby marked “advertisement” in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.
- Copyright © 2007 the American Physiological Society