Assessing the Metabolic Impact of Nitrogen Availability Using a Compartmentalized Maize Leaf

18 Pages·2016·1 MB·English

University of Nebraska - Lincoln
DigitalCommons@University of Nebraska - Lincoln
Chemical and Biomolecular Engineering -- All Chemical and Biomolecular Engineering,
Faculty Papers Department of
2014
Assessing the Metabolic Impact of Nitrogen
Availability Using a Compartmentalized Maize
Leaf Genome-Scale Model
Margaret Simmons
The Pennsylvania State University
Rajib Saha
National University of Singapore, [email protected]
Nardjis Amiour
Centre National de la Recherche Scientifique
Akhil Kumar
The Pennsylvania State University
Lenaïg Guillard
Centre National de la Recherche Scientifique
See next page for additional authors
Follow this and additional works at:http://digitalcommons.unl.edu/chemengall
Simmons, Margaret; Saha, Rajib; Amiour, Nardjis; Kumar, Akhil; Guillard, Lenaïg; Clément, Gilles; Miquel, Martine; Li, Zhenni;
Mouille, Gregory; Lea, Peter J.; Hirel, Bertrand; and Maranas, Costas D., "Assessing the Metabolic Impact of Nitrogen Availability
Using a Compartmentalized Maize Leaf Genome-Scale Model" (2014).Chemical and Biomolecular Engineering -- All Faculty Papers.
Paper 49.
http://digitalcommons.unl.edu/chemengall/49
This Article is brought to you for free and open access by the Chemical and Biomolecular Engineering, Department of at DigitalCommons@University
of Nebraska - Lincoln. It has been accepted for inclusion in Chemical and Biomolecular Engineering -- All Faculty Papers by an authorized
administrator of DigitalCommons@University of Nebraska - Lincoln.
Authors
Margaret Simmons, Rajib Saha, Nardjis Amiour, Akhil Kumar, Lenaïg Guillard, Gilles Clément, Martine
Miquel, Zhenni Li, Gregory Mouille, Peter J. Lea, Bertrand Hirel, and Costas D. Maranas
This article is available at DigitalCommons@University of Nebraska - Lincoln:http://digitalcommons.unl.edu/chemengall/49
Assessing the Metabolic Impact of Nitrogen
Availability Using a Compartmentalized Maize
Leaf Genome-Scale Model1[C][W][OPEN]
Margaret Simons2, Rajib Saha2, Nardjis Amiour, Akhil Kumar, Lenaïg Guillard, Gilles Clément,
Martine Miquel, Zhenni Li, Gregory Mouille, Peter J. Lea, Bertrand Hirel, and Costas D. Maranas*
Departments of Chemical Engineering (M.S., R.S., C.D.M.) and Bioinformatics and Genomics, Huck Institutes
of the Life Sciences (A.K.), Pennsylvania State University, University Park, Pennsylvania 16802; Institut Jean-Pierre
Bourgin, Institut National de la Recherche Agronomique, Centre de Versailles-Grignon, Unité Mixte de
Recherche 1318 Institut National de la Recherche Agronomique-Agro-ParisTech, Equipe de Recherce
Labellisée, Centre National de la Recherche Scientifique 3559, F–78026 Versailles cedex, France (N.A., L.G.,
G.C., M.M., Z.L., G.M., B.H.); and Lancaster Environment Centre, Lancaster University, Lancaster LA1 4YQ,
United Kingdom (P.J.L.)
Maize(Zeamays)isanimportantC plantduetoitswidespreaduseasacerealandenergycrop.Asecond-generationgenome-
scale metabolic model for the maiz4e leaf was created to capture C carbon fixation and investigate nitrogen (N) assimilation
4
by modeling the interactions between the bundle sheath and mesophyll cells. The model contains gene-protein-reaction
relationships, elemental and charge-balanced reactions, and incorporates experimental evidence pertaining to the biomass
composition, compartmentalization, and flux constraints. Condition-specific biomass descriptions were introduced that
account for amino acids, fatty acids, soluble sugars, proteins, chlorophyll, lignocellulose, and nucleic acids as experimentally
measuredbiomassconstituents. Compartmentalization ofthemodelisbasedonproteomic/transcriptomic dataandliterature
evidence.WiththeincorporationofinformationfromtheMetaCropandMaizeCycdatabases,thisupdatedmodelspans5,824
genes, 8,525 reactions, and 9,153 metabolites, an increase of approximately 4 times the size of the earlier iRS1563 model.
Transcriptomic and proteomic data have also been used to introduce regulatory constraints in the model to simulate an
N-limited condition and mutants deficient in glutamine synthetase, gln1-3 and gln1-4. Model-predicted results achieved 90%
accuracywhencomparingthewildtypegrownunderanN-completeconditionwiththewildtypegrownunderanN-deficient
condition.
Maize(Zeamays),alsoknownascorn,isanessential 2013 (Monaco et al., 2013), and MetaCrop 2.0 in 2012
dual-use food and energy crop. Maize production is (Schreiber et al., 2012), there is a need for an updated
increasing at the greatest rate among all cereals, with a genome-scale metabolic model (GSM; Saha et al., 2011)
worldwidetrendof0.06tonsha21year21(Leveauetal., that will integrate all newly available information from
2011)andarecord877milliontonsproducedinthe2011- diversesources.Theintegrationofthisinformationwith
2012 fiscal year (International Grains Council, 2013). experimental transcriptomic data, proteomic data, and
Withtherecentcompletionofthemaizegenomein2009 biomasscompositionmeasurementsobtainedwithwild-
alongwiththecreationandcurationofdatabasessuchas type plants grown under optimal nitrogen (N+ WT)
2
MaizeGDB in 2011 (Schaeffer et al., 2011), MaizeCyc in conditions and limited nitrogen (N WT) conditions
(Amiouretal.,2012),aswellastwoGlnsynthetase(GS)
mutants grown under optimal nitrogen (N), gln1-3 and
1This work was supported by the U.S. Department of Energy gln1-4(Martinetal.,2006),hasprovidedamoreaccurate
(grantno.DE–FG02–05ER25684)andbytheAgenceNationalpour assessment of N metabolism within the maize leaf.
la Recherche Genoplante program Maize and Yield (grant no. Moreover,sinceintegrationoftranscriptomic,proteomic,
GNP05015G). and metabolomic data appeared not to be straightfor-
2Theseauthorscontributedequallytothearticle. ward (Amiour et al., 2012, 2014), the development of a
*[email protected]. model could help to identify putative candidate genes,
Theauthorresponsiblefordistributionofmaterialsintegraltothe proteins, and metabolic pathways contributing to plant
findingspresented in this article in accordancewith the policy de-
growth and development.
scribed in the Instructions for Authors (www.plantphysiol.org) is: MaizeisaC plant thatovercomestheinefficiencies
CostasD.Maranas([email protected]). 4
[C]Somefiguresinthisarticlearedisplayedincoloronlinebutin of Rubisco, to capture oxygen over the preferred CO2,
blackandwhiteintheprintedition. by separating the photosynthetic carbon fixation process
[W]TheonlineversionofthisarticlecontainsWeb-onlydata. intotwocelltypes:thebundlesheathandmesophyllcells.
[OPEN]Articlescanbeviewedonlinewithoutasubscription. In comparison with C plants, this separation allows C
3 4
www.plantphysiol.org/cgi/doi/10.1104/pp.114.245787 plants to have a lower rate of photorespiration, a higher
Plant Physiology(cid:1), November 2014, Vol. 166, pp. 1659–1674, www.plantphysiol.org (cid:3)2014AmericanSociety ofPlant Biologists. All RightsReserved. 1659
Simons et al.
rate of photosynthesis at high light intensities (under 2010), rice (Goff et al., 2002; Yu et al., 2002), Populus
standard air and temperature conditions), and a higher trichocarpa(Tuskanetal.,2006),sorghum(Patersonetal.,
photosynthetic nitrogen use efficiency (NUE; Christin 2009), Theobroma cacao (Tuskan et al., 2006), and maize
and Osborne, 2013; Driever and Kromdijk, 2013; Peter- (Schnable et al., 2009). Gene annotations of the whole-
hansel et al., 2013; Sage, 2014; Wang et al., 2014). A genome sequences have been used to determine the re-
C -specific maize GSM could provide insight into actionswithinanorganismandtherefore build a GSM.
4
NmetabolismandprovidecuesforimprovingNUE(i.e. FBA calculates all reaction fluxes in a metabolic net-
thevegetative biomassor grain yieldproducedperunit work based on the optimization of an objective func-
of N present in the soil). Since N is the major limiting tion (typicallythemaximizationofthebiomass yield).
factor in agricultural production among mineral fertil- A quasi-steady state is assumed, and flux constraints
izers (Vitousek et al., 1997; Hirel et al., 2007; Andrews aresetbasedonthespecificmediumorthereversibilityof
and Lea, 2013; Andrews et al., 2013) and NUE is esti- reactionsderivedfromthermodynamics.Incorporationof
mated to be far below 50% in cereal grains (Raun and omics data into GSMs is achieved through appropriate
Johnson, 1999), improving NUE is essential for improv- constraints on fluxes that restrict metabolic flows to only
ingoverall productivity in maize (Hirel and Gallais, condition-relevantphenotypes.
2011). Amiour et al. (2012) experimentally determined During the last few years, multiple methods have been
150genetranscripts,40proteins,and89metabolitesthat developed to integrate omics data into GSMs. Proteomic
are significantly different between the N+ WT and andtranscriptomicdatahavebeenusedtoapplyfluxcon-
2
N WTconditionsduringthevegetativestageofgrowth. straints on corresponding reactions determined by gene-
N utilization is strongly linked to the GS enzyme, as protein-reaction (GPR) associations. The GIMME (Becker
all N, either in the form of nitrate or ammonium ions, andPalsson,2008),iMAT(Shlomietal.,2008),andMADE
is channeled through the reaction catalyzed by the GS (Jensen and Papin, 2011) algorithms use a switch ap-
enzyme(Martinetal.,2006;Cañasetal.,2010;Hireland proach to turn on/off reactions based on expression
Gallais, 2011; Andrews et al., 2013). The mesophyll cell- levels. The GIMME algorithm turns offreactions based
specific GS1-3 isozyme is involved in synthesizing Gln on a user-specified threshold of the expression level.
after nitrate reduction from the vegetative state until the TheiMATalgorithmturnsonaminimalsetofreactions
plant reaches maturity. Leaf aging induces the synthesis associatedwithlowexpressiondatainordertoachieve
of the bundle sheath-specific GS1-4 isozyme. Conse- a user-specified metabolic function. The MADE algo-
quently,Martinetal.(2006)hypothesizedthattheGS1-4 rithmincorporatesrelatedexperimentaldatasetsintothe
isoform is used in the reassimilation of ammonium dur- model to activate or repress reactions based on the pro-
ingprotein degradation in senescing leaves. During gressionoftheexperimentalconditions.Adifferentclass
vegetative growth in the leaf tissue, DNA microarray of algorithms, known as the valve approach, was de-
datarevealedthat243genetranscripts,46proteins,and velopedtoincorporateproteomicandtranscriptomic
48 metabolites exhibited significant differences in the data by constraining the allowable flux ranges of reac-
gln1-3 mutants and 107 gene transcripts, 14 proteins, tions. The E-Flux method incorporates a user-specified
and 18 metabolites displayed substantial differences in function to convert gene expression data to flux con-
thegln1-4mutants(Amiouretal.,2014).Inthissecond- straints(Colijnetal.,2009).Finally,thePROMalgorithm
generation maize model, we explore the effect of the (ChandrasekaranandPrice,2010)usesmultipledatasets
computational knockout of genes encoding for GS1-3 toconstrainfluxbounds(i.e.allowablefluxranges)based
and GS1-4 isozymes using flux balance analysis (FBA) ontheprobabilitiesassociatedwithgeneactivityamong
to elucidate the role of GS in N metabolism. all data sets. Lee et al. (2012) integrated gene expression
FBA of GSMs is used to model organism-specific databyminimizingthedifferencebetweenthepredicted
metabolism by simulating the internal flow of metab- flux levels and gene expression data over all reactions
olites. The number of GSMs for plants has increased with corresponding expression levels. Using the Yeast
rapidly, with models available for Arabidopsis thaliana 5model(Heavneretal.,2012)forSaccharomycescerevisiae,
(Poolman et al., 2009; de Oliveira Dal’Molin et al., Lee et al. (2012) compared the predicted fluxes with ex-
2010a), barley (Hordeum vulgare) seed (Grafahrend- perimentally determined exometabolome fluxes using
Belau et al., 2009), maize (de Oliveira Dal’Molin thecoefficientofdeterminationr2.Theauthorsachieved
et al., 2010b; Saha et al., 2011), sorghum (Sorghum bi- r2valuesof0.87and0.96at75%and85%ofthemaximal
color; de Oliveira Dal’Molin et al., 2010b), sugarcane biomass level, respectively. In comparison, the authors
(Saccharum officinarum; de Oliveira Dal’Molin et al., generatedabestFBAsolution,whichmaximizesr2over
2010b), rapeseed (Brassica napus; Pilalis et al., 2011), all feasible solutions generated for FBA, and achieved
and rice (Oryza sativa; Poolman et al., 2013). These r2valuesof0.2and0.58at75%and85%ofthemaximal
models rely on annotation information to assemble biomass level, respectively. These advancements per-
comprehensive compilations of all reactions and metab- taining to the integration of omics data with GSMs has
olites known to occur within the organism. Currently, enabled more accurate model predictions.
whole-genome sequencing has been completed for ap- In this work, we describe the reconstruction of a
proximately 40 vascular plants, including A. thaliana second-generation maize leaf model and the incorpo-
(ArabidopsisGenome Initiative, 2000), Arabidopsis lyrata ration of omics data into the model with the goal of
(Huetal.,2011),soybean(Glycinemax;Schmutzetal., improving the understanding of N metabolism. Both
1660 Plant Physiol. Vol. 166, 2014
A Maize Leaf Model Applied to Nitrogen Metabolism
the primary and secondary metabolic pathways of RESULTSANDDISCUSSION
maize are included, by combining information from
EffectofNConditions onBiomassComponents
MetaCrop (Schreiber et al., 2012), MaizeCyc (Monaco
et al., 2013), and the earlier iRS1563 (Saha et al., 2011) Biomass components were measured in the N+ WT
models. In comparison with the iRS1563 model, this condition as well as for each N background (N2 WT,
second-generation model spans an additional 4,261 gln1-3, and gln1-4). Table I and Figure 1 display the
genes and 6,540 reactions. The increased number of composition of the classes of biomass metabolites, and
genesandreactionsenablestheinclusionofadditional Supplemental Table S1 indicates the specific biomass
pathways such as fructan biosynthesis, siroheme bio- measurementsinallmodeledconditions.Asexpected,in
synthesis, and ubiquniol-9 biosynthesis. The model the majority of cases, the N2 WT condition produced a
accountsforthetwomajorcelltypesintheleaf(i.e.the smaller concentration of biomass components than the
bundlesheathandmesophyllcells).Thebundlesheath N+ WT, gln1-3, and gln1-4 conditions. However, the
cell contains seven compartments: the cytosol, mito- concentrationofaminoacidsproducedwasabout5times
chondrion, peroxisome, chloroplast stroma, plasma higher in the gln1-4 mutant than the gln1-3 mutant,
membrane, thylakoid membrane, and vacuole. The resulting in comparable amino acid concentrations be-
mesophyll cellcontainssixcompartments:thecytosol, tweenthegln1-4 mutantandN+WTaswellasbetween
mitochondrion, chloroplast stroma, plasma mem- the gln1-3 mutant and N2 WT. The similar amino acid
brane, thylakoid membrane, and vacuole. Compart- concentrationsbetweenthegln1-4mutantandtheN+WT
mentalization is based on maize-specific experimental conditioninthevegetativestagehelptoconfirmthatthe
proteomic and transcriptomic measurements (Majeran GS1-4 isozyme is essential in plant maturity and has a
et al., 2005; Friso et al., 2010; Li et al., 2010; Chang smaller effect compared with the GS1-3 isozyme at the
et al., 2012), as opposed to the A. thaliana-based vegetativestage.Asexpected,theconcentrationofstarch
compartmentalization adopted in the previous iRS1563 was higher in the N2 WT condition than in the N+ WT
maize model (Saha et al., 2011). Light reactions have condition.UndertheN2WTcondition,thebreakdownof
been expanded from an aggregate reaction (as de- starch is limited by the amount of N available (Tercé-
scribedintheiRS1563model)tomultiplereactionsfor Laforgue et al., 2004; Amiour et al., 2012). Due to the
each complex with the inclusion of a thylakoid mem- limited N available, the starch is stored rather than bro-
brane compartment. In contrast to the C GEM maize ken down to produce other biomass components. The
4
model (de Oliveira Dal’Molin et al., 2010b), which fo- stained micrographs depicting the starch visible in the
cuses exclusively on primary metabolism in maize, N+WT,gln1-3mutant,andgln1-4mutantconditionsare
the developed model also spans secondary metabo- available in Supplemental Figure S1. The condition-
lismbyincludingallreactionsknowntooccurwithin specific biomass concentrations have been incorporated
the maize leaf tissue. The model includes as many as in the maize leaf model to more accurately represent
763 secondary metabolism reactions (without includ- metabolism under each condition.
ingduplicatecountingduetocompartmentalization).
Through the incorporation of omics data, regulatory
restrictions are introducedin the model to switch-off/
on reactions under the N+ WT and N2 WT conditions DevelopmentoftheSecond-GenerationMaizeLeafModel
and two GS knockout mutants (gln1-3 and gln1-4) in The second-generation maize leaf model was devel-
the vegetative stage, during which the plant absorbs oped using a combination of gene, protein, and reaction
and assimilates N for root and leaf biomass produc- information from thepreviously developedmaizemodel
tion (Amiour et al., 2012, 2014). Reactions linked to iRS1563(Sahaetal.,2011),biologicaldatabasessuchasthe
genes or proteins with significantly different expres- Kyoto Encyclopedia of Genes and Genomes (Kanehisa
sion levels between the N+ WT and N2 WT condi- et al., 2014), MaizeCyc (Monaco et al., 2013), and Meta-
tions, as well as the gln1-3and gln1-4 mutants versus Crop (Schreiber et al., 2012), as well as published litera-
the N+ WT condition, are conditionally turned on or ture sources. The model contains 5,824 genes and 8,525
off accordingly. The metabolite pool is simulated by reactions, a significant increase from the iRS1563 model,
maximizing the total flux through a metabolite (i.e. which contained 1,563 genes and 1,985 reactions. The
fluxsum)asaproxyforthemetaboliteturnoverrate second-generationmaizemodelissplitintotwocelltypes
(ChungandLee,2009).Thedirectionalchangesofflux- (i.e. the bundle sheath and mesophyll cells). The bundle
sum levels between the N2 WT condition and the N+ sheath cell is further divided into seven compartments,
WTcondition,aswellastheGSmutantconditionsand while the mesophyll cell contains six compartments
theN+WT condition,arequalitativelycomparedwith (Fig.2).Ofthe8,525reactionsinthemodel,3,892reactions
the directional change in experimentally measured are unique, as duplicated counts due to compartmentali-
concentration levels. These analyses reveal similar zationhave been disregarded. Of these 3,892 unique
trends to the recently developed flux imbalance anal- reactions, 1,012 reactions were assigned localization in-
ysis (Reznik et al., 2013), which makes use of dual formation based on transcriptomic and proteomic data
variable values associated with metabolite balances to (Majeran et al., 2005; Friso et al., 2010; Li et al., 2010;
infer the effect of concentration changes on the objec- Chang et al., 2012). Light reactions were adjusted to
tive function value. modeltheflowofprotonsacrossthethylakoidmembrane
Plant Physiol. Vol. 166, 2014 1661
Simons et al.
TableI. Experimentalcontentofclassesofmetabolitesindifferentconditions
Thebiomasscomponentsweredeterminedexperimentallyforeachoftheconditions(N+WT,N2WT,gln1-3mutant,andgln1-4mutant).Values
aremeansofthreereplicatesunlessindicatedbytheasterisk,indicatingthattworeplicatemeasurementsweretaken.Biomassmeasurementsforthe
specificmetaboliteswithineachclassaredisplayedinSupplementalTableS1.
BiomassComponents N+WT N2WT gln1-3Mutant gln1-4Mutant
Biomassyield(gdrywt) 6163.5 1562 6464.5 6561.5
Solubleaminoacidcontent(mmolg21drywt) 0.073260.0170 0.026160.0040 0.0212460.00100 0.0930360.00640
Proteincontent(mgg21drywt) 132.660.7 5863 125.2561.7 140.3965.72
Fattyacidcontent(mgg21drywt) 43.364.4 16.662.2 45.1* 16.361.1
Starchcontent(mmolg21drywt) 0.15260.005 0.19960.007 0.08560.011 0.10760.006
RNAcontent(mgg21drywt) 3.7860.19 0.9260.10 1.0560.09 1.7760.11
DNAcontent(mgg21drywt) 8.31560.270 2.5360.10 9.6260.22 5.4860.13
Solublecarbohydratecontent(mmolg21drywt) 0.23560.012 0.11260.013 0.19860.041 0.19360.023
Cellwallcarbohydratecontent(mgg21drywt) 0.3260.03 0.2660.09 0.18760.017 0.2960.07
Chlorophyllcontent(mgg21drywt) 1.8760.16 0.6960.06 1.7160.14 1.8560.08
TotalN(%gdrywt) 4.3660.08 1.8060.15 4.2860.10 4.3160.12
tothechloroplaststroma,torepresentthepHdifferential lipids(SupplementalTableS2).Compilingtranscriptomic
betweencompartments,andtodescribetheconversionof andproteomiccompartmentalizationdatawithliterature-
light to ATP (Nelson and Cox, 2009). The mitochondrial based pathways yielded a model of 4,103 reactions,
electrontransportchainwassimilarlyupdatedtoinclude leaving 2,880 unique reactions, still with their locali-
the proton exchange of ATP synthase between the inter- zations unknown.
membrane space and the mitochondrial matrix (Taiz, Once reactions were compartmentalized based on
2010).Finally,303specificreactionswereaddedtomodel transcriptomic data, proteomic data, and published lit-
glycerolipid synthesis, as shown in Supplemental Figure erature,thereactionsweredividedintotwogroups.The
S2 and Supplemental Table S2 (Moore, 1982; Murata, first group (core set) includes reactions with known lo-
1983; Murata and Tasaka, 1997; Mekhedov et al., 2000; calizations, while the second group (noncore set) spans
Bachlavaetal.,2009;Li-Beissonetal.,2010;Rollandetal., reactionsknowntooccurwithinthemaizeleafbutwith
2012).Tothebestofourknowledge,thisisthefirstplant no localization evidence. Whenever possible, core reac-
model to include detailed glycerolipid synthesis. Aggre- tionswereunblockedbyfirstaddingreaction(s)fromthe
gate reactions were included to link specific two-tailed noncore set to one or multiple compartment(s) and sec-
glycerolipids to the experimentally measured single ond appending intercellular or intracellular transporters
Figure 1. Weight percentage of biomass compo-
nents. The weight percentage for each class of
metabolites experimentally measured contributing
tobiomasssynthesisisdisplayed.Thecomposition
isdisplayedfortheN+WT(A),N2WT(B),gln1-3
mutant(C),andgln1-4mutant(D)conditions.The
measurementsforspecificcomponentswithineach
class of metabolites are shown in Supplemental
TableS1.[Seeonlinearticleforcolorversionofthis
figure.]
1662 Plant Physiol. Vol. 166, 2014
A Maize Leaf Model Applied to Nitrogen Metabolism
Figure 2. Number of metabolic and transport
reactions distributed between compartments in
thebundlesheathandmesophyllcelltypes.The
numbersofmetabolicandtransportreactionsare
shownforeachcompartment.Integralmembrane
proteins are counted for the compartment in
which the main biotransformation occurs. For
example, the ATP synthase associated with the
mitochondrialelectrontransportchainiscounted
asametabolicreactioninthemitochondrion,not
the inner mitochondrial membrane (IMM). [See
onlinearticleforcolorversionofthisfigure.]
(see “Materials and Methods”). By following this ap- remainingblockedcorereactionsandbiomassformation
proach, 1,032 unique reactions with previously un- by adding reactions from similar organisms (Krumholz
known localizations were assigned to compartments et al., 2012) and model organisms (i.e. rice ssp. japonica,
and 729 transporters were added. The remaining 1,848 Brachypodium distachyon, sorghum, and A. thaliana). By
uniquereactionswereassignedtocompartmentsbased addingfiveuniquereactionsfromsimilarorganisms,the
on available pathway information or assigned to the flux through three additional reactions known to be in
cytosol of both the bundle sheath and mesophyll cells. maizewasresolved.Thesereactionswereallinvolvedin
With all the reactions assigned to specific compart- the formation of Glu from His through urocanic acid.
ments, thermodynamically infeasible cycles that were The model is provided in a Microsoft Excel format in
generated due to the overly permissive inclusion of re- Supplemental Table S3 and in Systems Biology Markup
actions in the model, as well as lack of reaction direc- Language format in Supplemental Table S4.
tionality information, were subsequently identified and
eliminated. By first restricting the directionality of reac-
tions and second removing reactions, it was possible to
IncorporationofTranscriptomic andProteomic Datain
eliminate all thermodynamically infeasible cycles in the
theModel
model.Bythisprocess,werestrictedthedirectionalityof
36reactionsandremoved2,055reactionsfromthemodel InordertomoreaccuratelymodeltheN+WT,N2WT,
(Table II). Upon the resolution of thermodynamically and GS mutant conditions in maize, GPR associations
infeasiblecycles,attemptsweremadetounblockthe mapped the gene transcripts and proteins that were
TableII. Numberofreactionsaftereachmodelcreationandcurationstep
Theoriginaltwodatasetsarethecoresetandthenoncoreset,whichcombinetoformthefinalmodelstatistics.Thetotalnumberofmetabolic,
transport,exchange,andbiomassreactionsaredisplayedaftereachprocessduringmodelcuration.Metabolicreactiontotalsincludeduplication
fromcompartmentalization.
DataProcessing MetabolicReactions TransportReactions ExchangeReactions BiomassReactions
Initialdata
Coreset 3,002 418 82 85
Coresetplusmanuallycreatedpathways 3,264 469 285 85
Noncoreset 18,951 0 0 0
Processesperformed
Compartmentalizationalgorithm 3,971 1,198 285 85
Manuallydeterminedcompartmentalization 9,005 1,198 285 85
Thermodynamicallyinfeasiblecycles 7,033 1,115 285 85
SimilarorganismGapFill 7,040 1,115 285 85
Finalmodel
Second-generationGSM 7,040 1,115 285 85
Plant Physiol. Vol. 166, 2014 1663
Simons et al.
statisticallyexpressedatalowleveltoreactionsthatwere catalyzed by catalase (Boamfa et al., 2005). These two
turned off in the model. However, no essential reactions reactionshavea very slight effecton biomassformation,
to the model, which are required for biomass formation, asbiomassyielddropsbylessthan1%.Asexpected,we
were altered. For example, the d-aminolevulinic acid de- findthat manyof thereactionsthat correspondto genes
hydratase reaction was experimentally determined to be that are significantly down-regulated in the N+ WT con-
higherintheN+WTcondition,suggestingthatitshould dition do not hinder biomass formation. In the N2 WT
2
berestrictedintheN WTcondition.However,whenthe condition, none of the reactions have an effect on the
flux through the d-aminolevulinic acid dehydratase re- biomassyield,suggesting,asexpected,thatthedecreased
actionisrestrictedtozero,biomasscannotbeformed,as amountofNisthemainlimitingfactorinbiomassyield.
this reaction produces porphobilinogen, a precursor to In the gln1-3 mutant condition, three of the 100 reac-
chlorophyll (Gupta et al., 2013). Due to the incomplete tions, which are switched off based on omics data,
information available in published literature or affect the biomass yield. These three reactions are the
databases regarding possible alternative routes of the glyceraldehyde-3-phosphate dehydrogenase, Fru-bisP
production and degradation of a specific metabolite, aldolase, and Fru-bisphosphatase reactions. The ca-
regulatingreactionsthatareessentialtothemodelwill pacity of glyceraldehyde-3-phosphate dehydrogenase
restrict biomass synthesis. Based on experimental ev- toformamultienzymecomplexinthechloroplastsfora
idence, the fluxes through 83 reactions in the N+ WT range of plants is regulated by environmental condi-
2
condition, 20 in the N WT condition, 100 in the gln1-3 tions such as the light/dark transitions (Howard et al.,
mutant, and nine in the gln1-4 mutant were restricted. 2011). Glyceraldehyde-3-phosphate is synthesized dur-
The reactions regulated in the N+ WT condition mainly ing carbon fixation in photosynthesis, and 1,3-bisphos-
correspondtoreactionsknowntooccuronlyunderstress pho-D-gycerate (i.e. 3-phospho-D-glyceroyl phosphate)
andareexpressedatalowlevelincomparisonwiththe can besynthesizedfrom 3-phospho-D-glycerate.ATPis
2
N WTandmutantconditions.Reactionsthathavebeen requiredfortheconversionof3-phospho-D-glycerateto
down-regulatedbasedonomicsdataareindicatedinthe 1,3-bisphospho-D-glycerate catalyzed by 3-phospho-D-
model file (Supplemental Table S3). N perturbations glycerate kinase in the bundle sheath chloroplast. This
within the leaf tissue were modeled by combining the reactionisanimportantenergy-requiringreactioninthe
incorporationoftranscriptomicandproteomicdatawith Calvin-Benson cycle, as it is essential that the enzyme
the unique biomass composition for each condition. immediately metabolizes 3-phospho-D-glycerate, the
The minimal set of reactions, whose elimination productoftheRubisco reaction.Thisconclusionisalso
causesadecreaseinbiomassyield,wasdeterminedfor consistentwiththe findings that3-phospho-D-glycerate
the N+ WT, N2 WT, gln1-3 mutant, and gln1-4 mutant 1-phosphotransferase is sensitive to changes in energy
conditions.Therearesixreactionsacrosstheconditions state (Nakamoto and Edwards, 1987). The Fru-bisP al-
that encompass the minimal set of reactions, as sum- dolasereaction,whichisinvolvedintheCalvin-Benson-
marizedinTableIII.Ofthe83reactionswithrestricted Bassham cycle and the glycolysis pathway, can be
flux in the N+ WT condition, only two reactions were bypassed using the sedoheptulose 1,7-bisphosphate/
identified to affect biomass yield. These two reactions D-glyceraldehyde-3-phosphate-lyase reaction, which
are the conversion of ethanol to acetaldehyde through catalyzes the synthesis of sedoheptulose 1,7-bisphos-
either ethanol oxidoreductase involving NAD+ or a phateusingdihydroxyacetonephosphate(i.e.glycerone
hydrogen peroxide-dependent oxidation of ethanol phosphate) and D-erythrose 4-phosphate (Lakshmanan
TableIII. Summaryofreactionsthataffectbiomasssynthesis
Theminimumsetofreactionsthataredown-regulatedasaresultoftheinclusionofproteomicandtranscriptomicdataandaffectbiomasssynthesis
isdisplayed.Thecorrespondingconditionisdisplayedforeachreactionaswellastheroleofthereaction.
Reaction ConditionAffected AffectoftheReaction
Ethanoloxidoreductaseandethanolcatalase N+WT Producesacetaldehyde,alleviatingfluxthrough
pyruvatedecarboxylase
Glyceraldehyde-3-phosphatedehydrogenase gln1-3mutant Participatesinglycolysisandcarbonfixationbut
isnotrequired,as3-phospho-D-glycerate
kinasecanrestorefluxto1,3-bisphospho-D-glycerate
Fru-bisPaldolase gln1-3mutant ParticipatesintheCalvin-Benson-Basshamcycle
butcanbebypassedthroughthesedoheptulose
1,7-bisphosphate/D-glyceraldehyde-3-phosphate
lyasereaction
Fru-bisphosphatase gln1-3mutant DecreasesATP-ADPratio,switchesmetabolism
fromSuctostarchsynthesis,andinhibits
photosynthesisathighCO levelsin
2
A.thaliana
Rib-5-Pisomerase gln1-4mutant AffectscellulosesynthesisinA.thaliana
1664 Plant Physiol. Vol. 166, 2014
A Maize Leaf Model Applied to Nitrogen Metabolism
et al., 2013). The decreased expression of the cytosolic 2010),Serbiosynthesis(HoandSaito,2001),andtheurea
Fru-bisphosphatasereactionhasbeenshowntodecrease cycle (Mérigout et al., 2008) must decrease compared
the ATP-ADP ratio, lead to the switch from Suc to withtheN+WTcondition.Cholinebiosynthesis(McNeil
2
starchsynthesis,andinhibitphotosynthesisathighCO et al., 2001) is decreased in the N WT condition, in-
2
levels in A. thaliana, resulting in the inhibition of plant creasedinthegln1-3mutant,anddecreasedinthegln1-4
growth (Strand et al., 2000). Finally, the regulatory mutantcondition.FluxthroughIleandLeubiosynthesis
2
restrictions for the gln1-4 mutant involve only nine (McCourt and Duggleby, 2006) is lower in the N WT
reactions, of which one affected the biomass drain condition, higher in the gln1-3 mutant condition, and
(i.e.Rib-5-Pisomerasereaction).ThelackoftheRib-5-P lowerinthegln1-4mutantconditioncomparedwiththe
isomerase reaction has been experimentally shown to N+WTcondition,asexpectedbytheproportionofthese
cause premature death and affect cellulose synthesis biomass componentsin the various conditions. The flux
inA.thaliana(Howlesetal.,2006;Xiongetal.,2009).A through the glyoxylate cycle (Schnarrenberger and
comparison of the number of reactions that affect the Martin, 2002), stearate biosynthesis (Li-Beisson et al.,
GS mutants suggests that at the vegetative stage, the 2010), and urate degradation (Ramazzina et al., 2006) is
impact of the gln1-4 mutation is less severe than that higher in the gln1-3 mutant condition compared with
occurring in the gln1-3 mutant. Such a finding is not the N+ WT condition. Val biosynthesis (McCourt and
surprising, since it has been shown that the gene en- Duggleby, 2006)is lowerin thegln1-3 mutant condition
coding the GS1-3 isozyme is constitutively expressed compared with the N+ WT condition. Flux through
irrespective of the leaf development stage and that the glutathione biosynthesis/degradation, Trp biosynthesis
expression of the gene encoding the GS1-4 isozyme is (TzinandGalili,2010),uracildegradation(Zrenneretal.,
much lower and only enhanced at later stages of leaf 2006),andXyldegradation(Pennaetal.,2002)ishigher
development(Hireletal.,2005).Althoughonlyasubset in the gln1-4 mutant compared with the N+ WT condi-
ofreactionsaffectthebiomassproductionintheN+WT, tion. Glu is converted to glutathione through two ATP-
gln1-3 mutant, and gln1-4 mutant conditions, the addi- dependent steps requiring the addition of Cys and then
tional regulation will have an effect on the flux pre- Gly. Glutathione is a vitally essential protectant against
dictions within the model. oxidative stress, heavy metals, and xenobiotics (Noctor
etal.,2012;Rahantaniainaetal., 2013). Severalroutesof
glutathione breakdown have been proposed, including
FluxRangeVariations amongConditions theformationofCysandGlythroughcysteinyl-Gly.The
Cys is then degraded to form pyruvate, helping to alle-
Thefluxrangeofeachreactionwasdeterminedinthe viatethegln1-4mutation.Theincreasedfluxesassociated
N+ WT, N2 WT, gln1-3 mutant, and gln1-4 mutant with Xyl (from 1,4-b-D-xylan) and uracil degradation
conditions under the assumption that biomass is maxi-
generate a larger pool of xylulose-5-phosphate and
mized.ThefluxrangeofareactionintheN2WT,gln1-3 b-Ala, respectively. Finally, phenylpropanoid biosyn-
mutant, and gln1-4 mutant conditions was compared
thesis (Vogt, 2010) is lower in the gln1-4 mutant condi-
withthefluxrangeintheN+WTreferenceconditionto tion compared with the N+ WT condition. The majority
determine reactions with flux ranges that must deviate
of the changes in these pathways are directly related to
from the N+ WT flux range. This indicates that the flux
differencesintheproportionofthebiomasscomponents
through the reaction must change as a result of the
between the modeled conditions.
limited N or mutation. Overall, the flux through 202
2
reactionsintheN WTconditionisnotcontainedwithin
the flux range of the N+ WT condition, 765 reaction
ComparisonofModelPredictionswithMetabolomicData
fluxesinthegln1-3mutantdivergefromtheN+WTflux
range, and 678 reaction fluxes in the gln1-4 mutant The metabolomic data were compared with flux pre-
mustchangefromtheN+WTfluxrange(Supplemental dictions within the modelin each of the variousN back-
2
Table S5). In all three N backgrounds (i.e. the N WT, ground conditions. The increasing or decreasing trend of
gln1-3 mutant, and gln1-4 mutant conditions), the flux the metabolite concentration, displayed in Figure 3, was
comparedwiththeN+WTreferenceconditiondecreases qualitatively compared with the change in the flux-sum
undermaximumbiomassthroughthechlorophyllcycle, rangedeterminedbythemodel,asdisplayedinFigure4.
chlorophyllide a biosynthesis, farnesyl diphosphate bio- Thefluxsumisameasureoftheamountofflowthrough
synthesis, methylerythritol phosphate pathway, and the reactions associated with either the production or
tetrapyrrole biosynthesis. Tetrapyrrole biosynthesis, consumption of the metabolite. A variability analysis of
chlorophyllide a biosynthesis, and the chlorophyll cycle the flux sum was performed, and flux-sum ranges, nor-
linktheproductionofchlorophyllfromGlu(Tanakaand malizedbythebiomassrate,thatdonotoverlapbetween
Tanaka, 2007; Kim et al., 2013). The methylerythritol the N background condition and the N+ WT condition
phosphate pathway and farnesyl diphosphate biosyn- wereanalyzed.Anincrease/decreaseinthefluxsum(i.e.
thesis lead to a reactant required for the production used as a proxy for the metabolite pool) of a metabolite
of chlorophyll a from chlorophyllide a (Lange and betweentheN2WTconditionandthe N+ WT condition
Ghassemian,2003).InbothoftheGSmutantconditions, and between the two GS mutants and the N+ WT con-
thefluxthroughchorismatebiosynthesis(TzinandGalili, dition was compared with the metabolite concentration
Plant Physiol. Vol. 166, 2014 1665
Simons et al.
when the flux ranges are similar to the wild-type condi-
2 2
tion,asintheN WTcondition.BetweentheN WTand
N+WTconditions,onlyapproximately7%ofthereactions
activeineitherconditionhavefluxrangesatthemaximum
biomass that do not overlap. In the gln1-3 and gln1-4
mutant conditions, the fluxes are significantly perturbed,
with 49% and 45% of the active reactions at maximum
biomass resulting in nonoverlapping ranges compared
with the N+ WT condition, respectively. The accuracy of
flux sum in the gln1-3 mutant and gln1-4 mutant condi-
tions with omics-based constraints incorporated reaches
53% and25%,with eight of 15 metabolitespredictedcor-
rectlyandoneoffourmetabolitespredictedcorrectlyinthe
gln1-3 and gln1-4 mutant conditions, respectively. This
Figure 3. Number of metabolites in each condition that statistically levelofpredictionaccuracyisfarbelowwhatwasseenfor
variedfromtheN+WTconditionatthevegetativestage.Thenumbersof 2
N WT, suggesting a tenuous connection between con-
metabolites that experimentally significantly increased (up arrows) or
centration changes and gene expression levels when the
decreased(downarrows)incomparisonwiththeN+WTconditionare
2 geneticbackgroundchanges.
displayedforeachoftheNconditionstested(i.e.N WT,gln1-3mutant,
We explored the efficacy of the flux-sum method un-
and gln1-4 mutant conditions). The metabolites are shaded based on
whethertheyareinvolvedincarbon(C),N,orothermetabolism.[See derdifferentgeneticbackgroundsforamuchmorewell-
onlinearticleforcolorversionofthisfigure.] studied and data-rich organism (i.e. Escherichia coli) to
explorewhetherthedissonancebetweengeneexpression
changes. Figure 4 demonstrates the importance of levels and concentrations was maize specific or applied
restricting fluxes based on transcriptomic and proteomic broadly.Weappliedflux-sumvariabilitytotheIshiietal.
data.IntheN2WTcondition,theaccuracychangesfrom (2007) fluxomic and metabolomic data using the
13%to90%whenthefluxconstraintsbasedonomicsdata iAF1260 (Feist et al., 2007) E. coli model. Two single-
are incorporated. Without the incorporation of these con- gene knockout mutants (i.e. ppsA and glk) were com-
straints, all flux-sum ranges normalized by the biomass pared with the wild-type condition, and predicting the
2
rate are predicted higher in the N WT condition. The directional change of the metabolite pool size was met
identified flux-sum levels are included in Supplemental with less than 50% accuracy in each condition. This
Table S6. The flux-sum variability approach is able to impliesthatchangesinthegeneticbackgroundseemto
predictthechangeinmetabolitepoolsizesmoreaccurately cause concentration changes thatare not predictable by
Figure4. Effectofomics-basedregulationontheflux-sumpredictioncomparedwiththeexperimentaltrendinmetaboliteconcen-
tration.Theaccuracyinpredictingtheincreasing(uparrows)ordecreasing(downarrows)trendinmetabolitechangebetweentheN
backgroundconditionandtheN+WTconditionisdisplayed.Byrestrictingthereactionfluxbasedonthetranscriptomicandproteomic
data,theaccuracyofthequalitativetrendinmetabolitepoolsizebetweentheN2WTandN+WTconditionsincreases.Beforeadding
omics-basedconstraints,themodelwasabletocorrectlypredictthedirectionofchangein13%ofthemetabolitesmeasuredinthe
N2WTconditioncomparedwiththeN+WTcondition.Theaccuracyincreasesto90%whenomics-basedconstraintsareincluded.
Theflux-summethodisnotabletoaccuratelyrepresentthegln1-3andgln1-4mutantconditions,suggestingthatthegeneticback-
groundaffectstheabilityoftheflux-summethodtopredictmetabolitechanges.[Seeonlinearticleforcolorversionofthisfigure.]
1666 Plant Physiol. Vol. 166, 2014