Team:Toronto/Model

https://static.igem.org/mediawiki/2019/f/f1/T--Toronto--Reaction_equations.png

Bioreactor




Goals


The processes of main interest are:
  1. The enzyme kinetics between PETase, MHETase, and PET (Tank #1).
  2. The dilution process (Tank #2).
  3. The chemical reactions between Ammonium Hydroxide (NH4OH) and Terephthalic acid (TPA) (Tank #3).
A continuous model (i.e. based on differential equations) was used for each process.

Models


Tank #1 - Enzyme Kinetics
The following model is based on Michaelis-Menten kinetics. For an introduction, refer to Ingalls [1].

The reaction between PET and PETase produces MHET, BHET, and TPA. However, as one can see from figure 2.B in Yoshida et al. [2], the amount of produced MHET vastly outnumbers the produced BHET and TPA. Since MHET is our main interest for the second reaction $MHET + MHETase$, we assume that the amount of BHET and TPA produced from the $PET + PETase$ reaction is negligible.

Figure 1. shows the chemical equations governing the reactions under the previous simplification and the Michaelis-Menten assumption.
where the $k$ values represent the reaction rates. (1) is the PET-PETase reaction. $C_1$ is the complex generated by PET and PEtase. (2) is the MHET-MHETase reaction. $C_2$ is the complex generated by MHET and MHETase.

We keep track of the seven chemical species in (1) and (2) (We are taking TPA and EG as a single one for simplicity). This leads to seven ODEs describing the concentrations of the species over time, namely:

$$S_1(t)' = -k_1S_1E_1 + k_{-1}C_1 \tag{3}$$ $$E_1(t)' = -k_1S_1P_1 + (k_{-1} + k_2)C_1 \tag{4}$$ $$C_1(t)' = k_1P_1E_1 - (k_{-1} + k_2)C_1 \tag{5}$$ $$S_2(t)' = k_2C_1 - k_3P_2E_2 + k_{-3}C_2 \tag{6}$$ $$E_2(t)' = -k_3P_2E_2 + (k_{-3} + k_4)C_2 \tag{7}$$ $$C_2(t)' = k_3P_2E_2 - (k_{-3} + k_4)C_2 \tag{8}$$ $$P(t)' = k_4C_2 \tag{9}$$
where $S_1$ is the concentration of PET, $E_1$ is the concentration of PETase, $C_1$ is the concentration of the complex PET+PETase, $S_2$ is the concentration of MHET, $E_2$ is the concentration of MHETase, $C_2$ is the concentration of the complex MHET+MHETase, and $P$ is the final products TPA and EG. The apostrophe denotes derivative respect to time t.

To simplify our system of equation (3)-(9), we assume a conservation of enzymes in reactions (1) and (2), namely:

$$E_1' = - C_1' \tag{10}$$ $$E_2' = - C_2' \tag{11}$$
which simplifies our model to five ODEs by not explicitly keeping track of the enzyme concentrations anymore. The final model for tank #1 is:

$$S_1(t)' = -k_1E_{1}(0)S_1 +C_1(k_1S_1 + k_{-1}) \tag{12}$$ $$C_1(t)' = k_1E_{1}(0)S_1 - C_1(k_1S_1 + k_{-1} + k_2) \tag{13}$$ $$S_2(t)' = k_2C_1 - K_3E_2(0)S_2 + C_2(K_3S_2 + k_{-3}) \tag{14}$$ $$C_2(t)' = k_3E_2(0)S_2 - C_2(k_3 + k_4 + k_3S_2) \tag{15}$$ $$P(t)'= k_4C_2 \tag{16}$$
where $E_1(0)$ and $E_2(0)$ are the initial concentration of PETase and MHETase, respectively.

Recall that the main goal of this model is to know the concentration of TPA and EG produced after time $t$ has passed. To answer the question, one would solve the equations (12)-(16) numerically to obtain the function $P(t)$.

Tank #2 - Dilution of growth medium
Recall that we are considering the process in tank #2 as a two-stages process.

In the first stage, tank #2 contains an initial volume of water $H_0$. A constant flux $L_{in}$ (mass/time) brings the growth medium with the products from tank#1 into tank #2. In practice, the growth medium would be transported by water. We denote the constant water flux carrying the medium into the tank as $Q_{in}$ (volume/time). Hence, $L_{in} = Q_{in}C_L$, where $C_L$ is the concentration of medium in the water flux.

We want to know how long we must wait for the concentration of the diluted solution to be a chosen number $c << 1$. In a real scenario, $c$ would depend on the design of the bioreactor.

Since $L_{in}$ is constant, the amount of mass brought into the tank after time $t$ is $tL_{in}$. By the same reasoning, the amount of water brought into the tank is $Q_{in}t$. At $t = t^*$, the solution reaches the desired concentration. Hence,

$$t^*L_{in}/(H_0 + Q_{in}t^*)= c \xrightarrow{} t^* = H_0c/(L_{in} -Q_{in}c) \tag{17}$$
Equation (17) gives the time one must wait to start pumping the solution out of tank #2 into tank #3, given initial volume $H_0$, desired concentration $c$, and fluxes $L_{in}$ and $Q_{in}$. We can also express $t^*$ in terms of the concentration $C_L$. The required time is:

$$t^* = \frac{H_0}{Q_{in}}\frac{1}{(C_L/c - 1)} \tag{18}$$
In practice, $C_L$ is determinated by the model of tank #1 (equations (12)-(16)) since $P(t) = C_L$. $H_0$, $Q_{in}$, and $c$ are determined by the design of the bioreactor.

In the second stage, when $t > t^*$, we start pumping out the solution inside the tank. We denote the out-flux by $F_{out}$ (volume/time). This carries a concentration $c$ with it into tank #3.

Since we have an out flux, the volume of water inside the tank changes over time while $L_{in}$ keeps bringing media into the tank. To keep the concentration constant over time, we must introduce a second constant in-flux $H_{in}$ carrying water into the tank. We want to know what $H_{in}$ is.

Since we have an outflux, the volume of water inside the tank changes over time while $L_{in}$ keeps bringing media into the tank. To keep the concentration constant over time, we must introduce a second constant in-flux $H_{in}$ carrying water into the tank. We want to know what $H_{in}$ is.

$$\frac{L(t)}{H(t)} = c \tag{19}$$
Since $c$ is constant, we differentiate (19) to get:

$$L(t)'H(t) = L(t)H(t)' \tag{20}$$
where $'$ denotes derivative respect to time. To solve (20), we derive equations for $H(t)'$ and $L(t)'$.

We start with $L(t)'$. The mass of media inside the tank changes by the in-flux $L_{in}$ coming from tank #1 and by the mass carried by the out-flux $cF_{out}$. Hence:

$$L(t)' = L_{in} - cF_{out}\tag{21}$$
The volume of water inside the tank changes by both in-fluxes $L_{in}$ and $H_{in}$, and the out-flux $F_{out}$. Hence:
$$H(t)' = Q_{in} + H_{in} - F_{out}\tag{22}$$
We substitute (21) and (22) back in (20) and use (19) to get:

$$L_{in} -cF_{out} = c(Q_{in} + H_{in} - F_{out})\tag{23}$$
After rearranging, we get:

$$ H_{in} = \frac{L_{in}}{c} - Q_{in}\tag{24}$$
Equation (24) returns the amount of water per unit time that we must bring into the tank to keep the concentration constant once we start pumping the solution out of tank #2 and into tank #3, given the $c$, $L_{in}$, and $Q_{in}$ (all constants depending on the design of the bioreactor). One can also express $H_{in}$ in terms of the concentration of media $C_L$ by using $L_{in} = C_LQ_{in}$.

Tank #3 - TPA Precipitation
As in tank #2, tank #3 has the same two-stages process.

In the first stage, the tank contains an initial volume $V_N$ of $NH_4OH$ and a constant in-flux $F_{in}$ (volume/time) coming from tank #2. $F_{in}$ carries a concentration $c$ of TPA and EG. We refer to the medium as TPA and EG since we assume it is diluted enough to neglect the medium itself. $c$ was defined in tank #2.

To not to add an extra tank to get rid of the $NH_4OH$ that is pumped out with EG, we must wait until the ratio $\frac{C_N}{C_{H_2O}} = c_2<< 1$, where $C_N$ is the concentration of $NH_4OH$, $C_{H_2O}$ is the concentration of water, and $c_2$ is a parameter smaller than one (similar to $c$ in tank #2). We can assume that no $NH_4OH$ flows to the next tank, provided $c_2$ is small enough. We are also assuming that $C_N$ doesn't change over time during this stage for simplicity, i.e. we do not take into account the chemical reactions going between TPA and $NH_4OH$ until the next stage. As in tank #2, we want to know the time $t^{**}$ we must wait to start pumping the solution out of tank #3.

We write the ratio of concentrations in terms of a ratio of volumes as follows:

$$c_2 = \frac{C_N}{C_{H_2O}} = \frac{M_N}{M_{H_2O}}\frac{V_{H_2O}}{V_N} \Longrightarrow \frac{V_{H_2O}}{V_N} = c_M\tag{25}$$
where $c_M = \frac{M_{H_2O}}{M_N}c_2$.

Hence, we use the same reasoning as with the first stage of tank #2 to calculate the required time $t^{**}$ and obtain:

$$\frac{F_{in}t^{**}}{V_N} = c_M \longrightarrow t^{**} = \frac{V_Nc_M}{F_{in}}\tag{26}$$
Equation (26) gives the required time one must wait to start pumping the solution out of the tank, given the initial volume $V_N$ of $NH_4OH$, the modified-desired concentration $c_M$, and the in-flux $F_{in}$. These are determined by the design of the bioreactor.

During the second stage, when $t > t^{**}$, we start pumping out the solution inside the tank. We denote the out-flux by $G_{out}$ (volume/time). We assume that all TPA reacts with $NH_4OH$ so that $G_{out}$ only carries EG into the next distillation column. By equation (25), we can assume that EG is transported only by water. As with tank #2, we must introduce an in-flux $N_{in}$ (volume/time) of $NH_4OH$ to keep precipitating TPA and to mantain equation (25). We want to know what $N_{in}$ is. Moreover, we also want to know how much di-ammonium salt will be generated after time $t$ by the chemical reaction between TPA and $NH_4OH$.

We start calculating the concentration $C_S$ of the salt (mass/volume) over time. The following model uses the common ODEs-based framework for chemical reactions. For an introduction, refer to Ingalls [1].

The chemical reaction is described by the following formula:


where $C_8H_{12}N_2O_4$ is the di-ammonium salt [4] and $k$ is the reaction rate.

We keep track the three species in equation (27). This leads to three ODEs describing their concentration over time, namely:

$$ C_T(t)' = \frac{c}{2}\frac{F_{in}}{V} - kC_TC_N\tag{28}$$
$$C_N(t)' = \frac{N_{in}[N]}{V} - kC_TC_N\tag{29}$$
$$C_S(t)' = kC_TC_N\tag{30}$$
where $V$ is the volume inside the tank (we take it as constant), $[N]$ is the concentration of $NH_4OH$ carried into the tank by $F_{in}$, $C_T$ is the concentration of TPA, and $C_N$ is the concentration of $NH_4OH$. Notice that in equation (28) we take $\frac{1}{2}\frac{cF_{in}}{V}$ instead of $\frac{cF_{in}}{V}$. According to Palm et al [3], the reaction in equation (2) generates TPA and EG in a 1:1 ratio. Hence, $\frac{1}{2}$ of the concentration $c$ is the concentration of TPA brought by $F_{in}$.

Constants $F_{in}$, $[N]$, and $V$ depend on the design of the bioreactor and $k$ is an experimental parameter. However, $N_{in}$ depends on the dynamics of the fluxes themselves.

We need one more equation for $N_{in}$ to close the system (28)-(30). Following the derivation of equation (24) for tank #2, we start with equation $\frac{C_N(t)}{C_{H_2O}(t)} = c_2$ and differentiate respect to time to obtain:

$$C_N(t)'C_{H_2O} = C_N(t)C_{H_2O}(t)'\tag{31}$$
Recall that the variables in equation (31) are concentrations and not volumes as in equation (20). We are forced to not use volumes since equations (28)-(30) follow from the *Law of Mass Action* [1]. $C_N(t)'$ is already determined by equation (29). To derive an equation for $C_{H_2O}$, we consider the amount of mass that is transported into the tank to obtain:

$$VC_{H_2O}(t)' = F_{in} + N_{in} - F_{out} \longrightarrow C_{H_2O}(t)' = \frac{F_{in} + N_{in} - F_{out}}{V}\tag{32}$$
Substituting (32) and rearranging terms we get:

$$N_{in} = \frac{c_M(F_{in} - F_{out}) + VkC_NC_T}{[N] - C_M}\tag{33}$$
Equation (33) gives the in-flux of water bringing $NH_4OH$ that we need to keep the ratio of concentrations constant.
Hence, we can obtain $N_{in}$ and $C_S(t)$ by solving the system of equations (28)-(30) with (33).

Results


Tank #1
To solve equations (12)-(16), we used MATLAB's ODEs solver *ode45*. Some of the kinetic constants $k$ could not be found in the literature, specifically the ones for the substrate-enzyme reaction in equation (1) and (2). Hence, we took parameters chosen at random and initial conditions $(S_1, C_1, S_2, C_2, P) = (100, 0, 0, 0, 0)$. Figure 1. shows the result of one calculation.


Figure 1. Enzyme kinetics over time. Parameters: $k_1 = 2$, $K_{-1} = 1$, $k_2 = 2$, $k_3 = 2$, $k_{-3} = 1$, $k_4 = 1$,$E_1 (0) = E_2(0) = 5$.

We can see in Figure 1. that the products EG and TPA saturate at 25 units of time at a concentration $C_L = 100$ units of concentration.

Although solutions with parameters chosen at random does not output data we can be entirely confident with, the behavior of the graphs over time gives an idea of how the interactions are taking place and at what timescale. For example, we can see in Figure 1. that the PET concentration $S_1$ and the PET-PETase complex $C_1$ stop mattering after $t \approx 10$ units of time, while the concentration of MHET $S_2$ and the concentration of the MHET-MHETase complex $C_2$ becomes irrelevant after $t \approx 22.5$. This allows one to define windows in which the dynamics of a certain part of the complete enzymatic reaction are negligible and to focus on specific interactions of interest.

Tank #2
To visualize equation (18), we chose $x = \frac{H_0}{Q_{in}}$ and $y = \frac{C_L}{c}$ as variables to obtain a three-dimensional plot. Figure 2. shows the result.


Figure 2. Plot of equation (18). x-axis is the ratio of the original volume inside the tank $H_0$ to the volume introduced by $Q_{in}$. y-axis is the ratio of the concentration of media $C_L$ to the chosen concentration $c$ for equation (19).

Figure 2. shows that as the ratio of concentrations approaches one, the required time increases. Hence, we would like to choose $c << C_L$ to decrease the time as much as possible. The condition is quite intuitive since we must dilute down the media enough to consider it negligible.

Figure 2. also shows that when $0 < \frac{H_0}{Q_{in}} < 5$, the time increases linearly, but for $\frac{H_0}{Q_{in}} > 5$ it stays approximately constant. Although one could argue that the result is unintuitive because we could expect a longer time if $\frac{H_0}{Q_{in}}$ is large, the situation is not so clear. We are looking at the ratio, not at the individual behavior of $H_0$ and $Q_{in}$. We could drastically change both of them without affecting their ratio and we would not see an effect. It suggests that it is a better strategy to study the ratio of volumes instead of the volumes themselves.

$\quad$ We also visualize equation (24) for the behavior of the in-flux $H_{in}$. Figure 3. shows a three-dimensional plot with $c= \frac{1}{50}$.


Figure 3. Plot of equation (24). x-axis is the mass in-flux $L_{in}$ coming from tank #1 with the products (EG, TPA). y-axis is the in-flux $Q_{in}$ of water that is transporting the products.

One can see directly from equation (24) that $H_{in}$ has a linear relation with $Q_{in}$ and $L_{in}$. This is confirmed by Figure 3. Moreover, the angle the plane makes with the x/$L_{in}$-axis is $50$. To calculate the angle, we differentiate equation (24) respect to $L_{in}$ to get:

$$\frac{\partial H_{in}}{\partial L_{in}} = 1/c$$
Equation (34) says that the smaller we pick $c$, the larger we need $H_{in}$ to be to keep equation (19) constant.

Tank #3
Figure 4. shows a plot of equation (26), the time that we must wait to start pumping the solution out of tank #3.

Figure 4. Required time to reach a negligible concentration of $NH_4OH$, given the flux $F_{in}$ and the initial volume in the tank $V_N$. We set $C_M = \frac{1}{50}$ for concreteness.

Figure 4. shows that the smaller the in-flux $F_{in}$ is, the longer we must wait. Similarly, the larger the initial volume $V_N$ is, the longer we must wait.

Instead of solving equations (28)-(30) and (33) simultaneously, we substitute equation (33) into equation (29) to obtain:

$$C_N(t)' = C_MkC_TC_N + \frac{[N]C_M(F_{in} - F{out})}{V([N]-C_M)}$$
Since we do not have access to many of the parameters in equation (35) because they depend on the design of the bioreactor (such as $F_{in}$ and $F_{out}$), we set them to arbitrary constants $n_1$ and $n_2$ to obtain:

$$C_N(t)' = n_2C_TC_N + n_1$$
where $n_1 = \frac{[N]C_M(F_{in} - F{out})}{V([N]-C_M)}$ and $n_2 = C_Mk$.

Similarly, in equation (28) we set $h_1 = \frac{c}{2}\frac{F_{in}}{V}$ to obtain:

$$C_T(t)' = h_1 - kC_TC_N$$
To obtain the concentration of di-ammonium salt $C_s(t)$, we solved equations (30),(36), and (37) together using MATLAB's ODEs solver *ode45*. We assume that $F_{in} > F_{out}$ and $[N] > C_M$ to get $n_1>0$. Figure 5. shows the result of the simulation.

Figure 5. Chemical reactions over time for $n_1 > 0$. Initial conditions $(C_N(0), C_T(0), C_S(0)) = (0,0,0)$.

Figure 5. shows that after $t \approx 3$ units of time, the concentration of TPA is negligible in this example. Moreover, one can see that while $C_T \ne 0$, the concentration $C_N$ of $NH_4OH$ and the concentration $C_S$ of the product grow exponentially. When $C_T = 0$, the grow is linear. An exponential growth is faster than a linear one. Hence, Figure 2. suggests prologing the amount of time $C_T \ne 0$ leads to a higher concentration of the di-ammonium salt in a shorter time.

To calculate the required in-flux $N_{in}$, we use the values calculated for $C_N$ and $C_T$ in the previous simulation. We set equation (33) to:

$$N_{in} = k_1 + k_2C_NC_T$$
where $k_1 = \frac{c_M(F_{in} - F_{out})}{[N] - C_M}$ and $k_2 = \frac{Vk}{[N] - C_M}$.

Then, we evaluate equation (38) with the values for $C_N$ and $C_T$. Hence, $k_1, k_2 > 0$. Figure 6. shows how the function $N_{in}$ looks like with randomly chosen $k_1$ and $k_2$.

Figure 6. In-flux $N_{in}$ over time. $k_1$ and $k_2$ are the parameters for equation (38). The rest of the parameters were used to solve equation (30), (36), and (37).

Figure 6. shows that $N_{in}$ starts increasing, reaches a maximum at $t \approx 1.8$ units of time, and then decreases. The plot suggests $N_{in}$ approaches an asymptote at $N_{in} \approx 6$ units of concentration per unit volume. It is promising to see that such asymptote exists. It means that although $N_{in}$ changes over time by the chemical reactions (equations (28)-(30)), it eventually stays constant like the in-flux described by equation (24) in tank #2. In fact, the asymptote exists even after changing $k_1$ and $k_2$ a bit. Refer to Appendix #1 to see other graphs.

Conclusion


In conclusion, we developed three models to study the enzyme kinetics in tank #1, the media dilution in tank #2, and the chemical reactions in tank #3. All were continuous models (ODEs framework) complemented with classical Michalis-Menten kinetics (tank #1) and chemical reactions (tank #3). The models were solved using MATLAB's ODEs solver *ode45*. Multiple parameters depend on the design of the bioreactor, which we are not concerned with, and some kinetic parameters are not available. Randomly chosen parameters were used. Our solutions and graphs do not show precise numerical values, but the possible behavior the different quantities of interest show given the parameters. We believe these behaviors will provide valuable information for a more accurate design.

Next Steps


After the precipitation of TPA in tank #3, EG is transported into the distillation column. Recall that the concentration of the substance in tank #3 is chosen to neglect the ammonium hydroxide that is used. Moreover, we are also assuming that all TPA reacts. Hence, the only chemical compound transported to the distillation column is EG. To extract EG, we can evaporate the volume of water inside the column. The boiling point of EG is 193 °C while water is 100 °C [5]. Hence, we can raise the temperature to a middle point between them and wait for the evaporation to take place. The next step is modelling the evaporation process. Questions of importance to answer are how long we must wait, what the purity of EG will be, and how the design of the bioreactor could improve these parameters.


References


[1] Ingalls, B. (2013). Mathematical Modeling in Systems Biology. Retrieved from https://www.math.uwaterloo.ca/~bingalls/MMSB/Notes.pdf.
[2] Yoshida, S., Hiraga, K., Takehana, T., Taniguchi, I., Yamaji, H., Maeda, Y., … Oda, K. (2016). A bacterium that degrades and assimilates poly(ethylene terephthalate). Science, 351(6278), 1196 LP – 1199. https://doi.org/10.1126/science.aad6359.
[3] Palm, G. J., Reisky, L., Böttcher, D., Müller, H., Michels, E. A. P., Walczak, M. C., … Weber, G. (2019). Structure of the plastic-degrading Ideonella sakaiensis MHETase bound to a substrate. Nature Communications, 10(1). https://doi.org/10.1038/s41467-019-09326-3.
[4] Rod, V., Bazant, V., Sir, Z. (1972). U.S. Patent No. US3694497A. Retrieved from: https://patents.google.com/patent/US3694497A/en. [5] msds ethylene


Appendices

Appendix #1: Tank 3, Additional Flux Graphs