Development of a kinetic metabolic model: application to Catharanthus roseus hairy root

A kinetic metabolic model describing Catharanthus roseus hairy root growth and nutrition was developed. The metabolic network includes glycolysis, pentose-phosphate pathway, TCA cycle and the catabolic reactions leading to cell building blocks such as amino acids, organic acids, organic phosphates, lipids and structural hexoses. The central primary metabolic network was taken at pseudo-steady state and metabolic flux analysis technique allowed reducing from 31 metabolic fluxes to 20 independent pathways. Hairy root specific growth rate was described as a function of intracellular concentration in cell building blocks. Intracellular transport and accumulation kinetics for major nutrients were included. The model uses intracellular nutrients as well as energy shuttles to describe metabolic regulation. Model calibration was performed using experimental data obtained from batch and medium exchange liquid cultures of C. roseus hairy root using a minimal medium in Petri dish. The model is efficient in estimating the growth rate.

network was taken at pseudo-steady state and metabolic flux analysis technique allowed reducing from 31 metabolic fluxes to 20 independent pathways. Hairy root specific growth rate was described as a function of intracellular concentration in cell building blocks. Intracellular transport and accumulation kinetics for major nutrients were included. The model uses intracellular nutrients as well as energy shuttles to describe metabolic regulation. Model calibration was performed using experimental data obtained from batch and medium exchange liquid cultures of C. roseus hairy root using a minimal medium in Petri dish. The model is efficient in estimating the growth rate. reaction rate for the flux x (mmol/d or mmol/g DW/d) m GLC partial reaction rate for glucose (mmol/ g DW/d) m max(x) maximum reaction rate for flux x (mmol/ g DW/d) m max GLC maximum partial reaction rate for glucose (mmol/g DW/d)

Introduction
In vitro plant biotechnology offers a controlled environment and has been widely studied for phytopharmaceuticals and recombinant proteins production. However, the low productivity and the poor reproducibility of the cultures are still limiting the economical feasibility of such in vitro bioprocesses. Moreover, the lack of reproducibility of the cultures significantly complicates process validation and acceptance by the regulatory agencies, and thus the potential to rapidly put a product to market. Different approaches have succeeded in decreasing the technological risk associated with in vitro plant biotechnology. The introduction of elicitors such as chitin [1] and jasmonic acid [2] has shown to enhance significantly the production level for many plant species cultured as cell suspensions and hairy roots. The use of an extractive phase allowed the simplification of harvesting procedures [3,4]. More recently, cell engineering studies have shown its potential in improving cell catalytic capacity towards the production of secondary metabolites [5][6][7] and in in vitro molecular farming for recombinant human proteins [8]. However, high variability levels in cell growth and production of biomolecules of interest are still observed. The genetic flexibility of plant cells may partly explain these phenomena [9] and inadequate culture management may also be involved.
A plant cell has the ability to accumulate nutrients and metabolites which are involved in the regulation of its metabolic pathways. The key for improving plant cell culture reproducibility may thus rely on a better understanding of the links that are exerting between a plant cell physiological state and its potential towards growth and production of a biomolecule of interest. This understanding (following that of Bailey [10]) could take the form of a descriptive and predictive metabolic model. Such structured model may then be either useful to enhance our understanding of cell behaviour, in identifying possible regulatory roles [11], as well as being a tool for defining adequate controlled culture conditions. Metabolic modelling has been applied to plant cells for studying specific metabolic sub-networks such as photosynthesis [12,13], respiration [14], cellulose biosynthesis [15] and lipid biosynthesis [16]. These studies have clearly showed the importance of the energy shuttles on the control of the metabolic pathways. In addition, some nutrients are known to be involved in the regulation of the cell metabolism. Intracellular inorganic phosphate (Pi) plays a central role in the regulation of enzymes activity through phosphorylation/dephosphorylation processes, ATP/ ADP concentration ratio, starch synthesis and storage, and in the flux distribution between the glycolysis and the pentose-phosphate pathway (PPP). Intracellular nitrate and ammonium are known to affect amino acid (AA) production.
The large capacity of plant cells for nutrient and metabolite accumulation plays a crucial role in cell growth and biomolecule production, as observed for Daucus carota hairy root [17] and Eschscholtzia californica suspension cells [18]. Recently, we have developed a kinetic model based on intracellular nutrients such as inorganic phosphate, nitrate and sugars which showed to be efficient in simulating carrot hairy root growth for different culture media composition [17]. Therefore, the aim of this work was to include the description of metabolic pathways to the nutritional model in order to describe plant cells behaviour from the estimation of the cells physiological state, including nutritional and metabolic states. Catharanthus roseus was studied as a model biological system. Cell nutritional state in Pi, nitrogen (NO 3 À and NH 4 + ) and carbohydrates (sucrose, fructose, glucose and starch) is described. The hypothesis of a central primary metabolism at steady state has been proposed based on literature [19,20]. Using the metabolic flux analysis (MFA) approach, a model reduction [21] was applied on the central primary metabolism network and resulted in independent pathways. A second network includes transient fluxes such as for nutrient uptake and storage, energy shuttles management and root cells growth. Metabolic regulation of the fluxes from energy shuttles and nutrients is included. The hairy root specific growth rate is described as a function of the content in cell building blocks such as amino acids (including proteins), lipids (LIP), organic acids (ORA), organic phosphates (OP) (including nucleic acids) and structural hexoses (STH). Batch and medium exchange cultures of C. roseus hairy root were performed and the experimental data were used for model calibration.

Model general structure
The model has been first developed by Tikhomiroff [22]. The cell metabolic network (Fig. 1) is divided into two interlinked sub-networks as the stationary (SPMP) (Fig. 2) and the transient (TPMP) primary metabolic pathways (Fig. 1). The SPMP includes glycolysis, PPP, the TCA cycle and the catabolic reactions leading to the cell building blocks. The cell building blocks are amino acids and peptides which were taken as a unique pool of AA, ORA, OP, LIP and STH. The TPMP network is linked to the SPMP network and describes cells growth and nutrient transport between medium and intracellular volumes. Compartmentalization of nutrients and metabolites among the cytosol, the vacuole and other organelles is not included in the  Table 3. Kinetic description of the resulting fluxes is presented in Table 4 Fig. 2 Fluxes in the stationary primary metabolic pathways (SPMP). Flux numbers refer to the stoichiometric biochemical reactions of Table 1 model and a single cell population was considered to describe the hairy root cells pools. This simplification already showed to be efficient to model hairy root growth and nutrition [17]. The secondary metabolism is simplified to fluxes leading to two pools, one accounting for the global pool in secondary metabolites derived from tryptamin (TRY SM ), and one accounting for the global pool in secondary metabolites derived from secologanin (SEC SM ). The model is thus composed of interlinked metabolic networks that are at steady (SPMP) and transient (TPMP) states, and is described by a mass balance using the stoichiometric matrix and the hairy root specific growth rate: where S is the stoichiometric matrix, U is a vector containing reaction fluxes, l is the hairy root specific growth rate and M is a vector containing the concentration in cellular metabolites and nutrients. Root mass with time can then be estimated both kinetically and from a mass balance on all cell constituents.

Stationary primary metabolic pathways
The pseudo-steady-state assumption for the central metabolism was based on observations from Rontein et al. [19] and Stitt and Fernie [20] and proposed to simplify the model development. The original SPMP metabolic network has 31 fluxes ( Fig. 2; Table 1), which were reduced to 20 independent pathways ( Fig. 3; Table 2) using the method proposed by Simpson et al. [21] and Stephanopoulos et al. [23] and the following simplifications. Briefly, the minimal number of independent feasible metabolic pathways is determined with a group of metabolites assumed at steady state: G6P, F6P, R5P, G3P, E4P, CHO, PEP, PYR, A-CoA, OXO and OAA. The ''N'' matrix (11·31) is first constructed from the 31 biochemical reactions and the 11 metabolites at steady state. The determination of the independent pathways then requires the knowledge of the kernel matrix ''K'' (31·number of independent pathways), which describe each pathway as a linear combination of the 31 biochemical reactions of the SPMP. This matrix is the non-trivial solution of the equation:  [39] b [24] c [40] d [27] e [41] f [26] Fig. 3 Reduced independent biochemical pathways from the SPMP. Pathway numbers refer to the stoichiometric biochemical reactions of Table 2. Kinetic description of the resulting fluxes is presented in Table 4. The 20 independent biochemical pathways were obtained from the 31 biochemical reactions of the SPMP described in Fig. 2 and Table 1, using the reduction method of Simpson et al. [21] and Stephanopoulos et al. [23], as described in Sect. 3 Since this equation has an infinite number of solutions, both ''N'' and ''K'' matrices are decomposed each into two matrices: with '' N 1 '' and '' K 2 '' selected as square invertible matrices. Starting using identity matrix for ''K 2 '', ''K 1 '' can then be determined by modifying ''K 2 '' such that possible independent pathways can be identified: Using the complete matrix K, we can build the independent pathways presented in Fig. 3. In this work, the size of this matrix is 31·20. The complete method to obtain the pathways can be found in Ref. [22]. The nutrient concentrations described are glucose (GLC), fructose (FRU), sucrose (SUC), ammonium (NH4) and inorganic phosphate (Pi). The anaplerotic pathways are simplified to the transformation of phosphoenolpyruvate (PEP) into oxaloacetate (OAA). The metabolite pools resulting from the SPMP are starch (STA), STH (cell wall and membrane constituents), OP (nucleotides, phospholipids, nucleic acids), total AA including that of peptides, ORA, pools from tryptamin (TRY SM ) and from secologanin (SEC SM ). Starch biosynthesis (flux v(r14) in Table 2) and catabolism (v(r13)) are included. AA (v(r5)) are assumed to be synthesized from oxoglutarate, since it is the site for ammonium fixation and thus the initiation step in the AA biosynthesis [24], even though there are AA synthesized from other pathways. Tryptamin (TRY SM ) biosynthesis (v(r1)) has been explicitly described because this AA is a precursor, with secologanin (SEC SM ) (v(r2) and v(r10)), to the formation of the secondary metabolites. OP biosynthesis (v(r11)) has been simplified as illustrated in Fig. 3 and is accounting for the four nucleotides. The biosynthesis of STH results from multiple biochemical reactions [25] but it is simplified to the condensation of G6P with R5P (v(r12)). The biosynthesis (v(r7)) and catabolism (v(r15)) of LIP is linked to acetyl-coenzyme A [26] as for ORA biosynthesis (v(r16)) and catabolism (v(r17)). The kinetic expressions of the fluxes are described in the next section.

Transient primary metabolic pathways
The TPMP network includes the biochemical reactions that cannot satisfy the pseudo-steady-state hypothesis. These include nutrients transport and accumulation, root cells growth and energy shuttles dynamics. The general structure of the TPMP and its interactions with the SPMP are illustrated in Figs. 1 and 2. The Table 2 Independent fluxes of the SPMPs after pathways reduction using MFA approach a The biochemical reactions are denoted with an ''r'' to indicate that they are obtained from pathways reduction a Stoichiometric coefficients were taken from Ref. [24] stoichiometric mass balances are presented in Table 3 and the kinetic terms used in fluxes' regulation are described in Table 4. Extracellular nutrients such as sucrose (ESUC), glucose (EGLC), fructose (EFRU), ammonium (ENH4), nitrate (ENO3) and Pi (EPi) are included. Sucrose is hydrolysed into fructose and glucose extracellularly by apoplastic invertases (Table 3, v(28)), or intracellularly as described in the SPMP ( Table 2, v(r19)). It should be noted that sucrose biosynthesis ( Table 2, v(r20)) is not included into the model because there are no evidence that hairy roots have active chloroplasts. v(r19) then represents the net flux of sucrose hydrolysis. Alternative glycolytic pathways [27] are described and the carbohydrate sources include sucrose, glucose, fructose and starch. However, since all the biosynthesis pathways of the SPMP are defined with glucose as precursor, all available intracellular carbon sources have first to be converted into glucose. Fructose conversion (v(r3)), starch catabolism (v(r13)) and sucrose hydrolysis (v(r19)) processes are thus feeding the intracellular glucose pool. Plant cell adaptation mechanism to Pi deficiency is described with the fluxes for the degradation of pyrophosphate (v(26)) and OP (v(23)) into Pi. The plant cell storage capacity for sugars, NO 3 À , NH 4 + , Pi and cell building blocks is described. Respiration is described (v(30)), however, oxygen is assumed to be non-limiting as discussed below. Energy loss associated with maintenance and other reactions that are not included into the SPMP is contained in flux v (29).

Metabolic regulation
Metabolic regulation of two kinds is integrated into the model: that associated with the energetic status of the cells and that related to the nutritional state of the cells, variable denoted as ''h i '' below. Each flux (m) is regulated as follows, according to a multiplicative kinetic of each The order of the regulation kinetics can also be adjusted from the term a i to account for multi-steps mechanisms.
Since we established the independent pathways for the SPMP with arbitrary conditions, it is not possible to compare the structure of the resulting kinetics with literature. However, the identified independent pathways have been selected for their feasibility, and their combination is mathematically equivalent to the whole metabolic system (31 fluxes), given that the pseudosteady-state hypothesis is acceptable. Furthermore, a sigmoid function is used for nutrients acting like switches to avoid large discontinuities when solving the differential equations, which would have been the case using a simple on/off switch type. Continuous functions are also more representative of biological behaviour [28]. The sigmoid function is also used for imposing maximum accumulation levels (see Table 4): The parameter ''a'' defines the steepness of the function, ''M i '' is the concentration of the nutrient involved in the regulation of a flux and ''M i,t '' is the nutrient concentration threshold.
The role of energy shuttles in the regulatory mechanisms is described in Table 4. NADH, NADPH and NTP (the sum of ATP, GTP, CTP and UTP) as well as their reduced forms (NAD, NADP and NDP) are used. The sum in energy shuttles per g DW of roots in both mðr20Þ is fixed at 0 and no sucrose synthesis is considered a Reaction number refers to pathways of SPMP in Fig. 3, and of TPMP in Fig. 1 b Not shown fluxes are related to exchange between the cytoplasm and the vacuole and were not included in this model c Use of NTP, NDP, NADH, NAD, NADPH and NADP in flux kinetics was taken from Ref. [24] d [42] e [43] f [44] oxidized and reduced forms (NADH + NAD; NDP + NTP; NADPH + NADP) are taken constant with time. However, the energetic status of the cells [NADH/ (NADH + NAD); NTP/(NDP + NTP); NADPH/ (NADPH + NADP)] was not assumed constant. ATPase proton pumps that are linked to H + /nutrient co-transport and involved in the control of intracellular pH are thus included. ATP consumption for transmembrane transport of glucose (v(21)), fructose (v(22)), sucrose (v(24)), nitrate (v (27)) and Pi (v(31)) is thus included. A unique value for the affinity constant for NTP (K mNTP ) is used. The same strategy of using unique affinity constants for NAD, NADH, NADP and NADPH is also applied.
Root cells nutritional state in Pi, NH4 and in sugars are involved in metabolic flux regulation. Intracellular Pi is involved in enzyme activation/deactivation processes through dephosphorylation/phosphorylation, in the biosynthesis of NTP (v(30)) and of many other metabolites as described below. Intracellular NH4 is involved in AA biosynthesis (v(r5)) ( Fig. 3; Table 4). Michaelis-Menten kinetics with NH4 and Pi are at power two since a second-order mechanism is involved. NH4 has first to be transformed into NH2 radical, then into an AA which is finally integrated into a protein structure. Management of carbon source is crucial for plant cells since glucose, fructose, sucrose and starch can be used. However, since the model is based only on glucose, biochemical reactions converting other sugars into glucose are described ( Fig. 3; Table 4). Starch biosynthesis (v(r14)) is controlled by the total intracellular available sugar concentration, and the sum in intracellular GLC, FRU and SUC has to be above 0.2 mmol/g DW for the STA biosynthesis to occur. A maximal STA storage capacity of 1.5 mmol/g DW was observed experimentally (see Sect. 7). Starch is degraded (v(r13)) when the GLC concentration is below 0.2 mmol/g DW. PPP (v(r4)) is regulated by NADP/ (NADP+NAPH) ratio, which has to be above 0.5 for the pathway to be active. Pyruvate kinase regulation (v(r6)) is under the control of intracellular Pi concentration with a flux increase at 1 mmol Pi/g DW and below. G3P conversion into PEP (v(r8)) is regulated by Pi concentration with a flux increase at 1 mmol Pi/ g DW and below. TCA cycle (v(r9)) is regulated by the NAD/(NAD+NADH) ratio, and a ratio that is above 0.5 will induce a flux increase. GLC (v(21)) and FRU (v(22)) uptake require a STA concentration that is below 0.07 and 0.02 mmol/g DW, respectively. OP synthesis (v(r11)) is regulated by intracellular Pi concentration with a flux increase at above 0.1 mmol Pi/g DW. OP degradation (v(23)) is controlled by intracellular Pi concentration with a flux increase at below 0.045 mmol Pi/g DW. LIP (v(r15)) and ORA (v(r17)) degradation into A-CoA is controlled by NADH concentration with fluxes increase at below 1 mmol NADH/g DW. Phosphofructokinase regulation (v(r18)) is controlled by intracellular Pi with a flux increase below 1 mmol Pi/g DW. Respiration rate (v (30)) is controlled by NDP/(NTP+NDP) ratio with a flux increase below a ratio of 0.5. Finally, regulation of root cells growth rate (v(32) is controlled by the intracellular concentration in the cell building blocks such as OP, LIP, amino acids and peptides (AA), ORA and STH. Monod model was used for LIP, AA and ORA. A hybrid Moser-Monod model was used for OP and STH with each kinetic term at power 4 and 1.25, respectively (Table 4). Since OP and STH are crucial to cell growth, the steepness of the affinity for both cell building blocks was increased.

Culture conditions
Liquid cultures of hairy roots were performed in Petri dish as described in Jolicoeur et al. [17]. The major problem in culturing hairy roots is the difficulty to obtain a representative sample of the root network since the roots develop a highly dense interlinked bed. Petri dish culture allowed distributing the roots inoculum in a series of dishes with each dish taken as a single sample. Root growth was also limited by the use of a low salt medium, thus preventing the occurrence of dense root network. Previous results in Petri dish [17] and in bioreactor [3] suggest that there was no oxygen limitation in the cultures of this work.
Hairy roots of C. roseus L. G. Don were established as described by Bhadra et al. [29], with Agrobacterium rhizogenes strain A 4 . Hairy roots were transferred every month in Petri dishes in 25 mL of minimum medium [30] supplemented with 3% (w/v) sucrose and with a threefold KH 2 PO 4 (0.352 mM). In the batch culture experiments, approximately 0.125 g fresh weight (FW) of hairy roots was inoculated in each Petri dish containing 25 mL of minimum medium. In the medium exchange cultures, approximately 0.125 g FW of hairy roots was inoculated in Petri dish containing 25 mL of minimum medium. The medium was renewed at 3-or 2-day intervals to avoid depletion of Pi in the medium. Whole Petri dishes were harvested in triplicates (n=3) periodically and taken as distinct samples. Liquid from each dish was filtered at 0.45 lm (Millipore, Billerica, Massachusetts) and stored at À20°C for further analysis. Roots were filtered under vacuum on a glass fiber filter (47 mm diameter Glass Microfiber filters GF/D, Whatman, #1823 047) and rinsed three times with 20 mL of de-ionized water. The filtered roots were removed from the filter and weighed for FW in a disposable aluminum dish (Fisher Scientific, # 08-732) on a precision balance (Sartorius). Fresh roots were immediately frozen into liquid N 2 and stored in liquid N 2 for further analysis. Root samples were then freeze dried (Duratop and Duradry, FTS Systems Inc., Stone Ridge, NY, USA), weighted for DW measurement and grinded (mortar and pestle) for further analysis. All further analyses were performed using freeze-dried roots.

Amino acids analysis
Approximately 2 mg of freeze-dried roots were extracted in 1 mL of 2% w/v 5-sulfosalicylic acid. The samples were sonicated for 15 min and then centrifuged at 12,000g for 5 min. The supernatant was analysed for AA by HPLC using a modified method from Gombert et al. [31] as described in Benslimane et al. [32]. The column is a high-efficiency Nova-Pak TM (C18, 4 lm). Precolumn derivatization of the AA was performed using AccQ.Fluor (6-aminoquinolyl-N-hydroxysuccinimidyl carbamate, or AQC), which is an N-hydroxysuccinimideactivated heterocyclic carbamate, and the derivatized AA were quantified via a fluorescence detector. a-Aminobutyric acid was used as an internal standard. The standard error of measurement using this method was routinely below 5%. AA concentration was then calculated since the peak for each individual AA accounted for the combination of free and that from peptides and proteins.   [17] Canada, Guelph, Ontario, Canada), coupled with a Prevail Carbohydrate ES All-Guard 4.6 mm·7.5 mm guard column (Alltech Canada, Guelph, Ontario, Canada), was used at a column temperature of 35°C. The injection volume was 20 lL. The mobile phase consisted of acetonitrile and water 75:25 (v/v) at 1.0 mL/min. The cell pellet obtained from the soluble sugars extraction was re-suspended in 1 mL de-ionized water and sterilized at 121°C for 15 min together with a 1 mL 6 g/L starch solution in de-ionized water. Samples were allowed to reach room temperature and their volume was readjusted to 1 mL with de-ionized water if necessary. Starch calibration standards from 0 to 6 g/L were prepared by serial dilutions of the starch solution in de-ionized water. Calibration standards and samples were diluted 1:1 with an amyloglucosidase solution (Sigma, St. Louis, Missouri, USA, Cat. #S9144) and incubated for 15 min in an ultrasound bath at 60°C. Samples and standards were centrifuged for 10 min at 16,000g. A 10 lL aliquot of the supernatant was directly transferred into a spectrometric cuvette with 500 lL of ''Glucose Infinity'' reagent (Sigma, Cat. #17-25). After incubation (15 min) at room temperature, 500 lL of 100 mM KH 2 PO 4 buffer pH 7.5 was added followed by spectrophotometric reading at 340 nm.  , NH 4 + , K + , Na + , Ca 2+ ) were analysed using a Dionex HPLC system (Dionex Canada Ltd., Oakville, Canada) equipped with an isocratic pump, an automated sampler AS-3500 and a pulsed electrochemical detector in the conductivity mode, controlled by the Dionex A1-450 software for cations and the Dionex Peaknet software for anions. Anions were separated using a 4·250 mm IONPAC AS14A-SC analytical column, an IONPACAG14A-SC guard column and a ASRS-1 anion self-regeneration suppressor to improve the signal-to-noise ratio. The mobile phase consisted of an aqueous buffer of 2 mM Na 2 CO 3 /1 mM NaHCO 3 solution flowing at a rate of 1.0 mL/min. Cations were separated using a 4·250 IONPAC CS-12 analytical column, a IONPAC CG-12 guard column and a CSRS-1 cation self-regenerating suppressor. The mobile phase was an aqueous 20 mM methanesulphonic acid solution flowing at a rate of 0.9 mL/min.

Model simulations
The model simulations were done using the Matlab software (The MathWorks Inc., Natick, MA, USA). The differential equations system was integrated through the Ordinary Differential Equation solver ode15s.m. The model error minimization was performed by means of manual and algorithm-based methods. The former was done to find good initial estimates based on literature (when available). The latter was then used to reduce the global error between the estimates and the experimental data. The algorithm used was the lsqcurvefit.m subroutine (Optimisation toolbox, Matlab) based on the Levenberg-Marquardt algorithm [33].

Results and discussion
Model calibration and determination of the kinetic parameters Model calibration was performed using experimental data from C. roseus hairy root batch and medium exchange cultures. Each experiment (one batch and one medium exchange) was performed in triplicate. A total of 14 samples were taken during the batch experiment and 17 during the medium exchange. Parameter values were obtained from experimental data (maximum uptake rates, maximum growth rate, sucrose maximum hydrolysis rate and some of the maximal accumulation levels) and from literature, taking values from other plant species when unavailable for C. roseus (Table. 5, 6, 7). Adjustment of the unknown parameters as well the values taken from literature was performed manually in parallel using a non-linear least-square algorithm (lsqcurvefit) from Matlab software. The least-squares criterion was applied on all the data points (triplicates taken as mean in (M i ) exp ). The sum of squares were weighted as described below.
where e i ¼ P Error terms were weighted using the maximum (max(M i ) exp ) and the minimum (min(M i ) exp ) experimental values (exp) from either batch or medium exchange cultures. This method was preferred to a method using experimental error because the experimental error was very high for some data points and varied with experimental data (Table 8). This method has also showed previously to perform adequately to calibrate a nutritional model on hairy root data [17]. All the experimental data points from batch and medium exchange culture were used for the non-linear leastsquare regression. The total number of experimental data was then of 411 with 13 different measurements for each of the 31 samples plus eight measurements for intracellular AA. Many set of parameters were tested as initial guesses to improve the fitting. A combination of parameter values (Table. 5, 6) showed to minimize the global error for combined data from batch and from medium exchange cultures, using a convergence criterion of 1·10 À6 on the variation of the global error. From the calculation of R 2 values ( Table 8) for each of the 28 state variables that were fitted (14 state variables in two experiments), it seems that the best fit is obtained for cells growth and nutrient transport and storage (NO3, Pi, NH4). However, the sensitivity analysis also suggests that some parameters can still be optimized (Fig. 4). A critical problem with such a descriptive model is a high number of parameters. The model has 35 maximum reaction rates and 26 affinity constants and the leastsquare regression was only performed on only 14 of the 31 state variables of the model. It is then possible that the identified solution corresponds to a local minimum. More measurements on metabolite concentration are thus necessary to improve the performance of the model  and for being able to perform an adequate identifiability analysis on the parameters. In this context, the proposed model is probably over-parameterized. However, the goal of this work was not to establish a minimal model but a descriptive model that can be used as a tool useful at improving our understanding of plant cell culture.

Simulation of the cells physiological state
Hairy root growth The model allows estimating root growth from flux m(32) as well as from the total estimated mass of the root cells constituents per dish (R M i on a mass basis), taking average molecular weights for the cell building blocks (Fig. 5). However, the total mass calculated at initial condition (t=0) was 25% lower than that of the inoculum with 0.0148 g DW per dish as compared to 0.0197 g DW per dish, respectively. The difference of À0.0049 g DW per dish can be explained from the molecular species that are present in the cells but which are not included in the model. Therefore, the estimation of root mass with time obtained from m(32) was compared with that calculated from total cell mass, but adding a constant correction factor of +25%. Model simulations then describe hairy root growth for both batch and medium exchange cultures. Both root mass estimates from m(32) and RM i are superimposed (Fig. 5). However, the sum of the molecular species is overestimating measured biomass at the end of batch culture because of cell accumulation in diverse compounds while growth slowed down. Surprisingly, both batch and medium exchange cultures behaved similarly with a maximum specific growth rate of 0.035 d À1 and a growth cessation at around 40 d. Periodic medium renewal has avoided nutrient depletion in macronutrients and sugars (Fig. 6). However, the plateau observed for the medium exchange culture may have resulted from a limitation in micronutrients, which were not measured. Therefore, only experimental data obtained before the occurrence of the plateau for the medium exchange culture were used for model calibration and other analyses. The strategy of describing cell growth as a function of intracellular content in building blocks thus showed to perform adequately. This result may support the global model structure as well as the hypothesis of pseudo-steady state for the SPMP. The central primary metabolic network, described here as SPMP, was observed at pseudo-steady state by Rontein et al. [19] for most of in vitro culture duration for tomato suspension cells [20]. However, this hypothesis has to be further investigated.

The cells building blocks
The measured total AA concentration decreased with time for the batch culture and reached a plateau for the medium exchange culture (Fig. 5a). Simulated AA concentrations seem to follow the trend for experimental data but were overestimated for both batch and medium exchange cultures. The deviations cannot be explained, however, by exuded AA and proteins since no AA were detected in the used culture media. Values of 0.63±0.10 to 0.82±0.10 mmol AA/g DW were measured for the batch culture which corresponds to 0.088±0.014 and 0.12±0.014 g AA/g DW, respectively (taking an average amino acid MW of 140 g/mol). Higher values were measured in the medium exchange culture with 0.63±0.15 to 1.1±0.14 mmol AA/g DW, which represents 0.12±0.021 and 0.16±0.020 g AA/g DW, respectively. Free AA content of 0.21 mmol AA/g DW (estimated at 0.029 g AA/g DW) is reported in literature for tobacco suspension cells using B5 medium [34] and 0.15 mmol AA/g DW (estimated at 0.021 g AA/g DW) in potato tuber [35,36]. In tomato cell culture, the total protein content was of 0.1-0.4 g proteins/g DW [19]. The total AA mass per root dry weight that have been measured, which represents the sum of free AA and that of peptides and proteins, is thus within the range reported in literature. The use of a minimal medium may explain a low value for total AA. The other cell building blocks were not measured but the model estimations were close to literature data. Total ORA concentration was estimated at 0.15 mmol/g DW as compared to 0.42 mmol/g DW reported by Farre´et al. [35]. Total STH were estimated at 0.14 mmol/g DW as compared to a range of 0.3-1.6 mmol/g DW (on the basis of glucose units). Maximal simulated rates of biosynthesis were compared with that from literature. The rate of AA biosynthesis (v(r5)) was estimated at 0.12 mmol/g DW/ d, which was within the range of 0.094-0.36 mmol/ g DW/d reported by Rontein et al. [19] for tomato suspension cells. STH biosynthesis (v(r12)) was slow with 0.009 mmol/g DW/d as compared to 0.36 mmol/ g DW/d [19]. Maximal net rate of starch (STA) biosynthesis (v(r14)) was 0.20 mmol/g DW/d as compared to 0.15 mmol/g DW/d [19]. Differences in metabolic fluxes may be explained from the differences in the specific growth rate. Rontein et al. [19] have reported a specific growth rate of 0.4 d À1 for tomato suspension cells in B5 medium as compared to 0.035 d À1 in this study for C. roseus using a minimal medium. The same minimal medium has been reported previously to induce lower growth rate as compared to standard media for carrot hairy root [17].

The cells nutritional state
Extracellular sucrose (ESUC) concentration profile was simulated by the model for both the batch and the medium exchange cultures (Fig. 6), thus suggesting that the sucrose hydrolysis modelling strategy was adequate.
In the case of glucose (EGLC) and fructose (EFRU), simulations for the batch culture followed the trend of experimental data but before growth cessation. Then, estimated glucose and fructose levels stayed high as compared to experimental data. This may suggest that the model is underestimating glucose and fructose  uptake at reduced growth. However, simulated concentrations in intracellular sucrose (SUC) and glucose (GLC) were following similar trends than experimental data (Fig. 7). However, model simulated constant starch (STA) concentration and overestimated fructose (FRU) concentration. In the case of the medium exchange culture the model simulated adequately experimental data on extracellular and intracellular sucrose and intracellular glucose. Extracellular glucose and fructose concentrations were underestimated and starch and fructose concentrations were overestimated. A maximal sucrose rate of hydrolysis (v(28)) of 0.72 mmol/g DW/d was estimated as compared to 1.7 mmol/g DW/d for D. carota hairy root using same minimal medium [17]. Contents of 0.25 mmol SUC/ g DW (batch) and 0.7 mmol SUC/g DW (medium exchange) were obtained as compared to 0.05 mmol SUC/g DW [19], 0.26 mmol SUC/g DW [36] and 0.4 mmol SUC/g DW [34]. For free glucose, cell contents of 0.3-0.05 mmol GLC/g DW were obtained as compared to 0.6 mmol GLC/g DW [19], 0.24 mmol GLC/g DW [36] and 0.18 mmol GLC/g DW [34]. For fructose, cell contents of 0.1-0.25 mmol FRU/g DW were observed as compared to 0.009 mmol FRU/ g DW [35], 0.21 mmol GLC/g DW [34] and 0.5 mmol GLC/g DW [19]. Starch accumulation reached 1.25 mmol STA (based on glucose)/g DW as compared to 0.055 [34] and 0.1-0.6 mmol STA (based on glucose)/g DW [19]. Model simulations of extracellular Pi and NO 3 À (ENO3) concentrations followed experimental data for the complete duration of the batch and the medium exchange cultures (Fig. 6). Hairy root growth (simulated and experimental) has ceased concurrently to ENO3 depletion suggesting this ion to be limiting in the batch culture. Intracellular Pi and NH4 were simulated (Fig. 7). It is interesting to note that the final intracellular Pi concentration in the batch culture seemed to reach a plateau at 0.05 mmol/g DW which was simulated. This plateau may be close to the Pi level that is essential to maintain the endogeneous metabolism [37]. A maximum Pi accumulation plateau was also observed for the medium exchange culture around 0.2 mmol/ g DW. This value is close to the maximum value of 0.23 mmol/g DW measured for carrot hairy root [17, and references therein]. Maximum intracellular nitrate of 0.6 mmol NO3/g DW was measured for both cultures. Intracellular nitrate was kept constant for the medium exchange and decreased to its initial value for the batch. Intracellular ammonium reached a value of 0.12 mmol NH4/g DW for the medium exchange culture and decreased to 0.09 mmol NH4/g DW after 30 d. For the batch culture, intracellular ammonium reached a maximum value (0.08 mmol NH4/g DW) at day 5 and then decreased closely to its value at inoculation. The model overestimated intracellular NO3 concentration and simulated experimental intracellular concentrations in NH4 for both the batch and medium exchange cultures. A more precise description of the secondary metabolism and of the AA biosynthesis may be required to improve model simulation of nitrogenous compounds.

The cells energetic state
Model has simulated energy rich cells with high levels in NTP and NADH at exponential growth (Fig. 8). However, estimated cell content in NADPH showed an opposite trend with lower values at exponential growth. This result is interesting since NADPH is produced through the pentose-phosphate pathway, which is mostly active at secondary metabolism. The simulated cell energetic level decreased after 20 d at the initial value (inoculation) in batch culture and was constant for the medium exchange culture. Since the medium was renewed every 3 days for the medium exchange experiment, a pseudo-steady state was expected. NTP/ (NTP+NDP) and NADH/(NAD+NADH) ratios of $75% have been reported for tobacco suspension cells [34]. Farre´et al. [35] measured ratios of $73% for NTP/ (NTP+NDP) in potato tubers. These authors have measured an NTP level of 0.00013 mmol/g DW as compared to 0.00143 and 0.00074 mmol/g DW that was obtained from simulations at exponential growth and stationary growth, respectively. We have recently measured by in vivo NMR ATP levels of 0.009 mmol/g DW at day 1 and 5, and of 0.0065 mmol/g DW at day 10 for suspension cells of E. californica cultured using B5 medium [38]. Energetic levels simulated by the model were then similar to that found in literature. However, variation of energy shuttles concentration has to be further investigated.

Conclusion
A metabolic model capable of describing hairy root growth from the estimation of the cells physiological state was developed. The model includes the central metabolism, the primary metabolic pathways (SPMP) assumed at steady state and a network for the TPMP. At this point, the model is simplified and intracellular compartmentalization processes into the different cell compartments and organelles are not included. Nevertheless, the metabolic model showed to perform efficiently in simulating hairy root growth and nutrition. The use of intracellular concentrations in nutrients and co-substrates as well as the cells energetic state seems an efficient strategy in describing regulation of the metabolic fluxes. However, more experimental is required for improving the model structure as well as parameter values which may be regulated with cell physiological state. Finally, the model will be applied to other plant species as well as cell suspensions, and it will be studied as a tool to describe transient processes such as metabolic regulation.