| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | TABLE OF CONTENTS |
Department of Bioengineering, University of California at San Diego, La Jolla, California
Correspondence: Address reprint requests to Bernhard O. Palsson, Dept. of Bioengineering, University of California, 9500 Gilman Dr., La Jolla, CA 92093-0412. E-mail:palsson{at}ucsd.edu.
| ABSTRACT |
|---|
|
|
|---|
| INTRODUCTION |
|---|
|
|
|---|
A recently developed dimension within the constraint-based modeling approach is uniform random sampling of the steady-state flux space (Almaas et al., 2004
; Wiback et al., 2004
). This approach is used to fully determine the range of possible steady-state fluxes allowed in the network under defined physicochemical constraints. A Monte Carlo sampling procedure was applied to the metabolic network of the human red blood cell, the modeling of which has reached an advanced state (Jamshidi et al., 2001
; Kauffman et al., 2002
; Lew and Bookchin, 1986
; Mulquiney et al., 1999
; Mulquiney and Kuchel, 1999
, 2003
; Price et al., 2003b
; Schuster and Holzhutter, 1995
; Schuster et al., 1988
; Wiback and Palsson, 2002
). Monte Carlo sampling has also proven very useful in analyzing the general properties of networks by testing their robustness to parameter variation (Alves and Savageau, 2000a
,b
,c
). The approach utilized herein has the utility to identify the selection of independent measurements to determine the state of a biochemical network and to predict systemic effects from the reduction in a maximal reaction rate to simulate enzymopathies.
| MATERIALS AND METHODS |
|---|
|
|
|---|
![]() | (1) |
![]() | (2) |
i, generated on each of the basis vectors, bi, as
![]() | (3) |
|
![]() | (4) |
The next sets of constraints used to define the steady-state solution space are minimum and maximum flux rates through each of the reactions,
![]() | (5) |
and the
for that reaction. These
and
constraints thus segment the null space defined by Eq. 4 (Fig. 1 B).
In this article, the elementary forward and reverse reactions are combined into a net reaction. Reactions can thus have a negative
to indicate that the reaction is being used in the direction opposite to that defined as positive in S. The
values are generally based on experimental measurements. For irreversible reactions, the
values are set to zero, and for reversible reactions the
is set to
.
Elimination of redundant constraints
Many reaction
levels cannot be reached in a steady state because the saturation of other reactions is more constraining upon the reaction flux, vi, than its own saturation state. Thus, many of the
constraints are redundant from a systems point of view, and do not affect the size of the solution space. Redundant
and
constraints that were not needed to define the steady-state flux space (i.e., these redundant constraints lay outside of more constraining
and
constraints) were eliminated and not needed for the generation of the sample points.
Choice of enclosing parallelepiped
Because each pair of
and
constraints form parallel hyperplanes, the shape of the null space leads naturally to the choice of a high-dimensional parallelepiped in which to enclose it. The set of possible parallelepipeds that could be used to enclose the steady-state flux space was chosen by forming the faces of the parallelepiped along the directions defined by these
and
constraints. Since each parallelepiped is defined by r planes which are chosen from the set of m
and
planes, the number of such parallelepipeds that could be used enclose the space is
![]() | (6) |
Minimizing the volume of the enclosing parallelepiped
Checking the volume of every possible parallelepiped to find the smallest can become prohibitive due to the large number of parallelepipeds (i.e., Eq. 6). Therefore, an alternate approach was used. This algorithm chooses its first direction based on the set of
and
constraints that are the closest together based on Euclidian distance. Then, the next direction is chosen by determining the smallest parallelogram that can be formed by choosing the next set of
and
constraints. The third direction is chosen as the set of constraints that forms the smallest parallelepiped formed using three sets of
and
constraints and so on until the parallelepiped fully encloses the solution space (Fig. 1 D).
Uniform random sampling of points
Uniform random points were generated within the solution space by randomly sampling within the enclosing parallelepiped. Randomly sampling within the parallelepiped was accomplished by uniformly sampling along each of the edges, bi, defining the parallelepiped and picking the point resulting from the sum of this set of weightings. Each point in the space is uniquely defined by weightings on the edges spanning the parallelepiped (Eq. 3). The weighting,
i, on each basis vector, bi, was uniformly selected by generating a pseudo-random number, f, between 0 and 1 and then assigning each weight as
![]() | (7) |
and
constraints to verify whether the point falls in the solution space or not. If the point satisfies all constraints, it is a valid solution and is kept in the set. If the randomly generated point does not satisfy all the necessary constraints, it is excluded. The fraction of total points generated that fall within the space, the "hit" fraction, was used to calculate the absolute volume of the steady-state flux space (see Fig. 1 D).
Volume calculation of steady-state flux space
The volume of the steady-state flux space can be calculated by multiplying the volume of the enclosing parallelepiped by the fraction of generated points that falls within the solution space
![]() | (8) |
2, is given as
![]() | (9) |
is the ratio of sampled points that fall inside the solution space (hit fraction) and n is the total number of sample points generated. The relative error of the hit fraction estimate,
, is calculated as the ratio of the standard deviation,
, to the mean, µ,
![]() | (10) |
, as would be expected.
Red blood cell metabolic network
The red blood cell metabolic network used in this study consists of 48 reactions utilizing 39 metabolites. Reversible reactions were not decoupled into forward and reverse reactions, and thus were allowed to take on negative values as discussed above. The dimension of the null space of this stoichiometric matrix was 11. This red blood cell network differs slightly from the one previously used to study extreme pathways (Wiback and Palsson, 2002
) in that the metabolic loads are represented as turnover reactions, rather than as exchange fluxes and also differs from the network studied in Wiback et al. (2004)
in that it contains reversible reactions. The
values utilized in this study were taken from Wiback and Palsson (2002)
, with the
of all reactions without a stated
in Wiback and Palsson (2002)
being set at an arbitrarily high value of 1000 so that none of these reactions was limiting on the system. Thus, all reaction fluxes were limited by either a measured
value or an uptake rate. The
values were set to
for reversible reactions and to zero for all irreversible reactions, unless there was a minimum physiological demand, as detailed in the next section.
Physiologic conditions
In addition to limitations on flux values due to
values, certain known physiologic demands were used to identify the range of potential flux values. For example, the red blood cell is obligated to produce a basal level of ATP to run the sodium potassium pump to balance against the natural diffusion rate of sodium and potassium that "leaks" through the membrane. Thus, the minimum ATP production rate in this study was set to 1 mM/h to reflect this fact. The maximum ATP production rate was set at 1.5, as a conservative bound, higher than the maximum value observed in a wide range of calculations based on a full-scale kinetic model (Jamshidi et al., 2001
). Also, a minimal amount of NADH production is necessary to convert Met-hemoglobin into its functional state. The minimum NADH production rate was set to 0.4 mM/h based on the oxidation rate of iron in hemoglobin. This minimal level of NADPH production was set to 0.05 mM/h for the turnover of GSSG to GSH to combat a minimal level of reactive oxygen species. Lastly, the minimum flux through DPGase was set to 0.3 mM/h, since this flux always operates near saturation levels because of the high intracellular concentration of 23DPG and the slow rate of this reaction (time constant on the order of half a day).
Convergence of statistics with increasing samples
To determine that sufficiently large samples were taken to accurately assess each of the statistical properties computed, samples were taken until the statistics being assessed ceased changing with increased sample size relative to the error deemed appropriate for the property being calculated. At the sample sizes used in this study, completely new sets of random samples were also taken with the results remaining unchanged for the significant figures related, and this resampling was done numerous times. Thus, as should be the case, all results presented herein are independent of the particular set of uniform random samples used, and are generic properties of the metabolic solution space. Specific details on the number of repetitions and the size of the samples are given along with each of the calculated results in the figure captions.
Computation and implementation
Computations for this study were performed on Dell Workstations (either a Dell Precision 340 or Dimension 8200) or on a Linux box (Dual Athlon MP 2400, 2 GB RAM). The program for finding the enclosing parallelepiped and for generating points was performed using MATLAB (The MathWorks, Natick, MA) and an interface with the linear programming package LINDO API (Lindo Systems, Chicago, IL). As an indication of the efficiency of the calculations, the computation of 250,000 uniform random samples within the red blood cell steady-state flux space using the Dell workstations was performed in <30 s. The calculation of 1,000,000 points inside the steady-state flux space was performed on the Linux box in <50 s. The approach herein was very fast compared to the previous sampling approach used in Wiback et al. (2004)
, which describes the calculation of 250,000 uniform points in the red blood cell metabolic network as taking "over a week" of calculation (Dell Dimension 8200). Thus, the approach presented herein represents an
20,000-fold increase in calculation efficiency over the approach used in Wiback et al. (2004)
. Increased sampling efficiency is important because it allows for more detailed and precise calculations, increased capacity to study improbable regions of the steady-state flux space, and better error assessment.
| RESULTS |
|---|
|
|
|---|
Distribution of flux values and space segmentation
The set of candidate steady-state flux distributions through each reaction in the red blood cell metabolic network were represented as a histogram of all possible flux values and displayed on the metabolic map (Fig. 2). Each histogram presents one-dimensional information on its x axis, in terms of the extent of possible values for that particular flux. The y axis represents the "size" of space in the other r1 dimensions resulting from slicing the metabolic solution space along a specific value of the flux through the indicated reaction. Thus the "height" of each histogram represents r1 dimensional data (see illustration in Fig. 1 D).
|
|
and
is the most probable with the upper and lower possible flux values approaching a zero probability, such as the fluxes through reactions PFK, PK, and LD. Histograms for reactions that do not contain a zero flux value are essential under the conditions examined. None of these histograms can have more than one peak due to the convexity of the steady-state flux space. The shape of these histograms is highly informative about how a reduction in the allowable range of each flux would affect the remaining size of the steady-state flux space. The percentage of the histogram area remaining in the reduced range is the same as the percentage of the remaining steady-state flux space.
Correlation of flux measurements: use in designing experiments
Uniform random sampling of the steady-state flux space allows for the calculation of the correlation coefficient (rij) between any two fluxes (vi and vj) in the network. Thus, sampling provides a straightforward means of not only calculating perfectly correlated subsets (
), but also of identifying well-correlated, but not perfectly correlated reaction sets. The matrix of squared pairwise correlation coefficients for all the RBC metabolic fluxes was computed. The fluxes can be ordered such that the "correlated reaction sets" (defined here as
) are listed in order of decreasing number of fluxes in each set (Table 1).
|
|
A second important factor to consider when choosing fluxes for measurement to best determine the flux state of the red blood cell is to measure a set of fluxes that are uncorrelated to each other (
values close to zero). Such a set ensures that the flux measurements are not providing duplicate information. Thus, by iteratively choosing to measure fluxes that are uncorrelated to each other and to what has already been measured, the correlation matrix guides the selection of independent flux measurements.
A third factor that determines how informative an experimental measurement will be is where in the solution space the measurement falls. For example, if a flux measurement falls in an improbable region of the histogram, it will be much more constraining for determining the flux state of the system than a measurement taken in a probable region. The numerical range of flux measurement is generally unknown before the measurement is taken, and thus the reduction in the size of the solution space is known only after the experimental measurement is taken. However, this information can be taken into account using previous experience and intuition to evaluate if a flux measurement is likely to occur in a region deemed improbable by, and thus highly informative to, the current status of the model.
Thus, the Monte Carlo sampling procedure provides three important criteria for designing a set of informative measurements to determine the state of the system:
Systemic effects of simulated enzymopathies
The sampling procedure can be used to track the networkwide changes that occur due to changes to a single, or a small number, of individual reactions. For example, a substantial number of single nucleotide polymorphisms (SNPs) have been found in genes of red blood cell enzymes, which significantly decrease the
values through many reactions in the red blood cell (Jacobasch and Rapoport, 1996
; Tanaka and Zerez, 1990
). Some of these SNPs have been correlated to chronic and nonchronic anemia (Grimes, 1980
; Thorburn and Kuchel, 1985
) and the systemic effects of the two most common SNPs in the red blood cell, pyruvate kinase (PK) and glucose-6-phosphate dehydrogenase (G6PDH), have recently been evaluated using a kinetic model of red cell metabolism (Jamshidi et al., 2002
).
Defects in glycolytic enzymes were simulated to demonstrate how a 75% reduction in their normal operating range (though a decreased
) would affect the network. For each reaction, an enzymopathy was simulated by setting a reduced maximum flux for the reaction,
as
![]() | (11) |
and
were the maximum and minimum flux values achievable in the unaltered system. The percentage of the steady-state flux space remaining after restricting each of these reactions to one-quarter their normal operating range is shown in Table 2. As discussed above, the shape of the flux histograms, shown in Fig. 2, are a means for immediately determining the impact of a reduction in the allowable range of any flux. Such impact can be evaluated by comparing the reduction in the size of the steady-state flux space seen in Table 2 with the histograms shown in Fig. 2. A reduction in the
of a reaction that has a right peak distribution will reduce the size of the steady-state flux space more than a reduction in the
of a reaction with a left peak distribution.
|
The
through PK was reduced such that the allowable flux range is decreased first to one-half and then to one-quarter of the original range. The systemic effects of the simulated enzymopathy through the PK reaction greatly affected the probability distributions through many reactions in the red blood cell metabolic network (Fig. 3). For example, the probability distribution for glucose uptake changed from being a right peak distribution to being a central peak distribution, meaning that its normal operating range in the red blood cell changed from being in its most probable region to being in an improbable operating region of the steady-state flux space. With the PK range constrained to one-quarter of the original range, the maximum possible value for glucose uptake also decreased, making the previous upper range impossible to the network.
Constraining the range of allowable PK values to one-quarter of the original range significantly changed the correlations between many metabolic fluxes (Table 3). For example, the correlation between glucose uptake and the flux through the pentose phosphate pathway increased dramatically, with the r2 value between HK and G6PDH increasing by 0.79. This change in correlation is due to the limitation on the steady-state flux allowable through glycolysis due to the decreased capacity of PK. The increase in the correlation between two fluxes can be seen in the example of G6PDH versus HK and GAPDH versus TPI (Fig. 4). A high correlation coefficient implies that the shape of the two-dimensional histogram will be narrow, whereas a low correlation coefficient implies that the two-dimensional histogram will be broad. In the case of G6PDH versus HK, the PK enzymopathy causes the r2 value to increase significantly. In the case of GAPDH versus TPI, the r2 value decreases significantly with the simulated PK enzymopathy. Also of interest, the direction of the ridge representing the most probable values with the simulated PK enzymopathy changes direction significantly compared to the dominant direction of the correlation without the PK enzymopathy.
|
| DISCUSSION |
|---|
|
|
|---|
will have on the size of the steady-state flux space, and thus on how constraining to the system a reduced flux capacity through a specific reaction would be.
The degree of dependence between all of the fluxes in the red blood cell metabolic network was determined by calculating correlation coefficients from the uniform random samples. Different methods have been used to calculate correlated subsets in a metabolic network, including extreme pathway analysis (Papin et al., 2002
; Price et al., 2002
), elementary mode analysis (Pfeiffer et al., 1999
), and the linear programming-based flux-coupling finder (Burgard et al., 2004
). The correlated reaction subsets calculated by extreme pathway and elementary modes correspond only to sets that are perfectly correlated. The sets calculated using the flux-coupling finder (Burgard et al., 2004
) are more informative and scalable to larger networks. In contrast to previous methods, the direct calculation of correlation coefficients from the Monte Carlo sampling procedure allowed for the stratification of pairwise correlations among all fluxes between 0 and 1 (
) or between 1 and 1 (rij). Thus, the unbiased degree of independence can be determined for any two fluxes under any set of conditions.
Experimental design can be guided from the results of uniform random sampling within a constrained metabolic solution space. Both the probability distributions and the calculated correlation coefficients are condition-dependent, and thus change with the addition of more experimental data. One implication of this fact is that it is best to include as much known information as possible into the model before sampling. Subsequently, the model predicts which combinations of measurements provide independent information and which fluxes are correlated to the highest number of other fluxes under the studied conditions. Since the choice of fluxes to measure providing independent information can change with the experimental conditions, it is important that the Monte Carlo sampling procedure be used iteratively to determine which flux measurements are the most independent at each step. Thus, following the results of one set of experimental measurements, the sampling procedure can be redone to reevaluate the most important measurements given what is already known. Of course, those fluxes that are perfectly correlated (rij = 1) will remain perfectly correlated, regardless of condition. The volume of the solution spaces also changes, giving a quantitative measure for how many potential steady-state fluxes are in agreement with the measured data. As the volume of the steady-state flux space approaches zero, the internal flux distribution inside a cell is determined. Other constraint-based methods also exist for guiding experimental design and could perhaps be used in conjunction with the method proposed herein. One such method enables the calculation of the optimal set of fluxes to measure to minimize the impact of experimental error on the prediction of the steady-state flux state (Savinell and Palsson, 1992a
,b
).
The effects of enzymopathies on the capabilities of a metabolic network can be studied using the Monte Carlo sampling procedure. Single nucleotide polymorphisms (SNPs) or other genetic defects can impair enzyme function. This impairment can be the result of such factors as an enzyme having a lowered
or by a lower rate constant which can effectively lower the allowable maximum flux based on the concentration of the substrates in the red blood cell. Importantly, since the probability distributions of individual fluxes and the correlations between them change under differing conditions, the Monte Carlo sampling procedure can be used to study both normal and pathological cases under a variety of different environments. Clinical outcomes of enzymopathies can obviously involve a great many factors outside of those modeled. However, model-driven studies describe the impact of an enzymopathy on the known metabolic network uncoupled from other considerations. Thus, model-driven assessments provide evidence for whether or not a metabolic explanation is sufficient to account for a clinical outcome or not. Indeed, the apparent lack of an explanation for a clinical outcome based on what is known about a metabolic network can lead to novel biological understanding and provide a basis for novel hypotheses. One example of this is the observation that the maximal enzymatic activity of G6PDH vastly exceeds what seems to be needed in the network (Salvador and Savageau, 2003
).
Constraining the flux range for the glycolytic enzymes showed a general trend seen that those simulated enzymopathies that severely restricted the steady-state flux space had been shown experimentally to cause anemia, whereas those that were less restrictive did not (Table 2). The exception to this trend was lactate dehydrogenase (LDH), the inhibition of which greatly restricted the steady-state flux space, but which has not been found to lead to anemia (Tanaka and Zerez, 1990
). This may be because restriction of LDH does not significantly restrict the glycolytic rate, since pyruvate can leave the cell as the endpoint of glycolysis. However, the restriction of LDH does affect the cell's ability to control NADH levels. Another observation was that, although G6PDH is by far the most common enzymopathy in the red blood cell, only a small fraction of G6PDH enzymopathies lead to anemia. Although fewer in number, a higher fraction of PK enzymopathies than G6PDH enzymopathies lead to anemia. This fact may be attributable, in part, to the fact that restricting PK is generally much more constraining on the metabolism as a whole than is G6PDH.
Strong similarities were seen when comparing the results of the simulated enzymopathies with results from an interesting earlier study using a full-scale kinetic model of glycolysis (Holzhütter et al., 1985
). In Holzhütter et al. (1985)
, a kinetic model of glycolysis was used to study the impact of lowering enzyme activity in a "metabolic homeostasis function." Of the reactions studied in Holzhütter et al. (1985)
, one of them, HK, was found using the Monte Carlo sampling procedure to have a right peak distribution, trailing off to zero probability as the flux level through this reaction decreased, and another of the reactions, AK, was left peak, trailing off to zero probability as the flux level through this reaction increased. The kinetic model predicted that constricting HK had the most impact on the homeostasis of the red blood cell and that restrictions of AK had the least, just as predicted with the distributions given with the Monte Carlo sampling procedure. The fluxes shown by the kinetic model to have effects on homeostasis between the extremes of HK and AK, all had either central peak distributions, or were constrained such that the difference in probabilities over the flux range were not highly variable, matching expectations from the sampling procedure. Thus, results of the Monte Carlo sampling procedure corresponded well to results from a kinetic model which requires measurement of a large number of kinetic enzyme data to formulate. The need for few kinetic parameters becomes essential as the Monte Carlo sampling procedure is applied to organisms where kinetic data is much sparser than in the well-characterized red blood cell.
The uniform distribution of points within the solution space used in this study is not meant to imply that the flux steady states of cells in a population are likewise distributed uniformly within the range of allowable steady states. Indeed, it seems highly unlikely that a population of cells would be distributed uniformly in this manner, since certain regions of the solution space are expected to be preferred or excluded based on additional, unmodeled demands and physicochemical constraints under which cells operate. Rather, uniformity is used to clearly grasp the implications of applying the physicochemical constraints to the reconstructed metabolic network. Although not utilized for the purposes of this study, the Monte Carlo sampling procedure could potentially be used to study distributions of populations in cells. The informative aspect of such an approach would be to see how the population density deviated from a uniform distribution within the range of allowable steady states for the in silico cell. By identifying which portions of the steady-state flux space were favored within the population of cells, hypotheses for the preferred states could be formed and the importance of additional physicochemical constraints could be assessed.
Although the present study focused on studying the steady-state flux space, the approach detailed herein can be equally well applied to study concentration or kinetic spaces associated with a reconstructed biochemical network. Also any additional physicochemical law can be used to further confine sets of candidate solutions. If the imposed constraints can be represented as linear equations, they can be easily implemented into this framework and, aside from computational limitations, uniform random sampling can easily be done within the space. Once the set of candidate solutions are generated, any type of nonlinear constraint, such as stemming from regulation or thermodynamics, can be applied as a postprocessing step (as long as dimensionality of the space is not reduced) by eliminating candidate solutions that do not satisfy the imposed constraint. Although "elimination" algorithms such as the one used in this study have certain advantages for sampling smaller networks, it is likely that any "elimination" approach to sampling will be inadequate for sampling genome-scale networks. This difficulty in generating samples in high-dimensional objects occurs because the ratio of the size of an enclosing object to the size of the enclosed object increases rapidly as dimension increases. However, genome-scale networks can be sampled (Almaas et al., 2004
) using alternate, but similar, methods for generating the set of candidate solutions, such as Monte Carlo Markov-chain samplers (Chen and Schmeiser, 1993
; Kaufman and Smith, 1998
; Lovasz, 1999
; Zabinsky et al., 1993
), and such methods have been well studied and continue to improve. Thus, the methods described herein for analyzing uniform random samples from solution spaces can be utilized for studying genome-scale networks.
Taken together, uniform random sampling provides a broadening of the constraint-based approach by allowing for the unbiased and detailed assessment of the capabilities of reconstructed networks subject to the imposition of physicochemical constraints. Uniform random sampling of solution spaces allows for a detailed framework in which to design experiments to determine the operating state of the cell, a context in which to study the resulting experimental data, as well as detailed insight into what network behaviors are allowed by the imposition of the stated constraints upon the reconstructed network. Clearly understanding model predictions is important because a cell-scale model provides the most compact and quantitative representation of a large-scale hypothesis of how a cell works. Thus, the Monte Carlo sampling procedure coupled with ongoing network reconstruction provides a powerful engine to drive experimental work and the iterative model-building process that is at the heart of systems biology.
| APPENDIX: ABBREVIATIONS |
|---|
|
|
|---|
|
TABLE 5 Reactions
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
| ACKNOWLEDGEMENTS |
|---|
|
|
|---|
Support for this research was generously provided by a grant from the National Science Foundation (NSF/BES-01-20363).
Submitted on March 18, 2004; accepted for publication June 29, 2004.
| REFERENCES |
|---|
|
|
|---|
Alves, R., and M. A. Savageau. 2000a. Comparing systemic properties of ensembles of biological networks by graphical and statistical methods. Bioinformatics. 16:527533.
Alves, R., and M. A. Savageau. 2000b. Extending the method of mathematically controlled comparison to include numerical comparisons. Bioinformatics. 16:786798.
Alves, R., and M. A. Savageau. 2000c. Systemic properties of ensembles of metabolic networks: application of graphical and statistical methods to simple unbranched pathways. Bioinformatics. 16:534547.
Bonarius, H. P. J., G. Schmid, and J. Tramper. 1997. Flux analysis of underdetermined metabolic networks: the quest for the missing constraints. Trends Biotechnol. 15:308314.[CrossRef]
Burgard, A. P., and C. D. Maranas. 2001. Probing the performance limits of the Escherichia coli metabolic network subject to gene additions or deletions. Biotechnol. Bioeng. 74:364375.[CrossRef][Medline]
Burgard, A. P., and C. D. Maranas. 2003. Optimization-based framework for inferring and testing hypothesized metabolic objective functions. Biotechnol. Bioeng. 82:670677.[CrossRef][Medline]
Burgard, A. P., E. V. Nikolaev, C. H. Schilling, and C. D. Maranas. 2004. Flux coupling analysis of genome-scale metabolic network reconstructions. Genome Res. 14:301312.
Burgard, A. P., S. Vaidyaraman, and C. D. Maranas. 2001. Minimal reaction sets for Escherichia coli metabolism under different growth requirements and uptake environments. Biotechnol. Prog. 17:791797.[CrossRef][Medline]
Chen, M. H., and B. Schmeiser. 1993. Performance of the Gibbs, hit-and-run and Metropolis samplers. J. Comput. Graph. Stat. 2:251272.[CrossRef]
Covert, M., I. Famili, and B. Palsson. 2003. Identifying the constraints that govern cell behavior: a key to converting conceptual to computational models in biology? Biotechnol. Bioeng. 84:763772.[CrossRef][Medline]
Edwards, J. S., and B. Palsson. O.2000. Metabolic flux balance analysis and the in silico analysis of Escherichia coli K-12 gene deletions. BMC Bioinformat. 1:1. Epub 2000 Jul 27.[CrossRef]
Edwards, J. S., M. Covert, and B. Palsson. 2002. Metabolic modelling of microbes: the flux-balance approach. Environ. Microbiol. 4:133140.[CrossRef][Medline]
Edwards, J. S., R. U. Ibarra, and B. O. Palsson. 2001. In silico predictions of Escherichia coli metabolic capabilities are consistent with experimental data. Nat. Biotechnol. 19:125130.[CrossRef][Medline]
Edwards, J. S., and B. O. Palsson. 1999. Systems properties of the Haemophilus influenzae Rd metabolic genotype. J. Biol. Chem. 274:1741017416.
Famili, I., J. Forster, J. Nielsen, and B. O. Palsson. 2003. Saccharomyces cerevisiae phenotypes can be predicted by using constraint-based analysis of a genome-scale reconstructed metabolic network. Proc. Natl. Acad. Sci. USA. 100:1313413139.
Forster, J., I. Famili, B. O. Palsson, and J. Nielsen. 2003. Large-scale evaluation of in silico gene deletions in Saccharomyces cerevisiae. Ohmics. 7:193202.
Grimes, A. J. 1980. Human Red Cell Metabolism. Blackwell Scientific Publications, Oxford, UK.
Hasty, J., D. McMillen, and J. J. Collins. 2002. Engineered gene circuits. Nature. 420:224230.[CrossRef][Medline]
Hasty, J., D. McMillen, F. Isaacs, and J. J. Collins. 2001. Computational studies of gene regulatory networks: in numero molecular biology. Nat. Rev. Genet. 2:268279.[Medline]
Holzhütter, H. G., G. Jacobasch, and A. Bisdorff. 1985. Mathematical modelling of metabolic pathways affected by an enzyme deficiency. A mathematical model of glycolysis in normal and pyruvate-kinase-deficient red blood cells. Eur. J. Biochem. 149:101111.[Medline]
Ibarra, R. U., J. S. Edwards, and B. O. Palsson. 2002. Escherichia coli K-12 undergoes adaptive evolution to achieve in silico predicted optimal growth. Nature. 420:186189.[CrossRef][Medline]
Jacobasch, G., and S. M. Rapoport. 1996. Hemolytic anemias due to erythrocyte enzyme deficiencies. Mol. Aspects Med. 17:143170.[CrossRef][Medline]
Jamshidi, N., J. S. Edwards, T. Fahland, G. M. Church, and B. O. Palsson. 2001. Dynamic simulation of the human red blood cell metabolic network. Bioinformatics. 17:286287.
Jamshidi, N., S. J. Wiback, and B. O. Palsson. 2002. In silico model-driven assessment of the effects of single nucleotide polymorphisms (SNPs) on human red blood cell metabolism. Genome Res. 12:16871692.
Kauffman, K. J., J. D. Pajerowski, N. Jamshidi, B. O. Palsson, and J. S. Edwards. 2002. Description and analysis of metabolic connectivity and dynamics in the human red blood cell. Biophys. J. 83:646662.
Kaufman, D. E., and R. L. Smith. 1998. Direction choice for accelerated convergence in hit-and-run sampling. Opt. Res. 46:8495.
Lew, V. L., and R. M. Bookchin. 1986. Volume, pH, and ion-content regulation in human red cells: analysis of transient behavior with an integrated model. J. Membr. Biol. 92:5774.[CrossRef][Medline]
Lovasz, L. 1999. Hit-and-run mixes fast. Math. Program. 86:443461.[CrossRef]
Meyer, C. D. 2000. Matrix Analysis and Applied Linear Algebra. Society for Industrial and Applied Mathematics, Philadelphia, PA.
Mulquiney, P. J., W. A. Bubb, and P. W. Kuchel. 1999. Model of 2,3-bisphosphoglycerate metabolism in the human erythrocyte based on detailed enzyme kinetic equations: in vivo kinetic characterization of 2,3-bisphosphoglycerate synthase/phosphatase using 13C and 31P NMR. Biochem. J. 342:567580.[CrossRef][Medline]
Mulquiney, P. J., and P. W. Kuchel. 1999. Model of 2,3-bisphosphoglycerate metabolism in the human erythrocyte based on detailed enzyme kinetic equations: computer simulation and metabolic control analysis. Biochem. J. 342:597604.[CrossRef][Medline]
Mulquiney, P. J., and P. W. Kuchel. 2003. Modeling Metabolism with Mathematica: Detailed Examples including Erythrocyte Metabolism. CRC Press, Boca Raton, Fl.
Papin, J. A., N. D. Price, and B. O. Palsson. 2002. Extreme pathway lengths and reaction participation in genome-scale metabolic networks. Genome Res. 12:18891900.
Papin, J. A., N. D. Price, S. J. Wiback, D. A. Fell, and B. O. Palsson. 2003. Metabolic pathways in the post-genome era. Trends Biochem. Sci. 28:250258.[CrossRef][Medline]
Pfeiffer, T., I. Sanchez-Valdenebro, J. C. Nuno, F. Montero, and S. Schuster. 1999. METATOOL: for studying metabolic networks. Bioinformatics. 15:251257.
Price, N. D., J. A. Papin, and B. O. Palsson. 2002. Determination of redundancy and systems properties of Helicobacter pylori's metabolic network using genome-scale extreme pathway analysis. Genome Res. 12:760769.
Price, N. D., J. A. Papin, C. H. Schilling, and B. O. Palsson. 2003a. Genome-scale microbial in silico models: the constraints-based approach. Trends Biotechnol. 21:162169