{% extends "layout.html" %} {% block title %}Model{% endblock %} {% block lead %}Our team modeled the decursin biomanufacturing system from micro to macro. We worked on multiple models and employed a stacking approach, where the output of one model was used as the input of another. Thus, we were able to optimize the system for large-scale production of our target metabolite - decursin.{% endblock %} {% block page_content %}
Modeling is an inherent part of the engineering cycle, even more so in the context of biomanufacturing, a discipline that relies on optimization when it comes to assessing feasibility. For our group, modeling was the driving force behind our work. Its importance was emphasized by the feedback we received during our meetings with stakeholders, which inquired regarding the production rates and costs. Thus, helping us identify which aspects needed to be modeled.
We worked on multiple models and employed a stacking approach, where the output of one model was used as the input of another. This was done to optimize the system for large-scale production of our target metabolite - decursin. Our biomanufacturing process will be the synthesis of decursinol from umbelliferone. Later on, decursinol can be converted chemically and cheaply to decursin.
Our models can be divided into four categories:
Figure 1: Schematic overview of Angel Root’s models. UMB precursor diffuses into a recombinant E. coli BL21 cell, expressing PT to convert UMB to DMS, and a combination of XimD and XimE to convert DMS into decursinol. Modeling simulations were run for the diffusion element (1), enzyme expression (2), in-silico directed evolution to improve XimE’s specificity (3), and the process’s enzymatic kinetics (4)
Before we set out to design the biomanufacturing process itself, our project started with a single basic assumption: the precursor UMB can diffuse into the cell, where the enzymatic reaction takes place. To examine this assumption, we wanted to model the diffusion process (1). To this end, we turned to the work done by iGEM CLS_CLSG_UK 2020 and expanded upon their efforts. We made the needed adjustments according to our needs and used the prediction of almost-immediate diffusion as output. We then wanted to establish a way to predict production rates, by focusing on three factors: bacterial population size, enzyme expression rates, and their kinetics. This led us to the development of a fermentation model (2), whose prediction of total protein expression was used as an output. Regarding the kinetics of the enzymes, previous work has shown that XimE and XimD are not substrate-specific[1]. To address this, and improve XimE’s catalytic efficiency in the decursin pathway, we employed an approach of using computational tools to model the substrate binding affinity of different mutations (3). The assumption of diffusion-based kinetics (1) and the protein concentration from (2) were used as input for an enzyme kinetics model (4). This model took into account the co-dependence of the enzymes in the pathway using statistical tools and was used to predict production rates. Furthermore, the insights regarding the limiting factor and byproducts resulted in a redesign of our plasmid system and were used to assess its feasibility. Afterward, the final model was fitted to the experimental data from the wet lab, to improve its predictions. All of the graphs displayed here are generated using our MATLAB code, which is open source and available here: Gitlab link.
In the anabolic pathway of DEC, the first metabolite is UMB which should be supplied externally, and be able to diffuse the membrane of E. coli cells to start the process. Therefore, we wanted to show the occurrence of diffusion using known equations and mathematical derivatives.
Table 1: Description of Model Symbols and Values
Parameters | Units | Description | Value |
---|---|---|---|
$ J $ | $ \frac{mol}{sec\cdot m^{2}} $ | Diffusion influx |
$-$ |
$ D $ | $ \frac{m^{2}}{sec} $ | Diffusion coefficient | $-$ |
$ C_{out} $ | $ \frac{mg}{L} $ | Umbelliferone concentration outside the cell | $ 100 $[1] |
$ C_{in} $ | $ \frac{mg}{L} $ | Umbelliferone concentration inside the cell | $-$ |
$ z $ | $ m $ | Thickness of E. coli membrane | $ 4\cdot10^{-9} $[5] |
$ A $ | $ m^{2} $ | Surface area of E. coli cell | $ 6\cdot10^{-12} $[4] |
$$ \beta $$ | Partition coefficient | $ 10^{1.896} $[6] | |
$ n $ | $ mol $ | Number of umbelliferone moles | $-$ |
$ V $ | $ m^{3} $ | Volume of E. coli cell | $ 10^{-18} $[4] |
$ t $ | $ sec $ | Time | |
$ k $ | $ \frac{m^{2}\cdot kg}{K\cdot sec^{2}} $ | Boltzmann constant | $1.380649\cdot10^{-23}$[7] |
$ T $ | $ K $ | Temperature | $ 305 $ |
$ R $ | $ m $ | Maximum radius of umbelliferone molecule | $ 5.06\cdot10^{-10} $ [6] |
$ \eta $ | $ Pa\cdot sec $ | Viscosity of medium | $ 796\cdot10^{-6} $[33] |
By the stated assumptions 1 and 2, the above equations can be rewritten as follows:
$$ (1.3) \ \frac{dn}{dt}=-DA\frac{(C_{in}-C_{out})}{z} $$To consider the membrane permeability it is necessary to consider the partition coefficient, which is the ratio of the solubility of the molecule in lipid and water and derived experimentally.
$$ (1.4) \ \frac{dn}{dt}=A\frac{\beta D}{z}(C_{out}-C_{in}) $$Dividing both sides by V:
$$ (1.5) \ \frac{dn}{dt}\cdot\frac{1}{V}=A\frac{\beta D}{Vz}(C_{out}-C_{in})\rightarrow\frac{dC}{dt}=A\frac{\beta D}{Vz}(C_{out}-C_{in}) $$Dividing both sides by $(C_{out}-C_{in})$:
$$ (1.6) \ \frac{dC}{dt}\cdot\frac{1}{(C_{out}-C_{in})}=A\frac{\beta D}{Vz} $$Solving first order equation:
$$(1.7) \ \int\frac{dC}{(C_{out}-C_{in})}=\int A\frac{\beta D}{Vz}dt $$Under the assumption that the only variable is $C_{in}$:
$$ (1.8) \ -ln(C_{out}-C_{in})=A\frac{\beta D}{Vz}\cdot t+a $$Raising both sides of the equation by the power of exp:
$$ (1.9) \ (C_{out}-C_{in})=be^{-A\frac{\beta D}{Vz}\cdot t} $$When $b=e^{a}$
At the time 0 sec concentration of umbelliferone inside the cell is 0 $\frac{mol}{m^{3}}$:
$$(1.11) \ 0=-b+C_{out}\rightarrow b=C_{out} $$Inputing equation (11) into equation (10):
$$ (1.12) \ C_{in}=C_{out}(-e^{-A\frac{\beta D}{Vz}\cdot t}+1) $$The value for diffusion coefficient wasn't found, therefore we decided to evaluate it using the Stokes-Einstein Law for diffusion[7]:
$$ (1.13) \ D=\frac{kT}{6\pi\eta R} $$Inputting equation (13) into equation (12), by using assumption (5):
$$ (1.14) \ C_{in}=C_{out}(-e^{-A\frac{\beta kT}{6\pi\eta RVz}\cdot t}+1) $$Inputting values into the final equation:
$$ (1.15) \ C_{in}=C_{out}(-e^{-12.4754\cdot t}+1) $$According to the work by Bu, X.L et.al [1], we assume an initial concentration of $100 \frac{mg}{L}$ outside of the cell. Thus $\text{C}_{\text{out}} = 100 \frac{mg}{L}$.
Using these values we get the following graph from our simulation:
Figure 2: Diffusion of UMB into E. coli BL21. UMB diffusion was simulated using Fick’s law. intracellular and extracellular UMB concentrations reach equilibrium in about 0.2 seconds. Simulation represents a final concentration of 100mg/L.
The graph shows the concentration of umbelliferone inside the cell in mg/L vs. time in sec. In less than one second the concentration inside the cell reaches the initial concentration outside the cell. Therefore, the model has shown that diffusion occurs at a time scale of less than a second under all our assumptions (fig.2). Consequently we can use this conclusion and neglect the diffusion time in future models and developments.
We aimed to model the growth behavior of our bacteria in a fermentor. We used this model to predict the maximum amount of protein that will be produced by the bacteria. This amount will be significant in our enzyme kinetics model. To model the fermentation process we use the Monod growth kinetics[9]. This theory links substrate (glucose), product (protein), and bacteria concentrations.
Table 2: Description of Model Symbol and Values as Inferred from Literature
Symbol | Units | Description | Value |
---|---|---|---|
$ \mu $ | $ \cfrac{1}{hour} $ | Specific growth rate |
$-$ |
$ X $ | $ \cfrac{g_{bacteria}}{liter} $ | Biomass concentration | $-$ |
$ t $ | $ hour $ | Time | $-$ |
$ P $ | $ \cfrac{g_{protein}}{liter} $ | Protein concentration | $-$ |
$ S $ | $ \cfrac{g_{glucose}}{liter} $ | Substrate (glucose) concentration | $-$ |
$ q_{p} $ | $ \cfrac{g_{protein}}{g_{bacteria}\cdot hr} $ | Specific protein formation |
$-$ |
$ m_{s} $ | $ \frac{g_{glucose}}{g_{bacteria}\cdot hr} $ | Maintenance coefficient | $ 0.07 $ [9] |
$ K $ | $ \frac{1}{hr} $ | Specific protein destruction | $ 0.01 $ [11] |
$ \mu_{max} $ | $ \frac{1}{hr} $ | Maximum growth rate | $ 0.27 $ [12] |
$ k_{s} $ | $ \frac{g}{liter} $ | Substrate concentration when $\mu=0.5\mu_{max}$ | $ 180\cdot10^{-7} $ |
$ Y_{x/s} $ | $ \cfrac{g_{bacteria}}{g_{glucose}} $ | Yield of bacteria per substrate | $ 0.37 $ [13] |
$ Y_{p/s} $ | $ \cfrac{g_{protein}}{g_{glucose}} $ | Yield of protein per substrate | $ 0.2 $ [14] |
$ Y_{p/x} $ | $ \cfrac{g_{protein}}{g_{bacteria}} $ | Yield of protein per bacteria | $ 0.46 $ [12] |
We describe our system with a set of ordinary differential equations (ODEs)[15]:
The following equation describes the change in biomass concentration over time:
$$ (2.1) \ \dfrac{dX}{dt}=\mu\cdot X$$We use the Monod equation, which describes the specific growth rate:
$$ (2.2) \ \mu=\mu_{max}\cdot\frac{S}{K_{s}+S}$$Then we define the following ODE that describes the change in protein concentration:
$$ (2.3) \ \dfrac{dP}{dt}=q_{p}\cdot X-K\cdot P$$Calculate a constant that describes the specific protein formation:
$$ (2.4) \ q_{p}=Y_{p/x}\cdot\mu$$The needed energy for glucose-consuming cellular activities can be described in terms of glucose (substrate) consumed per time. The energy required for cell growth:
$$ (2.5) \ \cfrac{\mu\cdot X}{Y_{x/s}}$$The energy required for maintenance of the bacteria:
$$ (2.6) \ m_{s}\cdot X$$The energy required for the synthesis of enzymes (product):
$$ (2.7) \ \cfrac{q_{p}\cdot X}{Y_{p/s}}$$Adding it all together, the rate of glucose consumption is:
$$ (2.8) \ -\cfrac{dS}{dt}=\cfrac{\mu\cdot X}{Y_{x/s}}+m_{s}\cdot X+\cfrac{q_{p}\cdot X}{Y_{p/s}}$$Figure 3: Bacterial concentration vs time. The input for the model is a time span of 10 hours with initial bacterial concentration of zero [g/L]. The plot shows exponential growth up until around four hours.
Figure 4: Glucose concentration vs time. The input for the model is a time span of ten hours with an initial glucose concentration of 10 [g/L]. Glucose is fully consumed at around four hours
Figure 5: Protein concentration vs time. The input for the model is a time span of 10 hours with initial protein concentration of zero [g/L]. The plot shows exponential increase in protein concentration up until around four hours. After four hours, the protein concentration decreases as no more protein is being produced, and natural degradation takes place.
As expected, in a batch process the glucose eventually runs out. (fig.4). Once there is no more glucose in the system, bacterial growth stops. That causes protein production to stop as well. Because there are no more new proteins, the relative effect of protein natural degradation increases and a reduction of total protein levels can be seen. The model predicts peaks of bacteria (fig.3) and protein (fig.5) concentrations at around 4 hours after induction. The predicted protein concentration is ~0.86g/L. The predicted protein concentration is on par with data from the literature for the BL21 strain, with an induction time similar to our experiment[12] .
In order to make the production of decursin cost-effective, throughout the modeling and biosynthetic design, we tried to focus on making the anabolic pathway more efficient. One of the ways used nowadays is directed protein mutagenesis, a robust method to design proteins with desirable functions, in our case greater production of decursin. Therefore, for this purpose computational tools were used to identify mutations that will improve our proteins and would help to shorten the time of future lab work.
The desirable target in the case of decursin production would be the remodeling of one of the key enzymes in the pathway, XimE. It is the final enzyme in our anabolic pathway, converting epoxidated 7-demethylsuberosin, (Fig.6), into decursinol. Using directed mutagenesis, we would like to achieve greater binding affinity of the substrate, which in turn, under our assumptions, leads to greater production of decursinol.
The computational approach of docking aims to predict the binding pose and affinity of a small molecule ligand and protein through scoring algorithms.
Figure 6: a proposed mechanism for the catalysis. First reaction from left to right: XimD converting 7-demethylsuberosin (DMS) into epoxidated 7-demethylsuberosin (Epoxidated DMS). Second reaction:XimE converting epoxidated 7-demethylsuberosin (Epoxidated DMS) into decursinol (DEC).
Our workflow consisted of performing structural prediction of a protein with a known single point mutation, conducting molecular docking of predicted mutant proteins of the same type with unchangeable ligand, and eventually analyzing and comparing results to reveal the desired mutation and its further impact on wet lab work.
First, we found the 3D structure and amino-acid sequence of XimE in the Protein Data Bank (PDB) and manually performed point mutations in desired positions to build the stock of mutant sequences destined for structural predictions. While accessing AlphaFold2 from ChimeraX [16], we retrieved the predicted structures of mutant proteins of the same type (XimE). AlphaFold[17] is a computational approach capable of predicting protein structures to near experimental accuracy, whereas ChimeraX[16] runs its script automatically on Google Colab servers. During our trials with docking between target ligand and XimE we noticed great differences in scores between WT protein retrieved from PDB, aka crystallographic structure, and WT protein predicted with AlphaFold, although both have the same amino acid sequence. Therefore, we decided to use the predicted wild-type (WT) protein from AlphaFold for comparison with predicted mutant proteins for our further investigations.
The next step was protein-ligand docking. It is a method in structural molecular biology, whose goal is to predict the optimal binding mode of a ligand with a protein of known three-dimensional structure. There are various programs that can be used for protein-ligand docking, however, more comfortable the fastest and the most dependable for us appeared to be CB-Dock2 used for blind docking because the protein of interest lacked the verified position of the active site. The program first automatically detects binding sites and customizes the docking box and then employs the latest version (1.2.0) of AutoDock Vina for the docking itself.[18][19] As a result, the software draws out the list of ranked binding sites of the ligand according to the Vina score, predicted binding affinity calculated in kcal/mol, which in turn can be translated into binding free energy ∆G.
Finally, we analyzed and compared all the scores and visualized the following interaction between mutant protein and target ligand using PyMol, a molecular visualization system on an open-source foundation.
You can find all of our files, and artifacts here: Supplements
XimE is a dimer made of two identical subunits, hence for the simplicity of the processes we ran the mutations on a single chain (chain a) as indicated in the literature[21][22][23].
Before the molecular docking, we needed to choose the suitable amino acids on which the mutation would be performed. To understand the region of the cavity where the ligand binds, we ran a wild-type (WT) simulation with our ligand and received an estimated location. We integrated the amino acids around the selected area with the literature and selected three amino acids that were most examined in mutagenesis experiments. In their paper, Bei-Bei H et al.[21] have identified Glutamic acid at position 136 (E136), Histidine at position 102 (H102), and Tyrosine at position 153 (Y153) as possible key residues for the pyran ring formation. Chuchu Jiang et al.[22] stated that mutagenesis in E136 or H102 didn’t eliminate the enzyme activity, reinforcing the choice.
For Y153 we also based on the work of Xu-Liang Bu et al. [23] that performed in vitro mutagenesis in Y153.
Moreover, we performed the initial docking of the WT protein with the target ligand and received the following list of residues participating in the bond as listed in Fig.6. It can be shown that the three selected amino acids are displayed.
Figure 7: A proposed mechanism for the catalysis. CB-Dock2 representation of XimE suggested the shown amino acids are involved in epoxidated 7-demethylsuberosin binding. List of residues participating in the bond: Y46, M47, L75, F93, A96, V97, H102, W118, V120, A122, R124, M130, A132, E136, N151, Y153, P157. In bold: selected positions for mutagenesis analysis.
Figure 8: Amino acids selected for mutagenesis analysis and their interaction with XimE's substrate, epoxidated 7-demethylsuberosin.(a) H102 residue making bond with ligand. (b) E136 residue making bond with ligand. (c) Y153 residue making bond with ligand.
The resulting scores obtained from CB-Dock2 were compared and summed in the supplements attached, to view them click here. The leading scores from each mutated amino acid are presented below in table 1.
It can be seen from the table that the initial score of the WT is -6.0 kcal/mol, and the favorable mutations’ scores are lower. This matches expectations as a more negative score indicate a higher affinity between the protein and the ligand[24].
The best score was obtained from the mutation H102P, with a score of -9.2 kcal/mol.
Table 3: Top Vina scores for the mutagenesis of each of the selected amino acid.
Mutation | Vina score [kcal/mol] |
---|---|
WT | -6.0 |
Y153R | -7.8 |
E136P | -8.4 |
H102P | -9.2 |
Figure 9: interaction of WT protein and H102P mutant protein with ligand. WT and H102P (mutant) protein aligned together interacting with 7-demethylsuberosin (in blue – WT protein with its ligand in the same color, in yellow – mutant protein with its ligand in the same color, in red – position of mutation, in pink – amino acids of WT protein not aligned, in grey - amino acids of mutant protein not aligned). (a) Rotating two proteins aligned together. (b, c) Closer look of areas not aligned betweeen two proteins. The ligand of mutant protein lays deeper in the cavity.
Proline is a much more rigid amino acid than others, and due to its characteristics, it often improves the stability of the protein as a whole[25], which may ensure a stronger and more stable bond between the XimE and the ligand.
Directed mutagenesis is a significant tool for improving the properties of proteins and the affinity between enzymes and their substrates. By using computational tools, it is possible to investigate a wide variety of mutations resulting in reducing time, manpower, and laboratory equipment. After targeting three amino acids- E136, Y153, and H102 we received different scores, with H102P being the lowest. As Proline is known for stabilizing tertiary structures of proteins it is a possible explanation as to why H102P is the optimal mutation.
The computational approach described here was desirable for us, because of the feasibility of testing the predictions and results in the wet lab. In our case, we were able to design primers that allow us to generate mutated enzymes based on the model results. By establishing the production rate of decursin as a criterion for comparison, we can assess the results of different mutations to the WT. The experiments we planned to validate our modeling work would be a good starting point for our future work.
We define the system as an expression system with three enzymes: PT, XimE, and XimD. PT catalyzes the reaction from UMB to DMS and XimE together with XimD catalyzes the reaction from DMS to DEC. To optimize the system for the production of DEC, a precursor of decursin, we needed to model the enzymatic reactions and characterize their dynamics.
Using our fermentation model we were able to estimate the total amount of protein produced in the system. According to our design, the bacterial population expresses all three enzymes. The next step was to understand the distribution between the different enzymes that catalyze the enzymatic reactions, this will allow us to calculate the concentration of each enzyme in the system. We used experimental data found by Bei-Bei et al. (2019), which held experiments on a two-enzyme system that expresses only XimE and XimD.
Now, we can model two different systems, plug-in experimental data from the two-enzyme system, and based on that draw conclusions regarding the three-enzyme system.
Table 4: Description of Model Symbol and Values as Inferred from Literature
Symbol | Units | Description | Value |
---|---|---|---|
$\text{XimD length}$ | $\text{No. of amino acids}$ | Amount of amino acids in a XimD protein | $474$ [27] |
$\text{XimE length}$ | $\text{No. of amino acids}$ | Amount of amino acids in a XimE protein | $124$[27] |
$\text{PT length}$ | $\text{No. of amino acids}$ | Amount of amino acids in a PT protein | $400$[27] |
$N_{A}$ | $ \frac{\text{molecules}}{\text{mol}} $ | Avogadro constant | $\sim6.022\cdot10^{23}$ |
$\text{[Total Protein]}$ | $ \frac{mg}{ml} $ | Concentration of protein produced in the cell according to presented fermentation model [1] | $0.5$ |
$m_{\text{aa}}$ | $Da$ | Average weight of an amino acid | $110$ [26] |
$MW_{\text{XimD}}$ | $ Da $ | Molar mass of a single XimD protein molecule | $51,010$[27] |
$MW_{\text{XimE}}$ | $ Da $ | Molar mass of a single XimE protein molecule | $13,637$[27] |
$m_{\text{PT}}$ | $ Da $ | Mass of a single PT protein molecule | $44,221$[27] |
$v_{\text{cell}}$ | $ \mu m^{3} $ | Volume of a bacterium cell | $1$[29] |
$m_{2,\text{Total XimD}}$ | $ mg $ | Total mass of XimD protein in the two enzyme system | Calculated in model |
$m_{2,\text{Total XimE}}$ | $ mg $ | Total mass of XimE protein in the two enzyme system | Calculated in model |
$m_{3,\text{Total XimD}} $ | $ mg $ | Total mass of XimD protein in the three enzyme system | Calculated in model |
$m_{3,\text{Total XimE}}$ | $ mg $ | Total mass of XimE protein in the three enzyme system | Calculated in model |
$[\text{XimD}]_{2} $ | $ \mu M $ | Enzyme concentration in an two enzyme system | Calculated in model |
$ [\text{XimE}]_{2}$ | $ \mu M $ | Enzyme concentration in an two enzyme system | Calculated in model |
$ [\text{XimD}]_{3}$ | $ \mu M $ | Enzyme concentration in an three enzyme system | Calculated in model |
$[\text{XimE}]_{3} $ | $ \mu M $ | Enzyme concentration in an three enzyme system | Calculated in model |
$ m_{\text{total protein}}$ | $ mg $ | Total mass of enzymes in a single bacterium | Calculated in model |
$N$ | $ - $ | Number of molecules in volume V | Calculated in model |
$V$ | $ \mu M $ | Total mass of enzymes in a single bacterium | Volume |
Because our expression system contains two enzymes, it holds that:
$$ (4.1)\ m_{\text{total protein}}=[\text{Total Protein}]\cdot v_{\text{cell}}=m_{\text{XimE}}+m_{\text{XimD}} $$By our assumption, the distribution of mass is equal to the ratio between enzyme lengths:
$$ (4.2)\ \frac{m_{2,\text{Total XimD}}}{m_{2,\text{Total XimE}}}=\frac{\text{XimDlength}}{\text{XimE length}}$$This yields:
$$ (4.3)\ m_{2,\text{Total XimD}}=\frac{\text{XimD length}}{\text{XimD length}+\text{XimE length}}\cdot m_{\text{total protein}} $$ $$ (4.4)\ m_{2,\text{Total XimE}}=\frac{\text{XimE length}}{\text{XimD length}+\text{XimE length}}\cdot m_{\text{total protein}}$$We calculate the number of molecules from each kind to get the concentration of each enzyme in a single bacterium.
For XimD:
For XimE:
$$ (4.6)\ [\text{Xime}]_{2}=\frac{m_{2,\text{Total XimE}}\cdot\frac{N_{A}}{MW_{\text{XimE}}}}{v_{\text{cell}}}=\frac{4400}{1\mu m^{3}}\underbrace{=}_{\text{(Assumption 5)}}4.4\,\mu M $$In this case:
$$ (4.7)\ m_{\text{total protein}}=m_{\text{XimE}}+m_{\text{XimD}}+m_{\text{PT}}=0.5mg $$In a similar fashion, considering PT Length:XimD Length:XimE Length=3.2:3.8:1, we find:
$$ [\text{XimD}]_{3}=2.7\,\mu M $$ $$ [\text{XimE}]_{3}=2.8\,\mu M $$ $$ [\text{PT}]_{3}=2.7\,\mu M $$After estimating the concentration of the enzymes in our system we continued our modeling work to answer two questions, the first being how much of the UMB is converted to DEC. The second would be the temporal aspect: understanding how long the enzymatic reaction will take. Both of these aspects together would allow us to quantify a production rate. This understanding of the key factors in the biomanufacturing process will allow us to create a more optimized biological circuit design, and accordingly provide us with predictions for the cost and time of the production.
Table 5: Description of Model Symbol and Values as Inferred from Literature
Symbol | Units | Description | Value |
---|---|---|---|
$\text{[XimE]}$ | $\mu M$ | XimE enzyme concentration | $2.8 \text{[previous model]}$ |
$\text{[XimD]}$ | $\mu M$ | XimD enzyme concentration | $2.7 \text{[previous model]}$ |
$\text{[PT]}$ | $\mu M$ | Prenyltrasnfarase enzyme concentration | $2.7 \text{[previous model]}$ |
$v_{\text{UMB}}$ | $\frac{\mu M}{sec}$ | Rate of formation of UMB | $-$ |
$v_{\text{DMS}}$ | $\frac{\mu M}{sec}$ | Rate of formation of DMS | $-$ |
$v_{\text{DEC}}$ | $\frac{\mu M}{sec}$ | Rate of formation of DEC | $-$ |
$v_{\text{MAR}}$ | $\frac{\mu M}{sec}$ | Rate of formation of marmesin | $-$ |
$k_{m,\text{PT}\cdot\text{UMB}}$ | $\mu M$ | Michaelis constant | $20$[30] |
$k_{1,\text{cat}}$ | $\frac{1}{sec}$ | Rate of formation DMS from UMB | $0.435$[30] |
$[\text{UMB}]$ | $\mu M$ | UMB concentration | $-$ |
$[\text{DMS}]$ | $\mu M$ | DMS concentration | $-$ |
$[\text{DEC}]$ | $\mu M$ | DEC concentration | $-$ |
$[\text{MAR}]$ | $\mu M$ | Marmesin concentration | $-$ |
We defined the following reaction:
$$ (4.8)\ \text{PT}+\text{UMB}\overset{k_{1}}{\underset{k_{-1}}{\rightleftarrows}}\text{PT}\cdot\text{UMB}\overset{k_{1cat}}{\rightarrow}\text{PT}+\text{DMS} $$based on the Michaelis-Menten equation[30] which states that the rate of a reaction is proportional to the concentrations of the enzyme and the substrate:
$$ (4.9)\ v=\frac{d[P]}{dt}=\frac{k_{\text{cat}}\cdot[S]\cdot[E]}{k_{m}+[S]} $$we reached the following:
$$ (4.10)\ v_{\text{UMB}}=-\frac{k_{1,\text{cat}}\cdot[\text{UMB}]\cdot[\text{PT}]}{k_{m,\text{PT}\cdot\text{UMB}}+[\text{UMB}]} $$The penultimate reaction in the biosynthesis transforms DMS to DEC. It is catalyzed by two enzymes co-dependently: flavin-dependent monooxygenase (XimD) and SnoaL-like cyclase (XimE). Together, they catalyze the formation of the pyran ring (6-endo-tet cyclization) necessary for the synthesis of pyranocoumarins such as decursinol.
In vitro enzymatic assays were used [1] to show that XimD alone (without the presence of XimE) catalyzes the formation of furan ring (5-exo-tet cyclization) which produces marmesin (MAR) instead of DEC.
It is only when both XimE and XimD are present, that the pyran ring is formed and decursinol is produced.
Microbial host models, expressing both XimE and XimD, clarified that the two reaction paths occur simultaneously. As we wanted to optimize our system for the production of DEC, we wanted the system to prefer the reaction that produces DEC
A suggested mechanism proposes an unstable epoxide intermediate, that can spontaneously decay to marmesin in the absence of XimE.
Figure 10: A proposed mechanism for the reactions of XimE and XimD. (Bei Bei et al.[1]). DMS is converted to an epoxide intermediate by XimD, which then in the presence of XimE is converted to DEC. in the absence of XimE, the intermediate is converted MAR.
Table 6: Description of Model Symbol and Values as Inferred from Literature
Symbol | Units | Description | Value |
---|---|---|---|
$N$ | $-$ | Number of molecules | $-$ |
$c$ | $\mu M$ | Concentration of reaction complex | $-$ |
$c_{0}$ | $\mu M$ | Initial concentration | $-$ |
$\beta$ | $\frac{1}{kJ}$ | Thermal influence | $\frac{1}{k_{B}T}$ |
$k_{B}$ | $\frac{kJ}{K}$ | Boltzmann constant | $1.380649\cdot10^{-26}$ |
$T$ | $K$ | Temperature | $-$ |
$\Delta\varepsilon$ | $kJ$ | Binding energy | $-$ |
$\varepsilon$ | $-$ | A constant expressing the dependency between XimE and XimD | In the magnitude of $K_{D,3}$ |
$Z$ | $-$ | Partition function | $-$ |
$k_{2,cat}$ | $\frac{1}{s}$ | Rate of formation marmesin from DMS | $0.1$ (Estimating model coefficients) |
$k_{3,cat}$ | $\frac{1}{s}$ | Rate of formation DEC from DMS | $0.44$ (Estimating model coefficients) |
$K_{D,1}$ | $\mu M$ | Dissociation constant of XimD with DMS | $0.1$ (Estimating model coefficients) |
$K_{D,2}$ | $\mu M$ | Dissociation constant of XimE,XimD system with DMS | $0.1$ (Estimating model coefficients) |
First, we defined the following reaction:
$$ (4.11)\ \text{XimD}+\text{DMS}\overset{k_{2cat}}{\rightarrow}\text{MAR} $$ $$ (4.12)\ \text{XimD}+\text{DMS}+\text{XimE}\overset{k_{3cat}}{\rightarrow}\text{DEC}+\text{MAR} $$We use a state variable description approach of protein-substrate binding. [31]
Each substrate molecule can be bound to each one of the two enzymes (XimD and XimE). This gives us four optional states:
Thus, we define a partition function. This function takes into account every state that our system could be in:
$$ (4.13)\ Z=1+\frac{[\text{XimD}]}{K_{D,1}}+\frac{[\text{XimD}][\text{XimE}]}{K_{D,1}\cdot K_{D,3}}\varepsilon+\frac{[\text{XimE}]}{K_{D,3}} $$While we perform the following reduction:$ \frac{\varepsilon}{K_{D,1}\cdot K_{D,3}}=\frac{1}{K_{D,2}}$
According to the Gibbs model[3], the average number of molecules in a bounded state is:
$$ (4.14)\ N=c\frac{\frac{1}{c_{0}}e^{-\beta\Delta\varepsilon}}{Z} $$To make our model depend on a single parameter, we define the following notation:
$$ (4.15)\ \frac{1}{\frac{1}{c_{0}}e^{-\beta\Delta\varepsilon}}=c_{0}e^{\beta\Delta\epsilon}=K_{D} $$We define the rate of a reaction as the number of molecules in bound state times the catalytic rate constant:
$$ (4.16)\ v=\frac{d[P]}{dt}=k_{cat}\cdot N $$And plug (4.6) and (4.7) into (4.8) to receive:
$$ (4.17)\ v=\frac{k_{cat}\cdot c}{K_{D}\cdot Z} $$Then, the change in the concentration of the metabolites follows:
$$ (4.18)\ v_{\text{DMS}}=\frac{k_{1,\text{cat}}\cdot[\text{UMB}]\cdot[\text{PT}]}{k_{m,\text{PT}\cdot\text{UMB}}+[\text{UMB}]}-\frac{k_{2,\text{cat}}}{K_{D,2}}\frac{[\text{DMS}]\cdot[\text{XimD}]\cdot[\text{XimE}]}{Z}-\frac{k_{3,\text{cat}}}{K_{D,1}}\frac{\cdot[\text{DMS}]\cdot[\text{XimD}]}{Z} $$ $$ (4.19)\ v_{\text{DEC}}=\frac{d[\text{DEC]}}{dt}=\frac{k_{2,\text{cat}}}{K_{D,2}}\frac{[\text{DMS}]\cdot[\text{XimD}]\cdot[\text{XimE}]}{Z} $$ $$ (4.20)\ v_{\text{MAR}}=\frac{d[\text{MAR]}}{dt}=\frac{k_{3,\text{cat}}}{K_{D,1}}\frac{\cdot[\text{DMS}]\cdot[\text{XimD}]}{Z} $$Where the change in DMS concentration is defined as:
$$ (4.21)\ v_{\text{DMS}}=-v_{\text{UMB}}-v_{\text{DEC}}-v_{\text{MAR}} $$To estimate the coefficients, we use the results found by Bei-Bei et al[30]: Cloning of XimD and XimE under two pET plasmids resulted in a DEC:MAR ratio of 1:1. Assuming the ratio of $\frac{dp}{dt}$ correlates with the final concentration, we get:
$$ (4.22)\ 1=\frac{[\text{DEC}]_{\text{final}}}{[\text{MAR}]_{\text{final}}}\implies1=\frac{v_{\text{DEC}}}{v_{\text{MAR}}}=\frac{\frac{k_{2},cat}{K_{D,2}}}{\frac{k_{3,cat}}{K_{D,1}}}\cdot[\text{XimE}]_{\text{basal}}\implies\frac{\frac{k_{2},cat}{K_{D,2}}}{\frac{k_{3,cat}}{K_{D,1}}}=\frac{1}{[\text{XimE}]_{2}} $$Define:
$$\frac{k_{2,cat}}{K_{D,2}}=k_{\text{eff,2}}$$ $$\frac{k_{3,cat}}{K_{D,1}}=k_{\text{eff,3}}$$As mentioned before:
$$ \frac{k_{3,cat}}{K_{D,1}}=k_{\text{eff},3}=4.4\frac{k_{2,cat}}{K_{D,2}}=4.4k_{\text{eff},2} $$Therefore:
$$\frac{k_{3,cat}}{K_{D,1}}=k_{\text{eff,3}}=4.4\frac{k_{2,cat}}{K_{D,2}}=4.4k_{\text{eff,2}}$$For that reason we chose arbitrary values while maintaining 2 conditions:
Using these values we get the following graph from our simulation:
Figure 11: Production of metabolites in the three-enzyme system. (a) Production of DEC vs time. UMB is converted to DMS, which is then immediately converted to DEC and MAR. The final DEC:MAR ratio is 0.64. The enzyme conrentration is as described in the three enzyme system. (b) Ratio of DEC:MAR vs XimE concentration. XimE molar fraction is calculated out of the total protein molecules in the system. The XimE percent stands for the amount of XimE molecules out of 100 protein molecules. As the molar fraction of XimE rises, the ratio increases in favor of decursin.
Figure 11a shows, that the entirety of our precursor, UMB is converted to the products, and results in saturation of DEC and MAR. The model also shows the feasibility of the process, the calculated production rate of DEC is 0.3 grams per liter per hour. The model predicts that the expression of the three enzymes under the same promoter will result in a higher concentration of the byproduct marmesin compared to the target product decursinol (fig.11a). Moreover, we observe that when the molar fraction of XimE increases, a higher concentration of decursinol is obtained compared to marmesin (fig.11b). This prediction prompted us to modify the design of our plasmid, we included a constitutive expression of XimE and an inducible expression of the other two enzymes. By controlling the expression levels of XimD & PT, we ensure a surplus of XimE. There is no correlation between [DEC] [MAR] and [XimD]. The ratio between the metabolites is proportional to the concentration of XimE.
Table 6: Description of this Section Symbols and Values as Taken from qPCR Results
Symbol | Units | Description | Value |
---|---|---|---|
$\text{XimE}_{\text{Copy number}}$ | $-$ | Number of XimE transcripts | $-$ |
$\text{XimD}_{\text{Copy number}}$ | $-$ | Number of XimD transcripts | $-$ |
$\text{XimE}_{\text{CT}}$ | $-$ | Cycle number in which XimE cDNA level reached threshold of 0.3 | $18.9810$ |
$\text{XimD}_{\text{CT}}$ | $-$ | Cycle number in which XimD cDNA level reached threshold of 0.3 | $16.2080$ |
According to our design, we have a plasmid containing XimD and XimE under T7 promoters (ximE_ximD). The expression of XimD is regulated by the Lac repressor while the expression of XimE is constitutive. We measured the expression levels in our two-enzyme system in the wet lab. This was done by running a qPCR of samples containing our engineered bacteria under the induction of IPTG. For detailed results please refer to our results. To calculate the ratio between the expression levels of XimD and XimE, we used the threshold defined for the qPCR process. For each gene, we calculated a cycle threshold (CT) value in which it reached the threshold of fluorescence. At these CT values, the cDNA level is equal.
Thus:
$$ (4.23)\ \text{XimE}_{\text{Copy number}}\cdot2^{\text{XimE}_{\text{CT}}}=\text{XimD}_{\text{Copy number}}\cdot2^{\text{XimD}_{\text{CT}}}$$ $$ (4.24)\ \frac{\text{XimD}_{\text{Copy number}}}{\text{XimE}_{\text{Copy number}}}=\frac{2^{\text{XimE}_{\text{CT}}}}{2^{\text{XimD}_{\text{CT}}}}=\frac{2^{18.9810}}{2^{16.2080}}\approx6.8$$6.8 is the molar ratio between XimD and XimE when IPTG is present. To get mass ratio we will multiply 6.8 with the ratio between the mass of a single XimD and the mass of a single XimE:
$$ (4.25)\ \text{XimD}:\text{XimD}=6.8\frac{\text{mol}}{\text{mol}}\cdot\frac{MW_{\text{XimD}}}{MW_{\text{XimE}}}=6.8\cdot\frac{51,010}{13,637}\approx25.4$$consequently, we deduce that XimD:XimE is 25.4:1 when IPTG is present. We redo the process detailed in the two-enzyme system, and get:
$$[XimE]_{2}=0.23\mu M$$ $$[XimD]_{2}=5.5\mu M$$We ran the model again with the new concentrations. However, now we are modeling a two enzymes system which means we remove equation 4.10 from our ODE system.
Figure 12: Production of DEC vs time with the concentrations calculated from qPCR results. The predicted production of DEC and MAR in the two-enzyme system (ximE_ximD) results in a DEC:MAR ratio of 0.05
Figure 12 shows a very high concentration of the unwanted product MAR. These wet lab results greatly emphasize to us the importance of creating a construct where the concentration of XimE will be as high as possible. We discuss ways to achieve this goal in the model implications.
After studying the results of the different models, we can proceed and test the predictions in the wet lab. To test our predictions regarding how different concentrations of XimE and XimD affect the end products of the manufacturing. We planned several plasmid designs, which can be compared according to the final concentrations of the end products. Each design allows for different expression rates of the XimE and XimD enzymes:
In the scope of iGEM, we focused on approach 3, which was the first one to be implemented in our design. {% endblock %} {% block references %}