Math modelling

The whole system can be thought to have 3 modules:

1. Arsenic importer: Independent module
2. Gene expression: Dependent on the output from 1st module
3. Chelation: Dependent on the output from both 1st and 2nd module

Assumptions:

1. One Arsenic Species: For simplicity, let’s assume import rate, chelation rate etc. are same for both “As(III)” and “As(v)”. Therefore, we will not consider “Aox” in this model, because this only lowers the toxicity level of “As” by changing its form and thus providing a less toxic environment for both the bacteria and gut cells.This conversion of As(III) to As(V) can be modeled easily by michales menten kinetics. But here in this model we are more focused with the dynamics of chelation.
2. Phytochelatin gene under Ars promoter: T7 promoter is used just to amplifying the amount of synthetic phytochelatin expression. The system is equivalent to the operon shown here, with a higher protein production rate of phytochelatin.
   

In this section we will model those 3 modules individually.

1. Modeling Arsenic Importer

Approach 1: Using Michales Menten kinetics

Importer is the Enzyme here. Substrate is Arsenic present outside the cell membrane and the Arsenic present inside the cell can be treated as Product of the reaction. Importer binding with Arsenic outside is considered to be a bidirectional reaction. And from Importer Arsenic complex, the formation of Arsenic inside is considered as unidirectional reaction.

\(As_{out} = \) Arsenic outside
\(As_{in} = \) Arsenic inside
\(Imp = \) Importer
\(As_{out}*Imp = \) Importer and Arsenic outside complex
\(K_f = \) Forward reaction constant for the bidirectional reaction
\(K_b = \) Backward reaction constant for the bidirectional reaction
\(K_2 = \) Forward reaction constant for the unirectional reaction

Kinetic equations:

\(As_{out} + Imp \rightarrow As_{out}*Imp\)  with reaction constant \(K_f\)
\(As_{out}*Imp \rightarrow As_{out} + Imp\)  with reaction constant \(K_b\)
\(As_{out}*Imp \rightarrow As_{in} + Imp\)  with reaction constant \(K_2\)

Ordinary Differential Equations (ODEs) to simulate over time:

$$\frac{d[As_{out}]}{dt} = -K_f [As_{out}] [Imp] + K_b [As_{out}*Imp]$$ $$\frac{d[Imp]}{dt} = -K_f [As_{out}] [Imp] + K_b [As_{out}*Imp] + K_2 [As_{out}*Imp]$$ $$\frac{d[As_{out}*Imp]}{dt} = K_f [As_{out}] [Imp] - K_b [As_{out}*Imp] - K_2 [As_{out}*Imp]$$ $$\frac{d[As_{in}]}{dt} = K_2 [As_{out}*Imp]$$

Solutions of those ODEs are the following:

Say an initial concentration of Arsenic outside is provided, and that is \(0.01 M\), and initial concentration of Arsenic inside is \(0.00 M\). The concentration outside in the gut is assumed to remain as initial, over time.

Drawbacks:

Now if we consider the Arsenic concentration inside the gut to be a constant value over time, then the cells will absorb Arsenic forever. But in reality, this is not possible. Arsenic is dangerous for cell, it can’t take Arsenic forever without having any preventive mechanism against arsenic.

Approach 2: Using Diffusion Equation

The flux of Arsenic through the cell membrane is considered here. As diffusion goes along a concentration gradient, it will not allow the cell to take more arsenic once Arsenic concentration inside the cell reaches to the Arsenic concentration outside the cell. That overcomes the previous drawback for Michales Menten model.

Flux of Arsenic through the membrane is \(J = As (Area)^{-1} (Time)^{-1}\) $$J = D \cdot \frac{As_{out}-As_{in}}{dt}$$ where \(D = \) diffusion coefficient, \(d = \) width of cell membrane.

$$\frac{dAs}{dt} = J \cdot A$$ where \(A = \) area covered by Importers on Cell Membrane.

By simplifying the Diffusion Equation, we get: $$\frac{dAs}{dt} = D \cdot \frac{As_{out}-As_{in}}{d} \cdot A$$ $$\frac{dAs}{dt} = R_{in} \cdot (As_{out}-As_{in})$$
where \(R_{in} = \frac{D \cdot A}{d}\)

Consider that \(R_{in}\) is the Arsenic intake rate, through importer. Further assume that \(As_{max}\) is the maximum Arsenic concentration one cell can have, if Arsenic level goes beyond that, the cell dies.
So, the modified equation would be: $$\frac{dAs}{dt} = R_{in} \cdot (As_{max}-As_{in})$$

The graph however has a similar trend like before.

2. Modeling Phytochelatin Expression

Assumptions:

1. We will assume that the phytochelatin production directly depends on the Arsenic concentration inside the cell. As Phytochelatin production is restricted by ArsR ( which represses expression of Ars operon ), we will set a threshold of Arsenic under which Phytochelatin expression would be zero. Once the cell overcomes the threshold Phytochelatin production linearly depends on increasing Arsenic concentration inside the cell. It will contribute to the Gain Term of the equation.
2. Considering the degradation of Phytochelatin, we will need to incorporate a loss term as well to get the exact dynamics. For now, assume there is no chelation. Once we will model the Chelation, we will incorporate another loss term for that.
   

Say the threshold has been set to \(As_{thr} = 0.1 M\)
Phytochelatin production rate \( = \gamma\)
Degradation rate of phytochelatin \( = \alpha\)
As concentration inside the cell \( = [As]\)

$$\frac{d[PC]}{dt} = 0 - \alpha[PC], \quad\quad if [As] \lt As_{thr}$$ $$\frac{d[PC]}{dt} = \gamma ([As]-As_{thr}) - \alpha[PC], \quad\quad if [As] \gt As_{thr}$$

3. Modeling the Chelation

Say one Phytochelatin chelates one Arsenic.

\(As = \) Arsenic
\(PC = \) Phytochelatin
\(PC*As = \) Phytochelatin and Arsenic complex
\(K_f = \) binding constant for forward reaction
\(K_b = \) binding constant for backward reaction

Kinetic equation:

\(As + PC \rightarrow PC*As\)  with binding constant \(K_f\)
\(PC*As \rightarrow As + PC\)  with binding constant \(K_b\)

ODEs to solve through time:

$$\frac{d[As]}{dt} = -K_f [As] [PC] + K_b [PC*As]$$ $$\frac{d[PC]}{dt} = -K_f [As] [PC] + K_b [PC*As]$$ $$\frac{d[PC*As]}{dt} = K_f [As] [PC] - K_b [PC*As]$$

Summarizing all the equations

Case 1: hmt1 was not considered

Change of Arsenic Concentration:

Gain term: Import due to importer + Release of Arsenic from Chelated Arsenic complex due to reverse reaction $$R_{in} (As_{max}-As_{in}) + K_b [PC*As]$$ Loss term: Chelation of Arsenic $$-K_f [As][PC]$$ So, the change of Arsenic Concentration inside the cell through time would be: $$\frac{d[As]}{dt} = R_{in} (As_{max}-As_{in}) + K_b [PC*As] - K_f [As][PC]$$

Change of Phytochelatin Concentration:

Gain term: Production of Phytochelatin when Arsenic concentration inside the cell is above the threshold value \(As_{thr}\) + Release of Phytochelatin from Chelated Arsenic complex due to reverse reaction cell is above the threshold value $$\gamma ([As]-As_{thr}) + K_b [PC*As]$$ Loss term: Loss due to degradation of Phytochelatin. + Loss due to chelation because the portion of Phytochelatin that has already chelated some amount of Arsenic is now unavailable as free Phytochelatin to chelate more free Arsenic $$-\alpha [PC] - K_f [As][PC]$$ So, the change of Phytochelatin Concentration inside the cell through time would be: $$\frac{d[PC]}{dt} = \gamma ([As]-As_{thr}) + K_b [PC*As] - \alpha [PC] - K_f [As][PC]$$

Change of Chelated Arsenic Concentration:

Gain term: Production of chelated arsenic due to chelation $$K_f [As][PC]$$ Loss term: Loss due to reverse reaction $$-K_b [PC*As]$$ So, the change of Chelated Arsenic Concentration inside the cell through time would be: $$\frac{d[PC*As]}{dt} = K_f [As][PC] - K_b [PC*As]$$

Simulating the ODEs we got the following results:

All the simulations were done with the set of following parameters and initial concentrations:

Parameters
\(R_{in} = 0.01 /s\)
\(\gamma = 0.7 /s\)
\(\alpha = 0.01 /s\)
\(K_f = 0.5 /s\)
\(K_b = 0.25 /s\)

Concentrations
\([As]_0 = 0.0 M\)
\([PC]_0 = 0.0 M\)
\([PC*As]_0 = 0.0 M\)
\([As_{thr}] = 0.1 M\)
\([As_{max}] = 0.9 M\)

Img-1

Img-2

The Images above are showing the amount of chelated Arsenic at different time ranges. From Img-1 it can be easily observed that the chelation starts at the very point from where the curve started to leave the Concentration=0.0 axis. At a larger time scale the concentration of Chelated Arsenic inside the cell reaches to an equilibrium (Img-2), when forward chelation reaction is balanced by the reverse reaction.

Img-3

Img-4

These two images are showing the same plots with different x ranges, i.e. with different intervals of time. The Red line shows the one when Arsenic chelation happens inside the cell, the black one is when there is no chelation mechanism inside cell.
In Img-3 one can easily see that Chelation starts from the time point, when the two curves starts to diverge. What happens for the red line, is basically due to chelation, the rate of increase of Arsenic concentration inside the cell, decreases. But as we have considered a reverse reaction in Arsenic Chelation, so once the Chelated Arsenic Concentration starts to reach an equilibrium at a further time point, no more effective chelation happens then and the Arsenic Concentration inside the cell reaches to the maximum value \(As_{max}\)(Img-4). After that no Arsenic is taken by those engineered bacterial cells.

Img-5

Img-6

These 2 images show the free Phytochelatin concentration inside the cell. From Img-5 it can be easily observed that Phytochelatin production starts at the particular time point previously noted. From this time point only chelation also starts. But, again at a larger time scale the Phytochelatin production saturates because its production due to arsenic is balanced by the degradation and chelation.

Discussion:

Well, the system is efficient to chelate Arsenic, but when the cell reaches at an equilibrium then it can chelate no more Arsenic. Now the question is – Can we make these cells more efficient to chelate more amount of Arsenic?

To answer this question we need to introduce one more Module to the previous system. And that is “hmt1”.

Case 2: Introducing hmt1 along with the previous system

What hmt1 does is basically transports already chelated Arsenic to the periplasmic space of the membrane.

Assumptions:

1. For a particular time-range consider the periplasmic space to be an infinite sink.
2. Further assume that the probability of hmt1 to transport the currently available chelated Arsenic from the cell to peripheral space, is 1. Therefore, all the chelated arsenic available inside the cell is readily transported to the periplasmic space without keeping any chance for the reverse reaction to happen.
   

All the ODEs would be same like the previous system, with only one modification that reverse reaction is not present.

Change of Arsenic Concentration:

$$\frac{d[As]}{dt} = R_{in} (As_{max}-As_{in}) - K_f [As][PC]$$

Change of Phytochelatin Concentration:

$$\frac{d[PC]}{dt} = \gamma ([As]-As_{thr}) + K_b [PC*As] - \alpha [PC] - K_f [As][PC]$$

Change of Chelated Arsenic Concentration:

$$\frac{d[PC*As]}{dt} = K_f [As][PC]$$

All the simulations were done with the same set of parameters and initial concentrations like the Case-1, only the reverse reaction is missing here:

Parameters
\(R_{in} = 0.01 /s\)
\(\gamma = 0.7 /s\)
\(\alpha = 0.01 /s\)
\(K_f = 0.5 /s\)
\(K_b = 0.0 /s\)

Concentrations
\([As]_0 = 0.0 M\)
\([PC]_0 = 0.0 M\)
\([PC*As]_0 = 0.0 M\)
\([As_{thr}] = 0.1 M\)
\([As_{max}] = 0.9 M\)

Img-7

Img-8

Both the images represent the same thing at different time ranges. Both the Free Arsenic concentration and the Free Phytochelatin concentration inside the cell reaches to an equilibrium after a couple of damped oscillations.

Reasons behind Oscillation:

From Img-8 it is clearly visible that from a time point greater than zero, the Phytochelatin production starts and thus the chelation. Actually, at this point Arsenic concentration inside the cell overcomes the threshold to switch on Ars promoter.
So after that time point the rate of chelation increases with increasing Phytochelatin concentration, and competes with the rate of import of Arsenic inside the cell. Eventually when these two rates cancel each other, the change of Arsenic inside the cell becomes zero. But the Phytochelatin production is still there.
After the time point when change of Arsenic concentration becomes zero, the chelation rate dominates and free Arsenic concentration decreases which limits the Phytochelatin production again and this time Phytochelatin production rate is balanced by degradation and chelation rates. So, net change of Phytochelatin inside cell becomes zero.
Then the Free Phytochelatin concentration inside the cell starts to decrease following the decrease of arsenic concentration. But due to the decrease of Phytochelatin concentration inside the cell leads to the decrease of chelation rate, and again the import rate of Arsenic starts to dominate. Change of Arsenic concentration inside the cell again becomes zero, followed by another increase. Phytochelatin concentration again starts to increase following the trend of Arsenic concentration.
Thus after 2-3 damping oscillation the system reaches to an equilibrium. But here the Arsenic concentration inside the cell remains very less than the \(As_{max}\) which is a safe situation for the cell, and also leading the cell to have a nonzero Arsenic absorption rate for a very large amount of time, until the whole periplasmic space and the cell become filled with chelated arsenic.

Img-9

Img-9 represents the change of free Phytochelatin concentration with the change of free Arsenic concentration inside the cell. And it converges to a point, that represents the equilibrium. And the oscillation is represented by the spiral curve.

Img-10

Img-10 represents the concentration of chelated Arsenic inside the cell. We can see this is a monotonically increasing plot. That tells that the chelation continues inside the cell through a larger time span, and it will continue until the cell becomes filled up with chelated Arsenic.

Thus, introducing hmt1 makes the system more efficient to chelate arsenic through a larger amount of time span. This thing we can theoretically guess also. Therefore, results of the mathematical model assure about our expectation about the system and also it provides an insight to increase chelation efficiency of the engineered cell by introducing hmt1 into the system.

Answer to some of the probable questions we expect from the system

1. How much number of Bacteria is needed to Absorb “X” amount of Arsenic, in given time “T”? (T is the time until the bacterial cell gets fully filled with chelated Arsenic)
   
   \(\rightarrow\) Our bacteria had a nonzero absorption rate of As, after the system gets an equilibrium. At equilibrium, the rate of Arsenic absorption \(= 0.00786774993215 \approx 0.008 M/s\).
   Let us assume number of Bacteria required is \(N\), \(N*0.008*T = X\)
2. How better is the Engineered Bacteria to Absorb As, in competition with the other Cells in Gut?
   
   \(\rightarrow\) Short Time Scale, before Equilibrium is reached.
   For Normal Cells, increase of Arsenic inside cell will eventually decrease the rate of Arsenic uptake. For the engineered bacteria the increase of Arsenic starts to balance by the chelation soon, so the decrease in uptake rate is not that much.
   \(\rightarrow\) Larger Time Scale, after Equilibrium is reached.
   Absorption rate of Arsenic for Normal Cells become 0 as Arsenic concentration inside them reaches to the maximum level. Absorption rate of the engineered bacteria is nonzero till a long time span.

Thus, in competition with the other cells in the gut the engineered bacteria are more efficient to absorb Arsenic.

Insights to design the Probiotic pill

While designing the Probiotic Pill we need to think about these following parameters:

1. Number of bacteria in the pill.
2. Time the pill will spend inside the body.
3. Number of pills one required.

Well, we can derive each of the parameters from the Mathematical Model.

Let us start with the 3rd parameter. To know how much pill is required per day, at first, we need to set a value for the amount of Arsenic we want to chelate from the gut, by one pill. Say this concentration of chelated bacteria by one pill is \(X\). And the concentration of arsenic in the gut is \(C\). Then the total number of pills required to chelate all the Arsenic from the gut would be \(N_p = C/X\) . If one takes 3 pills per day, then he/she needs to take those pills for \(N_p/3\) days.

Say \(T\) is the time required for the bacteria of one pill to become fully filled up with chelated Arsenic inside them. Coming to the 2nd parameter, the time the pill must spend inside the body should not be less than \(T\), (So that one pill can chelate the maximum amount it can) but the time it spends inside the body must not be greater than \(T\). Otherwise those cells might burst and spread Arsenic contamination in other parts of the excretory system or might be inside the intestine.

At last, we are going to set the first parameter, mentioned above, i.e. the number of bacteria required inside one pill. Assuming \(X\) is the amount of Arsenic chelated by each bacterium inside the pill. Further assume \(t \lt T\) is the time each bacterium spends inside the gut. We know that at equilibrium the bacterium has an absorption rate of \(0.008 M/s\). (This value we got from our simulations)
Let us assume that the number of bacteria in the pill is \(N\). $$N*0.008*T = X$$ $$N = \frac{X}{0.008*T}$$ Thus, our Math Model provides an very clear idea about the value of the parameters we need to consider while designing the probiotic pill.