Mathematical Modelling as a Tool for Optimized PHA Production

Poly(hydroxyalkanoates) (PHAs) are biodegradable intracellular polyesters synthesized by various eubacterial genera and some archaea1–5; further, the biosynthesis of PHAs in genetically modified yeasts6 and in modified plants was also reported7. PHAs occur in living systems as monopolymers (homopolymers) or copolymers (heteropolymers) with a variety of attributed pivotal functions in diverse ecosystems8, primarily serving as intracellular energy and carbon reservoirs9. They can be produced from renewable resources of the first generation (grains, sugar beet, sugarcane, starch, edible oils) and the second generation (molasses, lignocellulosic wastes, whey, biodiesel and biodiesel waste -i.e., glycerol and waste lipids)4,10,11. In the near future, it is expected that they will substitute some conventional plastics stemming from petrochemical sources, such as poly(ethylene) (PE) and poly(propylene) (PP)5. According to Gaden12, the rate patMathematical Modelling as a Tool for Optimized PHA Production


Introduction
Poly(hydroxyalkanoates) (PHAs) are biodegradable intracellular polyesters synthesized by various eubacterial genera and some archaea [1][2][3][4][5] ; further, the biosynthesis of PHAs in genetically modified yeasts 6 and in modified plants was also reported 7 .PHAs occur in living systems as monopolymers (homopolymers) or copolymers (heteropolymers) with a variety of attributed pivotal functions in di-verse ecosystems 8 , primarily serving as intracellular energy and carbon reservoirs 9 .They can be produced from renewable resources of the first generation (grains, sugar beet, sugarcane, starch, edible oils) and the second generation (molasses, lignocellulosic wastes, whey, biodiesel and biodiesel waste -i.e., glycerol and waste lipids) 4,10,11 .In the near future, it is expected that they will substitute some conventional plastics stemming from petrochemical sources, such as poly(ethylene) (PE) and poly(propylene) (PP) 5 .According to Gaden 12 , the rate pat-terns of product formation enclose three basic types in various bioprocesses: (1) growth-associated products arising directly from the energy metabolism of carbohydrates; (2) indirect products of carbohydrate metabolism; (3) products apparently unrelated to carbohydrate oxidation.
The aforementioned was studied for excretable products produced from carbohydrates.This was subsequently confirmed as valid for other types of C-sources (substrates related to the energy production and cell components synthesis).Although PHAs are intracellular (insoluble in cytosol) storage polymers, Gaden ' s types of kinetic behaviours can, under certain assumptions (circumstances) be applied to these compounds.
Regarding the kinetics of microbial growth and PHA production separately, three different types of microbial producers of PHAs can be distinguished, as follows: a) Strains displaying strict separation between the biomass growth phase and the PHA production phase normally provoked by N or P limitation (prototype organisms: Pseudomonas sp.2F, Methylomonas extorquens) 13 .b) Strains that accumulate certain quantities of PHA already under balanced nutritional conditions; also here, in the non-growth phase (usually provoked by N or P limitation), amplified PHA accumulation is observed (prototype organism: Cupriavidus necator).c) Strains displaying high PHA formation rates during growth phase without limitation of an essential growth component 2, [14][15][16] (prototype organisms: Azohydromonas lata DSM 1122, or Pseudomonas putida GPo1 ATTC 29347).
Based on the different kinetic properties of the production strains, different cultivation techniques have been applied in the past for their cultivation: batch, fed batch, repeated batch, repeated fed batch, or one-, two-, and multi-stage continuous processes 17 ; both pure and mixed microbial cultures (MMCs) where used.For example, microorganisms belonging to the mentioned groups (a) and (b) are well suited for two-stage continuous cultivation that requires extensive biomass growth in the first stage, followed by non-growth-associated PHA synthesis in the second stage [18][19][20][21][22][23][24] .As a novel approach, Atlić et al. 25 have recently tested a five-stage continuous process in a bioreactor cascade, specially designed to achieve different process conditions (concentrations of substrates and co-substrates, temperature, pH-value, dissolved oxygen [DO]) in each step of the cascade.From the point of view of such application, this cascade set-up was designed to act as a tool for the designing of novel biopolymers, e.g., blocky structured polymers with alternating "soft" and "hard" segments, or polymers displaying controlled molecular mass and polydispersity of molecular mass distribution.
Regarding the type of polymers, substrates, strains (wild types or genetically modified organisms), cultivation procedures and techniques, as well as the downstream processing systems for product recovery, numerous experiments were performed in the past in order to optimize PHA production both quantitatively and qualitatively.
Mathematical models are useful tools for optimizing and controlling microbial product formation and microbial metabolism, encompassing the modeling of cultivation techniques, the design of single cell metabolic models, or modelling of whole cell populations [26][27][28][29] .The segregated nature of biological systems and the complexity of cell reactions are cumbersome for mathematical treatment of processes in the bioengineering field.In addition, linking experimental data with mathematical modeling can reveal new aspects of microbial physiology, providing reasonable interpretations of results from experimental work.This way, the improvement of knowledge as well as the designing of new, more call-oriented experiments, can be achieved 30,31 .
Model theory classifies models into several classes and its antipodes: verbal/non-verbal, descriptive/explanatory, black-/grey-/white-box, unstructured/structured, non-mathematical /mathematical, deterministic/stochastic, discrete/continuous, and distributed (non-segregated)/segregated.The derived combinations of the mentioned model classes can be further grouped into logistic, kinetic, dynamic and cybernetic models.Usually, a specific model developed from a biological situation observed in living systems can belong to more than one class.For example, dynamic models can at the same time be structured, mathematical, continuous, and segregated.Kinetic models include two basic subcategories: They can be unstructured or structured.In literature, unstructured kinetic models are known as formal kinetic (or mechanistic) models; structured kinetic models can be divided into subclasses: morphological, compartmental, metabolic, and chemically related.Due to the plurality of microbial producers, the vast number of cultivation strategies, substrates, and products, and the existing broad spectrum of available models, a great number of interconnections (combinations) is conceivable.Therefore, very different model types have been applied to describe miscellaneous biotechnological situations regarding cultivation techniques, mass transfer, microbial growth kinetics, and metabolic network reaction kinetics.In the light of the aforementioned, it is practically impossible to discover a conceivable model type in the literature that has not yet been applied to describe PHA biosynthesis.Such highly extensive, comprehensive, and sophisticated material is difficult to systematize in order to achieve the desired and urgently needed overview for the relevant scientific community, and reach deeper insight into the essence of matter.Therefore, the aim of this work was to provide the respected readers with a comprehensive and critical synopsis of mathematical models that were applied for the microbial synthesis of PHA, together with the related cultivation strategies.Thereby, special attention is dedicated to the usefulness of applied models for improvement of processes, product yields, downstream processing, strains, and knowledge of the metabolic network.

Modeling approaches
Formal kinetic (unstructured and low-structured) models In the early phase of mathematical modeling of PHA biosynthesis, mainly formal kinetic unstructured, so-called "black box type", models were introduced [32][33][34][35][36] .These early works focused on the kinetic relationships between substrate (S), product (P) and biomass (X); hence, at the beginning, the optimization of the bioprocess itself was not the main purpose of modelling.At that time, microorganisms were cultivated by well-established standard microbiological techniques: Microbial broth usually contained the C-source as the growth-limiting substrate, whereas all other constituents were present in excess.For this type of cultivation, the standard differential mass balance equations for biomass (X), limiting substrate (S) and product (P; here: PHA) were established by applying the Monod equation for linking the substrate concentration to the specific growth rate.An application of appropriate yield and/or stoichiometric coefficients (biomass to substrate and product to substrate) was the common tool in these works.Such approaches were especially applicable for the growth-associated "(c)"type of PHA producers as described above.The reader can find an overview of these early works in the literature [37][38] .Soon after, these early mathematical models for PHA synthesis were corrected thanks to the works of Sonnleitner et al. 32 and Heinzle and Lafferty 34 , who developed the "two main compartments" strategy, which is applicable and still valid today.As the PHA mass fraction can account for the predominant part with mass fractions exceeding 0.9 of the whole biomass (X), the latter was structured by the mentioned authors into two components: (1) The catalytically active fraction consisting of membrane-and cytosolic proteins (structural proteins and catalytic enzymes), glycolipids, phospholipids, glycoproteins and nucleic acids (jointly known in literature as "residual biomass", i.e., the biologically active part Xr).
In addition, for some organisms it was observed that both the nitrogen (N) and phosphorus (P) sources act as limiting substrates, mainly triggering growth and product kinetics.So, N-and/or P-source limitation in the late phase of cultivation was postulated as the inductive triggering factor that provokes (or enhances) PHA synthesis; this is especially valid for the abovementioned producer types "(a)" and "(b)".By assuming that the rate of product formation is related both to the rate of growth and to the concentration of the biologically active part of biomass (Xr), a reasonably good correlation was obtained between experimental data and the model values.For this purpose, the Luedeking-Piret equation 39 was applied to describe product formation.Introducing growth-limiting substrates other than C-source (i.e., N-and/or P-sources) had certain effect on the mathematical expressions for specific growth and specific PHA production rate.In such cases, the specific growth rate of biomass (µ) was expressed using double-substrate Monod relations according to Megee et al. 40 .Regarding broth composition, this was extended to the multi-substrate Monod relation for more complicated mediums (containing complex substrates or mixture of different C-or N-sources).In some cases, the logistic equation was applied for the mathematical description of cell growth 41 .Unfortunately, these improvements were not the only changes needed in the kinetic modelling to achieve a valid model.
The introduction of new microbial species and strains constructed by genetic manipulation with specific kinetic performances in PHA biosynthesis, as well as the application of abundant inexpensive waste substrates (and by-products) in order to reduce production costs, required more advanced approaches in formal kinetic and low-structured modelling.Substrate inhibition was the second strong grip that was introduced in kinetic equations in order to adapt mathematical models to such biological/technological situations.Several different substrate inhibition terms reported by Hill 42 (sigmoidal), Aiba et al. 43 (exponential decreasing), as well as those reported by Yanno et al., Webb, Andrews and Noack (all reviewed by Moser 44 ) were applied.Two general application fields of substrate inhibition terms could be identified among published models in several papers: -The substrate inhibition of growth by C-sources caused by the impact of organic com-pounds (i.e., glycerol), pH-value effect (i.e., organic acids), osmotic pressure effect (sugars), and katabolic repression effect (typical in multi-substrate cultivations).
-The triggering of PHA production in the late phase of cultivation by limitation on N-and/or P-sources; high concentrations of the mentioned sources acted as an inhibiting factor on the specific PHA production rate, whereas a balanced supply of N and P resulted in predominant biomass formation at the expense of PHA formation.
The third important model adaptation was the introduction of the specific carbon-substrate consumption rate for generation of maintenance energy.For such purposes, the Guthke relation defined by simple Monod equation 44 was usually applied.Thus, maintenance energy was related to the concentration of C-sources.Furthermore, in experiments under aerobic conditions, it was recognized that the dissolved oxygen level (DO) influences the specific growth rate, but also acts on the intracellular NADH, NADPH and FADH 2 pools that are essential for the reduction of PHA precursors (e.g., Acetoacetyl-CoA), thus also affecting the PHA production rate.Hence, some authors 45 accomplished furt her improvement of formal kinetic and low-structured models for PHA synthesis by introducing DO in their models.
Finally yet importantly, product inhibition terms were introduced in mathematical modelling of PHA biosynthesis.As PHA granules can account for a significant part of whole cell mass, they may occupy a large share of the intracellular volume.In addition, some strains were found to stop PHA biosynthesis when the average mass fraction of PHA exceeded a certain value.Bearing this in mind, there were possible steric hindrances in PHA biosynthesis; therefore, some authors 35 implemented product inhibition terms in equations for the specific growth rate and/or in those for the specific PHA production rate.Usually used for such purposes were logistic equations, such as the equation according to Aiba et al. 43 , and the Luong equation 46,47 .
Simultaneously with the progress in kinetic modeling, the application of mathematical models in reactor designs, process optimization, substrate feeding strategies, copolymers production, one/two stage continuous cultivations and mixed cultures had increased.The reader can find all the situations described above, along with appropriate comments attached, in the review paper presented by Patnaik 48 who divided the described models into three groups: -Mechanistic models with pure culture, -Mechanistic models with mixed cultures, and -Cybernetic models After the year 2005 (the year when the mentioned paper was published), a number of more complex models was established.For example, Koller et al. 49 published a formal kinetic model for batch and fed-batch processes of copolymer production with Haloferax mediterranei, as well as a low-structured model for poly(3-hydroxybutyrate) (PHB) homopolymer synthesis from hydrolyzed whey by Pseudomonas hydrogenovora (today: Burkholderia fungorum).In the first case, the following assumptions were applied: i) Residual biomass (non-PHA biomass) is synthesized from all three main carbon sources i.e., glucose, galactose, γ-butyrolactone and from yeast extract as an obligate complex nitrogen source (independent growth on each substrate according to Monod relation was introduced).
ii) Direct influence of one C-source on the consumption rate of the two other C-sources does not exist (simultaneous independent consumption of C-sources).Note: this is valid for the special case described in this publication; in other cases, preferences of the organism for one C-source that is consumed prior to other, less easily convertible C-sources have been observed.
iv) PHA synthesis is inhibited by increasing intracellular mass fraction of PHA (known as steric disturbing effect).The equation according to Luong 46,47 was applied for this purpose.
Except for a few basic metabolic routes, metabolic pathways for the applied microorganism were not known in detail.This is especially valid for the galactose degradation pathway in archaea like the investigated organism H. mediterranei.Insufficient metabolic knowledge was the reason why authors had chosen formal kinetic modelling (instead of low-structured or high-structured metabolic flux models).In the second case (batch cultivation of Pseudomonas hydrogenovora), a low-structured metabolic mathematical model for fed-batch cultivation was established.In low-structured models, several metabolic reactions of particular metabolic pathways are lumped in one kinetically defined reaction (usually with kinetic properties of the slowest reaction or those that are the most sensitive to the regulating factors).In such models for PHA synthesis, except native biomass (residual biomass Xr + PHA), metabolic pools of intracellular metabolites and biocatalyst pools were structured.Here, the production of PHB from glucose and galactose as C-sources using NH 4 + and casein hydrolyzate as N-sources was studied.P. hydrogenovora (today: Burkholderia funghorum) is a typical strain having 3HB synthesis provoked by nitrogen limitation of growth; the total uncoupling of growth phase and 3HB synthesis phase (''non-growth-associated PHA production'') was evidenced.The following assumptions were necessary to establish a useful low-structured mathematical model: -Ac-CoA (originated from the central metabolic pathways -meaning from the metabolism of both sugars) is chosen as the ubiquitous precursor for PHB (because all sugar metabolizing (degrading) pathways of microbial metabolism i.e., Entner-Doudoroff, pentose-phosphate, Leloir, DeLey-Doudoroff, and EMP pathway) lead to EMP substrates (i.e., to pyruvate consecutively decarboxylated to Ac-CoA!).
-Ac-CoA constitutes the substrate for two competitive reactions, the thiolase reaction (first step in the PHB synthesis), and the citrate synthase reaction (substrates for synthesis of cell substances).
-The breakdown of the complex nitrogen source (deamination and degradation of casamino acids) leads to Ac-CoA.An independent Monod kinetic of Ac-CoA synthesis from both sugars was assumed for the modelling strategy.
-Negative feedback control mechanism of Ac-CoA synthesis was built-in as strategy for its own regulation (Luong type of inhibition pattern).
-Ac-CoA was foreseen for consumption toward biomass formation, energy supply (including maintenance energy), NADPH generation in TCA cycle, 3HB accumulation (thiolase substrate), a-ketoglutarate excretion, for the production of one chemically unknown (excreted) compound and for the PHB-polymerase catalytic unit synthesis.
-Ac-CoA synthesis from sugars is inhibited (but not stopped) by complex nitrogen source (application of the Jerusalimsky type of inhibitory influence).
-A small quantity of biomass is assumed to be synthesized directly from the complex nitrogen source (casamino acids were used as a separate C/N source for protein synthesis).
-PHB synthesis rate was assumed to be proportional to the intracellular biocatalyst concentration (PHB-polymerase complex).Polymerase synthesis rate from Ac-CoA and nitrogen sources (multiple Michaelis-Menten kinetics) with particular degradation rate (protein turnover) was assumed.
-The PHB-polymerase complex synthesis was assumed to be inhibited by the presence of high levels of the complex nitrogen source (Jerusalimsky type of kinetic equation); its activity was assumed to be started after the complex nitrogen source was almost completely depleted.
-Inorganic nitrogen source consumption was foreseen to be inhibited by complex nitrogen source (no consumption of NH 4 + in the first part of the growth phase when casamino acids were provided in sufficient amounts).
-The excreted metabolite a-ketoglutarate was assumed to be influenced by a negative feedback controller of its own synthesis.The same was applied for an additional unknown metabolite.NH 4 + and PHB polymerase were adopted as inhibitory agents of production of these two metabolites (Jerusalimsky type of equation).
The aforementioned approach is more advanced and more demanding if compared to the models published earlier, but it is relatively simple if compared to high-structured, metabolic, genome-scale models.This simplified approach in modelling appears to be a promising tool in the case of insufficiently known metabolic pathways and/or rate control mechanisms making the building of metabolic flux models impossible.
An additional example of new application of formal kinetic modeling for process optimization is the work published by Horvat et al. 50.In this work, the partially growth-associated production of PHA under nitrogen-limited growth was chosen as the modelling strategy, thus the Luedeking-Piret's model of partial growth-associated product synthesis was used as a working tool.Specific growth rate relations according to Megee et al. 40 and the Mankad-Bungay relation 51 , adapted to a double substrate (C and N source) limited growth, were tested.The first stage of an investigated bioreactor cascade was modelled according to the principle of a nutritionally balanced, continuous biomass production system; the second as a two substrate controlled process, while the three subsequent reactors were adjusted to produce PHB under continuous C-source fed, but under strict nitrogen deficiency (limitation).The simulated results of production, obtained by applied mathematical models with computational optimization, indicate that PHB productivity of the whole system could be increased significantly if certain experimental conditions are provided regarding the overall dilution rate, C-and N-source feed concentration.Additionally, a supplemental feeding strategy for switching from batch to continuous mode of cultivation was proposed to avoid the substrate inhibition.
Recently, Mozumder et al. 52 have contributed a mechanistic model describing the production of PHB by, similar to Horvat et al. 50, a pure culture of C. necator DSM 545, which was calibrated and validated for the two different substrates, glucose and surplus glycerol, derived from biodiesel production.Respectively to glucose and glycerol, non-growth-associated PHB production was triggered by low ni-trogen source concentration.Besides non-growth-associated PHB production, some growth-associated PHB synthesis was observed, although inhibited in the presence of nitrogen source.Biomass growth, PHB and non-linear product inhibition of PHB were included in the model.The established models have taken into account four main processes: biomass growth on carbon substrate, biomass growth by utilizing intracellular PHB, PHB production and maintenance energy.Three models were used for further model calibration and model selection: (A) represents biomass growth on substrate only; (B) takes into account biomass growth on intracellular PHB; (C) was based on model (B), but cell density limitation for biomass growth on both substrates and on intracellular PHB was additionally included.The mass balances of different components in the bioreactor comprised two crucial parts: macroscopic transport and biochemical conversion.All models containing well-known formal-kinetic terms for specific growth rate, specific product synthesis rate, product-, substrate-, and biomass inhibition were subjected to sensitivity analysis, model calibration and model validation.The developed mathematical model was successful in the simulation of PHB production, and in the prediction of dynamic behaviors related to heterotrophic biomass growth for the twophase pure culture system 52 .
In the field of formal-kinetic, mechanistic and low-structured modelling, in addition to the aforementioned, very different models can be found in the scientific database.They consider different strains, substrates, mutants, kinetics, and cultivation properties.Raje and Srivastava 53 have investigated continuous production of PHAs by Alcaligenes eutrophus (today: C. necator) B4383 in a nitrogen-limited culture.To define the specific growth rate µ, these authors have used a linear combination of Monod kinetic term and sigmoidal Hill term, additionally corrected with the inhibition term according to Luong containing N/C ratio and with PHA inhibition.The rate of PHA formation was arranged to be growth-and non-growth-associated.A few years later, Pathwardan and Srivastava 54 applied the identical model for the strain Ralstonia eutropha (today: C. necator) NNRL 14690.
Yu et al. 55 and Wang et al. 56 have worked with C. necator ATTC 17699 in studies of PHAs biosynthesis from volatile fatty acids (VFAs), particularly acetic, propionic and butyric acid, respectively.They applied Monod kinetics for growth rate, which was corrected by the cell activity impacted by the toxic activity of acids.Leudeking-Piret term was integrated into the model in order to describe PHA production rate.The biological system was tested on all three acids separately and on mixtures of acids.
Shahhosseini et al. 57 were active in the investigation of fed-batch cultivation of Ralstonia eutropha (today: C. necator) ATTC 17697 on fructose and NH 4 + .These authors have introduced the Monod term for the definition of specific growth rate with the ratio of nitrogen to fructose (N/C) as the function argument.The equation for specific growth rate µ was corrected by N/C ratio incorporated in the Luong inhibition term.As mentioned previously, Yu et al. 55 also incorporated the dependence of production rate on biomass concentration and growth rate.
R. eutropha (today: C. necator) NRRL B14690 strain was the subject of batch/fed-batch/continuous cultivation and mathematical modeling by Khanna and Srivastava [58][59][60][61] .For the microbial growth rate, these authors applied the double substrate (fructose and nitrogen) Hill (sigmoidal) kinetic expression corrected by the double substrate inhibition term according to Luong.The PHA production term was modeled in accordance with Leudeking-Pirets equation, and consumption of both limiting substrates (fructose and nitrogen/or phosphorus) was assumed for cell maintenance.These authors did not consider product inhibition.
Monod growth kinetics corrected with an inhibition term was combined by Patnaik 62 with the dispersion model in modelling of growth of R. eutropha (today: C. necator) ACM 1296 on fructose in a 2-liter bioreactor.Similar to Shahhosseini et al. 57 , Patnaik 62 applied the N/C ratio as an argument in the specific growth rate equation, thus the dependence on two substrates was included.Imperfect mixing of bioreactor content was represented by a Pecklet number with 20 as the optimal value.Shang et al. 63 conducted high cell density fedbatch culturing to produce PHB by Ralstonia eutropha (today: C. necator) NCIMB 11599 under phosphate limitation.High glucose concentration was found to be inhibiting for PHB synthesis.An unstructured model was proposed for predicting the cell growth, PHB synthesis, glucose, and P-source consumption with the phosphate concentration as the key factor triggering both accumulation of PHB and cell growth.Specific growth rate was expressed in accordance to double substrate limitation (C and P) Monod kinetics, and specific PHB production rate was characterized by combining Andrews kinetic term (employed glucose inhibition at high concentration), Jerusalimsky term (triggering by phosphate limitation) and logistic term (product inhibition caused by steric hindering).
Faccin et al. 64 tested four different models intended to be useful for PHAs synthesis by Bacillus megaterium on sucrose as the main carbon source.These authors used the models proposed by Mulchandani 36 , Raje and Srivastava 53 , as well as Khanna and Srivastava [58][59][60][61] as the starting point to con-struct the basic formulation of four models.They modified these models by: (i) changing the expression for the specific growth rate (μ); (ii) introducing a correction factor (Φ) for products formation rates; and (iii) introducing a term for cell death in the expression for residual biomass.
The first model was a modified Monod double substrate kinetic with sucrose and nitrogen as limiting substrates, without substrate inhibition and without Φ correction factor, the correction of product synthesis rate.In the second type, a Monod model equation for the specific growth rate with nitrogen as the only limiting substrate (sucrose limitation was not included), but with Φ as correction factor.The next two models were based on Monod double substrate kinetics with sucrose and nitrogen as limiting substrates; these models additionally take into account the influence of the pH-value (one according to Ghose and Tyagi 65 and the other according to Åkerberg et al. 66 ).Both were equipped with Φ factor.When included, the correction Φ factor was a hyperbolic time dependent variable based on the experimental observation that the production of PHB by B. megaterium displays a two-stage process regarding the carbon flux towards residual biomass or polymer formation.
Gahlawat and Srivastava 67 studied batch cultivation of Azohydromonas australica DSM 1124 for growth-associated PHB production.Based on experimentally achieved production kinetics, the mathematical model was established.Thereafter, by using continuous fermentation modeling principles, the model was extended to fed-batch production in order to identify nutrient feeding regimes to improve PHB accumulation.The authors used double substrate Monod kinetics for sucrose and nitrogen, and corrected by inhibition terms for substrates (Luong type for nitrogen and Jerusalimsky type for sucrose).Direct growth-associated PHA production dependence was incorporated.
Recently, Špoljarić et al. 68,69 have published two modeling articles dealing with PHA production by C. necator DSM 545.Glycerol (GLY) and fatty acid methylesters (FAME) from biodiesel production were the main substrates in these works.Mathematical models were established in order to optimize these processes.For glycerol, five different expressions for double substrate (C, N sources) dependent growth rate with appropriate substrate inhibition terms were tested: double substrate Monod relation was combined with different substrate inhibition terms known as Yanno et al., Webb, Andrews and Noack, Aiba et al. relations as reviewed by Moser 44 .As proposed by Guthke, the specific car-bon-substrate consumption rate for the maintenance of cells was defined according to simple Monod relation.Leudeking-Piret relation was applied for the production rate; the depletion of nitrogen source was the trigger for PHA synthesis.
In the second work, two low-structured mathematical models for fed-batch production of PHB and poly(hydroxybutyrate-co-hydroxyvalerate) (PHBV) on renewable substrates (GLY and FAME) were established.GLY combined with glucose and FAMEs enriched with valeric acid (VA) were the main C-substrates.The models were used for development/optimization of feeding strategies of carbon and nitrogen sources regarding PHA content and polymer/copolymer composition.
The following modeling principles were applied when glucose and GLY were used as C-sources: (a) Residual biomass (Xr1) is synthesized only from glucose when nitrogen source is available; the growth on GLY in the presence of glucose is approached as practically negligible, evidencing that the consumption of GLY severely affected the availability of glucose (b) The model reflected also the experimental lag-phase.
(c) PHB is synthesized from two carbon sources, glucose and GLY.In the exponential phase of growth, it was assumed that PHB synthesis from glucose is growth-associated with an experimentally determined time delay (biosynthesis after beginning of the exponential growth phase).In the fedbatch part of fermentation and in the stationary phase of growth, PHB synthesis from glucose switches to the non-growth-associated type, and nitrogen source depletion was assumed to be an activator of the non-growth-associated type of synthesis.The synthesis of PHB from GLY was assumed to be non-growth-associated, to occur at maximal rate, but inhibited if glucose is present in the broth.This was reflected by introducing an inhibition term.
(d) Maintenance is generated by metabolism of glucose.
(e) Nitrogen source was added for pH correction; the quantity of consumed nitrogen source is assumed to be proportional to the growth rate in the exponential phase.
The following modeling principles were applied for PHBV production under N-limited conditions when FAME and VA were used as carbon sources: i) Total biomass is divided into two ''compartments'': residual biomass (Xr2) and PHBV.
ii) Residual biomass (non-PHA part of biomass) is synthesized only from biodiesel (FAME) and NH 4 + (as C and N sources) because valeric acid was added only during the non-growth-phase.Growth was assumed to be limited by depletion of nitrogen source.At the same time, FAME and VA were assumed as sources for cell maintenance.Furthermore, the experimental lag phase was incorporated.
iii) PHBV is divided into two ''sub-compartments'': 3HB and 3HV.iv) 3HB is synthesized from biodiesel particularly during the growth phase (growth-associated 3HB synthesis) and dominantly during the nongrowth phase (non-growth associated 3HB synthesis).
v) 3HV is synthesized only from VA during the non-growth phase (non-growth associated 3HV synthesis).Conversion of FAME in 3HV was not detected during the growth phase.
vi) Part of the carbon sources (FAME and VA) are lost by synthesis of CO 2 and other minor metabolites.
Some general comments to all the models described above: the formal kinetic modelling of microbial growth was applied first as the simplest among the group of kinetic models.Because of their simplicity, such models are still widely used.When such mathematical systems are insufficient to represent the biological system, low-structured models are introduced that are still based on formal kinetic principles.It is a challenging task to strictly distinguish low-structured and high-structured models of microbial metabolism (here, the corresponding number of included reactions is not conventionally fixed).This partition line often depends on the level of our knowledge about the objects/subjects of modelling.It is clear that genome-scale metabolic models tend to be white-box models; therefore, they certainly belong to the high-structured ones, as they tend to contain the true number of reactions as present in the natural system.In this review, in accordance with a non-written rule, we attempted to separate the models into categories, as follows: I) No metabolic reactions, extracellular molecular species (substrates, products) and biomass itself are interesting variables and therefore the object of mathematical modeling → formal kinetic models; II) Extracellular molecular species, biomass and a few intracellular metabolic reactions more or less lumped, without genetic regulations → low-structured models; III) Complete, interconnected, highly networked metabolic pathways present, and/or genetic regulation of pathways present → high-structured models.
To provide the reader with a more complete picture of formal kinetic and low-structured kinetic models, an overview of such mathematical models is provided in Table 1.

Dynamic models in PHA biosynthesis by pure and mixed cultures
The growth of microbial biomass is strictly related to environmental conditions (e.g., substrate concentrations, temperature, pH-value, DO, etc.).These factors are not the only ones that act on the properties of a growing microbial population.Genetic factors and individualism of each cell inside the whole population are additional variables that determine microbial processes.Unstructured formal kinetic models do not take into account such "individual" factors.They are usually of the steady-state type (regardless of whether is it in stabile, exponential, nutrient balanced, in equilibrium with nutrient supplying, or in stationary phase of growth).The traditional (early) approach for modeling of PHA-producing microbial cultures was based on the figure of an unstructured, non-segregated biophase during exponential growth, substrates that are sufficient for most of the cultivation period and a gradientless environment, i.e., unrestrained growth, without the substrate(s) concentration as the limiting parameter(s).Unfortunately, such models are unable to predict the dynamic behavior of the microbial culture under the real range of conditions.In fact, they do not consider cell segregation and ignore the inhibitory effect of PHB on the dynamic behavior of the cell culture.Dynamic models account for time-dependent changes in the state of the system; they are typically represented by differential equations and, in the case of PHA biosynthesis, they are usually structured.This model type is appropriate for describing different transient states (e.g., from lag-phase to exponential phase of growth, exponential to stationary phase of growth).Concerning structuring, bioreactor compartments, biomass itself (producers/non-producers, young/old cells, plasmids harboring/non-harboring, PHA steric hindered/no hindered, gene expressed /non-expressed, biologically active/non-active), metabolites (intracellular/extracellular), and products (homopolymer/heteropolymer, intracellular/extracellular) can be found in literature as objects of structuring.The dynamic responses of genetically structured models were published earlier 30,31 .Detailed kinetic information for basic cellular processes, e.g., enzymatic reactions, protein-DNA and protein-protein interactions, are the prerequisite to describing the dynamics.By combining the known stoichiometry of metabolic pathway reactions with kinetics, such dynamic models can be established.This approach is far from being an easy endeavour.For example, Rizzi et al. 29 applied 22 ma-terial balances for metabolites and 23 enzyme rate equations to predict the intracellular and extracellular metabolite levels in transient phase caused by glucose pulse in a continuous culture of S. cerevisiae.This model was successful only in a timescale of 120 seconds, but enzyme synthesis/degradation of pivotal importance in larger timescales have not been included.Regulating activities in the aforementioned work were included at the individual enzymatic level, not on the genetic level.Further improvements of dynamic and metabolic models have solved some of the aforementioned problems.
Regarding dynamic models for PHA synthesis, the basic type could perhaps be the one proposed by van Aalst-van Leeuwen et al. 70 , which deals with activated sludge microorganisms.It is well known that activated sludge organisms from the wastewater treatment plants respond to feast/famine regimes by the production of the storage polymers glycogen and PHA.PHA seems to be the more common when cells are exposed to excess carbon source.A pure culture of Paracoccus pantotrophus LMD 94.21 was the object of the case study in this work, where the steady-state C-limited chemostat culture was switched to batch mode with a pulse of added acetate as C-source.As long as the external substrate (acetic acid) was present, growth and accumulation of PHB took place.After depletion of the external substrate, intracellular stored PHB was used as growth substrate.The accumulation of PHB was found to be strongly dependent on the growth rate of the organism before the acetate pulse addition.Mass-balanced PHB was correlated to the difference in acetate uptake rate and the acetate quantity required for growth.Using the C-mole convention, a metabolically structured model that adequately describes the observed kinetics of PHB formation/ consumption, has been set up.
Katoh and coworkers 71 published a work dealing with the dynamics and modeling of the fermentative production of PHB by a mixed culture of Lactobacillus delbrueckii and A. eutrophus (today: C. necator) in one fermenter ("one-pot-reaction"). Sugars were transformed by L. delbrueckii in lactate, which was the substrate for PHB synthesis in A. eutrophus.Metabolic flux distributions were computed for two cultivation phases: the cell growth and PHB production periods.It was found that when NH 4 + was abundant, the NADPH generated through isocitrate dehydrogenase in TCA cycle was predominantly utilized for the a-ketoglutarate-to-glutamate reaction.On the contrary, when NH 4 + concentration decreased under a certain level or completely depleted, NADPH tended to be utilized for the PHB production (acetoacetyl-CoA reductase reaction).Several mixed culture experiments, conducted to see the dynamics of the system, were fundamental for the development of mathematical models able to describe the dynamic behavior of mixed cultures.According to the authors, the model may be used for control strategy, process optimization, and process dynamic investigations.
Roussos and Kiparissides 72 developed a bivariate population balance model (PBM) to describe the dynamic evolution of a PHB-producing microbial batch culture.The population balance model was solved by three different numerical methods: continuous finite element method, discrete-continuous finite element method, and a discretized method.Different high-resolution finite difference schemes for the calculation of the cell flux term were applied.The applied numerical methods were evaluated respectively to the accuracy and computational requirements.For this purpose, the authors directly compared the ongoings at pivotal instants of time of the developed bivariate distribution to the respective quantities calculated that were obtained by a simple homogeneous model for the non-product inhibited case.Results obtained from the solution of the 2-D PBM were compared with those achieved by the 1-D PBM (considering the product inhibition reflected on the residual biomass growth rate).It was shown that a bivariate population balance modeling approach is required to describe the dynamic behavior of the microbial culture when the intracellular accumulated PHB inhibits the biomass growth rate.The 1-D PBM was not accurate enough to describe the dynamics of the microbial system.PBMs may consider that the total biophase is distributed into a number of individualities (different cells).Due to the high degree of detailing, PBMs represent very accurate tools to describe the complicated phenomena associated with cell growth, nutrient uptake, and product formation.They allow the modeling of cell division with partitioning of cell space and intracellular PHB into new cells.Furthermore, instead of average cell properties predicted by formal kinetic models, PBMs can even predict the distribution of a characteristic property over the whole cell population.Before the Roussos / Kiparissides work was published, among PBMs only the article by Mantzaris et al. 73 was dedicated to the PHB accumulation and its inhibition effect on the behavior and dynamic of microbial cell cultures.
As the content of PHA in microbial cells is an important factor for the development of activated sludge properties (e.g., sludge volume index, sedimentation velocity), PHA storage kinetics and simultaneous growth of cells were investigated for aerobic heterotrophic biomass under unsteady feast conditions 74 .An additional reason for this effort was the idea that activated sludge could act as a cheap source for PHA production.Insel and coworkers 74 investigated the short-term variations in growth and storage kinetics of activated sludge under disturbed feeding conditions.To achieve this target, a multi-component biodegradation model was used, supported by a multi-response modeling procedure and identifiability analysis.Under the condition of ample availability of external substrate, the heterotrophic biomass was able to increase the direct growth activity with simultaneous reduction of substrate storage capability.Reduced sludge age (SRT) from 10 to 2 days had increased the maximum specific growth rate from 3.9 to 7.0 per day without affecting the maximum storage rate.The alteration of sludge age provoked the elevation of value of half-saturation constant for growth from 5 to 25 mg COD L -1 .This was interpreted as an indication of a change in species distribution among the activated sludge consortium.The increase in primary growth metabolism in parallel with the reduction of PHA storage rate was validated by modelling two different sludge ages, both at availability of external substrate.The alteration of feeding conditions was found to differently influence storage and growth kinetics.Modifications in model structure were highlighted in order to reduce the tedious task of frequent recalibration efforts caused in variations of process conditions.
Except wastewater treatment-related PHA synthesis, its production from defined substrates was also the target of dynamic modeling.Wu et al. 75 investigated an optimal adaptive control mechanism of fed-batch processes of R. eutropha (today: C. necator).A dynamic model of the fed-batch PHB production by R. eutropha was prepared.Because of the need to detect "to uncertainties very sensitive kinetic parameters", the authors developed a robust and feasible feeding strategy.A specific optimal adaptive control strategy was established by solving the constrained discrete-time optimization algorithm using the genetic algorithm solver from Matlab software package.Based on simplified kinetics, the proposed control methods were implemented to the feed flow manipulation with stepwise changes.
A Lyapunov stability analysis, encompassing state estimation errors, and investigation of parametric uncertainties was performed.In addition, by using simulations, a satisfactory output tracking performance was achieved by a two-input control configuration.
In this work, a set of simplifying hypotheses was applied: -Nitrogen is the limiting substrate affecting growth kinetics.
-All nutrients are in excess except C-source (acetate) and N-source (NH 4
-The broth is perfectly mixed without diffusion and with negligible flock formation.
-DO is not a limiting factor.
-The present organisms have an average metabolic behavior that does not change significantly in time.
-A simple two-compartment cell model, consisting of PHB and active biomass more or less mixed cultures.
-The rate of volume loss due to evaporation is neglected.
The main results were as follows: -The optimal adaptive control action is obtained by solving the constrained discrete-time optimization algorithm using the genetic algorithm solver from Matlab.
-The proposed control scheme is based on the optimization of kinetic parameter estimation.
-The stability analysis of the closed-loop system with respect to specific constraints and specified learning rates was provided.
-It was demonstrated that a satisfactory output tracking performance is achieved by using a two-input control configuration.
An overview of dynamic oriented modeling approaches that consider PHA production is given in Table 2.

General
The main characteristic of metabolic mathematical models is that they should reflect the real biochemical metabolic situation of cells.From the complexity point of view, metabolic models refer to a very broad range: from simple cases, where the metabolic pathway with at least two or three enzymatic reactions is investigated, to more or less complicated metabolic networks that represents the main catabolic and anabolic routes of the cellular metabolism.Within very complex networks that consist of hundreds or more reactions grouped in networked, interconnected pathways, the metabolic control and regulatory mechanisms specific for living cells must be included in the model.These regulatory mechanisms are usually not completely elucidated; therefore, the establishment of fully mechanistic models that can describe the whole cellular metabolism and/or the cell life cycle is not possible at the present state of science.That is why all metabolic models contain some simplifications (e.g., lumped reactions, simplified regulating patterns, simplified kinetics, etc.).There is an intention in science to lift such "grey box" model types to the level of "white box" models.Development of advanced analytical techniques dealing with the measurement of intracellular metabolite concentrations, intracellular/membrane-bound enzyme activities, regulating factors and signaling molecules lead to the formulation of higher structured and more improved models that have more possibilities to "reflect" real cell physiology.The improvement of speed and processing power of computers as well as the development of metabolic software has opened the possibility to solve the very complex non-linear systems of differential equations, necessary for successful application of such models.Metabolic mathematical models that describe the cellular metabolism of PHA producers play a central role in elucidation of metabolism and substrate transport properties of existing strains; further, they are unavoidable in the development of new strains by methods of metabolic engineering.Among metabolic models, kinetic models (usually based on both enzyme and/or microbial kinetics as well as on the stoichiometry of metabolic reactions), stoichiometric models (based on the time invariant characteristics of metabolic networks), and cybernetic models (based on optimal nature of microbial processes and on the metabolic regulation by the cybernetic framework which has to be discriminated from the pure kinetic ongoings) 76 can be found in the scientific and professional literature.Furthermore, hybrid metabolic models that combine at least two of the aforementioned types of model properties are in use.Kinetic metabolic models function on the principles of biochemical reaction kinetics.If information on all the enzymatic reactions (of the whole metabolic framework) for an organism would be available, it would be possible to build and apply a detailed model in order to interpret experimental results, to predict the dynamic changes in the cells, and to predict the changes in the cultivation broth.Such models should reflect the shifts in the cultivation conditions and the eventual genetic changes.Usually, the lack of information on general regulatory mechanisms and the predominant availability of in vitro kinetic parameters instead of the desired availability of in vivo data are the limiting factors for such models.That is why established models are often limited to a certain biological situation, e.g., steady-state metabolic fluxes during exponential phase of growth, PHA degradation under C-source limitation, or resting cell activities.
Stoichiometric models are characterized by application of strict stoichiometric relationship considering extracellular substrates, products, and intracellular metabolites.To each molecular (or macroscopic) species (intracellular or extracellular), at least one "pool" is related as an object of mass balancing.Metabolic flux analysis (MFA) is widely performed to achieve the quantification of the intracellular fluxes of central metabolism.In MFA, mass balances are used to calculate the fluxes through the different branches of the metabolic network.The establishment of a stoichiometric matrix is the first step that follows after the metabolic fluxes; these fluxes can either be measured, or calculated.In the second step, the stoichiometric matrix and vectors of reaction rates must be multiplied to define the mass balance equations.After that, an algebraic manipulation of the stoichiometric matrix should be performed.According to different authors (depending on the objective of the analysis), the further work can be done in five different ways: -One way is to impose some constraints (e.g., by measuring some fluxes) so that the determined system can be solved by principles of linear algebra (flux balancing!).
-If the labelling state of some metabolites is measured, further constraints may be imposed.This leads to the combination of flux balancing with labelling balancing, so the system has to be solved numerically.
-In the case when a system cannot be constrained to a "determined" type, the linear optimization is a suitable tool for finding the maximum or minimum of a chosen objective function (OF).
-Applying the convex analysis in order to determine elementary flux modes (EFMs) 77,78 .
-By calculating the biochemically meaningful base vectors for the null space of the stoichiometric matrix.
By applying the procedures described above, it is possible to obtain a snapshot of the metabolic situation under a particular condition.The main drawback of stoichiometric models is their limited predictive power usually caused by the lack of regulatory information.Some basic principles of metabolic modeling can be found in the article published by Giersch 79 .
The cybernetic models as a third main group among metabolic models will be treated in a separate chapter.
Metabolic models related to PHA biosynthesis in activated sludge cultures Kinetic and stoichiometric models were extensively used in works dealing with MMCs from the wastewater treatment field.Yu and Wang 80 developed a metabolic cell model with five fluxes in order to depict the detoxification mechanism, mass transfer and acetyl-CoA formation related to acetic acid.The model encompassed the formation of PHB, active biomass, and CO 2 as three final metabolic products.Fluxes were measured under different conditions.Based on achieved results, it was concluded that, under the present high cell mass concentration, elevated extracellular acetate concentrations occur by increased metabolic fluxes resulting from acetate detoxification.It was concluded that the magnitude of detoxification fluxes is highly dependent on the intracellular acetate pool.
Veldhuizen and associate authors published an article about modelling of biological phosphorus and nitrogen removal that is closely connected to PHA synthesis and degradation in a full-scale activated sludge plant 81 .An activated sludge model 82 (ASM2) foreseen for the removal of COD, N and P was the basis for model development.The entire model was implemented by means of a SIMBA 3.0 computer software package (based on MATLAB/ SIMULINK).Performance and validity of the model was successfully tested with the help of experimental data elaborated both under aerobic and anoxic conditions in the full-scale wastewater treatment plant.
Yagci and associates 83,84 published a new metabolic model based on an activated sludge model (ASM2d) 85 dealing with PHA synthesis and degradation in a sludge microbial community under anaerobic conditions and acetate uptake performed by phosphate-and glycogen-accumulating organisms (PAOs and GAOs, respectively).The variable overall stoichiometry, based on the assumption that PAOs use the glyoxylate pathway to produce the reducing agents for PHA synthesis, was basically applied.A strong agreement was observed between experimental values and model predictions for PHA production, PHA composition, glycogen utilization, and P release.
The mathematical modeling of similar processes was also performed by Dias and coworkers [86][87][88] .Primarily, these authors have developed a simple two-compartmental metabolic model for PHA production by mixed microbial culture 86 .After that, they decomposed the metabolic network into separate EFMs to develop a dynamic/hybrid/semi-parametric model 87 that was able to simulate the application of propionic acid as C-source.Thereafter, the metabolic model for copolymer production from mixtures of acetic and propionic acid was established 88 .Material and energy balances were established on the basis of previously elucidated metabolic pathways.Equations were derived for the theoretical yields (related to cell growth and PHA production) as functions of the oxidative phosphorylation efficiency (i.e., P/O ratio).Two different feast and famine culture-enrichment strategies were studied: either with acetate or propionate.Metabolic flux analysis (MFA) and flux balance analysis (FBA) were performed to achieve the optimal feeding conditions, culture enrichment, and high quality PHA.
Johnson and associated scientists 89 used the basic model published before 70 for the model-based data evaluation of PHB-producing MMCs in aerobic, sequencing batch and fed-batch reactors.A two-step process was reported: (I) for biomass, and (II) for PHB enrichment.These authors implemented the following procedure: (1) Measurement of a sufficiently large set of parameters including off-gas concentrations.
(2) Corrections of measurements for effects of sampling and addition of liquids (pH control, substrate control).
(3) Calculation of oxygen uptake and CO 2 evolution rates.
(4) Balancing of the measured conversions.
(5) Evaluation of measurements by means of a metabolic model.
The feast phase was characterized by: (i) acetate uptake affected by the PHB inhibition, growth of biomass, maintenance, (ii) PHB production affected by PHB inhibition, CO 2 evolution, oxygen and NH 4 + uptakes.Additionally, the famine phase was characterized by PHB degradation, maintenance, biomass growth, CO 2 evolution, and oxygen and NH 4 + uptakes.Conversion of acetate, NH 4 + and the storage polymer (PHB) were described very accurately by the developed model.
Jiang and coworkers 90 investigated the production of PHBV copolyester from acetate and propionate by P. acidivorans-dominated mixed microbial cultures.The used model was based on the work performed earlier by Dias et al. 88 with the following improvements that were taken into account: 1. Metabolic reactions active during the famine phase were included into the model structure; 2. Biomass synthesis from acetyl-CoA and propionyl-CoA were separately specified; 3. Maintenance energy (in the form of ATP requirement) was used as an estimated parameter in the kinetics, whereas the P/O ratio was fixed; 4. The TCA cycle was assumed inactive when propionate was the sole carbon source.
The regulation for simultaneous acetate and propionate uptake was included in the model as a function of the acetate-propionate substrate mixture composition.
Simultaneous growth and storage kinetics variation (for the aerobic heterotrophic biomass under unsteady feast conditions) was also an object of modeling 74 that belongs to the topics of MMCs and PHA production.Tajparast and Frigon 91 studied the storage metabolism during feast/famine cycles of activated sludge treatment systems together with the foaming, bulking, and process optimization for the production of value added by-products (e.g., bio-plastics).PHB, glycogen, and triacylglycerols (TAGs) were examined during feast-famine cycles using two genome-scale metabolic models: (i) for Rhodococcus jostii RHA1 (iMT1174) and (ii) for Escherichia coli K-12 (iAF1260), both foreseen to simulate growth on glucose, acetate, and succinate.The target was to develop the proper objective function (OF) suitable for the model prediction of the storages compounds production under feast/famine cycle conditions.Using the flux balance analysis, combinations of three OFs were tested.The main OF was established to maximize growth rate, but for the time interval between feast and famine phases, two additional sub-OFs were introduced: minimization of biochemical fluxes, and minimization of metabolic adjustments (MoMA).When glucose and acetate were set as sole carbon sources, all (sub-)OFs predicted identical substrate-storage "associations" (i.e., glucose-glycogen and acetate-PHB).These predictions were in good agreement with experimental observations.Interestingly, in case of succinate as substrate, the predictions depended on the network structure: the metabolic model for E. coli predicted glycogen accumulation and the R. jostii model predicted PHB accumulation.The authors concluded that new modeling insights between metabolic predictions and popula-tion ecology would be necessary to predict metabolism properly.
Tamis et al. 92 provided an overview on currently described models, and developed concepts for a generalized, or, as they said, "more predictive" model.The focus was devoted to the feast/famine parts of the process.Model improvements were introduced for modeling of mixed substrates uptake, microbial growth in the feast phase, switching between feast and famine phase, PHA degradation, and modeling of the accumulation phase.The authors postulated that the scientific community would welcomed the establishment of an improved generalized model.
The mathematical modeling of PHA production by microbial cultures from activated sludge is tightly connected to modeling of the entire wastewater treatment processe.A short overview of published models dealing with the wastewater treatment is given in a review work published by Gernaey et al. 93 .Metabolic models related to PHA biosynthesis in activated sludge cultures and their characteristic properties are summarized in Table 3.

Metabolic models targeted for industrial PHA biosynthesis
One of the early published metabolic models dealing with PHAs was that by Daae and coworkers 94 targeted for PHA biosynthesis in plants.The mathematical simulation was programmed in SCAMP using the metabolic pathway simulation package established by Sauro 95 .This mathematical model was based on kinetics of the key enzymes responsible for PHA synthesis: ketothiolase, acetoacetyl-CoA reductase, and PHB synthase.For this purpose, the complex enzyme kinetic was used: a reversible ping-pong BiBi mechanism, reversible sequential ordered BiBi mechanism, and Michaelis-Menten equation for nonreversible system corrected by incorporated competitive inhibition term.Steady-state production rates and steady-state concentrations were achieved by dynamic simulation.Flux control coefficients and concentration control coefficients were calculated at the steady states according to rules of metabolic control analysis (MCA) 96,97 .This paper outlined the theoretical prospects of producing the copolymer PHBV in plant plastids.The model indicated that both the 3HV/3HB ratio and the copolymer production rate can vary considerably between dark and light conditions.Furthermore, it was illustrated that natural variations in substrates and cofactor levels have a significant influence on the production rate and copolymer ratio.
Dealing with biological copolyesters, Xu and coworkers 98 used a mathematical model for regulat-ing monomer composition of microbially synthesized medium-chain-length copolymers (mcl-PHA).The authors simplified the biochemical reactions and the metabolic network.To provide a general model, both reversible and irreversible reactions were considered.The pseudo-steady-state assumption was taken into account; consequently, the concentration of each intermediate did not vary in time.In addition, it was also assumed that the type II PHA synthase, encoded on the PhaC1 and PhaC2 operons, from P. oleovorans and P. aeruginosa did not differ too much in its substrate specificity.Because both prefer 3-hydroxyoctanoyl-CoA and 3-hydroxydecanoyl-CoA as substrates, the simplification to "one synthase kinetic" in the model was very reasonable.The model was segregated into the following parts: -single carbon source, -effect of chain length of the related carbon source, -quantitative selectivity of the enzymatic system, -effect of the enzyme levels on the copolymer composition, -mixed carbon source, and -genetically engineered pathways.The authors reveal that if a sole carbon source had been used, the monomer compositions would have depended on the enzyme selectivity, while the cultivation time, concentration of the carbon source and enzyme levels would have limited effects, provided that the enzyme levels did not shift in time after PHA synthesis had started.Quantitative selectivity of the enzyme system towards different intermediates was determined, so the highest selectivity towards the intermediates with 8-10 carbon atoms was evidenced.Furthermore, the model was able to predict the effect of mixed substrates and genetic engineering on copolymer composition.In addition, the model was further applied for the simulation of short-chain-length PHA (scl-PHA) copolymer synthesis, as well as the copolymers of 3HB and mcl-hydroxyalkanoates.The conclusion was that manipulation of two simultaneous and independent metabolic pathways responsible for the precursor formation could effectively regulate PHA copolymer compositions.
Copolymer synthesis was also the aim of the work performed by Iadevaia and Mantzaris 99 .These authors investigated a genetic network driven control of PHBV copolymer composition.The strategy was based on the modification of an artificial genetic toggle that controls the ratio of monomer concentrations.For this purpose, the authors applied a mathematical model consisting of three coupled parts for describing the dynamics of three modules: -Artificial genetic network, which should control the expression levels of the enzymes catalyzing the formation of the monomers, -Monomer formation kinetics, and -Monomer polymerization dynamics (PHA synthesis).
The monomer formation module was coupled to the driver genetic toggle module through the dependence on expression levels, whereas the polymerization was coupled to monomer formation due to the dependence of the initiation, elongation and transition rates on the monomer concentrations.The presented simulation results showed that the steadystate mass fraction of 3HB and molecular mass distribution of the copolymer PHBV were strongly affected by the ratio and absolute values of the chain elongation rates.At the same time, the dynamics of the polymerization process was related to and determined by the specific growth rate.According to the authors, the presented model can be used as a predictive tool for the rigorous design of bioprocesses in order to produce PHBV chains of desirable compositions.
Ao and associates 100 worked on a systematic method for constructing large-scale kinetic metabolic models.Initially, the authors had applied this method for modeling of the metabolism of Methylobacterium extorquens AM1.Its central metabolic network included formaldehyde metabolism, serine cycle, citric acid cycle, pentose phosphate pathway, gluconeogenesis, PHB synthesis and acetyl-CoA conversion pathway, respiration and energy metabolism.The authors reported that they had overcome "an outstanding difficulty in large-scale kinetic modeling: the requirement for a massive number of enzymatic reaction parameters".The latter was done through a systematic and consistent procedure of parameter finding.They were able to construct a kinetic model containing 80 metabolic (computational) species that were based on general biological considerations and incomplete experimental kinetic parameters.The following steps were performed: -The number of enzymatic parameters for the generic enzymatic rate equation was reduced to a minimum set (still preserving crucial characteristics!): -A set of steady-state fluxes and metabolite concentrations in the physiological range were applied to restrict the parametric space of enzymatic reactions; -If experimental values for enzyme saturation constants (and equilibrium constants) are unknown, the values under physiological concentrations should be chosen; -"Designing a dynamical exchange for the coupling between the metabolism represented in the model without including the rest", valid for models which do not cover the entire metabolic network of the organism.
The results of the work described above was a robust model with respect to fluctuations of the carbon source.Although the main target of the aforementioned work was not PHA itself, the applied method can be a suitable tool for the construction of a number of PHA-related metabolic models.
P. putida is a well-described PHA producer.This microorganism (strain KT2440) was studied on the level of genome-scale reconstruction and metabolic network analysis performed by Puchałka et al. 101 .A genome-scale constraint-based model of the metabolism was developed where network reconstruction and flux balance analysis (FBA) have enabled the pin-pointing of essential metabolic functions, structure of the metabolic network, identification of knowledge gaps, and refinement of gene annotations.FBA and flux variability analysis (FVA) have been used for analyzing the properties, potential, and model limits, allowing, under various conditions, the identification of metabolic key features such as: gene essentiality, model and network robustness, growth and yield.The validation of the model was performed using data from continuous cell cultures, C 13 -measurement of internal flux distributions, high-throughput phenotyping data, and gene knockout studies.Interestingly, it was revealed that biomass composition had negligible influence on the accuracy of predictions, whereas the structure of metabolic network was the main factor determining this category.The model network accounts for 877 reactions with 886 metabolites (824 intracellular and 62 extracellular), where 6 % of the reactions in the network are non-gene-associated, and 821 (94 %) reactions have at least one assigned gene (delineated in the gene-protein-reaction (GPR) relationships).The mentioned model was used to improve metabolic engineering strategies in production of PHAs.One of the strategies frequently exploited in the past was the increasing of intracellular availability of the central PHA precursor AcCoA.According to the authors, two methods were tested by modeling to increase the cellular AcCoA pool.The first was the maximization of AcCoA production by enhanced action of pyruvate dehydrogenase (PDH).In the second "an auxiliary reaction was introduced that consumed AcCoA (concurrently producing CoA, to avoid cofactor cycling artifacts) and that would represent the pooling of AcCoA".The value of 'AcCoA production' predicted by the first method includes AcCoA that is consumed subsequently in other reactions (e.g., biomass production).Contrarily, "the value of 'AcCoA pooling' predicted by the second method includes only Ac-CoA that is taken completely out of the system, and is therefore made available for PHA production but unusable for growth or other purposes".For further information, the reader should consult the cited reference 101 .
The same strain was also investigated by Sohn et al. 102 .Various characteristics of P. putida KT2440, such as capacity for PHA synthesis and degradation of aromatics were investigated using the constructed genome scale metabolic model (PpuMBEL1071).Although characterized as a strict aerobic bacterium, the physiological conditions required to ensure anaerobic survival were investigated.An extended survival under anaerobic stress was achieved by introducing the ackA gene from Pseudomonas aeruginosa and E. coli.In the model network, 1071 reactions can be found (among them 958 metabolic and 113 transport reactions) represented by 900 genes (percentage of genes represented =16.6), and 1044 metabolites.Metabolic reactions in the PpuMBEL1071 model were categorized into 11 functional groups, where the largest group of reactions (ca. 25 %) belonged to the synthesis and degradation of amino acids, the next largest group of reactions belonged to the carbohydrate metabolism (16 %), and 10 % related to more than 300 proteins (transporters).A detailed description of the process of the genome-scale metabolic modeling and in silico flux analysis can be found elsewhere 103 .Linear programming, mass conservation, and biochemical thermodynamics were needed to determine the fluxes.The fluxes of irreversible reactions were considered positive, and the negative flux was used for the opposite direction of the reactions.
A series of articles [104][105][106][107][108] published by Penloglou and coworkers, and Chadzidoukas et al. was dedicated to the modeling of PHA production by Alcaligenes latus (today: Azohydromonas lata).A mathematical model that incorporates the metabolic regulation was used to predict the simultaneous bacterial growth and polymer accumulation, respecting the polymerization/kinetic mechanism for build-up of the polymer chains.The applied framework was able to predict the dynamic rate of carbon source consumption, the concentration and the molar mass of the produced PHB, and the molar mass distribution throughout the culture time.Additionally, the model can predict the effect of different fermentation strategies (e.g., batch or fed-batch conditions) on productivity, product quality, polymer concentration and the PHB molar mass distribution.This model can be characterized as a dynamic metabolic/polymerization kinetic model that integrates two different sub-models: a length/time scales polymerization kinetic sub-model, and a metabolic sub-model.The connecting point between these two sub-models is the concentration of the active substrate of PHA synthase, namely, 3-hydroxybutyryl-CoA.The following assumptions were used: -Total biomass is considered to be structured in three compartments: the residual biomass (X), the monomer (3HB), and PHB (P); -The cultivation medium is chemically defined, all nutrients are available in excess except for carbon and nitrogen sources; -The effects of population heterogeneity and mass-transfer phenomena are negligible.
Furthermore, the macroscopic sub-model for bioreactor space was connected to the metabolic and polymerization model 106 , and used for the investigation of tailor-made production of PHAs.
PHAs are produced by a great number of natural microbial species.Some of them are genetically changed in order to produce more PHA or in order to be able to use other substrates.By gene manipulations and plasmid transfers, some host strains that do not naturally produce PHAs (or they produce insignificant quantities) have been constructed to be powerful PHAs producers.An example of this is recombinant E. coli DH5a used by Carlson and coworkers 109 .By transforming DH5a strain with two plasmids, two different strains were constructed.Plasmid pPT500 was a carrier of the native, threegene R. eutropha (today: C. necator) PHB operon from pAeT41 ligated into the pCR Blunt vector.The operon was expressed efficiently in E. coli by the native R. eutropha operon promoter, which resembles the E. coli (σ 70 ) system.This E. coli strain was referred to as PHB (+).The second plasmid, pCR-KT, was constructed by ligating the R. eutropha b-ketothiolase gene, without a promoter, into the pCR Blunt vector.The control strain harboring this plasmid expressed no PHB genes and was referred to as PHB(-).A theoretical biochemical network model was used to interpret the experimental results of the product secretion profiles, PHAs synthesis, as well as for study of E. coli network capabilities under anaerobic conditions.The latter were analyzed by the elementary mode model described earlier 110,111 ; for details, the reader is referred to these works.To conceive the main idea of this work, some modifications were added to the model, concerning acetate, ethanol, lactate, succinate, and glucose as potential substrates, the designation of some transport reactions from irreversible to reversible, and the designation of interconversion of acetate and acetyl coenzyme A (CoA) to be reversible.
The elementary mode analysis was run by the program METATOOL, and the Euclidean distances were used as the objective tool.Elementary mode analysis was used to determine the network model capability of predicting the anaerobic PHB production.Two hundred and two anaerobic PHB-producing modes were identified, among which ninety-eight modes made PHB without coproducing biomass.The anaerobic synthesis of PHB was demonstrated under both growth and non-growth conditions.One of the examined recombinant pathways had a significant effect on by-product secretion patterns.
One of the best examined PHA producers is R. eutropha (today: C. necator).Its native (wild-type) strain (H16) is able to use fructose and N-acetyl-glucosamine (NAG), but not glucose 112 .Some mutants developed from this wild-type strain (i.e., G +1 and H1G +3 ) have been found capable of utilizing glucose 113 .In the past, this wild-type strain and selected mutants were the target of many biological and technological investigations, so the related metabolism was an object of mathematical modeling.Based on the annotated genome, Park et al. 114 developed a genome-scale metabolic model of R. eutropha H16.The genome of R. eutropha H16 consists of two chromosomes and one plasmid, fully sequenced and annotated: chromosome1 (4,052,032 bp), chromosome 2 (2,912,490 bp), and mega plasmid pHG1 (452,156 bp).The stoichiometric part of the model, RehMBEL1391, related to 1391 reactions (inclusive 229 transport reactions) and 1171 metabolites.To refine the genome-scale metabolic model and perform its validation under environmental and genetic perturbations, the constraints-based flux analysis was applied.According to the authors, "the metabolic flux distribution and the changes of metabolic fluxes under several perturbed conditions were examined by flux variability analysis (FVA) that calculates minimal and maximal flux values of each reaction for an objective function of maximum cell growth rate".These calculated fluxes were compared with the control values in order to examine the changes of flux solution space for each reaction under genetic and environmental perturbations.In order to identify the gene knockouts suitable for the enhanced production of 2-methylcitric acid, the method of minimization of metabolic adjustment (MOMA) was applied using quadratic programming (QP).This way, the lithoautotrophic growth characteristics were investigated by varying the feeding ratios of gas (CO 2 and H 2 ) mixture.Additionally, the developed model was used for the research of strategies for PHB production under different cultivation conditions (pH-value and carbon/nitrogen source uptake ratios).Furthermore, the metabolic characteristics of the strain were analyzed under phosphofructokinase gene expression.Finally, "in silico" gene knockout simulations were used "to identify targets for metabolic engineering essential for the production of 2-methylcitric acid in R. eutropha (today: C. necator) H16.
Lopar and associated colleagues analyzed the five-step continuous production of PHB on glucose 115 as well as glycerol metabolism pathways by C. necator DSM 545 116 , respectively.Elementary flux modes, yield space analysis and a high-structured metabolic model were used as tools.Established was a high-structured, metabolic, kinetic type model consisting of 43 mass balance equations related to the intracellular compounds 115 .The metabolic state of cells cultivated in a continuously operated five-stage bioreactor cascade was analyzed with the help of elementary flux modes and two-dimensional yield space.Stoichiometric matrix and elementary modes were calculated using the program Metatool, version 5.1 originally established by Pfeiffer et al. 117 and further developed by Kamp and Schuster 118 .The metabolic yield analysis and weighting factors were obtained using Matlab (inclusive function "fmincon").Two different C-source feeding strategies were performed.Regarding PHB and biomass yields, the values of the more efficient feeding strategy were used as the data source for elementary modes and metabolic flux calculations, respectively.Metabolic fluxes were calculated from experimental yield data using a combination of elementary modes.For that purpose, the quadratic programming approach 119 was applied, where the sum of squared weighting factors was minimized.Excellent agreement of metabolic model with experiments was achieved for the growth-associated PHA synthesis phase of cultivation.Regarding biomass and PHA-related experimental yields, the "in silico" calculated data presented in yield space were in satisfactory agreement with the real metabolic situation along the 5-step reactor cascade.According to the authors, this was the first time that a multi-step continuous process of PHA biosynthesis was analyzed combining elementary flux modes, yield space analysis and high-structured metabolic model tools.
For analyzing metabolic pathways of C. necator DSM 545 when growing on glycerol, a model with a metabolic network consisting of 48 reactions was established in order to describe intracellular processes and PHB production 116 .Special attention was devoted to the branching points in sugar metabolism pathways, direction of reversible reactions, fluxes in gluconeogenesis route, and to the possible coenzymes (NAD or FAD) involved in the glycerol-3-P-dehydrogenase (GLY-3-P DH) reaction.Four sets of elementary modes were obtained, depending on whether the pair NAD/NADH or FAD/FADH 2 contributes to the reaction of GLY-3-P DH, and whether 6-phosphogluconate dehydrogenase (6-PG DH) is present or not.The established metabolic network and the related system of equations provided multiple solutions for the simultaneous synthesis of PHB and biomass.The number of solutions was increased further if either NAD/NADH or FAD/ FADH 2 were assumed to contribute in the reaction of GLY-3-P DH.As a major outcome, it was demon-strated that experimentally determined yields for biomass and PHB (with respect to glycerol) fit well to the values obtained "in silico" if assuming that the Entner-Doudoroff pathway (ED) dominates over the glycolytic pathway.Interestingly, this was also the case if the Embden-Meyerhof-Parnas pathway dominates over the ED.Such results suggest that the producing strain possessed multiple metabolic possibilities for glycerol uptake; the crucial points for further investigation of metabolic pathways were pointed out by the authors.
Apart from pure cultures, mixed microbial cultures have also been applied for PHA biosynthesis.For this purpose, Pardelha and coworkers 120 cultivated microbial consortia in broths based on cane molasses and volatile fatty acids (VFA).Metabolic modeling principles described for the MMC (originated from wastewater activated sludge) were also used in this case.The authors presented a dynamic metabolic model describing the uptake of complex mixtures of VFA.This model was outsourced from a previously published flux balance analysis model 121 that identified the minimization of TCA cycle activity as the key metabolic objective to predict PHA storage capacity.The model was calibrated by experimental data obtained for synthetic mixtures of volatile fatty acids or by data obtained from sugarcane molasses experiments.The MMC was selected using sugarcane molasses under feast and famine cultivation regime.This model gave an excellent fit between experimental and computed values (regression coefficients above 0.92).The introduced C-source consumption regulatory factor reflected well the decrease in acetyl-and propionyl-CoA, confirming the hypothesis that the minimization of TCA cycle is a key objective for the maximization of PHA production by MMC in feast and famine regimen.
Grousseau et al. 122 investigated the role of NA-DPH i.e., how growth controls the NADPH generation and availability, as well as the resulting impact on PHB when C. necator DSM 545 (H1 G+3 mutant) grows on butyric acid.Well-controlled phosphorus limited fed-batch cultures were carried out under non-inhibitory conditions with butyric acid as substrate.Kinetic and stoichiometric models were constructed based on known metabolic pathways, bibliographic studies and available databases (KEGGS and others).The metabolic descriptor was composed of two equation systems: (1) An anabolic network (53 reactions) for the synthesis of macromolecular cell components from 28 intracellular metabolites used to calculate the global molar stoichiometry of biomass and the mass fraction of its components; (2) A catabolic network (45 reactions, 46 intermediates) for carbon assimilation, inclusive the cen-tral metabolic pathways i.e., glyoxylate shunt, TCA cycle, gluconeogenesis, PHB synthesis.
The experiments and modeling were used to explore how microbial growth (controlled by phosphate availability) influences NADPH generation and availability, as well as the impacts on PHB yields.
It was assumed that butyric acid in its non-ionized form penetrates into the cells, where it is activated by the actions of butyrate kinase and phosphate butyryltransferase, after which butyryl-CoA is converted into acetoacetyl-CoA via oxidation, and then decomposed into acetyl-CoA (by ketothiolase) to produce biomass and energy or to be incorporated into PHB .
Modeling (i.e., a mathematical solution) was performed by solving the complete descriptor system of linear equations on the basis of the mass balance under pseudo-steady-state assumption.The model was used in two ways: (a) to simulate the theoretical yield that could be reached, (b) to calculate the intracellular flux based on experimental results.The results indicated that the anabolic demand allowed the NADPH production by Entner-Doudoroff (ED) pathway with maximal theoretical PHB production yield of 0.89 C mol C mol -1 ; further, without residual biomass production, NADPH regeneration is only possible via the isocitrate dehydrogenase from TCA cycle with a theoretical yield of 0.67 C mol C mol -1 .Furthermore, it was found that the maximum specific rate of NADPH production at maximal growth rate (for biomass synthesis requirements) was maximal in all conditions.This by consequence determines the maximal PHB production rate.
These results lead to the conclusion that "sustaining a controlled residual growth improves the PHB specific production rate without altering production yield".
A very interesting modeling paper was recently published by Mavaddat and coworkers 123 .A 3D simulation was made by using a computational fluid dynamics software package (FLUENT 6.3.26) that was combined with a metabolic model.The target was to investigate hydrodynamics and production of PHB in an airlift bioreactor 124 .An Eulerian approach was applied to model the gas-liquid interactions.This system, characterized by the combined effect of bubble breakup and coalescence, was solved by a population balance model implemented in the software.Biosynthesis of PHB was determined by the maximal forward reaction rates for thiolase, reductase, and synthase.[106][107] Azohydromonas lata DSM 1123 (metabolic and polymerisation model) See Ref. [102][103][104][105][106] See Ref. [102][103][104][105][106] See Ref. Glucose PHB LEGEND: X, biomass (residual), (MCA) 94,95 metabolic control analysis; (FBA) flux balance analysis, (FVA) flux variability analysis, (IFD) internal flux distribution, (MILP) mixed integer linear programming, (FA) flux analysis, (EFM) elementary flux modes, (YSA) yield space analysis, (QP) quadratic programing, (MMC) mixed microbial culture, (CIM) consistency index method, (NLWL) nonlinear weighted least squares, (CFD) computational fluid dynamic Also, molar concentration profiles of PHB and glucose within the bioreactor were obtained.
From the presented models, the reader can easily conclude that very different modeling approaches and a very broad spectrum of methods has been applied to achieve the desired target.Some important characteristic properties of the models described in this section are listed in Table 4.

Cybernetic models in PHA biosynthesis
As previously stated, with the present state of knowledge it is practically impossible in kinetic modeling to have all the data necessary to compile a fully developed, and, from a regulatory point of view, exact, perfect, mechanistic model.That is why various kinetic and mechanistic simplifications are introduced in models.Another way to involve regulatory mechanisms related to cell metabolism is to apply cybernetic principles of modeling [125][126][127][128][129] .Thus, cybernetic variables introduced into a kinetic model substitute the unknown mechanistic factors of the cell regulatory system; they represent the different cellular processes: induction, repression, inhibition, and activation.To that end, an OF must be created, assuming that the cell metabolism operates with a physiologically specific goal (optimal growth or maximal product yield, for instance).The optimality hypothesis describes that cells regulate the synthesis and activity of enzymes in such a way that a nutritional objective (the goal) is achieved in an optimal manner.During model development, these goal-functions can be optimized.These facts are summarized comprehensively in a review by Patnaik 130 .Hatzimanikatis et al. 131 presented yet another way of including regulatory aspects in a mathematical model: discrete variables were used to represent the presence or absence of possible regulatory loops.Thus, the model is formulated as a "mixed integer linear programming (MILP) optimization problem".Using a (log)linear approximation for a de facto "non-linear model", the authors verified a qualitative agreement between the model predictions and experimental results.In cases when mechanistic details are not attainable, fuzzy logic-based models 132,133 or neural networks 134,135 can be used to achieve simulation of metabolic system.However, to generate these models, such an approach requires an extremely high quantity of experimental data.
In this section, the cybernetic models dealing with PHAs are discussed.Neural networks will be a separate part in this review.
One of the first cybernetic works dealing with PHAs was the investigation of Yoo and Kim 136 .These authors assumed that the pathway of PHB synthesis exhibits the transcriptional control induced by environmental stress, i.e., nitrogen limitation.As the objective for this control, an optimization of acetyl-CoA utilization was chosen, thus ensuring the necessary degree of metabolic flexibility in cells catabolism.The state equation related to key protein synthesis was assumed dependent on the nonlinear control variable.Yoo and Kim 136 built their model based on two assumptions, as follows: (i) Cells may be divided into two components/ compartments, i.e., the residual cell mass (non-PHB part) and PHB (storage material); (ii) C-source is allocated to the enzyme synthesis system respecting the principle of considerable flexibility under nitrogen limitation.
Although the application of economic laws in biotechnology is a rare case, very surprisingly, Herrnstein's law from macroeconomics field 137 has found its place in cybernetic modeling.That means that the conventional cybernetic approach is fully adopted applying the rule: "when many resources are to be allocated to many activities (not necessarily an equal number), the fractional allocations are proportional to the fractional returns".Yoo and Kim 136 did not use the above approach, because the singular optimality criteria applied by these authors "make the model too stiff to integrate numerically".They proposed a nonsingular strategy and compared their predictions with those of Asenjo and Suk 35 and Mulchandani et al. 36 as well as with own experimental data.Such organization of the model was able to predict the "mixed-growth-associated synthesis of PHB" for the overall fermentation range.The described model was compared with earlier reported unstructured models by statistical tests.
It is also worthwhile mentioning the model proposed by Varner and Ramkrishna 129 , developed for the synthesis of storage polymers.This is a generic model with assumed storage pathways comprising two substrates and seven enzymes.The authors concluded that the nutrient environment is a determining factor for the pathway performance.To improve productivity, they investigated genetic alterations suggested by the model simulations.When the cybernetic criteria of over-expressing a key enzyme that catalyzes PHB synthesis and blocking of two unproductive pathways were used, a PHB yield of 90 % of the total biomass was predicted.This is significantly higher than the experimental values of 65-75 % obtained by conventional fed-batch and continuous cultivations.It is unclear if this discrepancy is the consequence of model weakness or if it is due to a "steric inhibition factor" in the cells that was not considered in the model.
Ferraz and coworkers 138 established a cybernetic model for PHA synthesis by A. eutrophus (today: C. necator) DSM 545, containing 53 parame-ters (among them 14 fixed and 39 adjustable) and 25 equations.A preliminary estimation of adjustable model parameters was achieved using model linearizations and simplifications as well as the treated experimental data.Thereupon, the model parameters were estimated by applying the flexible geometric simplex method proposed by Nelder and Mead 139 , known as the flexible polyhedron search.For a more effective application of this technique, the results of the preliminary search were used as the initial guess.The model considers the feedings, sample withdrawals and evaporation influences taking into account the phenomena related to the kinetics of active biomass growth and PHBV production.According to the authors, the proposed model adequately represents the characteristic phenomenology of the two-phase cultivation process; nevertheless, this model was unable to reflect well all the data and transient states achieved under tested conditions.
Gadkar et al. 140 used a similar modelling approach to develop their cybernetic model for R. eutropha (today: C. necator).They established the metabolic model that contains four main metabolic parts (pathways): glycolytic and pentose phosphate pathways for sugar metabolism, the PHB synthesis, the PHB degradation (depolymerisation) and the TCA cycle (biomass precursor's synthesis).It was proposed that the pathway of PHB degradation (i.e., depolymerisation) becomes active when the C-source (glucose) concentration drops and starts to be growth-limiting.Gadkar et al. 140 validated their model with those proposed by Yoo and Kim 136 .The targets of modelling were the multi-rate model predictive control algorithm, a continuous flow bioreactor with cell recycling, PHB productivity and the permeate flux through a micro-filter acting as cells separator.
The formerly described model by Ferraz et al. 138 was further elaborated and improved by Riascos and Pinto 141 .The authors proposed the use of orthogonal collocation to a simultaneous optimization approach.Initially, this methodology was implemented to the model of penicillin biosynthesis that was simplified to some extent; thereafter it was applied to a cybernetic structured model for PHA synthesis 138 .The model was intended for optimizing PHB and PHBV production by A. eutrophus (today: C. necator).The model is characterized by 11 states, 2 -4 control profiles, 3 C-sources (glucose, fructose, propionic acid) and one nitrogen source (NH 4 OH), with all feed rates selected as control variables.This time, the optimization problem was treated as the optimal control problem (OCP).In order to solve the singularity problem, the authors discretized the optimization problem on finite elements by orthogonal collocation.The results of this inves-tigation indicated that "the discretization of differential-algebraic equation systems (DAE) by orthogonal collocation in finite elements efficiently trans forms dynamic optimization problems into non linear programming (NLP) problems".Thus, the solution of complex problems with several control variables satisfied the approximation error tolerance 141 .According to Patnaik 48 ,there is evidence that the results achieved by the model described above are "inferior" to those presented by Gadkar et al. 140 .If so, the model complexity did not sufficiently point out the expected positive contributions.Contrarily, it seems that a counter-effect was generated.
The published cybernetic model of continuous PHB production by Alcaligenes eutrophus 141 (today C. necator) was extended to an industrial two-stage continuous process 142 .This model incorporated three main metabolic pathways -the glycolysis, pentose phosphate, and the TCA cycle pathwaysupposing that glucose and NH 4 + are converted by the cells into precursors required for cell growth and storage materials.The formation of PHB and its degradation pathways were also involved.Regarding cybernetic variables, two levels were definedthe "local" and the "global" level.Regarding the local cybernetic variables, they considered "the individual strands of the metabolic pathways (production of growth precursors, production of PHB etc.)".Furthermore, on the global level the competition between the different directions of metabolic pathways was related, such as between the production of precursors for the cell growth or for the PHB synthesis.Such a solution resulted in diverting a larger amount of precursor resources towards PHB synthesis, when excess glucose was present.Similar, for the C-limiting growth phase, the competition between the production of metabolites from glucose on the one side, and PHB degradation on the other side was taken into account.The technological process in this study was performed in two reactors in series, followed by a separation stage.In the first reactor, the bacteria were forced to grow on glucose and NH 4 + as the nitrogen source under non-limiting conditions.The output from this reactor was connected to the second reactor where only glucose (without nitrogen source) was fed.Any NH 4 + in the second reactor originated only from the unassimilated part that had left the first reactor, resulting in the exposure of bacteria to excess glucose, preferentially channeled into PHB.The output from the second reactor was subjected to a separation procedure where PHB was extracted from biomass.Furthermore, the remaining biomass residuals underwent degradation to convertible substrates that were subsequently supplied as feed streams to both reactors.A bifurcation analysis was performed with the help of the described cybernetic model as a useful tool to study the dynamics of the process.This allowed insight into the system features (i.e., steady states; limit cycles), and was very useful for an efficient reactor design and process control.
Cybernetic modeling was combined with elementary flux modes resulting in the hybrid model described by Kim et al. 143 .Such method was also applied for modeling of the metabolism of R. eutropha (today: C. necator) 144 concerning PHAs synthesis on carbohydrates.The latter is characterized by the strong participation of cybernetic part in the whole modeling system.As this model is a typical hybrid model, it will be considered later under the separate section "hybrid models".Cybernetic models are based on concepts of natural evolution.They are more feasible, and give more physiological fidelity than mechanistic models.Nevertheless, cybernetic models display at least two disadvantages.One is certainly their complexity that can be considered as the consequence of representation of included intracellular regulatory systems and across the cell walls active transport systems.
Additional problems are, for a given situation, possible multiple cybernetic goals.It is often unclear which of these goals is the most relevant for the given situation.When mixed microbial cultures or/and multiple substrates are involved, a conflict of goals is possible, and consequently, a conflict of results can occur.

Neural networks and hybrid models in modelling of PHAs biosynthesis
Neural networks have been applied in modelling of PHAs biosynthesis in order to overcome limitations of mechanistic and cybernetic models.For such models, a mathematical description of the biological/technological system is not obligatory.This is a certain benefit for industrial level of applications, but such models must be carefully designed and "trained" with full attention 145,146 .It is not a rare case that more than one neural network is pointed out as the best solution.To choose correctly the truly best one, an additional effort in the form of compiling hybrid models can be helpful.In order to achieve maximal effect, neural networks and partial (or fully developed) mathematical models of mechanistic or/and cybernetic type are combined in such models.
Neural networks and hybrid models in modelling of PHA biosynthesis by pure microbial cultures One of the first neural network/hybrid models dealing with PHA biosynthesis was developed by Peres et al. 147 .The method combines "first principles models" with modular artificial neural networks trained by the EM (Expectation-Maximization) algorithm (for details the reader should consult the works of Jordan and Jacobs 148 and Xu et a1. 149).The modular networks were used in order to treat the 'cells system' as a "highly complex network of metabolic reactions organized in modular pathways".After validation of the model by the experimental, laboratory scale, production process, the results indicated that the modular network, if trained with the EM algorithm, was able to organize itself in modules related to the basic biological pathways.In particular, the network was able to learn to discriminate between acetate-and internal reserves respiration.That meant that the network had detected the switch between the presence and absence of acetate considering the transient phase from normal extracellular acetate respiration (PHB storage) and internal storage consumption, i.e., PHB degradation.
Patnaik 150 performed the neural network modelling to enhance the production of PHB in an ideally mixed fed-batch bioreactor.The process was represented by a neural network and controlled by another network.A feed-forward artificial neural network (ANN) was used as a feedback controller, and seven different neural networks were used and studied for the process.The pair of networks (feed controller/ bioreactor) was applied with the intention to ensure the stepwise control of feed rates in time (i.e., solutions for the carbon-glucose and NH 4 + ), in order to maximize the concentration of PHB at the end of each time interval.Dispersion, in all cases represented by a Pecklet (Bodenstein) number ≈ 20, was found to be optimal.All seven ANN that represent the fermentation can be classified into one of three types: feed-forward (FF) trained by the back-propagation (BP) algorithm; radial basis (RB) i.e., generalized regression type and recurrent (Elman and Hopfield type).The achieved results suggested that enhanced final PHB concentration could be reached (ranging from 16 to 93 %).Regardless of the choice, the application of neural networks significantly improves the productivity under simulated industrial conditions.A further work was performed by the same author 151 where two neural networks had been employed to enhance PHB production in fed-batch system by R. eutropha (today: C. necator).One network was adapted to filter the noise in the two feed streams (C, N), whereas the second was foreseen to control their flow rates at the optimum dispersion.This arrangement had more than doubled the PHB concentration per reactor volume if compared to an ideal fermentation.While the unfiltered noise reduced biomass growth (20-25 %) and PHB formation (25-30 %), the used neural filter reversed this trend.In addition, in simulated experiments 152 , Gaussian noise was added to the flow rates of both C-and N-source.Lyapunov exponents were used for comparison of different noise filters influence (auto-associative neural, cumulative sum control chart-CUSUM, extended Kalman, Butterworth filter, respectively) versus noise-free result.Negative exponents indicated a stable fermentation.An auto-associative neural filter performed best, followed by the combination of a CUSUM filter and an extended Kalman filter, whereas Butterworth filters were inferior and inadequate.The same author tested three fundamental hybrid designs by comparison of the neural and mechanistic optimizations done for a fed-batch bioreactor 153 .In his study, the author proposed the general neural hybrid model.The Gaussian noise was used to challenge the non-ideal features in the two feed streams for glucose and NH 4 + ; the optimal finite dispersion of substrates was assumed and incorporated in the models.The one purely neural model and three hybrid designs were compared with mechanistic kinetic models.All three hybrid designs were superior to the neural and mechanistic approaches.Models that contained both mechanistic and neural components, generated higher biomass and PHB concentrations if compared with "pure" neural system.Among these, the model in which concentrations of C and N were represented by neural networks were found less suitable, whereas the model containing neural networks for the biomass and PHB performed superior to the first.Finally, the third, involving the hybrid representation with an optimally weighted combination of mathematical and neural forms for the kinetics, was found to outperform the other modelling approaches the best; it generated 140 % more biomass and 330 % more PHB than obtained after conventional mechanistic optimization.
Franz and associates 144 published the hybrid cybernetic modeling approach to describe PHB biosynthesis and its depolymerisation with simultaneous consumption in R. eutropha (today: C. necator).This mathematical model allows a systematic derivation of the model equations from elementary mode analysis.It was assumed that internal metabolites are in quasi-steady state, and PHB was considered as internal metabolite with slow dynamics.The vector of fluxes through EFMs was set to be liable to the control of vector of cybernetic variables.All relevant enzymes were divided into two sub-groups: inducible and constitutive.The vector of inducible enzyme synthesis rates was assumed to be controlled by the vector of separate cybernetic variables.In addition, it was adapted by an OF, covering the assumption of maximum carbon source uptake by the organisms (again, the model was assumed to be two-compartmental (PHB and "residu-al" non-PHB biomass)).Metabolic yield analysis was applied in order to reduce the number of EFMs.Due to the stoichiometry of the network, the concentration of the residual biomass (non-PHB part of cells) was linearly correlated, inter alia, with the metabolized nitrogen source; therefore, the yield space was reduced to two-dimensional.This approach was extended for non-quasi-stationary metabolites, i.e., PHB.The model was validated by separate experiments.Five different metabolic model configurations were tested concerning ratio of metabolic flux to biomass or PHB synthesis.Nonlinear analysis of continuous cultures was performed with the intention to identify nonlinear phenomena: oscillations 154 and multiple steady states 155 .The model showed good agreement with experimental data concerning PHB synthesis and its intracellular degradation.Multiple steady states were identified in the continuous bioreactor for a certain range of dilution rates/substrate loading.It was underlined that the predicted multiplicity region is of relatively small range, so that multiple steady states are realistically not probable in the practical work.
Recently, Zafar and coworkers published a series of articles [156][157][158] dealing with optimization of PHB and PHBV production by A. lata MTCC 2311 on cane molasses and on a mixture of cane molasses with propionic acid.In the first step, maximizing biomass and PHB production was targeted.During optimization of broth components in shaking flask experiments, among three carbon sources (sucrose, fructose, and glucose) and four nitrogen sources [(NH 4 ) 2 SO 4 , NH 4 Cl, urea, and NH 4 NO 3 ], sucrose and urea were found to be the best carbon and nitrogen sources, respectively.Further, the effects of sucrose, urea, and trace elements solution on biomass and PHB concentrations were investigated.Response surface methodology (RSM) and artificial neural network models (ANN) were applied to direct the experimental data obtained in accordance with the central composite design.The relationships among the screened variables were expressed in the form of quadratic polynomial function.All experiments were conducted in the central composite design (CCD) procedure, fitted by an empirical fully developed second-order polynomial model, representing in the form of surface response over the variables range.Each variable in the CCD was studied at three different levels; hence, a 2 3 -factorial design was applied.With three variables, 20 experimental sets were spawned and used for optimization of sucrose, urea, and trace elements solution concentrations.Two ANN models were used for modelling the influence of the mentioned substrate concentrations on biomass and PHB concentrations.A "feed-forward" architecture of ANN model (multilayer perception; MLP) was combined with propagation algorithm (BP) to reach the predictive models.Three medium concentrations acted as input variables, in contrary, biomass and PHB concentrations were outputs.Root mean square error (RMSE), model predictive error (MPE) and standard error of prediction (SEP), were calculated (respecting the bias (B f ) and accuracy (A f ) factors) to reach the desired fitting/prediction accuracy of models.After good prediction accuracy was reached by RSM and ANN based models, the genetic algorithm (GA) was applied to optimize the input space with an objective to maximize the process performances (in this case biomass and PHB concentrations).RSM and ANN models (combined with CCD) were successful in modelling/studying of the input variables interactions.The modeling ability and the optimization power of hybrid ANN-GA model provided higher accuracy (considering optimal substrate concentrations) than hybrid RSM-GA model.The experience collected by the work described above was applied by the same authors for the optimization of PHBV production through the simulations performed with artificial neural network (ANN), response surface methodology (RSM) and (GA).The predictions obtained by ANN were in good agreement with experimental findings and better than those achieved by RSM.Furthermore, the cited authors applied the same methodology for modeling and optimization of PHBV production using the same microorganism and cane molasses supplemented with propionic acid.The experiments were carried out in a stirred-tank reactor, conducted in a three-level factorial design by varying the impeller and aeration rates to investigate the effect of agitation and aeration regimes.Further, the data were fitted to quadratic polynomial equation and ANN; process variables were optimized by GA-coupled models.ANN and hybrid ANN-GA were found superior for modeling and optimization of process variables, respectively.

Neural networks and hybrid models in modelling of PHA biosynthesis by mixed microbial cultures
In the previously mentioned work, Dias et al. 87 studied the optimization of PHA production by mixed microbial cultures.This work was based on a detailed hybrid metabolic model.The metabolic network from this model was, by applying elementary flux modes (EFM), first decomposed into its fundamental pathways.The EFMs are the minimal set of reactions, able to define the metabolism of cells obtained using the FluxAnalyzer software.A step after that, the dynamical hybrid semi-parametric model was formulated.The EFM fluxes were interpreted in terms of metabolic consistency.The final model allowing characterization of the metab-olism dynamics and metabolic pathway analysis, showed that the propionate metabolism is more flexible (redundant) than the acetate metabolism.Using the classical macroscopic approach, four reactions were obtained for the "feast" phase: 3HB and 3HV production, cell growth, and maintenance on propionate.The EFM technique identified two additional pathways: simultaneous growth and 3HB production, and simultaneous 3HB and 3HV production.The EFMs obtained in the "famine" phase were in good agreement with the macroscopic reactions.
A neural network based methodology was carried out for the modelling of a sequencing batch reactor (SBR) for PHB production by MMCs 159 .In order to achieve successful process control and optimization, the data driven modelling method through operating regime decomposition, with the bootstrap aggregated neural networks, were used in this study.The target was to enhance model accuracy and its reliability by capitalization of empirical models, developed earlier from process data.When PHB is produced in SBR using MMC, two substrates (acetate and NH 4 + ) in the feed stream have an important, dominant influence on PHB production time course.Different process operation regimes must be applied depending on the concentrations of these substrates.Using bootstrap aggregated neural networks, the authors proposed a method for the classification of such operation regimes, and for building of neural network models in accordance with operational regimes.Nine different regime types were tested.Mathematical equations for the "border curves" that separate the "feast" from the "famine" operation regime area were established with the intention that they will display an important role in process optimization.The neural network model inputs were initial biomass concentration, acetate and NH 4 + feeding concentrations; the model output was the time when PHB production reaches its maximal value.The bootstrap aggregated neural networks were applied to predict the batch completion time.Bootstrap re-sampling was used to generate 40 replications of the model development data, and each replication was partitioned randomly into a "training" and a "testing" data set.These pairs (training / testing data set) were applied to build the neural network model, further, the optimal number of hidden neurons was determined for each neural network (using cross validation, i.e., the lowest sum squared errors, SSE).The combination of individual networks generated an aggregated neural network model (ANNM).The bootstrap aggregated neural network models provide a more robust prediction capability.According to the authors, the optimal SBR cycles should operate near border curve, differentiating the "feast" and "famine" operation re-gime areas.An optimal duration of the "famine" phase requires further investigation because a long "famine" phase positively influenced the stability of the cell culture, however resulting in the loss of some quantity of PHB.Some general comments on the neural networks and related hybrid models: Similar to neural network designs, to date the scientific community has not found a unique procedure (solution) for the targeted, exact, final design of hybrid neural networks models.It seems that bidirectional exchange of solutions and experience between experimental science and modelling performers will be the right way to reach a "satisfactory point".With the present state of knowledge, we are "doomed" to improve the models by experimental results, and contrarily, the simulation results achieved by the models should be "first aid" in the planning and interpreting of experiments.As all other model types, these models are a useful tool for detecting possible methods in investigation procedures.

Conclusion
It is not possible to generally state which model type is the best for the mathematical modelling of different problems connected with PHA production (bioreactor performance, metabolism, genetic engineering, cultivation conditions and kinetics).All the models presented in this review possess some range of validity, and are characterized by some range of adjustability (indicated by so-called "interpolation/ extrapolation power of model"; related to the applied range of cultivation conditions and inevitably connected to the used strain).It is sure that certain types of models are meant to be direct mathematical pictures of the real biological systems, hence "as similar as possible" (e.g., genome-scale metabolic kinetic models).Contrarily, for example, neural network models are constructed without the intention of being "direct pictures" of biological systems, but they are proposed to give simulation results as exact and compliant "as possible" if compared with experimental data.In this review, the presented mathematical models fully meet different functions.They are used for the modeling of processes in bioreactors, mixing, metabolic pathways, cell cycles, genetic purposes, and metabolic engineering, respectively.It seems that for "standard bacterial cultivations and everyday practices", the formal kinetic models (for simple cases) and "low-structured" models will, in most situations, be accurate enough and of great benefit.This is due to their simplicity and relatively low computational demand.In contrast, it seems that complex multi-substrate cultivation media require a higher degree of model struc-turing.In addition, for scientific purposes and advanced development of industrial equipment, real systems (technical or biological) will in the future be modelled by highly organized hybrid models.For example, computational fluid dynamics (CFD) in bioreactors offers many new tools for the computational characterization of fluid streaming and concentration/temperature/momentum fields.Otherwise, genome-scaled models including enzyme kinetic properties are coming ever closer to the real biological "picture" of cell metabolism.Combination of such models in hybrid-type models by implementation of the necessary degree of structuring could be the right way to obtaining a model capable of reflecting the broad spectrum of biological and physical diversities.This is certainly only an example.Cybernetic or neural network models, hybridized with kinetic or genome-scale metabolic models (as well as with CFD models), also have the potential to solve complex biotechnological situations.In general, there is still a lot of work to be done in the field of mathematical modeling dealing with biotechnological processes; in this context, mathematical modeling in the field of PHA only constitutes a prime-example for the variety of bioprocesses.
Ta b l e 3 -Characteristic properties of metabolic models dealing with PHA production by mixed microbial cultures.
Ta b l e 4 -Main characteristics of metabolic models targeted for industrial PHA biosynthesis and metabolism of producers