Flux balance analysis of metabolism during growth by osmotic cell expansion and its application to tomato fruits

Cell expansion is a signiﬁcant contributor to organ growth and is driven by the accumulation of osmolytes to increase cell turgor pressure. Metabolic modelling has the potential to provide insights into the processes that underpin osmolyte synthesis and transport, but the main computational approach for predicting metabolic network ﬂuxes, ﬂux balance analysis, often uses biomass composition as the main output constraint and ignores potential changes in cell volume. Here we present growth-by-osmotic-expansion ﬂux balance analysis (GrOE-FBA), a framework that accounts for both the metabolic and ionic contributions to the osmotica that drive cell expansion, as well as the synthesis of protein, cell wall and cell membrane components required for cell enlargement. Using GrOE-FBA, the metabolic ﬂuxes in dividing and expanding cells were analysed, and the energetic costs for metabolite biosynthesis and accumulation in the two scenarios were found to be surprisingly similar. The expansion phase of tomato fruit growth was also modelled using a multiphase single-optimization GrOE-FBA model and this approach gave accurate predictions of the major metabolite levels throughout fruit development, as well as revealing a role for transitory starch accumulation in ensuring optimal fruit development.


INTRODUCTION
Flux balance analysis (FBA), a method for predicting and analysing steady-state metabolic fluxes, has been widely applied in the study of plant metabolism (Sweetlove and Ratcliffe, 2011;Nikoloski et al., 2015;Gomes de Oliveira Dal'Molin and Nielsen, 2018). The approach requires a matrix of metabolic reaction stoichiometries and an objective function that represents the optimization target of the biological system (Feist and Palsson, 2010). In simple unicellular microorganisms growing in a nutrient-rich environment, where life-cycle events involve only growth and reproduction, the maximization of flux representing the accumulation of biomass components has proved to be a reasonable objective function (Varma and Palsson, 1994;Feist et al., 2007;Feist and Palsson, 2010). In complex biological systems, however, the accumulation of biomass may not be the primary purpose of every cell type in the organism. For example, the principal metabolic objective of fully expanded source leaves is the biosynthesis of sucrose and amino acids for the rest of the plant (Cheung et al., 2014). Moreover, organ development in plants typically involves phases of cell differentiation, cell division and cell expansion (Gonzalez et al., 2012), and these are not temporally synchronous, meaning that at different stages of organ development different mechanisms of growth are dominant. Metabolism, being closely related to the demands of the cell, is thus likely to vary between these stages of development (Sweetlove and Ratcliffe, 2011;Nikoloski et al., 2015). The cell expansion stage is responsible for the main increase in organ size and metabolic content in plants, and it is the dominant mechanism when growth is measured during plant phenotyping (Fahlgren et al., 2015;Tardieu et al., 2017). This makes understanding metabolism in expanding cells of great interest and relevance to breeding and crop engineering. In its conventional form FBA does not take account of the changing volume of the cell, indicating the need for a new FBA formulation if this objective is to be achieved.
To provide a biological context, we focused on tomato, Solanum lycopersicum, which has an extensive phase of cell expansion during fruit development. This process has been described in detail at the molecular-biochemical level (Valle et al., 1998;Carrari and Fernie, 2006;Legland et al., 2012;Biais et al., 2014), and the available biochemical data can be used to guide the new FBA formulation. Tomato is also an important model for fleshy fruit development and ripening, owing to its agronomic value, ease of cultivation and diploid genetics. FBA has been previously used to model steady-state snapshots of tomato fruit metabolism at different stages of development (Colombi e et al., 2015). This was achieved by using timeseries metabolomic data sets to constrain metabolite accumulation and degradation rates in a model of primary metabolism. This snapshot approach is informative and captures the accumulation of osmolytes for cell expansion in the form of experimental constraints, but it has some disadvantages. First, accurate predictions of metabolic state, including the accumulation of solutes during cell expansion, can only be made by imposing a large number of experimentally measured constraints. Second, the majority of the flux predictions are a direct consequence of the constraints imposed. For example, the rate of starch degradation observed in the model was the result of a constraint that dictates that starch is consumed at a set rate. Third, because each developmental point is modelled without reference to future needs, the context of the metabolic events at a given time point within the overall development programme is lost. In our previous work on fully expanded source leaves we highlighted a similar issue with modelling day and night metabolism separately, and demonstrated that modelling day and night phases simultaneously as a single-diel FBA problem led to better predictions of leaf metabolism (Cheung et al., 2014). Applying a similar multiphase FBA approach to metabolism during tomato fruit development should therefore improve the usefulness of the computational approach.
In this study, we developed a framework for modelling metabolism in expanding cells using FBA. This approach, growth-by-osmotic-expansion FBA (GrOE-FBA), accounts for the metabolic and ionic contributions to the osmotica that drive cell expansion, as well as the synthesis of protein, cell wall and cell membrane components required for cell enlargement. We show how GrOE-FBA can be used to identify the metabolic network fluxes in expanding cells, highlighting the major metabolic differences between dividing and expanding cells. We also show how GrOE-FBA can be combined with a multiphase single-optimization FBA approach to study plant metabolism during organ development. This approach provided evidence that transitory starch stores are necessary for optimal fruit development, and that prevention of starch accumulation during the early stages of fruit development may result in smaller fruits owing to the reduced phloem uptake rate during later developmental stages, where fruit expansion has been observed to be maximal.

Modelling cell expansion using GrOE-FBA
Cell expansion is driven by the accumulation of osmolytes and the resulting influx of water (Boyer et al., 1985). There are two key equations (see Appendix S1 for the derivation of all the equations) that describe the relationship between the osmotic content of a cell and its volume, and which ultimately allow FBA to be used to model metabolism in expanding cells. First, assuming for convenience that the intracellular solutions behave as ideal solutions, the total osmotic content of a cell at steady state is equal to the product of its osmolarity (C cell ) and volume (V cell ): where n i and m i are the van't Hoff factor and the number of moles of metabolite i, respectively, in the cell. This equation relates changes in cell volume to changes in osmotic content, and for predicting fluxes it needs to be used in tandem with a second equation that considers the distribution of solutes between the cytosol and the vacuole: where metabolites j and k are vacuolar and cytosolic, respectively. This equation states that at steady state, the ratio of the osmotic content of the vacuole and cytosol is equal to the ratio of their volumes, V v and V c , respectively. Equations 1 and 2 provide the link between volume and metabolic content at cellular and subcellular levels, and they imply that including a representation of osmotic content in an FBA model should be sufficient to account for volume changes. This was achieved by introducing two pseudo-metabolites in the metabolic network to represent the osmoles associated with the accumulation of osmotically active metabolites/ions in the vacuole and cytosol (O vac and O cyt , respectively; Figure 1a). The pseudometabolites take account of the expected differences in the van't Hoff factor by constraining the model so that there can be no net change in the total charge of the vacuole or cytosol. To implement these constraints, two new pseudoreactions, 'aggregator reactions', were included in the model. The first aggregator reaction: As well as accumulating osmolytes, expanding cells synthesize extra cell wall, cell membrane and proteins to maintain cell functions. Estimates of the additional biomass were obtained by creating a simple geometric representation of a cell, assuming: (i) cells are cubic; (ii) the cell wall is uniformly thick and composed of cellulose only; (iii) the cell membrane is uniformly thick; (iv) the fraction of the total protein content of the cell in the vacuole is negligible; and (v) the cytoplasmic protein concentration is maintained during cell expansion. Based on these assumptions, if V is the volume of the cube-shaped cell, b is the thickness of the cell wall and c is the thickness of the cell membrane (Figure 1b), then the amount of cellulose, phospholipid and protein in the cell can be estimated from the following equations (Appendix S1): where moles cellulose is the amount of cellulose in moles, q cellulose is the density of cellulose, MW glucose is the molecular weight of glucose, moles phospholipid is the amount of cell membrane phospholipids in moles, f PM is the fraction of the total lipid content found in the plasma membrane, q phospholipid is the density of cell membrane phospholipids, MW phospholipid is the molecular weight of membrane phospholipid, moles protein is the amount of protein in moles (based on how a unit protein is represented in metabolic models), C protein is the molar concentration of protein in the cytoplasm (based on the molar mass of a unit protein as represented in the metabolic model) and V cytoplasm is the volume of the cytoplasm. The demand for osmolytes and biomass components required to support a change in cell volume can be estimated from the difference in the amount of these metabolites for the initial and final cell volumes. For example, the cellulose demand for when a cell changes its volume from V 1 to V 2 can be calculated as follows: where moles cellulose; V1 and moles cellulose; V2 are moles of cellulose when the cell has a volume of V 1 and V 2 , respectively, and time ; V 1!V 2 is the time that the cell takes to change its volume from V 1 to V 2 . In this manner, by combining previously described osmotic constraints (based on Equations 1 and 2) with biomass constraints (based on Equations 3-5), it is possible to perform FBA while accounting for cell volume.
Validation of the equations used to estimate the change in biomass during cell expansion The cell wall, lipid and protein contents of tomato pericarp cells during fruit development were determined (Appendix S2) and the results were compared with predictions based on Equations 3-5 ( Figure 2). Pericarp cell membrane phospholipid composition was assumed to be composed of phosphatidyl ethanolamine (PE), phosphatidyl choline (PC) and phosphatidic acid (PA), based on published data in cherry tomato (Guclu et al., 1989), and the relative quantities of PE, PC and PA were used to estimate values for q phospholipid and MW phospholipid . The value of f PM was estimated from data published in Nicotiana tabacum (tobacco) leaves (Cacas et al., 2016) and was assumed to be constant throughout fruit development. Values of other parameters required in the equations were collected from published data or were measured experimentally (Table S1). Figure 2 shows that Equations 3-5 generated a good match to experimentally measured values for tomato fruit. The slight underestimation of fruit cellulose content is likely to arise from the assumption of the cells being geometric cubes with a constant cell wall thickness, which does not account for increased cell wall thickness around the corners of the cell (Legland et al., 2012). Thus, a simple geometric model can be used to predict increases in cellulose, phospholipid and protein content during cell expansion, and this allows these major biomass outflows from the metabolic network to be estimated, avoiding the need for direct measurements.

Application of GrOE-FBA to expanding tomato cells
Due to the energy cost of synthesising cellular macromolecules, such as proteins, lipids and carbohydrates (Schwender and Hay, 2012), it is generally assumed that the synthesis of new cells is more expensive than expanding existing cells (Taiz, 1992;Lynch, 2019). Even during cell expansion, however, new macromolecules need to be produced and the synthesis and accumulation of osmolytes to drive cell expansion also represents a significant cost. The GrOE-FBA framework allows the costs of cell expansion to be computed and the result can be compared with the costs of cell division computed using conventional FBA with a biomass objective. To do this, we modelled rapidly dividing tomato pericarp cells in culture using existing data (Rontein et al., 2002) and compared the resulting model with one of rapidly expanding pericarp cells in tomato fruit. Both models were constructed using a minor updated version (PlantCoreMetabolism_v1_2) of our previously published core stoichiometric model of primary plant metabolism (Shameer et al., 2018). Given the focus on tomato, the model contains tomato gene-to-reaction associations, although these were not used in this study. The dividing cell model was provided with the nutrients present in the cell culture medium: glucose as the sole carbon source, as well as ammonium, oxygen, phosphate and sulphate. The expanding cell model was provided with sugars and amino acids in the same relative proportions as reported for tomato plant phloem sap (Walker and Ho, 1977a;Valle et al., 1998), as well as oxygen, phosphate and sulphate. The same objective function of minimization of the sum of fluxes was applied to both models. Experimental data were used to constrain the rate of biomass production in the dividing cell model to 2 mg dry weight (DW) ml -1 day -1 , which was the fastest rate reported in the Rontein et al. (2002) study. The rate of cell expansion in the expanding cell model was set to that of fruit at 26 days post anthesis (DPA), where the fastest rate of fruit growth was observed . The initial exploration of the expanding cell model included an additional component to the objective function: to maximize the organic solute content while satisfying the osmotic constraint (see subsequent explanations and discussions). Figure 3 summarizes the results of this comparison for equivalent numbers of cells, and the complete list of predicted fluxes is provided in Appendix S3. Comparison of the metabolic fluxes in the two models revealed that dividing cells have higher fluxes of cellulose, protein, starch and lipid biosynthesis, as expected, whereas the expanding cells accumulated significantly more organic solutes, both in the vacuole and cytosol (Figure 3a,b). This accumulation of organic solutes in the expanding cell model increased the carbon demand of the model, and somewhat surprisingly, the overall rate of glucose consumption was higher in the expanding cell model (2.26 mg ml -1 day -1 ) than in the dividing cell model (1.43 mg ml -1 day -1 ). The flux maps (Figure 3a,b) reveal that a substantial proportion of this additional carbon (68% of carbon taken up by the cell) was transported into the vacuole to satisfy the osmotic requirement for cell expansion. But it is also apparent that there were high fluxes of glycolysis, and an even higher tricarboxylic acid (TCA) cycle flux, relative to the dividing cell model, suggesting that the energy demand of cell expansion is comparable with that of cell division. To provide a more precise quantitative comparison of energy costs, we constructed an energy budget for the two systems by collating all the fluxes in the models that required ATP consumption or led to ATP production ( Figure 3d). The total ATP demand of the expanding cells was 0.13 mmol ATP ml -1 day -1 and that of the dividing cells was 0.14 mmol ATP ml -1 day -1 . This showed that indeed both models have very similar energy demands. Tomato fruit development from 13.1 to 43.7 DPA (with mature green fruits observed at 42 DPA) is mainly the result of the expansion of existing cells, and from the breakdown of energy expenditure (Figure 3d) it is apparent that although the expanding cell model devoted considerably less ATP to the synthesis of protein, lipid and cell wall (12% compared with 44% in dividing cells), this was offset by a large increase in the cost of solute biosynthesis and accumulation, which represented 39% of the ATP budget in the expanding cell model compared with only 11% in dividing cells ( Figure 3d). Given that the maximization of organic solute was a component of the objective function of the expanding cell model, the very high rate of accumulation of glucose to satisfy the osmotic constraint may be a feature of how we set up the model. Although it is known that fruit cells contain high quantities of hexose sugars (Biais et al., 2014), it is likely that inorganic solutes contribute to the osmotic balance and this would affect the total energy cost of osmolyte production and accumulation. To explore this we removed the objective for the maximization of organic solute and just used the minimization of the sum of fluxes as the sole objective, as in the dividing cell model. This led to a significantly lower rate of carbon consumption in the expanding cells (0.63 mg C ml -1 day -1 ), but the overall ATP demand remained comparable with that of the previous model (0.11 mmol ATP ml -1 day -1 ). The predicted flux distribution ( Figure 3c) revealed that the model was exclusively using inorganic ions to satisfy the osmotic constraints.
In reality expanding cells are expected to accumulate both organic and inorganic solutes to facilitate cell expansion, so the predicted energy demand of the expanding and dividing cells are comparable. These results demonstrate that cell expansion is not necessarily a cheap option for plant growth, as previously argued (Taiz, 1992).
Many similarities between the three systems (dividing cells, expanding cells accumulating organic solute and expanding cells accumulating solely inorganic solute) were apparent. For example, in each system, more than 75% of the ATP demand was met by mitochondrial ATP synthesis (75, 73 and 78%. respectively) and non-growthassociated maintenance (NGAM) was found to be a significant energetic drain (45, 49 and 55% of the total ATP demands, respectively). The second largest energetic demand was the ATP cost of maintaining the plasma membrane proton gradient (PM-ATPase) during nutrient uptake (18, 23 and 29% of the total ATP demands, respectively). Although the basis for NGAM was identical in the three models, the PM-ATPase flux was dependent on the carbon demands of the model. As a result, the energy demand of the system is significantly influenced by the carbon demand for biosynthesis. Note that NGAM was calculated based on an assumed 3:1 ratio of ATP:NADPH expenditure on maintenance (Cheung et al., 2013), but the general conclusion about the comparative overall ATP expenditure in dividing versus expanding cells is not altered if different maintenance ATP:NADPH ratios are used ( Figure S1).

A multiphase metabolic model of primary metabolism in developing tomato fruit
To investigate the changing metabolic requirements for growth by cell expansion during organ development, pericarp metabolism was modelled using GrOE-FBA in the expansion and ripening stages of developing tomato fruit by dividing the time course from 8 DPA to red ripe stage (59 DPA) into 10 developmental substages. Ten copies of PlantCoreMetabolism_v1_2, each representing one of these substages, were combined to generate a 10-phase constraint-based model of primary metabolism in tomato pericarp cells. To allow for the accumulation and use of sugars (glucose, fructose and sucrose), amino acids (including c-aminobutyric acid, GABA) and organic acids (malate and citrate) in each phase, reactions transferring metabolites and ions from one phase to the next, called 'linker reactions', were added to the cytosol and the vacuole of the model. Similarly, the transfer of starch from one phase to the next (allowing for the accumulation and use of starch) was enabled via linker reactions in the plastid. The linker reactions and the multiphase nature of the model ensure that metabolism in each stage of fruit development depends not only on metabolism in the previous developmental stages but is also influenced by that of future developmental stages. Each phase in the multiphase model also had access to the xylem and the phloem, and was subjected to biomass requirements. A full schematic representation of the model is presented in Figure 4.
The initial metabolite content of the fruit pericarp at 8 DPA was established by adding source reactions (reactions that represent the import of metabolites from outside the modelled system) for sucrose, fructose, glucose, amino acids, GABA, organic acids and starch to the first phase of the model. The flux through these source reactions was constrained using the experimentally determined soluble metabolite content of 8-DPA fruits (Appendix S2). The soluble metabolite content of red ripe tomato fruits (59 DPA; Appendix S2) was used to constrain the relative proportions of the soluble metabolites in the final model phase of the model by introducing appropriately constrained sink reactions. Previously published data on vacuolar pH during fruit development in cherry tomatoes (Rolin et al., 2000) was used to set the vacuolar pH for the 10 phases and thus the abundance of the different metabolite charge states for each phase. A data-based constraint was also applied for the maximal rate of uptake of nutrients into the fruit from the phloem (see Experimental procedures). The multiphase developing fruit model was solved as a single parsimonious FBA (pFBA) problem, an approach similar to previously published multiphase (Cheung et al., Shaw and Cheung, 2018) constraintbased models, with the maximization of fruit soluble metabolite content as the objective. Fruit content and composition is key in increasing yield in tomato cultivation, and hence the maximization of flux through the reaction representing the accumulation of soluble metabolites in ripe tomatoes was deemed an apt objective function for the model. NGAM, nutrient uptake from the phloem and the osmotic constraints based on Equations 1 and 2 were implemented as described in the Experimental procedures.
A multiphase developing tomato fruit model predicts realistic final fruit metabolite content Figure 5 depicts the influence of the various constraints on the predictions of the multiphase developing fruit model. Optimal flux distributions of the model without the osmotic and biomass constraints to describe cell expansion and with no upper limit on phloem uptake rates predicted unrealistically high fruit organic solute content (approx. 5000 mOmol per fruit; Figure 5 upper left). The introduction of osmotic constraints on the system greatly reduced the final fruit organic solute content to 24.5 mOsmol per fruit ( Figure 5 upper right). This value was just 7.5% lower than the measured value (26.5 mOsmol per fruit). The imposition of biomass demand fluxes did not affect the predicted ripe fruit organic content, although it did cause small variations in predicted fruit content during the early stages of fruit expansion ( Figure 5 lower right). Limiting the phloem uptake rate had no effect on the predicted organic content of ripe fruits, but it did reveal a requirement for a transitory carbohydrate store in expanding fruits ( Figure 5 lower left). The maximum level of starch that was predicted to accumulate was 8.3 mmol hexose equivalents per fruit, corresponding to 1495 mg per fruit. This is 2.8 times the value of 530 mmol hexose equivalents per fruit reported in an experimental study (Petreikov et al., 2009).
In order to understand the factors driving metabolism during fruit development, predicted fruit metabolic content during development was analysed (Table 1). From Table 1 it can be seen that the soluble organic content contributes the most quantitatively significant fraction of fruit content, particularly during the later stages of fruit development. The optimal flux distribution and a fully constrained version of the multiphase model are available in Appendixes S3 and S4, respectively.
Comparing model predictions with experimental measurements Figure 6 compares the fruit content of selected metabolites predicted by GrOE-FBA with the experimentally determined values. The complete data set for 22 solutes is provided in Figure S2. The model was capable of making accurate predictions for glucose, fructose, glutamine and glutamate, despite being free to choose which metabolites it accumulated to satisfy the osmotic constraint. These four metabolites formed 88.6% of organic osmolytes in the pericarp cells and their accurate prediction by the osmotically constrained model suggests that their pattern of accumulation in developing tomato fruit represents an efficient way to drive osmotic expansion. Some other metabolite levels, including those of malate and citrate, were predicted with reasonable accuracy during the later stages of fruit development, whereas for others, including sucrose and GABA, it was not possible to predict the observed pattern of accumulation. This suggests that the accumulation of these metabolites is not primarily driven osmotically and that other aspects of fruit physiology may drive their accumulation during fruit development. For example, factors such as insect resistance, texture and flavour all have an impact on determining fruit composition during development (Tohge et al., 2014;Cohen et al., 2014;Takayama and Ezura, 2015). Translating these factors into constraints for a metabolic model would be challenging, as the measured content of some of the poorly predicted amino acids (aspartate and serine) fell within the FVA ranges of the model ( Figure S2), but applying extra constraints to the system may improve the predictive performance of the model.

A transitory carbohydrate store is required when phloem nutrient uptake is limited
If nutrient influx from the phloem was unlimited, then the nutrient requirement for fruit development was predicted to peak in phases 4 and 5 (23.3-33.5 DPA; Figure 7a). In contrast, when an upper bound on phloem uptake (Figure 7b) was imposed, based on experimental data, the model predicted the need for transitory carbohydrate storage (Figure 7d), which is a well-documented behaviour of tomato fruits (Ho and Hewitt, 1986;Gillaspy et al., 1993). When phloem uptake rate limits were removed from the Figure 5. Impact of model constraints on the predicted metabolite contents of developing Solanum lycopersicum (tomato) fruit. The charts show the sequential effect on the predicted metabolite contents of adding osmotic, biomass and phloem uptake constraints to the model. In the absence of these constraints the predicted total metabolite content greatly exceeded the expected value, whereas after applying the constraints there was close agreement between the measured and predicted values. AA, amino acids; OA, organic acids; SUG, sugars (glucose, fructose and sucrose).
© 2020 The Authors.  (Figure 7c), suggesting that the accumulation of starch was associated with phloem influx limits and metabolic demand. Moreover, preventing the accumulation of starch caused a 40% drop in the organic content of the fruit, demonstrating the significance of transitory starch storage in developing tomato fruits.

GrOE-FBA: proof of concept in the comparison of expanding and dividing cells
GrOE-FBA was developed to predict metabolic network fluxes in cells undergoing osmotically driven cell expansion. As a proof of concept, GrOE-FBA was used to highlight the differences between metabolism in dividing and expanding cells, and the analysis led to the conclusion that cell expansion does not necessarily put a smaller requirement on cellular resources than cell division. Although it is not surprising that differences existed between metabolism in dividing and expanding cells, quantitative network-level comparisons have not been made before. The two systems chosen here for comparison were dividing cells in heterotrophic cell suspension cultures and expanding cells in developing tomato fruits. A more appropriate comparison would be to model metabolism during the very early phase of tomato fruit development, when cell division dominates, and then again for the later phases, when expansion dominates. There is very little biochemical data to constrain such a cell division model, however, because of the very small size of the young plant organs in the cell division phase. Even so, the marked difference in the predicted metabolic fluxes of dividing and expanding cells (Figure 3) implies that conventional FBA studies on plant tissues, such as growing leaves (Yuan et al., 2016;Shaw and Cheung, 2018) and roots (Grafahrend-Belau et al., 2013;Shaw and Cheung, 2018), need careful interpretation, even if they include experimental metabolite constraints that would indirectly account for the generation of osmotic pressure for cell expansion.

GrOE-FBA: application to tomato fruit development
The GrOE-FBA approach was combined with multiphase modelling to analyse metabolism in developing tomato fruit. The approach required minimal constraints: a knowledge of cell volume changes during expansion, an estimate of the volume ratio of vacuole to cytosol, and an estimate of protein concentration (C protein , Equation 5). Using these constraints, it was possible to make surprisingly accurate predictions about the accumulation of biomass components and the principal osmotic metabolites that drive cell expansion in tomato fruit (Figures 5 and 6). The same approach could be applied to the growth of other organs such as leaves and roots.
Previous modelling studies taking into account cell expansion and metabolism were based on an analysis of enzyme kinetics and were capable of exploring only limited aspects of metabolism because of the difficulty of gathering or estimating large numbers of kinetic parameters . Previous work using FBA modelling to analyse metabolism in developing tomato fruits used separate models for each stage of development (Colombi e et al., 2015). Such highly constrained models can provide answers to the question 'What is happening?' but cannot answer the question 'Why is this happening?'. For example, with this approach it has been shown that the breakdown of transitory starch is the underlying cause of the respiratory climacteric in tomato fruit (Colombi e et al., 2017), but it was not possible to explain the purpose of the transitory accumulation of starch. This is because starch accumulation and degradation were constraints of the model and not predictions. In contrast, starch accumulation and degradation in the multiphase GrOE-FBA model are system-level predictions that depend on the metabolic demand on the system and the availability of nutrients from the phloem, throughout development.

GrOE-FBA: insights into tomato fruit development
Two major biological insights emerged from the analysis of the GrOE-FBA model of tomato fruit development. First, we demonstrate that transitory carbohydrate stores are only predicted by the model if a constraint is added that gradually decreases the rate of influx of nutrients from the phloem. Our model demonstrates that because of the decline in phloem influx during development, the accumulation of starch (or other carbohydrates such as sucrose) in the earlier stages of fruit development is required to meet the carbon demand in later stages. The question then arises as to why the phloem influx into a fruit declines below the level required to maintain its metabolic requirements. One possibility is that because of the staggered initiation of fruits along a truss, it is necessary to progressively reduce the share of phloem nutrients taken by maturing fruit to allow younger fruits, which have lower sink strength, to develop successfully. Temporal extension of starch synthesis has been reported to result in increased fruit size (Petreikov et al., 2006(Petreikov et al., , 2009. This is in agreement with our hypothesis. Increased sucrose demand owing to starch synthesis will increase the sink strength of smaller fruits. Once the fruits are bigger, the greater transitory starch stores will reduce their sink strength, resulting in more phloem constituents being available to the smaller fruits. The second biological insight is the relatively minor impact biomass demand had on developing fruit metabolism ( Figure 5). It is known that cell expansion results in an increase in cell wall and cell membrane content. The protein content of cells is also thought to increase in order to maintain the optimal concentration required for unhampered metabolism. In most FBA models of growing plant tissues published to date, the synthesis of these biomass components is the main drain on the metabolic system. The GrOE-FBA model demonstrates that for tissues growing by the expansion of existing cells, the demand for these biomass elements imposes only a small fraction of the metabolic cost. Instead, the dominant demand on metabolism was solute accumulation for osmotically driven cell expansion. The accumulation of solutes was responsible for most of the gain in fruit DW from 18.0 to 43.7 DPA (Table 1). From this we can conclude that the drivers for metabolism in expanding cells are significantly different from those in dividing cells, and that osmolaritybased constraints such as GrOE-FBA applied in this study are necessary to model metabolism in tissues growing by cell expansion.

GrOE-FBA: limitations of the approach
GrOE-FBA is a modified form of FBA and hence inherits its limitations (Sweetlove and Ratcliffe, 2011), including not being able to predict futile cycles. A substantial fraction of the energy consumption in a tomato cell suspension has been attributed to the operation of such cycles (Rontein et al., 2002) and the FBA model is unable to capture this energy expenditure. Although this might explain the discrepancy between the glucose consumption rate of the model (3.57 mg ml -1 day -1 ) and the value reported in the literature (glucose influx rate for day 4 is 5.34 mg ml -1 day -1 ), it should also be noted that the method used to estimate the futile cycles has been questioned because of the complicating effects of the subcellular structure of the metabolic network (Kruger et al., 2007).
Another limitation of GrOE-FBA is that it focuses on the osmotic significance of the accumulating solutes and makes no attempt to tackle the challenges involved in imposing constraints on metabolism related to insect resistance, texture and flavour. Nevertheless, introducing such constraints into the model might improve its ability to predict sucrose and organic acid levels. Finally, the GrOE-FBA model predicts metabolic fluxes by flux optimization and does not account for signalling and the regulation of metabolism, although these could be included in the form of experimental constraints on maximal enzyme activity. As a result, the model we present does not activate flux through energy dissipation mechanisms (such as alternative NADH dehydrogenases and alternative oxidase), and hence the GrOE-FBA model is unable to predict the respiratory climacteric (Colombi e et al., 2017). A computationally expensive FBA approach termed cost-weighted flux minimization has been shown to be capable of predicting alternative pathways in source leaves (Cheung et al., 2015). Such an approach might allow the multiphase GrOE-FBA tomato model to predict metabolically inefficient processes, such as uncoupled respiration and futile cycles.

CONCLUSION
GrOE-FBA, an FBA approach that uses volume-based osmotic constraints, allows metabolism to be modelled during cell expansion, and thus extends the scope of constraints-based metabolic modelling to growing tissues in which the expansion of existing cells is the dominant driver of growth. Understanding metabolism in expanding cells could help develop engineering strategies to optimize metabolism during the cell expansion phase of leaf and root development. The approach could also be used to study metabolism in guard cells where turgor-pressure-driven changes in cell size regulate stomatal aperture, in turn determining CO 2 assimilation and transpiration rates. GrOE-FBA could therefore help complement current kinetic models of guard cells (Hills et al., 2012) with a more holistic modelling approach. Recently, there has been increased interest in combining metabolic models with whole-plant developmental models (Marshall-Colon et al., 2017). Using GrOE-FBA to account for the metabolism of expanding cells could help to improve the relevance of the metabolic models in the context of whole-plant growth and potentially improve the predictive power of integrated metabolic whole-plant growth models. systematically removed. Each flower anthesis was recorded, and trusses were pruned at five developed fruits to limit fruit size heterogeneity. Tomato fruits were harvested at nine different developmental stages corresponding to 8,15,22,28,34,42,50,52 and 59 DPA, with the last four corresponding to mature green, turning, orange and red fruit developmental stages, respectively. Fruit samples were collected from the first two levels of trusses in the plant. All materials were frozen in liquid nitrogen and storage at À80°C until use.

Fruit biomass measurements and morphology analysis
Cell wall was extracted according to an established protocol (Ruprecht et al., 2011). The remaining insoluble material was washed with water and ethanol, air dried and weighed. Lipids were extracted from a known mass of ground tissue using the chloroform/methanol protocol (Bligh and Dyer, 1959). Protein, extracted with 6 M urea, 2 M thiourea buffer (Salem et al., 2016), was quantified using the Bradford assay (Bradford, 1976). The amino acid content of protein hydrolysates (50% w/v trichloroacetic acid, 6 M HCl, 100°C/24 h) (Antoniewicz et al., 2007) was determined by an established gas chromatography (GC)/electron impact time-of-flight mass spectroscopy (EI-TOF-MS) protocol (Luedemann et al., 2008;Fernie et al., 2011;Osorio et al., 2011). Starch content was determined by enzymatic digestion and spectrophotometric assay of the resultant glucose (Hendriks et al., 2003). Fruit height and diameter were measured using a calibrated digital vernier. The fresh weight to dry weight ratio (FW/DW) was determined by weighing fresh tissues and reweighing them after drying to constant weight in a forced-air oven at 65°C.

Gas chromatography mass spectrometry of fruit content
Metabolite analysis of fruit samples was performed using GC/EI-TOF-MS (Lisec et al., 2006). Plant material was extracted using a method described elsewhere (Osorio et al., 2011). Data analysis was performed using CHROMATOF 1.0 (Leco, https://www.leco.com) and TAGFINDER 4.0 (Luedemann et al., 2008). Cross-referencing of mass spectra was performed with the Golm Metabolome database (Kopka et al., 2005). Documentation of metabolite profiling data acquisition is reported following recommended guidelines . Curve fitting (Colombi e et al., 2015) was used to generate curves for metabolite content and fluxes based on experimental measurements.

Flux balance analysis and flux variability analysis
Metabolic reactions associated with phospholipid biosynthesis (PC, PE and PA), amino acid catabolism, nucleotide biosynthesis, b oxidation, lycopene biosynthesis and phytol metabolism were added to a previously published mass and charge-balanced model of primary metabolism in plant cells (Shameer et al., 2018) to generate the model PlantCoreMetabolism_v1_2. A complete log of all model curations and associated literature is presented in Table S2. Parsimonious FBA (pFBA) and flux variability analysis (FVA) functions available in COBRAPY 0.13.4 (Ebrahim et al., 2013) were updated to perform weighted pFBA and weighted FVA, respectively. Aggregator reactions used to measure fruit metabolite content were given a zero weighting, whereas all other reactions were given a weight of 1. COBRAPY in JUPYTER NOTEBOOK (PYTHON 2.7) was used to run and document all scripts. BINDER was used to generate interactive versions of the JUPYTER NOTEBOOKS used in the study. All scripts, data, models, associated files and the link to launch BINDER have been made available in a Github repository (https://github.com/ljs1002/Shameer-et-al-Predicting-metabolismduring-growth-by-osmotic-cell-expansion).

Estimating maintenance costs
Non-growth associated maintenance (NGAM) costs in all models were represented by flux through an ATPase and three NADPH oxidase (cytosolic, mitochondrial and plastidic) pseudo-reactions, constrained to a 3:1 (ATP hydrolysis:NADPH oxidase) ratio based on the results of a previously published study (Cheung et al., 2013). Dividing cells were assumed to have a carbon conversion efficiency (CCE) of 70% based on measurements for a range of plant tissues and species (Chen and Shachar-Hill, 2012). Fluxes through the NGAM reactions in the dividing cell model were gradually increased until a CCE of approximately 70% was achieved. When the NGAM ATPase flux was constrained to 0.062 mmol ml -1 day -1 , the CCE of the system was observed to be 69.92%. This NGAM ATPase flux was then used to constrain maintenance in the dividing cell model. The same NGAM constraints were imposed on the expanding cell model.
Each phase of the multiphase GrOE-FBA model was 5.1 days long. The mean respiratory cost of NGAM based on data reported by Walker and Thornley (1977) was assumed to be 1.7 mmol CO 2 per fruit per day throughout fruit development. Note that these results indicated that the maintenance costs per unit mass of fruit decreased during fruit growth, leading to an approximately constant maintenance cost per fruit. This is equivalent to 8.5 mmol CO 2 per fruit per 5.1 days. Using a model of primary metabolism in pericarp cell constrained so as to prevent the net synthesis of any metabolite, an NGAM ATPase flux of 26.2 mmol CO 2 per fruit per 5.1 days was found to result in a respiration rate of 8.5 mmol CO 2 per fruit per 5.1 days. This NGAM ATPase flux was used to constrain NGAM in all phases of fruit development.

Modelling primary metabolism in a dividing heterotrophic cell suspension
The PlantCoreMetabolism_v1_2 model was used to generate representations of metabolism in dividing tomato cells. Previously published data from tomato heterotrophic cell suspension cultures report a maximum growth rate of 2 mg DW ml -1 day -1 (Rontein et al., 2002). Data on the biomass content (cell wall, sugars, organic acids and proteins) of cells demonstrating the maximum growth rate, from the same study, was used to generate a biomass equation for the dividing cell system as in conventional FBA modelling. The lipid:cellulose content ratio reported in Arabidopsis cell cultures (Williams et al., 2010) was used to estimate the lipid content of the cells and the biomass equation was updated accordingly. The biomass accumulation rate of the model was constrained to 2 mg DW ml -1 day -1 . For the sake of simplification, the cell wall and cell membrane were assumed to be composed of only cellulose and PA, respectively. Constraints based on Equation 2 were introduced to promote solute partitioning between the cytosol and the vacuole. With the lack of data available for heterotrophic tomato cells, the value for V v /V c required for Equation 2 was calculated from data published on tomato pericarp subcellular volume fractions . NGAM was represented in the model as described earlier. Glucose was set as the sole carbon source and FBA, with the minimization of the sum of fluxes as the objective, was used to predict the optimal flux distribution.
Modelling primary metabolism in a system of expanding 26-DPA pericarp cells pericarp cells. The rate of cell expansion was estimated to be at its highest in 26-DPA tomato pericarp cells (0.174 nl cell -1 day -1 ; estimated from data in Beauvoit et al., 2014). For the sake of consistency, the number of cells in the expanding cell system was kept equal to that in the dividing cell suspension when the growth rate attained the maximum rate of 2 mg DW ml -1 day -1 (1.45 9 10 6 cells ml -1 ). The composition of soluble metabolites in 26-DPA pericarp cells were estimated from third-degree polynomial curves fitted to published experimental data (Colombi e et al., 2015). Equations 1 and 2 were used to introduce osmotic constraints in the model and the values for V cell and C cell were estimated from published data (Almeida and Huber, 1999;Beauvoit et al., 2014). The difference in osmotic content between DPA 26 and 27 was used to set the demand for osmolytes. Data on subcellular volume fractions published by Beauvoit et al. were used to calculate V v /V c at 26 DPA. The difference between cellulose, lipid and protein contents between DPA 26 and 27, as calculated from Equations 3-5, was used to set the demand for these biomass elements. NGAM was represented in the model as described earlier. Uptake of sucrose and amino acids from the phloem was permitted and FBA was used to predict metabolic fluxes that minimized the sum of fluxes while maximizing the organic content.
Implementation of GrOE-FBA constraints in the multiphase developing fruit model Tomato fruit development from 13.1 to 43.7 DPA (mature green fruits observed at 42 DPA) is mainly the result of the expansion of existing cells (Gillaspy et al., 1993) and GrOE-FBA was used to model metabolism during these stages. The model extended to 54 DPA to allow the final metabolite composition of ripe fruit to be incorporated as a constraint. Assuming the tomato fruit is composed of uniformly packed pericarp cells, Equation 1 can be transformed to the following form for a whole fruit: where i is a metabolite/ion in the fruit and V fruit is the volume of the fruit.
Constraints based on Equations 2 and 6 were imposed on each phase of the multiphase developing fruit model. Previously published data on subcellular volume fractions  were used to calculate V v /V c during the different stages of fruit development and C cell was estimated based on osmolalities reported in the literature (Almeida and Huber, 1999). Values for V fruit for all 10 phases of fruit development modelled were determined experimentally (Table S1).
The pericarp cell volume for each phase of the multiphase fruit model was calculated from published data , and cellulose, phospholipid and protein demand fluxes were imposed on the system using Equations 3-5, scaling them from µmol per cell to mmol per fruit.
The carbon influx rate into developing tomato fruits has been reported to be inversely proportional to fruit carbon content and size (Walker and Ho, 1977a,b). A hyperbolic curve was fitted to capture the relationship between carbon influx rates and fruit carbon content more accurately. Fruit carbon content has been reported to be linearly related to fruit volume (Walker and Ho, 1997b). Hence, the initial carbon content for each developmental phase was calculated based on fruit volume at the respective phases. These values were then used to predict the upper bounds for the flux of nutrients from the phloem in the respective phases of the developing tomato model.

AUTHOR CONTRIBUTIONS
SS, RGR and LJS developed the GrOE-FBA framework; ARF, RGR and LJS designed the experiments and supervised the study; JGV collected and analysed all biological data; SS performed the modelling and analysed the results; SS, JGV, RGR and LJS wrote the article; all authors approved the final version of the article for publication.

CONFLICT OF INTEREST
There are no conflicts of interest.

DATA AVAILABILITY STATEMENT
All computational code required to reproduce this study is available from the following Github repository: https:// github.com/ljs1002/Shameer-et-al-Predicting-metabolismduring-growth-by-osmotic-cell-expansion. All other data generated and used in this study are available from the supporting information.

SUPPORTING INFORMATION
Additional Supporting Information may be found in the online version of this article. Figure S1. Expanding cells show higher ATP demand associated with solute accumulation when compared with dividing cells, irrespective of the non-growth associated maintenance (NGAM) ATP hydrolysis and NADPH oxidase flux ratio constraints. Figure S2. Comparison of predicted and measured metabolite contents in developing tomato fruit. Table S1. List of parameters used in the study. Table S2. Model curation log for the PlantCoreMetabolism_v1_2 model. Appendix S1. Derivations for Equations 1-5.
Appendix S2. Tomato fruit fresh weight, fresh weight/dry weight ratio, organic solute content, starch content, protein content, cell wall content and lipid content during different stages of fruit development. Appendix S3. Predicted metabolic fluxes in dividing cell, expanding cell and developing fruit models. Appendix S4. Systems Biology Markup Language (SBML) version of the fully constrained multiphase developing tomato model.