Overview

modeling overview

Microbial metabolism can be very complex, which makes it difficult to predict the effect and viability of genetic modifications without conducting experiments in the lab. Optimization of microbial production systems, including our use of Methylobacterium extorquens AM1 as a polyhydroxyalkanoate (PHA) factory, can thus be very challenging. To increase the success rate of our interventions, we first improved our understanding of the metabolism of M. extorquens and predicted the most promising optimization strategies for our goal using computational modeling. Using a model that describes the entire metabolism of M. extorquens, we identified the carotenoid biosynthesis pathway as an interesting knockout target, which is responsible for the production of pink pigments. Moreover, we discovered 19 promising gene knockout targets in the carotenoid biosynthesis pathway of which we selected three to be tested in the lab. Ultimately, we were able to obtain a viable knockout strain, with the help of our computational design.

Goals

M. extorquens naturally produces PHA for energy storage. This is advantageous as it allows the organism to adapt to the nutrient availability in its environment1. However, for industrial-scale biotechnological applications, the production of PHA by M. extorquens needs to be optimized to become more cost-efficient. Compared to model organisms like Escherichia coli, M. extorquens is relatively unknown, although it is the most studied methylotrophic bacterium2. When using a non-model organism, computational approaches can provide valuable additive information on the efficacy and viability of genetic modifications or metabolic engineering strategies. Therefore, we decided to use a computational approach for designing a promising metabolic approach to enhance PHA production. Predicting these effects in silico, before testing them experimentally in the lab, can save time and resources and may help elucidate counter-intuitive metabolic engineering approaches.

To this aim, we used the genome-scale metabolic network model of M. extorquens, which describes its entire metabolism3. Using this model, we had two main goals. Firstly, we aimed to identify metabolic targets that need to be up- or downregulated to maximize PHA production. This showed us that the carotenoid biosynthesis pathway is an interesting knockout target, as downregulation of this pathway increased PHA production and because carotenoid production is most likely redundant in our PHA overproducing strain. Based on the results of this first iteration, we then aimed to predict the effect of specific carotenoid biosynthesis pathway gene knockouts on the carbon flux to the PHA production pathway. This was done to select promising gene knockout targets to be implemented into our experimental design .

Expert advice

Dr. Erika Tsingos, an expert in the field of metabolic network modeling, helped us set up our modeling approach. She advised us on implementing biological assumptions into the modeling constraints to mimic culturing conditions. For the knockout simulation, she recommended considering the effect on both PHA production and the growth rate to determine the viability of the intervention. Furthermore, she advised us to compare the fold change increase or decrease in reaction fluxes of the optimized strain to the wild type, instead of the ‘raw’ flux data.

Genome-scale metabolic model and data set

    Genome-scale metabolic models (GEMs) are computational mass-balanced reconstructions of the biochemical reactions present in the metabolism of an organism. They contain information on the compartments, reactions, metabolites, enzymes, and genes of a cell, based on genome annotations and experimental data. The integration of multi-omics data, including genomics, metabolomics, transcriptomics, proteomics, and epigenomics, has advanced the development of GEMs, making them a powerful tool in metabolic and biochemical engineering4. GEMs allow us to gain an understanding of cellular metabolism, predict the response of an organism to different environments, and guide the design of mutant strains with optimal features which do not (yet) exist in nature. The COnstraint-Based Reconstruction and Analysis (COBRA) toolbox can be used in MATLAB to analyze GEMs to provide insight into cellular metabolism by simulating the flow of metabolites through the network.

    The GEM of M. extorquens, named iRP911, consists of 1139 reactions and 977 metabolites3. The model is represented as a stoichiometric matrix (S) describing the mass-balanced reactions (n) in the metabolism and the metabolites (m) involved in each reaction (Figure S1). The flow of the metabolites through each reaction, also known as flux, is represented by vector v. The objective function (Z) defines one or multiple reactions that are to be maximized or minimized when solving a constraint-based analysis. It is mathematically represented by the desired flux (v) multiplied by vector c. The solution space is confined by the steady-state assumption (Sv = 0), which means that at a steady-state, the rate of production and consumption of a metabolite must be zero. As this allows for a large number of solutions, the solution space can be limited by applying upper and lower flux bounds (vub and vlb) and objective functions, based on biological assumptions5. The COBRA toolbox relies on this framework and can be used to predict various outcomes, including strain design and knockout simulation.

Fig. S1 | A genome-scale metabolic model (GEM) contains mass-balanced reactions. The mathematical formulation of a GEM consists of stoichiometric matrix S, which describes the reaction (n) and metabolites (m) involved in each reaction. The flux and objective reactions are described by vectors v and c, respectively. Figure adapted from Ng et al5.

Identification of metabolic targets to optimize PHA production

Aim

Identify up- and downregulation targets that improve the production of PHA by M. extorquens.

How

An OptForce analysis was performed on the genome-scale metabolic network model of M. extorquens.

Why

Predicting the effect of different metabolic engineering strategies before testing them experimentally in the lab can save time and may help elucidate counter-intuitive interventions.


PHA is produced by M. extorquens from C1 compounds such as methanol, via acetyl-CoA, which is subsequently used in the poly-β-hydroxybutyrate (PHB) cycle to produce PHB, a type of PHA (Figure 1). Acetyl-CoA and pyruvate, one of the precursors of acetyl-CoA, are key intermediates in the central metabolism of M. extorquens and are converted into a variety of metabolites, including PHA, fatty acids, amino acids, carotenoids, etc. We hypothesized that knocking out a downstream metabolic pathway that competes with the PHB cycle for the use of pyruvate and acetyl-CoA could alleviate the metabolic burden on M. extorquens and increase the availability of these key metabolites to be used in the PHB cycle. This strategy has previously been proposed for other bacteria to improve the production of pyruvate- or acetyl-CoA-derived bioproducts6.

Fig. 1 | Central carbon metabolism of M. extorquens. PHB cycle = Poly(β-HydroxyButyrate) cycle, EMCP = EthylMalonyl-CoA Pathway, TCA cycle = TriCarboxylic Acid cycle. Figure adapted from Schada von Borzyskowski et al7.

Based on this hypothesis, we have computationally identified such competing metabolic pathways that could be beneficial for our overarching aim: improving PHA production. Using the GEM of M. extorquens, we attempted to identify various gene targets which could improve PHA production by being down- or upregulated.

Model set-up and assumptions

Because we want to simulate the flow of metabolites through the metabolic network of M. extorquens as realistically as possible, we constrained the model based on the following biological assumptions:

  1. The bounds of the exchange reactions, the reactions describing the uptake and efflux of compounds, are set according to the medium components used for culturing the bacteria in the lab. In our case, we use a minimal medium, which contains methanol as its sole carbon source, such that M. extorquens feeds on methanol instead of sugar. Other medium components are salts. If an exchange reaction in the model describes the exchange of a compound that is not present in the medium we use, such as glucose or succinate, the bounds resembling the uptake of the compound are set to zero.
  2. To represent the growth of the bacterium, biomass synthesis is selected as the objective function. To simulate the targeted overproduction of PHA, the objective function of the ideal mutant is changed to biomass in combination with the PHA polymerization reaction.

To identify up- and downregulation targets, the tool OptForce was used, which compares the ideal PHA overproducing model with the wild type model.


Learn more about OptForce

    OptForce is an optimization tool from the COBRA toolbox that can be used to optimize microbial metabolic pathways to overproduce a targeted metabolite8. It compares the reaction fluxes of the wild type metabolic network with the reaction fluxes of the ideal mutant strain which targets the production of a desired compound, in our case PHA. Figure S2 shows a schematic representation of the output of the OptForce algorithm. If the flux of a reaction in the wild type model and overproduction model overlap, it means that it is possible to achieve the ideal overproducing strain without any intervention. However, if the flux ranges do not overlap, it means that an intervention is required for the overproduction mutant. In this case, two options are possible: if the range in the overproducing model is greater than the range in the wild type model, this specific reaction should be upregulated. If the range in the overproducing model is smaller than the range in the wild type model, the reaction should be downregulated or knocked out, depending on whether the optimal flux range reaches zero. The OptForce algorithm goes over all reactions in the wild type and mutant model and classifies the reactions in two subsets of reactions that must be upregulated (mustU) or downregulated (mustL). Besides identifying single upregulation or downregulation targets, OptForce determines combinations of two reactions that must both be upregulated (mustUU) or downregulated (mustLL) or one upregulated and one downregulated (mustUL).

Fig. S2 | Overview of must sets, describing the difference between the wild type model and overproducing model. Adapted from Ranganathan et al8.

Results and discussion

The output of the OptForce optimization tool consists of sets of reactions that should either be upregulated, downregulated or knocked out to achieve the desired overproduction of PHA. For our PHA overproduction simulation, a large number of possible downregulation targets was identified. These targets were genes involved in a variety of metabolic pathways, including the starch and sucrose metabolism and the biosynthesis pathways for a variety of amino acids (valine, leucine, isoleucine, phenylalanine, tyrosine, tryptophan, lysine, and histidine), nucleic acids (purine and pyrimidine), glycerophospholipids, and fatty acids. These pathways all branch off of the central carbon metabolism and thereby use carbon that could otherwise be used for other branches, including the PHB cycle. Furthermore, previous research has shown that fatty acid biosynthesis and PHA biosynthesis pathways are competitive in Pseudomonas aeruginosa, which could also be the case in M. extorquens9. Depending on gene essentiality, these reactions could be interesting metabolic engineering targets for further optimization of M. extorquens.

Besides single targets, combinations of two reaction targets were also identified by the OptForce algorithm. One of the interesting findings of this analysis suggests that a combined downregulation of reactions in the carotenoid biosynthesis pathway and reactions lowering 2-oxoglutarate levels may be beneficial to enhance the flux in the PHB cycle. Carotenoids are produced by M. extorquens from pyruvate, as shown in Figure 1. They are responsible for the characteristic pink color of M. extorquens and protect the bacterium from UV damage. The use of M. extorquens as a microbial cell factory for PHA requires a lab environment, where the risk of UV damage is low, which means the production of carotenoids is likely to be redundant and therefore a promising knockout target. 2-oxoglutarate is an important regulator in the tricarboxylic acid (TCA) cycle and its concentrations are known to fluctuate according to nitrogen and carbon availability10. At lower nitrogen concentrations, the intracellular 2-oxoglutarate concentration tends to be higher10. Interestingly, previous research has shown that nitrogen deficiency increases PHA production, which suggests that 2-oxoglutarate levels may indeed influence PHA production11,12. Based on these findings, combining nitrogen-deficient conditions with a carotenoid biosynthesis knockout may favor the production of PHA to an even higher extent than the two interventions can achieve separately.

Implementation

In this first iteration, we discovered that the carotenoid biosynthesis pathway is an interesting knockout target, as it was shown that downregulation of the carotenoid pathway enhanced production. Furthermore, the production of carotenoids is most likely an unnecessary characteristic for our PHA overproducing strain of M. extorquens. We therefore favored this knockout target over other identified targets. However, the carotenoid biosynthesis pathway in M. extorquens comprises over thirty reactions (Figure 2), and this computational approach did not tell us which specific gene knockout is most optimal for eliminating carotenoid production and enhancing PHA production. We therefore included a second modeling approach to specifically simulate knockouts in the carotenoid biosynthesis pathway and predict the effect on PHA production.

Fig. 2 | Carotenoid pathway of M. extorquens.

Effect of different carotenoid gene knockouts on PHA production

Aim

Simulate the effect of various gene knockouts in the carotenoid biosynthesis pathway on the production of PHA by M. extorquens.

How

Gene knockouts were simulated in the GEM of M. extorquens. Flux Balance Analysis was used to determine the effect of the knockout simulations on the reactions in the PHB cycle.

Why

Predicting the effect of knocking out different genes before testing them experimentally in the lab can guide the choice of the best knockout targets, increasing the chance of success.


Our second iteration aims to predict which gene knockout in the carotenoid biosynthesis is most optimal for enhancing PHA production. The genes involved in the carotenoid biosynthesis pathway were knocked out one by one within the genome-scale metabolic network model of M. extorquens to create 34 knockout models. Subsequently, the effect on the reaction fluxes involved in the PHB cycle was predicted using Flux Balance Analysis (FBA). FBA is a method from the COBRA toolbox that mathematically simulates the flow of metabolites throughout a metabolic network. It predicts the most optimal or likely flux at a steady-state (Sv = 0), based on the flux constraints and the objective function. FBA is used to predict growth rate, the internal flow of metabolites, and the exchange rate of compounds with the environment.

Model set-up and assumptions

Similar to the first modeling approach, the model needs to be constrained based on biological assumptions, to simulate the knockouts as realistically as possible.

  1. The exchange fluxes of the knockout simulation model are also constrained based on the medium components.
  2. Knockouts are simulated by constraining both the lower and upper bounds of the knocked-out reaction to zero.
  3. The objective functions of the wild type and knockout simulations combine growth and PHA production.

The fold-change effects of the knocked-out reactions on the biomass synthesis flux and PHB cycle fluxes compared to the wild type were assessed to predict the effect on growth and PHA production, respectively. We looked at five reactions involved in PHA production: the conversion of pyruvate into acetyl-CoA, the three reactions involved in the conversion of acetyl-CoA into PHB, via acetoacetyl-CoA and 3-hydroxybutyryl-CoA, and the depolymerase reaction responsible for breaking down PHB (Figure 1). The genes encoding the enzymes catalyzing the last four reactions are phaA, phaB, phaC and phaZ, respectively.

Results and discussion

Figure 3 shows the predicted fold-change effect of all carotenoid biosynthesis pathway knockouts compared to wild type M. extorquens. The most notable effect was observed for the phaA reaction with a maximum fold-change increase of 76.45, followed by 42.51. An increase in pyruvate to acetyl-CoA conversion was also observed for the same knocked-out reactions, but to a lower extent than the phaA reaction, with the most optimal increase of 1.61 times higher than the wild type flux, followed by 1.44. The phaB reaction was unaffected by all carotenoid biosynthesis pathway knockouts compared to the wild type. Unfortunately, the phaC reaction was negatively affected by the knockouts, although only a small decrease was observed, with the lowest fold-change being 0.91. A similar pattern was observed for the phaZ reaction, but since this gene encodes for the depolymerase of PHB which breaks down the desired product, this decrease in flux is advantageous. Interestingly, 19 knockout simulations showed the same effect. Combining the flux changes, these knockouts are predicted to result in a ~110-fold increase in total PHA flux compared to the wild type, based on multiplying the fold-changes of pyruvate to acetyl-CoA conversion, phaA, phaB and phaC and divided by the fold-change of phaZ. The effect on biomass synthesis was also taken into account to determine if the growth of the bacterium would be affected. Only a slight, but insignificant decrease was observed, 0.99 compared to the wildtype, showing that the potential knockouts should be viable. The carotenoid production was effectively decreased for most of the simulated knockout strains, except for reactions in the lower part of the pathway, branching off to produce different types of carotenoids. This, together with the fact that the most desirable effects are present in the upper part of the carotenoid biosynthesis pathway, suggests that gene knockouts in the upper half of the pathway are the most promising gene targets to increase PHA production and eliminate carotenoid production.

Fig. 3 | Heat map with fold-change difference of reaction fluxes related to the PHB cycle and growth of simulated knockouts compared to wild type M. extorquens based on FBA results using the GEM of M. extorquens. Desired effects are marked green, undesired effects are marked orange and no effect is marked yellow.

A visualization of the flux in the predicted optimal knockout strains is shown in Figure 4. As the phaC reaction is the only reaction negatively affected by the carotenoid biosynthesis pathway knockouts, further engineering of M. extorquens should be attempted to overcome this decrease. This suggests that combining a phaC overexpression strain with a carotenoid knockout strain could further improve M. extorquens for optimal PHA production.

Fig. 4 | Visualization of the effect of the simulated optimal knockout on reaction fluxes involved with PHA production (marked pink), determined by FBA.

Conclusion

This knockout simulation identified 19 promising gene targets which are predicted to result in a ~110 times increase in total PHA flux and elimination of carotenoid production. These genes encode enzymes that are mainly involved in the first half of the carotenoid biosynthesis pathway. We used these results to select three promising gene targets to create carotenoid biosynthesis pathway knockout strains in the lab.

Experimental implementation & Future prospects

From the 19 promising gene targets that showed a similar effect, three genes encoding enzymes that are involved in different parts of the upper half of the carotenoid biosynthesis pathway were selected to create a carotenoid knockout strain of M. extorquens: dxr, ispA and crtI (Figure 5). These genes encode for the enzymes 1-deoxy-D-xylulose-5-P (DXP) reductioisomerase, geranylgeranyl-diphosphate (GGPP) synthase, and phytoene dehydrogenase, respectively. This selection was included in the design of cycle C, which investigates the effect of the carotenoid biosynthesis pathway knockouts on PHA production in the wet lab. Read more about the implementation of these predictions in our experimental design on the Engineering page.

Fig. 5 | Selected knockout targets dxr, ispA and crtI which are predicted to effectively eliminate carotenoid production and enhance PHA production in M. extorquens.

As can be read on our Results page, a viable ΔcrtI strain was created and validated. Although the model suggested that the carotenoid production in this knockout would be reduced, this mutant strain was still pink, which means that the carotenoid production is not completely eliminated. In contrast to the ΔcrtI strain, the Δdxr mutant was unviable, although this was not predicted by our knockout simulation. The knockout simulation only showed a minor, unsignificant effect on biomass production (Figure 3). We thus hypothesize that dxr is an essential gene, and that this gene deletion is lethal, as previous research also suggested13. Lastly, an ΔispA strain could not be obtained due to cloning failure.

To validate the predictions of our gene knockout models, the PHA production of the knockout strains should be determined and compared to wild type M. extorquens to evaluate whether our predicted fold-change carbon flux has indeed led to an increase in PHA production. However, due to limited time availability, we could not quantify the PHA content and growth rate of our ΔcrtI strain. Nonetheless, we propose the Nile Red staining as a simple fluorescence-based quantification protocol, which we used for the PHA quantification in our other modified M. extorquens strains. Furthermore, the growth rate should be determined to see if this gene knockouts impairs the growth of M. extorquens. Based on these results, the model could be optimized to be able to improve the validity of the genome-scale metabolic network model and its predictions. For a more profound validation, multi-omics data could be used, which allows a more in depth investigation of the modeled metabolism.

Furthermore, future optimization attempts of the PHA overproducing strain could focus on redesigning our modified M. extorquens based on our modeled predictions combined with our experimental data. We have created a phaC overexpression strain, which showed a significant increase in PHA production. Combining this with a carotenoid biosynthesis pathway knockout could overcome the only limitation of the carotenoid pathway knockout, as it showed a predicted decrease in phaC flux, although the total PHA flux increased (Figure 3 and Figure 4). This conclusion was implemented in the redesign of cycle C.

References
  1. Tan, G.-Y., Chen, C.-L., Li, L., Ge, L., Wang, L., Razaad, I., Li, Y., et al. (2014). Start a Research on Biopolymer Polyhydroxyalkanoate (PHA): A Review. Polymers, 6(3), 706–754. MDPI AG. http://dx.doi.org/10.3390/polym6030706
  2. Chistoserdova, L., Chen, S. W., Lapidus, A., & Lidstrom, M. E. (2003). Methylotrophy in Methylobacterium extorquens AM1 from a genomic point of view. Journal of bacteriology, 185(10), 2980–2987. https://doi.org/10.1128/JB.185.10.2980-2987.2003
  3. Peyraud, R., Schneider, K., Kiefer, P., Massou, S., Vorholt, J. A., & Portais, J. C. (2011). Genome-scale reconstruction and system level investigation of the metabolic network of Methylobacterium extorquens AM1. BMC systems biology, 5, 189. https://doi.org/10.1186/1752-0509-5-189
  4. Ramon, C., Gollub, M. G., & Stelling, J. (2018). Integrating -omics data into genome-scale metabolic network models: principles and challenges. Essays in biochemistry, 62(4), 563–574. https://doi.org/10.1042/EBC20180011
  5. Ng, R.H., Lee, J.W., Baloni, P., Diener, C., Heath, J.R. and Su, Y. (2022). Constraint-Based Reconstruction and Analyses of Metabolic Models: Open-Source Python Tools and Applications to Cancer. Front. Oncol. 12:914594. doi: 10.3389/fonc.2022.914594
  6. Luo, Q., Ding, N., Liu, Y., Zhang, H., Fang, Y., & Yin, L. (2023). Metabolic Engineering of Microorganisms to Produce Pyruvate and Derived Compounds. Molecules, 28(3), 1418. MDPI AG. Retrieved from http://dx.doi.org/10.3390/molecules28031418
  7. Schada von Borzyskowski, L. , Carrillo, M., Leupold, S., Glatter, T., Kiefer, P., Weishaupt, R., et al. (2018). An engineered Calvin-Benson-Bassham cycle for carbon dioxide fixation in Methylobacterium extorquens AM1, Metab Eng, 47, pp. 423-433, 10.1016/j.ymben.2018.04.003
  8. Ranganathan, S., Suthers, P.F., Maranas, C.D., (2010). OptForce: An Optimization Procedure for Identifying All Genetic Manipulations Leading to Targeted Overproductions. PLoS Comput Biol 6(4): e1000744. https://doi.org/10.1371/journal.pcbi.1000744
  9. Pham, T.H., Webb, J.S., & Rehm, B.H. (2004). The role of polyhydroxyalkanoate biosynthesis by Pseudomonas aeruginosa in rhamnolipid and alginate production as well as stress tolerance and biofilm formation. Microbiology (Reading, England), 150 (Pt 10), 3405–3413. https://doi.org/10.1099/mic.0.27357-0
  10. Huergo, L.F., & Dixon, R. (2015). The Emergence of 2-Oxoglutarate as a Master Regulator Metabolite. Microbiology and molecular biology reviews : MMBR, 79(4), 419–435. https://doi.org/10.1128/MMBR.00038-15
  11. Orita, I., Nishikawa, K., Nakamura, S. et al. (2014). Biosynthesis of polyhydroxyalkanoate copolymers from methanol by Methylobacterium extorquens AM1 and the engineered strains under cobalt-deficient conditions. Appl Microbiol Biotechnol 98, 3715–3725. https://doi.org/10.1007/s00253-013-5490-9
  12. Chang, W., Yoon, J. & Oh, M.K. (2022). Production of Polyhydroxyalkanoates with the Fermentation of Methylorubrum extorquens Using Formate as a Carbon Substrate. Biotechnol Bioproc E 27, 268–275. https://doi.org/10.1007/s12257-021-0218-7
  13. Mo, X.H., Zhang, H., Wang, T.M. et al. (2020). Establishment of CRISPR interference in Methylorubrum extorquens and application of rapidly mining a new phytoene desaturase involved in carotenoid biosynthesis. Appl Microbiol Biotechnol 104, 4515–4532. https://doi.org/10.1007/s00253-020-10543-w