- © 2009 American Society of Plant Biologists
Abstract
The accumulation of storage compounds is an important aspect of cereal seed metabolism. Due to the agronomical importance of the storage reserves of starch, protein, and oil, the understanding of storage metabolism is of scientific interest, with practical applications in agronomy and plant breeding. To get insight into storage patterning in developing cereal seed in response to environmental and genetic perturbation, a computational analysis of seed metabolism was performed. A metabolic network of primary metabolism in the developing endosperm of barley (Hordeum vulgare), a model plant for temperate cereals, was constructed that includes 257 biochemical and transport reactions across four different compartments. The model was subjected to flux balance analysis to study grain yield and metabolic flux distributions in response to oxygen depletion and enzyme deletion. In general, the simulation results were found to be in good agreement with the main biochemical properties of barley seed storage metabolism. The predicted growth rate and the active metabolic pathway patterns under anoxic, hypoxic, and aerobic conditions predicted by the model were in accordance with published experimental results. In addition, the model predictions gave insight into the potential role of inorganic pyrophosphate metabolism to maintain seed metabolism under oxygen deprivation.
Cereal seeds are worldwide a major source for food and feed. The main storage compounds of cereal seeds are starch and storage proteins (Crocker and Barton, 1957; van Dongen et al., 2004), which are mainly derived from Suc and amino acids and synthesized in the endosperm (Bewley and Black, 1994). Due to the agronomical importance of seeds, seed metabolism has been studied intensively with respect to biochemistry, physiology, and genetics (Finnie et al., 2004; Lai et al., 2004; Sreenivasulu et al., 2004; Ho et al., 2005; McIntosh et al., 2007) and has also been a target for metabolic engineering to improve nutritional value (Mazur et al., 1999; Ye et al., 2000; Houmard et al., 2007; Wang et al., 2007).
To identify suitable targets for metabolic engineering, mathematical modeling of metabolism is particularly promising as it offers systems approaches to analyze the structure, dynamics, and behavior of complex metabolic networks. In plant research, the issue of modeling metabolism is constantly gaining attention, and several mathematical modeling approaches applied to plant metabolism exist (for reviews, see Giersch, 2000; Morgan and Rhodes, 2002; Poolman et al., 2004; Rios-Estepa and Lange, 2007).
Flux balance analysis (FBA) is a constraint-based modeling approach that allows the prediction of metabolic steady-state fluxes by applying mass balance constraints to a stoichiometric model (Edwards et al., 1999). FBA has the advantage of not requiring the knowledge of kinetic parameters; instead, only the stoichiometry of the metabolic network has to be known, and an objective function is needed to identify the optimal flux distribution among all possible steady-state flux distributions. The flux balance approach has been applied to a variety of different biological systems, such as bacteria (Schilling et al., 2002; Reed and Palsson, 2003), fungi (David et al., 2003; Famili et al., 2003), algae (Shastri and Morgan, 2005), and animals (Fell and Small, 1986; Ramakrishna et al., 2001; Cakir et al., 2007; Duarte et al., 2007), to study different aspects of metabolism, including the prediction of optimal metabolic yields and flux distributions (Varma et al., 1993a, 1993b), gene deletion lethality (Edwards and Palsson, 2000; Förster et al., 2003), and pathway redundancies (Van Dien and Lidstrom, 2002). Although metabolic flux determination is acknowledged to be an important part of plant metabolic engineering and will help to improve the understanding of plant metabolic networks (Schwender et al., 2004; Fernie et al., 2005), to the best of our knowledge no detailed FBA model for the prediction of metabolic fluxes in plant metabolic systems has been published as yet.
In this study, a computational analysis of cereal seed storage metabolism was performed with the major aim of getting insight into storage patterning in response to environmental and genetic perturbations. To provide a framework for in silico analysis, a stoichiometric model of primary metabolism in developing endosperm of barley (Hordeum vulgare) during starch accumulation was constructed. Barley has been studied extensively with respect to seed storage metabolism (Weschke et al., 2000; Rolletschek et al., 2004; Sreenivasulu et al., 2004; Wobus et al., 2005) and serves as a model plant for temperate cereals (Wobus et al., 2005), which makes it particularly suitable for modeling cereal seed metabolism. FBA was applied to the model to study (1) the effect of oxygen depletion (hypoxia) on grain yield and metabolic flux distribution, and (2) grain growth in response to enzyme deletion. To validate the model, it is shown that the simulation results are consistent with experimental results described in the literature.
RESULTS
A Model of Cereal Seed Storage Metabolism
With the aim of getting a systemic understanding of cereal seed storage metabolism, we constructed a metabolic network of primary metabolism in the developing endosperm of barley seed during starch accumulation by integrating biochemical, physiological, proteomic, and genomic data derived from the literature and databases (Supplemental Data Sets S1 and S2). The resulting stoichiometric model includes 234 metabolites and 257 reactions, of which 192 represent biochemical conversions and 65 represent transport processes. The reactions are compartmentalized between the extracellular medium and the intracellular compartments cytosol, mitochondria, and amyloplast. Table I lists the network characteristics of the constructed model. Key features of the metabolic network are described below.
Network characteristics
Number of metabolites and reactions included in the model, and their locations.
In developing barley seeds, the primary substrate for the accumulation of storage compounds (starch and storage proteins) is Suc (Sreenivasulu et al., 2004), which is mostly cleaved by Suc synthase during the main storage phase (Weschke et al., 2000). The products of Suc cleavage are either metabolized through glycolysis and/or incorporated into starch. The key intermediate of the starch biosynthesis pathway, ADP-Glc, is synthesized from hexose phosphates by the cytosolic isoform of ADP-Glc pyrophosphorylase and transported into the amyloplast for starch synthesis via an ADP-Glc transporter (Thorbjørnsen et al., 1996; James et al., 2003; Patron et al., 2004).
Seed-specific amino acid synthesis, a key determinant for storage protein accumulation, is based on the import of Asn and Gln, which provide the nitrogenous component for de novo amino acid biosynthesis (Bewley and Black, 1994). Storage product accumulation is accompanied by high transcriptional activity of genes involved in energy metabolism (i.e. glycolysis, fermentation, and citrate cycle; Sreenivasulu et al., 2004). In addition to providing energy and reducing equivalents, these pathways yield precursors for Glc and amino acid metabolism (Owen et al., 2002).
Further information about the definition of the metabolic model can be found in “Materials and Methods.” A metabolic map of the network (Supplemental Fig. S1), the set of reactions included in the network (Supplemental Data Set S1), and the terminology used for the network reactions and metabolites (Supplemental Data Set S3) are given in the online supplemental material.
Seed Storage Metabolism in Response to Oxygen Limitation
Barley seeds are known to develop under hypoxic conditions during intermediate and storage phases (Rolletschek et al., 2004). To elucidate the possible role of oxygen and Suc supply for storage patterning in developing barley seeds, a phenotypic phase plane (PhPP; Edwards et al., 2002) was computed that depicts the metabolic behavior of seed metabolism at various levels of oxygen and Suc availability (Fig. 1 ). PhPP analysis is a method used to obtain the range of optimal flux distribution (i.e. optimal phenotypic behavior) in response to changes in two environmental parameters such as substrate and oxygen uptake rates. By defining the environmental parameters under investigation as two axes on an (x,y) plane, the optimal flux distribution is computed for all points in the plane and lines are drawn to demarcate regions of constant flux distributions (i.e. phenotypic phases). Based on this procedure, the optimal flux distributions in the phase plane are classified into a finite number of regions, each with a distinct metabolic pathway utilization pattern, corresponding to a different optimal phenotypic behavior (Edwards et al., 2001, 2002). Thus, PhPP allows a quick and comprehensive overview of the optimal use of metabolism in the growth environment studied (see “Materials and Methods” for details).
PhPP analysis. The PhPP revealing the dependence of seed storage metabolism on Suc and oxygen supply. To obtain the range of optimal flux distributions (i.e. optimal phenotypic behavior) in response to changes in the Suc and oxygen supply, PhPP analysis was performed by representing the Suc and oxygen uptake rates as two axes on the (x,y) plane and by computing the optimal flux distributions for all points in the plane. Based on the demarcation of regions of constant flux distributions, the optimal flux distributions in the PhPP are classified into a finite number of regions, each with a distinct metabolic pathway utilization pattern, corresponding to a different optimal phenotypic behavior. In the given example, the PhPP is divided into five phases (P1–P5), each characterized by a distinct optimal use of seed metabolism. The LO, shown as a dashed line, defines the optimal ratio of the substrate uptake rates for maximal biomass yield. The dotted lines depict the growth conditions used to perform simulation scenarios 1 (Sim 1) and 2 (Sim 2). Dots on these lines correspond to the substrate uptake rates (Suc uptake rate and oxygen uptake rate) used to simulate the metabolic flux distribution characterizing a given phase (see “Materials and Methods” for details).
As shown in Figure 1, the two-dimensional PhPP is divided into five phases (P1–P5), each one characterized by a distinct metabolic pathway utilization pattern (Edwards et al., 2001). Phases 1 to 4 are characterized by dual substrate limitation. In these phases, an increase in the supply of either oxygen or Suc will increase biomass production. Phase 5 is defined as a futile phase (Edwards et al., 2001), with the metabolic behavior in this region being characterized by futile cycles. The physiological relevance of these carbon cycles, which appear to waste energy, is still controversial, and many studies have proposed various physiological roles, as reviewed by Portais and Delort (2002). In the futile phase, an increase in Suc supply will increase biomass production, while an increase in oxygen supply will decrease the growth rate. The line of optimality (LO), which lies on the boundary between P4 and P5, defines the optimal ratio of the substrate uptake rates for maximal biomass production. For aerobic growth, LO represents the optimal oxygen uptake required for the complete oxidation of Suc to maximize biomass production.
Grain Growth
To study grain growth with respect to oxygen limitation, the Suc uptake rate was fixed at 8 μmol g−1 dry weight h−1, a value derived from experimental reports (Felker et al., 1984), and the oxygen uptake rate was varied from completely anaerobic to fully aerobic conditions (Fig. 2 ). The in silico analysis shows that the growth rate linearly increases with increasing oxygen uptake rate (Fig. 2, P1–P4). The maximal growth rate of 0.003 h−1 is reached under optimal growth conditions (Fig. 2, LO). Further increasing the oxygen uptake rate results in a decrease of the growth rate due to futile cycles (Fig. 2, P5). The predicted growth rates are in the scope of published experimental results (Felker et al., 1983; Quarrie et al., 1988), which range from 0.003 to 0.007 h−1 at 14 d after anthesis.
Grain growth depending on oxygen supply varying from completely anaerobic to aerobic conditions. Simulations were performed with the Suc uptake rate fixed at 8 μmol g−1 dry weight h−1. The vertical lines demarcate the phenotypic phases P1 to P5. The LO, shown as a dotted line, defines the optimal ratio of the substrate uptake rates for maximal biomass yield.
Metabolic Flux Maps
To illustrate the metabolic behavior characterizing the different phases, simulations were run with the substrate uptake rates fixed at representative values for each phase as described in “Materials and Methods” (simulation scenario 1). The corresponding metabolic flux distribution patterns, represented as carbon flux, are shown in Figure 3 . A table listing the simulated flux values for each phase is given in Supplemental Data Set S7.
Carbon flux maps depicting the key uptake/excretion rates and fluxes within central metabolism of the anoxic phase, phase 1 (A), the hypoxic phase, phase 3 (B), and the aerobic phase, LO (C). Simulations were performed using the growth conditions outlined for simulation scenario 1. Cytosolic and plastidic glycolysis have been collapsed to obtain a condensed representation of the metabolic network (demarcation, gray dashed line). Key metabolites taken up by the model (Suc, Asn, and Gln) are highlighted in gray circles. Metabolites excreted by the model (CO2, ethanol, and lactate) or incorporated into biomass (starch, storage proteins, β-glucan, arabinoxylan, and cellulose) are highlighted in gray rectangles. Reactions are as follows: 1, Suc synthase; 2, cytosolic ADP-Glc pyrophosphorylase; 3, plastidic ADP-Glc pyrophosphorylase; 4, phosphoenolpyruvate carboxylase. AcCoA, Acetyl-CoA; ADPglc, ADP-Glc; AraXyl, arabinoxylan; B-glucan, β-glucan; Cel, cellulose; Cit, citrate; Eth, ethanol; Frc, Fru; F6P, Fru-6-P; Fum, fumarate; G1P, Glc-1-P; G6P, Glc-6-P; Lac, lactate; Mal, malate; OAA, oxaloacetate; 2OG, 2-oxoglutarate; PEP, phosphoenolpyruvate; 3PG, 3-phosphoglycerate; Pyr, pyruvate; UDPglc, UDP-Glc.
To further characterize each phase, a shadow price analysis (Edwards et al., 2002) of the key metabolites was performed (Table II ). Mathematically, the shadow price measures the increase in the value of the objective function due to a small increase in the availability of a given metabolite. In this study, in which the objective function is the maximization of growth, the shadow price denotes the usefulness of a metabolite toward increasing the growth rate, with a negative shadow price indicating that a metabolite is limiting (i.e. the increase of the metabolite will increase the growth rate) and a positive shadow price indicating that a metabolite is available in excess. Thus, shadow price analysis supports the interpretation of a given metabolic behavior and therefore helps to characterize a given simulation scenario (see “Materials and Methods” for details).
Shadow prices
Shadow prices calculated for the phenotypic phases and the LO. c, Cytosolic; m, mitochondrial; n.c., not computed, as LO defines the optimal ratio of substrate uptake rates for maximal biomass production; p, plastidic.
In the following, each phase of the PhPP is described individually and the changes in the metabolic behavior are specified.
Anoxia: Phase 1
Under anoxic conditions, the metabolic flux distribution is characterized by the lack of respiration and high flux through fermentation (Fig. 3A). Flux through ATP-consuming biosynthetic processes such as starch and amino acid synthesis is low, resulting in low biomass production. Amino acid synthesis is based on the import of Gln and Asn, and generation of substrates for starch synthesis is solely mediated by cytosolic AGPase, as flux through the plastidic AGPase is directed toward glycolysis. Glycolytic flux is high and mainly directed through the ATP-utilizing reactions ATP:Fru-6-P 1-phosphotransferase (PFK) and pyruvate kinase (PK; Fig. 4A ). Pyrophosphate:Fru-6-P 1-phosphotransferase (PFP) and PFK form a futile cycle in which ATP is converted to inorganic pyrophosphate (PPi). There is no cyclic flux through the citrate cycle. Except for plastidic NADH, the redox carriers NADH and NADPH are available in excess, whereas energy in the form of ATP is limited, as indicated by the respective positive shadow prices (Table II).
Carbon flux maps depicting the key fluxes within the cytosolic glycolysis and Suc-to-starch pathways of the anoxic phase, phase 1 (A), the hypoxic phase, phase 3 (B), and the aerobic phase, LO (C). Simulations were performed using the growth conditions outlined for simulation scenario 1. The representation of glycolysis is restricted to the ATP-utilizing reactions PFK and PK and the PPi-utilizing bypass PFP and PPDK to focus on differences between ATP- and PPi-dependent glycolysis in response to oxygen supply. Reactions are as follows: AGPase, ADP-Glc pyrophosphorylase; FK, fructokinase; PFK, ATP:Fru-6-P 1-phosphotransferase; PFP, pyrophosphate:Fru-6-P 1-phosphotransferase; PK, pyruvate kinase; PPDK, pyruvate orthophosphate dikinase; SuSy, Suc synthase; UGPase, UDP-Glc pyrophosphorylase. Metabolites are as follows: ADPglc, ADP-Glc; Frc, Fru; F6P, Fru-6-P; F1,6BP, Fru-1,6-bisP; G1P: Glc-1-P; PEP, phosphoenolpyruvate; Pyr, pyruvate; UDPglc, UDP-Glc.
Hypoxia: Phase 2
The onset of oxygen exposure is characterized by the induction of respiration and decline of glycolytic flux. Metabolism is still highly oxygen limited; therefore, there are no major changes in the flux distribution.
Hypoxia: Phase 3
Further increasing the oxygen supply results in the utilization of the complete citrate cycle, accompanied by a decrease of fermentative fluxes and an increase of flux through ATP-consuming processes such as starch synthesis (Fig. 3B). In contrast to phases 1 and 2, the futile cycle formed by PFK and PFP is no longer active and glycolytic flux is additionally carried by the PPi-utilizing bypass of PFP and pyruvate orthophosphate dikinase (PPDK; Fig. 4B). The Suc availability decreases with increasing aerobioses (Table II), indicating that Suc has a growing importance as a limiting factor with increasing oxygen supply.
Hypoxia: Phase 4
The metabolic flux distribution pattern denoted by this phase is characterized by the cessation of fermentation: ethanol and lactate are no longer secreted, so more carbon can be utilized for growth. In this phase, the citrate cycle is the major source of energy. In contrast to phases 1 to 3, glycolytic flux is restricted to the PPi-utilizing bypass of PFP and PPDK (Fig. 4C). In addition, plastidic AGPase now operates in the direction of starch synthesis, resulting in higher biomass production. In contrast to phases 1 to 3, redox equivalents are limiting (Table II).
Aerobiosis: LO
The metabolic flux map denoted by LO corresponds to the optimal metabolic flux distribution pattern without limitations on the availability of oxygen, thus allowing optimized biomass production. In contrast to phases 1 to 4, Gln is the only form of nitrogen delivered to the endosperm and ATP is imported into the amyloplast to provide energy for ATP-consuming biosynthetic processes.
Hyperoxia: Phase 5
In this futile phase, oxygen excess is inhibitory toward obtaining maximal biomass production (Table II). Phase 5 is excluded from further analysis, as the physiological condition of oxygen excess does not occur under in vivo conditions.
To take into account the fact that carbohydrate import into seeds is dramatically decreased under hypoxic or anoxic conditions (van Dongen et al., 2004), a second simulation scenario was performed by decreasing the Suc input concomitant to the decreasing oxygen consumption rate, as described in “Materials and Methods” (simulation scenario 2). The corresponding flux values are given in Supplemental Data Set S7. As depicted in Figure 1 and Supplemental Data Set S7, the resulting flux distribution pattern between the various pathways in the different oxygen conditions corresponds to the flux distribution described for simulation scenario 1, while the overall flux is reduced.
Reaction Essentiality
For each of the metabolic networks specific for the phases described above, the effect of in silico knockouts of enzymatic reactions involved in central metabolism (i.e. glycolysis, citrate cycle, pentose phosphate pathway, fermentation, oxidative phosphorylation, Suc-to-starch pathway) was studied to determine the importance of the reactions under different growth conditions. It should be mentioned that an enzyme deletion in one of the mostly linear and nonredundant pathways for amino acid or cofactor synthesis would logically result in a loss of the ability to grow. Therefore, this analysis was restricted to central parts of the metabolic network. Furthermore, it should be noted that the in silico knockout of enzymatic reactions as described in this paper is not analogous to gene deletion, as a gene family might contain more than one member and the deletion of only one member of the family would not affect cell growth.
Table III shows that for each of these networks only a small number of reactions are essential for growth, reflecting the flexibility of the metabolic networks to provide biosynthetic precursors required for biomass production. Three enzymatic reactions belonging to three different pathways were determined to be essential under all conditions: (1) fumarate hydratase (citrate cycle), (2) phosphoenolpyruvate carboxylase (anaplerosis), and (3) fructokinase (Suc breakdown). In addition, all reactions involved in the nonoxidative phase of the pentose phosphate pathway and the Rubisco bypass, and one reaction involved in oxidative phosphorylation (H+-exporting ATPase), were determined to be essential for the anoxic phase (phase 1).
Enzyme deletion
Results of the in silico deletion study performed on reactions of central metabolism in the metabolic networks of phases 1 to 4 and the LO. The given values correspond to the percentage of optimal biomass production (complete enzyme set) achieved after deleting a particular enzyme. The essential in silico predictions (percentage of optimal biomass production = 0) are in boldface type, and growth-inhibiting reactions (percentage of optimal biomass production < 95%) are marked with an asterisk. Ana, Anaplerosis; Glyc, glycolysis; oxP, oxidative phosphorylation; PPP, oxidative pentose phosphate pathway; Rub, Rubisco bypass; Suc, Suc-to-starch pathway; TCA, citrate cycle.
Furthermore, the model predicts that the deletion of certain reactions involved in energy metabolism, namely oxidative phosphorylation, citrate cycle, and glycolysis, has a growth-inhibiting effect under most conditions (phases 2–4 and LO). Eliminating reactions involved in the oxidative phosphorylation pathway (H+-exporting ATPase, NADH dehydrogenase, cytochrome c oxidase) leads to a major decrease (60%–91%) in growth efficiency, while deleting the citrate cycle enzymes aconitate hydratase, isocitrate dehydrogenase (NAD+), pyruvate dehydrogenase, and citrate synthase has only a minor growth-inhibiting effect (93%–99%) under these conditions. With the exception of the cytosolic phosphoglucose isomerase (94%–99%), none of the enzymes involved in glycolysis was determined to have a growth-inhibiting effect.
With respect to reactions involved in the synthesis of biomass components, the deletion of the cytosolic AGPase, UDP-Glc pyrophosphorylase, and Suc synthase leads to a major decrease (79%–89%) in the growth efficiency under all growth conditions. The remaining reactions were determined to be unessential, since the phase-specific metabolic networks maintained the capability of achieving a growth rate of 95% to 100% compared with simulations with the complete enzyme set.
DISCUSSION
Modeling Seed Storage Metabolism
With the aim of getting a systemic understanding of cereal seed storage metabolism, we constructed a stoichiometric model of primary metabolism in the developing endosperm of barley seeds during starch accumulation. The constructed model includes central metabolism (glycolysis, pentose phosphate pathway, citrate cycle), amino acid metabolism, starch synthesis, and some minor pathways. To our knowledge, this is the first detailed attempt at stoichiometric modeling of seed metabolism to date, although a few smaller models of storage metabolism in other plant species exist, including sugarcane (Saccharum spontaneum; Uys et al., 2007), potato (Solanum tuberosum; Poolman et al., 2004), and canola (Brassica napus; Schwender et al., 2003).
To evaluate the predictive capabilities of the model, published experimental results were compared with the simulation results. With respect to prediction of grain yield, the computed growth rates were in good agreement with the published literature, indicating that the objective function used for model simulation is appropriate for biologically meaningful predictions and that the constructed model has the potential to simulate cereal seed metabolism. By providing an initial framework for studying seed storage metabolism in silico, different aspects of storage patterning were computationally analyzed, which are discussed in more detail in the following sections.
Cereal Seed Storage Metabolism in Response to Oxygen Depletion
Oxygen depletion appears to be a common feature within dense or metabolically active plant tissues, such as developing seeds and grains (Rolletschek et al., 2002, 2004; Geigenberger, 2003; van Dongen et al., 2004). To study the effect of oxygen depletion on storage patterning in developing barley seeds, growth simulations under varying oxygen conditions were performed. To validate the model, the simulation results were compared with published experimental results. In addition, the model predictions were used to test controversial hypotheses on the role of PPi metabolism in oxygen-depleted tissues published in the literature in order to benefit from the potential of mathematical modeling to verify and extend the understanding of controversially discussed biological processes.
Seed Metabolism under Anoxic Conditions
Under fully anaerobic conditions, the model showed characteristic anaerobic metabolic behavior (Kennedy et al., 1992; Geigenberger, 2003): (1) inhibition of respiration, (2) induction of fermentation, and (3) stimulation of glycolysis (Pasteur effect), resulting in a decrease in the cellular energy state and an increase in the redox state [i.e. the NAD(P)H-NAD(P)+ ratio]. In response to the surplus of reducing equivalents, the operation of pathways generating NAD(P)H is limited, as indicated by an incomplete citrate cycle and the inactivity of the oxidative part of the pentose phosphate pathway.
The role of citrate cycle enzymes in anaerobic metabolism in plants has been examined very little, and in most cases with respect to anoxia-tolerant species (Fox and Kennedy, 1991; Kennedy et al., 1992). However, mitochondria of most plants have been reported to degenerate or fail to develop properly without oxygen (Ueda and Tsuji, 1971; Oliveira, 1977), indicating inoperative mitochondrial activity. In agreement with the incomplete citrate cycle predicted by the model, (1) in canola embryos it has been shown that cyclic flux from the citrate cycle is absent under most conditions (Schwender et al., 2006; Junker et al., 2007), and (2) invertebrates have been shown to commonly utilize a partial citrate cycle during anoxia (Hochachka, 1986).
With respect to glycolysis, the model predicts a PFK-PFP substrate cycle (Stitt, 1990) in which energy is directed from ATP to PPi and, thus, to PPi-consuming reactions. PFP catalyzes a near equilibrium reaction (Stitt, 1990) and is usually characterized by high activity in young developing tissues and starch-storing tissues (Xu et al., 1989). The precise metabolic function of PFP is still unclear and controversially discussed (Stitt, 1990; Hajirezaei et al., 1994). Several metabolic functions have been proposed, including a role in glycolysis (Duff et al., 1989; Hatzfeld et al., 1989), gluconeogenesis (Fahrendorf et al., 1987; Paul et al., 1995), and operation in a cycle with PFK to produce PPi required for PPi-consuming reactions, such as Suc mobilization via Suc synthase and UGPase (Stitt, 1990). Several investigations support the existence of such a cycle in cereal seed metabolism. Studies about the role of PFP and PFK in the endosperm of developing seeds of wheat (Triticum aestivum; Mahajan and Singh, 1990, 1992) indicate that PFP is primarily involved in the generation of PPi. A similar mechanism was suggested in anoxic rice (Oryza sativa) coleoptiles (Gibbs et al., 2000). By providing evidence for the gluconeogenetic direction of PFP in anoxic tissues, the authors support the view of PFK being a significant control point for glycolytic flux under anoxia. Studying the role of PPi-dependent glycolytic enzymes during anoxia in anoxia-tolerant rice and anoxia-intolerant Arabidopsis (Arabidopsis thaliana), Huang et al. (2008) proposed PFP in rice coleoptiles to function (1) in the glycolytic direction in the early stage of anoxia in order to accelerate glycolysis in response to the anoxic energy crises and (2) in the gluconeogenetic direction during long-term anoxia in order to slow net glycolysis to conserve carbohydrates. In agreement with the model predictions, the authors suggest that PFP operates in a cycle with PFK to produce PPi required for PPi-consuming processes in anoxia-intolerant Arabidopsis. Nevertheless, as there is still no clear evidence about the precise metabolic function of PFP under anoxia, more studies are needed to provide further evidence for the model predictions.
Seed Metabolism under Hypoxic Conditions
In general, the simulation results under hypoxic conditions are consistent with the main qualitative physiological characteristics reported for cereal seed storage metabolism under hypoxic conditions (Rolletschek et al., 2004, 2005; van Dongen et al., 2004). With respect to starch metabolism, the model correctly predicts the Suc-to-starch pathway reported from barley seed metabolism (Thorbjørnsen et al., 1996; Weschke et al., 2000; Beckles et al., 2001; Emes et al., 2003) by predicting that (1) Suc degradation is restricted to the Suc synthase pathway and (2) synthesis of ADP-Glc, which is the main precursor for starch synthesis, is predominantly catalyzed by the cytosolic isoform of AGPase. The prominent role of cytosolic AGPase with respect to starch synthesis is furthermore confirmed by the in silico deletion studies, which agree with experimental results in which a low starch content was observed in a barley mutant lacking the cytosolic small subunit of AGPase (Johnson et al., 2003).
With increasing oxygen supply, the model predicts a shift in direction of the PFP reaction, which now operates in the glycolytic direction, as proposed by several studies with heterotrophic cereal tissues (Doehlert et al., 1988; Roscher et al., 1998). Thus, the PPi-utilizing bypass of PFP and PPDK now acts in parallel to the ATP-utilizing reactions PFK and PK. Although PPDK is found to be highly abundant in developing seeds of graminaceous cereals such as barley, wheat, and rice (Meyer et al., 1982; Aoyagi et al., 1984; Nomura et al., 2000; Chastain et al., 2006), its metabolic function in these organs is not yet known. Different putative metabolic functions for cytosolic PPDK operating in the pyruvate-to-phosphoenolpyruvate direction have been proposed for cereal seeds, suggesting that the enzyme is involved either in providing phosphoenolpyruvate for the refixation of respiratory CO2 or in amino acid interconversions (Meyer et al., 1982; Aoyagi et al., 1984). Contrary to these hypotheses, Chastain et al. (2006) postulated a putative function for the enzyme working in the phosphoenolpyruvate-to-pyruvate-forming direction. Studying endosperm-localized cytosolic PPDK in developing seeds of rice, the authors proposed PPDK to function as an efficient mechanism for glycolytic ATP synthesis in oxygen-depleted tissues by converting AMP to ATP, resulting in a net formation of two ATPs compared with the single ATP generated by PK (Plaxton, 2005). With respect to the hypoxic nature of cereal seeds and the consequent low ATP-ADP ratio (van Dongen et al., 2004; Rolletschek et al., 2004), a similar role of PPDK in energy metabolism of hypoxic rice seeds has been postulated by Kang et al. (2005). In contrast, Huang et al. (2008) suggest PPDK to operate in a cycle with PK in anoxia-tolerant rice coleoptiles to provide PPi required for PPi-consuming reactions. Nevertheless, in agreement with the model predictions, the authors propose PPDK to operate in the glycolytic direction in anoxia-intolerant plant species such as Arabidopsis.
Assuming a PPDK- and PFP-using glycolysis as predicted by the model, the bioenergetic efficiency is furthermore enhanced due to the yield of five ATP molecules per Glc molecule, instead of the two produced by conventional glycolysis (Mertens, 1993; Plaxton, 2005). The role of PPi as an alternative energy donor in the cytosol of plants has been addressed in several reports (Geigenberger et al., 1998; Stitt, 1998; Plaxton, 2005). Although different studies indicate an important role of PPi metabolism in young growing tissues and in stress conditions including anaerobiosis (for review, see Stitt, 1998), more studies are needed to provide further evidence for the model predictions.
A different response of plant metabolism to unfavorable conditions such as hypoxia is reported to be the accumulation of γ-aminobutyric acid (GABA; Kinnersley and Turano, 2000; Bouché and Fromm, 2004). The associated GABA shunt is proposed to function as an alternative NAD+-independent bypass for Glu entry into the citrate cycle (Plaxton, 2005). Further roles as well as aspects of its regulation have been discussed recently (Fait et al., 2008). In contrast to the observations that the GABA shunt in developing soybean (Glycine max) seed is associated with hypoxia (Shelp et al., 1995), this pathway was not active under the simulated hypoxic conditions. However, in agreement with the model prediction, Inatomi and Slaughter (1975) showed that in barley embryo under water imbibition, GABA accumulated markedly only when oxygen was available.
Seed Metabolism under Aerobic Conditions
Under fully aerobic conditions, the model showed characteristic aerobic metabolic behavior of cereal seeds provided with ambient oxygen (van Dongen et al., 2004): (1) up-regulation of respiratory energy production, resulting in an increase in the cellular energy state; (2) subsequent increase of storage metabolism, leading to an increase of phloem transport toward the seed; and (3) increase of seed dry weight due to extensive storage metabolism.
With respect to glycolysis, the model predicts PPi-dependent glycolysis to be the only form of glycolysis operating under fully aerobic conditions. Although no experimental data on this point have been reported from barley seeds, the bioenergetic advantages of PPi-dependent glycolysis are a possible explanation for this model prediction.
Altogether, the simulation results support the idea of PPi availability determining the net flux through PFP and PPDK, as proposed by previous reports (Stitt, 1990; Gibbs et al., 2000): (1) under anoxic conditions, a substantial requirement for PPi results in PFP operating in the direction of PPi synthesis to compensate for PPi deficiency; (2) under hypoxic conditions, increasing PPi availability results in a net flux of PFP and PPDK operating in the direction of PPi consumption; and (3) under aerobic conditions, sufficient PPi availability allows the energy-saving mode of PPi-dependent glycolysis to be the only form of glycolysis. This interpretation is supported by simulation studies, showing that in silico-generated PPi deficiency leads to PFP and PPDK operating in the direction of PPi synthesis (data not shown).
Although no experimental data on PPi deficiency are available from oxygen-depleted barley seeds, low flux through the major source of PPi production (cytosolic AGPase) in conjunction with high flux through the major source of PPi consumption (UGPase) suggest low PPi availability in anoxic and hypoxic tissues. Although further experiments are required to evaluate the proposed hypothesis, the given examples show the potential of the model to verify or extend the understanding of controversially discussed biological processes by looking at systemic stoichiometric constraints only. Furthermore, it reveals the advantage of systems-oriented modeling by giving insight into complex biological processes provided by the integration of multiple data, which in this form cannot be obtained by studies restricted to the analysis of single enzymes.
Effect of Enzyme Deletion on Seed Storage Metabolism
In silico knockout analysis provides an efficient method to study the importance of the reactions in the metabolic network and to gain insight into metabolic changes caused by enzyme deletion. The in silico deletion study demonstrated the flexibility of the metabolic network of barley to compensate for enzymatic perturbations. For all of the growth conditions examined, only three enzymatic reactions were determined to be essential: (1) fumarate hydratase (citrate cycle); (2) phosphoenolpyruvate carboxylase (anaplerosis); and (3) fructokinase (Suc breakdown).
Plant fumarate hydratase has been documented to be an important control point in the citrate cycle (Behal and Oliver, 1997). Unlike the situation in microbial systems, such as yeast, for which a variety of knockout studies exist to date (Przybyla-Zawislak et al., 1999; McCammon et al., 2003; Kokko et al., 2006), only a single transgenic study, analyzing the effect of fumarase inhibition, has been carried out in plants (Nunes-Nesi et al., 2007). That study revealed that antisense inhibition of mitochondrial fumarase activity in tomato (Solanum lycopersicum) plants has a negative impact on photosynthesis and leads to a consequent reduction in total plant biomass. However, biochemical analysis revealed only minor alterations in leaf metabolism, indicating that the residual activity of fumarase in the transgenic plants is still sufficient to ensure relatively normal metabolic function. This is supported by the discovery that most enzymes in central metabolism are present at activities far beyond the actually mediated flux (Junker et al., 2007). Thus, knockout studies are needed to support further evidence for the model prediction. Nevertheless, knockout analysis in yeast suggested an important role of fumarase in central metabolism, with cells lacking fumarase being essentially unable to grow on any nonfermentable carbon source (Przybyla-Zawislak et al., 1999; McCammon et al., 2003; Kokko et al., 2006).
In heterotrophic tissues (i.e. developing seeds and fruits) of C3 plants, phosphoenolpyruvate carboxylase plays an essential role in the anaplerotic replenishment of citrate cycle intermediates by providing precursors for several biosynthetic pathways, including the biosynthesis of amino acids (Chollet et al., 1996). Several groups have attempted to manipulate assimilate partitioning in heterotrophic tissues of C3 plants (Lebouteiller et al., 2007; Radchuk et al., 2007) or to improve CO2 fixation in autotrophic tissues of C3 plants (Gehlen et al., 1996; Ku et al., 1999) by overexpressing phosphoenolpyruvate carboxylase. In contrast, only a few reports exist on plants with decreased phosphoenolpyruvate carboxylase activity, and studies of knockout mutants have not been carried out to date. In the inhibition studies, transgenic potato plants with reduced phosphoenolpyruvate carboxylase activity (30% of wild-type activity [Rademacher et al., 2002] and 50% to 70% of wild-type activity [Gehlen et al., 1996]) revealed no significant effect on whole plant and tuber growth. However, metabolite analysis showing only little alteration in leaf metabolism suggests an insufficient reduction of phosphoenolpyruvate carboxylase activity to affect metabolism. Nevertheless, the severe reduction of growth rate in ppc knockout Escherichia coli strains, despite of the presence of other anaplerotic reactions not present in plants, points to an essential role of phosphoenolpyruvate carboxylase in central metabolism (Peng et al., 2004; Fong et al., 2006).
With respect to the identified essential enzymatic reaction involved in Suc breakdown, antisense inhibition studies document the predominant role of fructokinase in starch metabolism (Davies et al., 2005), thus providing experimental evidence for the model predictions. Antisense inhibition of potato fructokinase resulted in a reduced tuber yield with a substantial shift in tuber metabolism, suggesting an important role for fructokinase in maintaining a balance between Suc synthesis and degradation (Davies et al., 2005).
Taken together, the in silico knockout results demonstrate an important property of the cereal seed metabolic network, namely that there are only few essential enzymatic reactions in central metabolism. Spatial compartmentation of biochemical pathways, which exists in all higher organisms, and isoforms expressed in different compartments could be an explanation for the observed network flexibility. For example, none of the enzymes of the glycolytic pathway, which is present in cytosol and plastid, was predicted to be essential. In contrast, in microorganisms lacking compartmentation, such as E. coli, several glycolytic enzymes are considered to be essential based on in silico deletion studies and experimental observations (Fraenkel, 1996; Edwards and Palsson, 2000). These results indicate that redundancy properties of the plant metabolic network play a major role in network robustness (Barkai and Leibler, 1997; Barabási and Oltvai, 2004; Kitano, 2004). By accomplishing the same step in a metabolic pathway in different ways, this metabolic flexibility supports plants growing and/or surviving under suboptimal environmental conditions (Plaxton, 1996, 2005).
PERSPECTIVES
The potential of stoichiometric modeling has been shown in a variety of organisms (Fell and Small, 1986; Varma et al., 1993a, 1993b; Edwards and Palsson, 2000; Ramakrishna et al., 2001; Schilling et al., 2002; David et al., 2003; Famili et al., 2003; Förster et al., 2003; Reed and Palsson, 2003; Shastri and Morgan, 2005; Cakir et al., 2007; Duarte et al., 2007). Here, FBA was applied to a stoichiometric model of cereal seed metabolism to gain insight into storage patterning in developing barley seeds. Although the predictions of the computational analysis were found to be in good agreement with published experimental results, the current model may be improved and refined. For example, additional experimental results, such as in vivo measurements of metabolic fluxes (Schwender et al., 2004), are needed for experimental verification to allow us to further improve and update the metabolic construction.
Despite the known limitations of FBA, including the constraints of not incorporating regulatory events and of predicting optimal behavior only, which may not reflect suboptimal growth in vivo (Edwards and Palsson, 2000), the results presented here indicate that the constructed model has the potential to simulate cereal seed metabolism. Thus, by providing an initial framework for studying cereal seed storage metabolism in silico, in future applications the model can be used to verify and extend the understanding of complex processes, to generate and test hypotheses, and to explore in silico scenarios, which eventually should allow us to find suitable targets for the improvement of grain and yield quality.
MATERIALS AND METHODS
Metabolic Network Construction
The metabolic network of primary metabolism in the developing endosperm of barley (Hordeum vulgare) seed was constructed in a stepwise manner by integrating biochemical, physiological, and genomic data through manual curation of an extensive survey of the scientific literature and online databases. Only pathways of primary seed metabolism necessary to generate major biomass components (>1% of total dry weight) were included in the model using the biomass composition of barley grains reported in OECD (2004) as a basis (see Supplemental Data Set S4 for more details). The literature references supporting each reaction (Supplemental Data Set S1) and a list of the databases inquired during the construction process (Supplemental Data Set S2) are given in the online supplemental material. In case no data were available to support the presence of a reaction in developing barley seeds (e.g. not all of the enzymatic steps in a biosynthetic pathway have been characterized in barley), data referring to closely related monocotyledon species such as wheat (Triticum aestivum), rice (Oryza sativa), and maize (Zea mays) were considered to substantiate the reaction. Unless data on the irreversibility of a reaction were available, all reactions were considered to be reversible. Multienzyme complexes, such as pyruvate dehydrogenase, were modeled as a single reaction by merging the reactions of the respective subunits. With respect to compartmentation, reactions located in compartments other than the extracellular medium, cytosol, mitochondria, and amyloplast, as well as reactions without information on their localization, were considered to be cytosolic due to the fact that the number of reactions annotated for other compartments is very small. The set of reactions included in the network (Supplemental Data Set S1) and a systems biology markup language (Hucka et al., 2003) version of the model (Supplemental Data Set S6) are given in the online supplemental material. Model construction was done using the information system Meta-All (Weise et al., 2006). Model and flux visualization was done using the VANTED system (Junker et al., 2006).
FBA
FBA is a constraint-based modeling technique developed to characterize the production capabilities and systemic properties of a metabolic network based on the mass balance constraints (Edwards et al., 1999). Assuming metabolic steady state, the system of mass balance equations derived from a metabolic network of n reactions and m metabolites can be represented as follows:
with
where S is the stoichiometric matrix (m × n) and v is a flux vector of n metabolic fluxes, with αi as lower and βi as upper bounds for each vi, respectively. FBA uses the principle of linear programming to solve the system of mass balance equations by defining an objective function and searching the allowable solution space for an optimal flux distribution that maximizes or minimizes the objective.
Objective Function
The identification of the appropriate objective function is one of the main subjects with FBA. In most of the application studies, the maximization of growth (i.e. biomass yield) has been assumed to be the main objective of metabolism (Edwards and Palsson, 2000; Schilling et al., 2002; Price et al., 2004) and has proven useful in predicting in vivo cellular behavior (Edwards et al., 2001; Ibarra et al., 2002).
For the analysis presented here, the maximization of growth per flux unit, which is a combination of the maximization of growth and the minimization of overall flux, was used as the objective function. According to the principle of flux minimization (Holzhütter, 2004), the minimization of overall flux (i.e. the squared sum of all fluxes) ensures maximal enzymatic efficiency, resulting in an efficient metabolic flux distribution. The hypothesis underlying the minimization principle postulates that cells and whole organisms gain functional fitness by fulfilling their functions with minimal effort (i.e. the economic usage of the available sources; Holzhütter, 2004). Applications of this approach and its derivatives have been shown to be in good agreement with the main in vivo cellular behavior of prokaryotic and eukaryotic systems, thus leading to flux values of biological relevance (Bonarius et al., 1996; Holzhütter, 2004; Cakir et al., 2007; Schuetz et al., 2007).
The objective function was computed as a two-step optimization process, where the first step is to maximize growth (linear optimization) and the second step is to minimize the overall intracellular flux (nonlinear optimization). The computation was performed by adding the objective value (growth rate) of the first optimization as an additional constraint to the second optimization. The second (nonlinear) optimization step not only achieves an efficient channeling of metabolites but also allows handling of the fundamental FBA problem of alternative optima (i.e. multiple flux distributions with the same optimal value identified by the objective). Due to its underlying quadratic optimization, the nonlinear optimization step does not produce alternative optima, thus eliminating those possibly produced by the linear optimization step (Schuetz et al., 2007).
The growth objective was mathematically defined as a flux drain composed of all biosynthetic precursors and cofactors (e.g. amino acids, ATP, etc.) required for biomass production. Growth was modeled as a single reaction (i.e. the biomass reaction) converting all of the precursors into biomass. To estimate the demands of each of the biosynthetic precursors, the biomass composition of barley grains reported in OECD (2004) was used as a reference. A table listing the biomass components and their respective stoichiometric coefficients in the biomass equation of the model is given in Supplemental Data Set S4. Energy requirements (i.e. growth- and nongrowth-associated ATP maintenance requirements) for growth were determined by reviewing the relevant literature. Due to the lack of experimental data, growth-associated maintenance was estimated to be 5.36 mmol ATP g−1 dry weight based on the framework proposed by Amthor (2000), as detailed in Supplemental Data Set S5. Experimentally determined nongrowth-associated maintenance in barley seeds was found to be as low as 7 × 10−6 to 14 × 10−6 mmol ATP g−1 dry weight h−1 (Penning de Vries, 1975) and thus was ignored.
Simulation Conditions
Simulations were performed in MATLAB using the COBRA toolbox (Becker et al., 2007). To allow for the computation of the nonlinear optimization, the toolbox was extended by the CLP solver (COIN-OR Linear Program Solver; http://www.coin-or.org/Clp/index.html) and a MATLAB routine comprising the two-step optimization process.
All simulations were performed using the conditions outlined in this section. The exchange flux for Suc was constrained between 0 and 8 μmol g−1 dry weight h−1, a value taken from experimental reports (Felker et al., 1984). In these [14C]Suc-labeling experiments, the Suc uptake in developing grain of barley was measured to be 5 and 8 μmol g−1 dry weight h−1 for 17 and 5 d after anthesis, respectively. Relating to maximum substrate uptake, the higher value was used to constrain the model. Unless otherwise stated, the exchange flux for metabolites allowed to enter (i.e. Asn, Gln, oxygen, H2S) or leave (i.e. lactate, ethanol, carbon dioxide) the system was always unconstrained in the net inward or outward direction, respectively (Supplemental Data Set S1). Fluxes through all internal reactions were always unconstrained.
Growth under oxygen deprivation was simulated by varying the oxygen uptake rate between 0 and 30 μmol g−1 dry weight h−1.
To illustrate the metabolic behavior characterizing the different phenotypic phases (see below), two simulation scenarios were run. Simulation scenario 1 was performed by fixing the Suc uptake rate (SUR) at 8 μmol Suc g−1 dry weight h−1 and the oxygen uptake rate (OUR) at representative values for each phenotypic phase (units of μmol g−1 dry weight h−1): P1, OUR = 0; P2, OUR = 1; P3, OUR = 4; P4, OUR = 8.5; LO, OUR = 8.9. Simulation scenario 2 was performed by fixing SUR and OUR as follows (units of μmol g−1 dry weight h−1): P1, SUR = 5, OUR = 0; P2, SUR = 5.3, OUR = 0.8; P3, SUR = 6.1, OUR = 3.5; P4, SUR = 7.7, OUR = 8.4; LO, SUR = 8, OUR = 8.9.
Simulations for the shadow price analysis and the deletion studies (see below) were performed using the conditions outlined for simulation scenario 1.
PhPP Analysis
PhPP analysis is a method used to obtain a global perspective on the optimal flux distributions and how they are affected by changes in two environmental parameters of specific interest (Edwards et al., 2001, 2002). To define a PhPP, the respective values of the environmental conditions under investigation are represented on the axes of an (x,y) plane and the optimal flux distribution is computed for all points in the plane. To demarcate regions of constant flux distributions (i.e. phenotypic phases), shadow prices are computed for the two-parameter space and demarcation lines separating the phenotypic phases in the phase plane are drawn based on changes in the shadow prices (Edwards et al., 2002). Mathematically, the shadow price (γi) depicts the sensitivity of the objective function (Z) to changes in the availability of a metabolite (bi):
thus indicating the usefulness of the metabolite toward increasing the objective function. The shadow price can be negative, zero, or positive, with a negative shadow price indicating that a metabolite is limiting (i.e. the increase of the metabolite will increase the objective function) and a positive shadow price indicating that a metabolite is available in excess. Thus, shadow price analysis allows us to identify qualitative shifts from one optimal flux distribution to another due to the fact that shadow prices are constant within each phenotypic phase and will be different in the other phases.
In this study, the PhPP was computed by calculating the shadow price for each point (interval, 1 μmol g−1 dry weight h−1) in the Suc-oxygen plane. In addition, to support the interpretation of the phenotypic phases derived from the PhPP analysis, a shadow price analysis of the key metabolites (Suc, oxygen, ATP, H+, NADH, and NADPH) was performed.
Deletion Studies
To predict essential reactions in the phase-specific metabolic networks, enzymatic reactions involved in core metabolism (i.e. glycolysis, citrate cycle, pentose phosphate pathway, fermentation, oxidative phosphorylation, Suc-to-starch pathway) were subjected to deletion. To simulate deletion, the flux through the corresponding reaction was set to zero. An enzymatic reaction was considered to be essential if its deletion led to the complete loss of the ability to grow.
Supplemental Data
The following materials are available in the online version of this article.
Supplemental Figure S1. Metabolic map of the model.
Supplemental Data Set S1. Set of reactions included in the constructed network and the respective references supporting each reaction.
Supplemental Data Set S2. Online databases inquired during the construction process.
Supplemental Data Set S3. Terminology used for network reactions and metabolites.
Supplemental Data Set S4. Biomass composition considered for the determination of the stoichiometric coefficients in the biomass equation in the metabolic model.
Supplemental Data Set S5. Description of the computation of growth-associated maintenance.
Supplemental Data Set S6. Systems biology markup language file of the model.
Supplemental Data Set S7. Flux distributions for all reactions computed for simulation scenarios 1 and 2.
Acknowledgments
We thank Christian Klukas for valuable help with the VANTED system and the reviewers for constructive and helpful comments.
Footnotes
-
The author responsible for distribution of materials integral to the findings presented in this article in accordance with the policy described in the Instructions for Authors (www.plantphysiol.org) is: Eva Grafahrend-Belau (grafahr{at}ipk-gatersleben.de).
-
↵1 This work was supported by the German Ministry of Education and Research (grant nos. 031270–6A and 031504–4A).
-
↵[W] The online version of this article contains Web-only data.
- Received September 9, 2008.
- Accepted October 30, 2008.
- Published November 5, 2008.
LITERATURE CITED
