Team:Grenoble/Modeling/Amplification/ODE

From 2012.igem.org

(Difference between revisions)
(Created page with "{{:Team:Grenoble/Templates/Modeling}} <html xmlns="http://www.w3.org/1999/xhtml" xml:lang="en" lang="en"> <body id="Modeling"> <div id="cadre"> <section> </section> </div> </bod...")
 
(76 intermediate revisions not shown)
Line 4: Line 4:
<body id="Modeling">
<body id="Modeling">
<div id="cadre">
<div id="cadre">
 +
<section>
<section>
 +
<h1>Preliminary</h1>
 +
We will use the quasi steady state approximation (QSSA) then. The idea is that there are quick reactions, such as enzymatic ones, complexations, etc… And there are slow reactions such as protein production. We assume that the evolution speed of an element that is created only by quick reaction is null.
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/c/c6/Qssa_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
Indeed, we have these types of evolution for the biological elements. The ones involved only in quick reactions are most of the time in a steady state, and there jump from one steady state to an other has an infinite speed, which does not interest us.
 +
</section>
 +
<section>
 +
<h1>Goal</h1>
 +
In this part, we want to answer to three questions:
 +
</br>
 +
</br>
 +
 +
<div style="margin-left: 30px;">
 +
<a href="#1"><img src="https://static.igem.org/mediawiki/2012/4/49/1_mod.png" alt="" /> What is the sensitivity of our system?<br/>
 +
<a href="#2"><img src="https://static.igem.org/mediawiki/2012/1/1e/2_mod.png" alt="" /> What is the time response?
 +
<br/>
 +
<a href="#3"><img src="https://static.igem.org/mediawiki/2012/5/57/3_mod.png" alt="" /> What steady states will our system always reach?
 +
</div>
 +
 +
</br>
 +
</section>
 +
<section>
 +
<h1>The system</h1>
 +
</br>
 +
</br>
 +
<b>We give the explanations on all the parameters used here in the <a href="https://2012.igem.org/wiki/index.php?title=Team:Grenoble/Modeling/Amplification/Sensitivity">parameters part</a>.</b>
 +
</br>
 +
</br>
 +
Here is the schema of the real system, in orange are the reactions which did not appear in the simplified system of the overview:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/8/82/Schema_system_grenoble.png" alt="" /></center>
 +
</br>
 +
cAMP is the quorum sensing molecule. When we put some cAMP out of the system, it enters into the system. Then, it complexes with CRP to create (CRP-cAMP), which is the transcription factor of the gene arac. When some  AraC is created, it will complex with arabinose to create AraC*. AraC* is the active form of Arac. Arac* with (CRP-cAMP) are the transcription factors of the gene cya. Then when some protein of adenylate cyclase is produced, it will catalyze the production of  cAMP.
 +
</br>
 +
</br>
 +
We first have the complexation of cAMP with CRP to get (CRP-CAMP). It is modeled by this biological reaction:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/9/96/Complexation_crp_camp_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
Thus, we get the evolution speed of (CRP-CAMP), r<SUB>(CRP-cAMP)</SUB>:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/0/0e/Eq2_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
Then, we have the conservation equations:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/0/03/Eq3_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
Thus, we get:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/d/de/Eq4_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
Eventually, we get:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/d/d6/Eq5_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
where <img src="https://static.igem.org/mediawiki/2012/2/22/Eq71_grenoble.png" alt="" />
 +
</br>
 +
We have : <br/>
 +
<br/>
 +
<center><img src="https://static.igem.org/mediawiki/2012/f/fc/Eq7_grenoble.png" alt="" /></center>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/b/b7/Eq8_grenoble.png" alt="" /></center>
 +
</br>
 +
Eventually, we have to choose between the two solutions
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/e/e7/Eq9_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
These two solutions are positive, because:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/d/d4/Eq10_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
To chose which solution is the good one, we know that:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/1/13/Eq11_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
<img src="https://static.igem.org/mediawiki/2012/e/eb/Eq12_grenoble.png" alt="" />
 +
</br>
 +
</br>
 +
Eventually, we get:
 +
</br>
 +
</br>
 +
<div class="encadre" id="red">
 +
<center><img src="https://static.igem.org/mediawiki/2012/8/85/Eq13_grenoble.png" alt="" /></center>
 +
</div>
 +
</br>
 +
</br>
 +
(CRP-cAMP) is the transcription factor of the gene <i>araC</i>. When it appears in the network, it activates the production of the protein AraC. This is modeled by a Hill function. In addition, there is some outflow linked to the promoter paraBAD, which is the promoter regulating <i>araC</i>, thus there is a basal production of AraC. We take into account this basal production, because we need to know if because of them our system will always be turned on, thus useless. AraC is also naturally degraded by the bacterium. Thus, we get as the equation of evolution of AraC concentration:
 +
</br>
 +
</br>
 +
<div class="encadre" id="red">
 +
<center><img src="https://static.igem.org/mediawiki/2012/8/84/Eq14_grenoble.png" alt="" /></center>
 +
</div>
 +
</br>
 +
</br>
 +
Then, the protein AraC complexes with arabinose to create AraC active, written AraC*. It is modeled by the following chemical equation:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/a/a8/Complexation_arac_arabinose_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
We get the evolution of AraC* r<SUB>AraC*</SUB>:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/8/88/Eq15_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
With the QSSA, we get:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/e/e2/Eq16_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
In addition we have the conservation equation of AraC:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/d/d1/Eq17_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
We assume that we have <img src="https://static.igem.org/mediawiki/2012/8/82/Eq18_grenoble.png" alt="" />
 +
</br>
 +
</br>
 +
<div class="encadre" id="red">
 +
<center><img src="https://static.igem.org/mediawiki/2012/3/32/Eq19_grenoble.png" alt="" /></center>
 +
</div>
 +
</br>
 +
</br>
 +
<img src="https://static.igem.org/mediawiki/2012/e/eb/Eq20_grenoble.png" alt="" />
 +
</br>
 +
</br>
 +
Then, AraC* with (CRP-cAMP) are the transcription factors of  of the gene ca. When they appear in the network the protein Ca is produced. The product of two hill functions models this. For the same reasons as for AraC we take into account the basal production of the adenylate cyclase. In addition it is degraded by the bacterium.
 +
Thus, we get the equation:
 +
</br>
 +
</br>
 +
<div class="encadre" id="red">
 +
<center><img src="https://static.igem.org/mediawiki/2012/5/55/Eq21_grenoble.png" alt="" /></center>
 +
</div>
 +
</br>
 +
</br>
 +
Eventually, Ca catalyzes the production of cAMP with ATP. We have the following enzymatic reaction:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/7/73/Eq22_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
We use the Michaelis Menten model for an enzymatic reaction, thus we get the evolution of cAMP r<SUB>cAMP</SUB>:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/6/6d/Eq23_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
In addition, we assume that [ATP]>>K<SUB>M</SUB>. Eventually we get:
 +
</br>
 +
</br>
 +
<div class="encadre" id="red">
 +
<center><img src="https://static.igem.org/mediawiki/2012/a/a7/Eq24_grenoble.png" alt="" /></center>
 +
</div>
 +
</br>
 +
</br>
 +
In addition, the bacterium naturally degrades cAMP. Finally there is a quorum sensing term. However, we are not modeling the quorum sensing here, thus we assume that [cAMP<SUB>out</SUB>] is constant. Thus, we get:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/6/60/Eq25_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
The QSSA enables us to have r<SUB>cAMP</SUB> = 0. Then, we have:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/3/33/Eq26_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
<span style="text-decoration:underline;">Remark:</span>
 +
</br>
 +
</br>
 +
We do not take into account in the equations the growth of the bacterium. Indeed, the bacterium grows as long as it has some glucose. And, as long as there is some glucose, the bacterium will not use the arabinose. However, if the bacterium does not use the arabinose, the protein AraC can not be actived, and thus no adenylate cyclase is produced. The bacterium begins to use the arabinose when the whole glucose has disappeared. But it does not grow with the arabinose.
 +
Indeed, the biologists in order to check the “AND gate” behavior, the biologists built, see protocol <a href="https://2012.igem.org/Team:Grenoble/Biology/Protocols/AND_test">protocol "AND gate test"</a> . Here we give the biological graphs of the absorbance and the graph of the RFU  in function of the time for arabinose and cAMP maximum:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/2/23/Graphe1_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/a/ad/Graphe2_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
We can observe that until 100 minutes the bacterium is growing (the absorbance increases), and we do not get any RFU signal. Then, the bacterium almost stops growing, and thus we begin to get a signal.
 +
</section>
 +
<section>
 +
<h1 id="1"><img src="https://static.igem.org/mediawiki/2012/4/49/1_mod.png" alt="" /> Sensitivity </h1>
 +
</br>
 +
</br>
 +
The first question was: what is the sensitivity of our system?
 +
</br>
 +
</br>
 +
Even if we do not know the exact value of all the parameters, we have enough information on them to be able to have a good evaluation of the sensitivity of our system.
 +
</br>
 +
</br>
 +
To answer this question, we plot the evolution of the steady state of adenylate cyclase in function of [cAMP<SUB>out</SUB>], the initial concentration of cAMP, which is assumed to be constant. <a href="https://static.igem.org/mediawiki/2012/f/f1/Resol_Eq_cAMP_grenoble.zip">Here</a> are the scripts that enable us to plot the graphs. We solve the differential equations to get the steady state, because if we wanted to solve a set of equations we would have had to give an initial point. If we had given 0, matlab would have stayed at this point, and we could not give another initial point without solving the equations.
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/f/f3/Graphe3_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
<span style="text-decoration:bold;">Conclusion:</span>
 +
</br>
 +
</br>
 +
The sensitivity of our system is 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span> of initial cyclic AMP. When we introduced this quantity in the system, the bacteria will turn on.
 +
Then, the next question is to know when we are under this value in how much time we will be able to observe that one bacteria turned on.
 +
</section>
 +
<section>
 +
<h1 id="2"><img src="https://static.igem.org/mediawiki/2012/1/1e/2_mod.png" alt="" /> Temporal evolution</h1>
 +
</br>
 +
</br>
 +
To evaluate the time it will take to be able to detect a signal, we need to plot the evolution of the adenylate cyclase in the time for an initial concentration of cAMP<SUB>out</SUB> ≥ 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span>. We first give the graph with cAMP<SUB>out</SUB> = 10<SUP>-3</SUP> mol.L<span class="exposant">-1</span> :
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/e/e6/Graphe4_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
We notice that after 200 minutes, we will be able to begin to detect a signal.
 +
</br>
 +
</br>
 +
Then, we want to see the behavior of the system around the threshold. We give the evolution of the adenylate cyclase in the time in function with cAMP<SUB>out</SUB> = 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span> :
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/8/8c/Graphe5_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
Here, to be able to begin to detect a signal, we should wait around 1300 minutes. So even if our system can detect cAMP<SUB>out</SUB> = 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span>, we may not be able to say if it is a real detection or a false positive. We will be able to answer this question with the stochastic part.
 +
</br>
 +
</br>
 +
Then we give one temporal evolution of the adenylate cyclase for cAMP<SUB>out</SUB> = 10<SUP>-8</SUP> mol.L<span class="exposant">-1</span>. It is bellow the threshold, but because of the basal values, we want to see exactly what happens.
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/1/13/Graphe6_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
We notice that Ca is expressed, and that half of the steady state is quickly reached (around 100 minutes). However the value of the concentration of the adenylate cyclase at steady state is too low to enable us to get a signal.
 +
</br>
 +
</br>
 +
<b>Conclusion:</b>
 +
The more cAMP<SUB>out</SUB> is introduced in the system, the fastest the answer is. The quickest answer would be 200 minutes to reach half of the maximal expression of ca enabling us to get a signal.
 +
Because of the basal values, the adenylate cyclase is always expressed. Thus, we will make a steady state study of the system. This is useful seeing that with the sensitivity graph we could not see the low expression of adenylate cyclase and it is only in the temporal part that we could notice it, so we need a real study.
 +
Then, the stochastic part will be really important to be sure that there will not be too many false positives because of these basal values.
 +
</section>
 +
<section>
 +
<h1 id="3"><img src="https://static.igem.org/mediawiki/2012/5/57/3_mod.png" alt="" />Steady state study</h1>
 +
</br>
 +
</br>
 +
Here, we want to find the steady states of the system, to understand why our system works. At the steady states, the derivatives of our system are null:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/9/92/Eq27_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
If we solve this system of equations, we will get the steady states of the system. Because of the non linearity they are not obvious. However, we need to give an initial point to the matlab solver, and it will return only the closest steady state. Because we want to know all the steady states, we will use a different method.
 +
</br>
 +
</br>
 +
Equation (1) gives us AraC in function of Ca. If we plot this function, we plot all the steady states of AraC, while Ca is varying. Equation (2) gives us Ca in function of AraC. If we plot this function, we plot all the steady states of Ca while AraC is varying. These functions are called the isoclines of the system. The intersections of these two graphs represent the steady states of the system.
 +
</br>
 +
</br>
 +
<a href="https://static.igem.org/mediawiki/2012/9/9e/Steadys_state_study.zip">Here</a> you can find the scripts we worked with in this part. First, I give the isoclines with cAMP<SUB>init</SUB> = 10<SUP>-5</SUP> mol.L<span class="exposant">-1</span>.
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/d/dc/Graphe7_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
Isoclines of Ca and AraC with cAMP<SUB>init</SUB> = 10<SUP>-5</SUP> mol.L<span class="exposant">-1</span>.
 +
</br>
 +
</br>
 +
</br>
 +
In black we have Ca isocline, in purple, AraC one, it will always be the same color code. We notice that there is only one steady state. It is high enough to enable us to get a signal.
 +
</br>
 +
</br>
 +
Then I give the isoclines with cAMP<SUB>init</SUB> = 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span>.
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/4/49/Graphe8_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
Ca and AraC isoclines with cAMP<SUB>init</SUB> = 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span>.
 +
</br>
 +
</br>
 +
It seems that there are two steady state. However, with a zoom around 0, we have:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/0/08/Graphe9_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
Ca and Arac isoclines with cAMP<SUB>init</SUB> = 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span>, zoom around 0.
 +
</br>
 +
</br>
 +
We notice that there is no low steady state, but only a high one.
 +
</br>
 +
</br>
 +
Eventually, I give the isoclines with cAMP<SUB>init</SUB> = 10<SUP>-7</SUP> mol.L<span class="exposant">-1</span>:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/8/85/Graphe10_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
Ca and AraC isoclines with cAMP<SUB>init</SUB> = 10<SUP>-7</SUP> mol.L<span class="exposant">-1</span>.
 +
</br>
 +
</br>
 +
We have the same graph as previously. We zoom around 0 again:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/6/60/Graphe11_ampli_grenoble.png" alt="" /></center>
 +
</br>
 +
Ca and AraC isoclines with cAMP<SUB>init</SUB> = 10<SUP>-7</SUP> mol.L<span class="exposant">-1</span>, zoom around 0.
 +
</br>
 +
</br>
 +
We can still see a high steady state, but this time there also is a low one!
 +
</br>
 +
</br>
 +
To actually understand why our system works, we need to study the stability of our steady states. If they are stable when our system stops to one of these steady states, it will always stay there. It works like this:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/7/79/Stab_point_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
To answer the question of the stability, we follow the classic dynamical system method. We first linearize the system around the steady states. We have:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/8/81/Eq28_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
The equilibrium points are <img src="https://static.igem.org/mediawiki/2012/5/50/Eq29_grenoble.png" alt="" />
 +
</br>
 +
We compute the matrix <img src="https://static.igem.org/mediawiki/2012/b/be/Eq30_grenoble.png" alt="" />
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/4/4e/Eq31_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
Thus, we have
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/2/2e/Eq32_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
Let us call λ<SUB>1</SUB>  and λ<SUB>2</SUB> the eigenvalues of the matrix. Then, we get the solution of the system around the equilibrium points:
 +
</br>
 +
</br>
 +
<center><img src="https://static.igem.org/mediawiki/2012/7/7d/Eq33_grenoble.png" alt="" /></center>
 +
</br>
 +
</br>
 +
We give the computed values of the steady states and of the eigenvalues for:
 +
</br>
 +
</br>
 +
- cAMP<SUB>init</SUB> = 10<SUP>-5</SUP> mol.L<span class="exposant">-1</span> :
 +
</br>
 +
</br>
 +
AraC steady state = 0.167058129527727 10<SUP>-4</SUP> mol.L<span class="exposant">-1</span>
 +
</br>
 +
</br>
 +
Ca seady states = 0.1837444563636 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span>
 +
</br>
 +
</br>
 +
λ<SUB>1</SUB> =  -0.006000000912526 λ<SUB>2</SUB> = -0.005763188664176
 +
</br>
 +
</br>
 +
</br>
 +
- cAMP<SUB>init</SUB> = 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span> :
 +
</br>
 +
</br>
 +
AraC steady state = 0.166879570344986 10<SUP>-4</SUP> mol.L<span class="exposant">-1</span>
 +
</br>
 +
</br>
 +
Ca seady states = 0.1832826298080 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span>
 +
</br>
 +
</br>
 +
λ<SUB>1</SUB>=  -0.006000000910603  &emsp; λ<SUB>2</SUB>= -0.005745344108236
 +
</br>
 +
</br>
 +
</br>
 +
- cAMP<SUB>init</SUB> = 10<SUP>-7</SUP> mol.L<span class="exposant">-1</span> :
 +
</br>
 +
</br>
 +
AraC steady state = 0.182361098919416 10<SUP>-6</SUP> mol.L<span class="exposant">-1</span>
 +
</br>
 +
</br>
 +
Ca seady states = 0.249177541683 10<SUP>-9</SUP> mol.L<span class="exposant">-1</span>
 +
</br>
 +
</br>
 +
λ<SUB>1</SUB>=  -0.006000006994365  &emsp; λ<SUB>2 </SUB>= -0.002117175391388
 +
</br>
 +
</br>
 +
</br>
 +
The steady states are globally stable if and only if the real part of the proper values is < 0 .
 +
</br>
 +
</br>
 +
Their real part are always < 0 ! The steady states  are stable.
 +
</br>
 +
</br>
 +
<b>Conclusion:</b>
 +
Now, we can be sure that our system will not always be turned on. In function of the quantity of initial cAMP, our system will stay at a low or a high steady state.
 +
</section>
 +
<section>
 +
<h1>Conclusion</h1>
 +
The sensitivity of our system is around 10<SUP>-7</SUP> mol.L<span class="exposant">-1</span>. To be able to know if it is satisfying, we need to link it with the signaling part.
 +
In addition, to know if our system is fast we need to link this part with the signaling too. That is what we are going to do in the next part.
</section>
</section>
</div>
</div>

Latest revision as of 03:47, 27 September 2012

iGEM Grenoble 2012

Project

Preliminary

We will use the quasi steady state approximation (QSSA) then. The idea is that there are quick reactions, such as enzymatic ones, complexations, etc… And there are slow reactions such as protein production. We assume that the evolution speed of an element that is created only by quick reaction is null.



Indeed, we have these types of evolution for the biological elements. The ones involved only in quick reactions are most of the time in a steady state, and there jump from one steady state to an other has an infinite speed, which does not interest us.

Goal

In this part, we want to answer to three questions:


The system



We give the explanations on all the parameters used here in the
parameters part.

Here is the schema of the real system, in orange are the reactions which did not appear in the simplified system of the overview:


cAMP is the quorum sensing molecule. When we put some cAMP out of the system, it enters into the system. Then, it complexes with CRP to create (CRP-cAMP), which is the transcription factor of the gene arac. When some AraC is created, it will complex with arabinose to create AraC*. AraC* is the active form of Arac. Arac* with (CRP-cAMP) are the transcription factors of the gene cya. Then when some protein of adenylate cyclase is produced, it will catalyze the production of cAMP.

We first have the complexation of cAMP with CRP to get (CRP-CAMP). It is modeled by this biological reaction:



Thus, we get the evolution speed of (CRP-CAMP), r(CRP-cAMP):



Then, we have the conservation equations:



Thus, we get:



Eventually, we get:



where
We have :



Eventually, we have to choose between the two solutions



These two solutions are positive, because:



To chose which solution is the good one, we know that:





Eventually, we get:



(CRP-cAMP) is the transcription factor of the gene araC. When it appears in the network, it activates the production of the protein AraC. This is modeled by a Hill function. In addition, there is some outflow linked to the promoter paraBAD, which is the promoter regulating araC, thus there is a basal production of AraC. We take into account this basal production, because we need to know if because of them our system will always be turned on, thus useless. AraC is also naturally degraded by the bacterium. Thus, we get as the equation of evolution of AraC concentration:



Then, the protein AraC complexes with arabinose to create AraC active, written AraC*. It is modeled by the following chemical equation:



We get the evolution of AraC* rAraC*:



With the QSSA, we get:



In addition we have the conservation equation of AraC:



We assume that we have





Then, AraC* with (CRP-cAMP) are the transcription factors of of the gene ca. When they appear in the network the protein Ca is produced. The product of two hill functions models this. For the same reasons as for AraC we take into account the basal production of the adenylate cyclase. In addition it is degraded by the bacterium. Thus, we get the equation:



Eventually, Ca catalyzes the production of cAMP with ATP. We have the following enzymatic reaction:



We use the Michaelis Menten model for an enzymatic reaction, thus we get the evolution of cAMP rcAMP:



In addition, we assume that [ATP]>>KM. Eventually we get:



In addition, the bacterium naturally degrades cAMP. Finally there is a quorum sensing term. However, we are not modeling the quorum sensing here, thus we assume that [cAMPout] is constant. Thus, we get:



The QSSA enables us to have rcAMP = 0. Then, we have:



Remark:

We do not take into account in the equations the growth of the bacterium. Indeed, the bacterium grows as long as it has some glucose. And, as long as there is some glucose, the bacterium will not use the arabinose. However, if the bacterium does not use the arabinose, the protein AraC can not be actived, and thus no adenylate cyclase is produced. The bacterium begins to use the arabinose when the whole glucose has disappeared. But it does not grow with the arabinose. Indeed, the biologists in order to check the “AND gate” behavior, the biologists built, see protocol protocol "AND gate test" . Here we give the biological graphs of the absorbance and the graph of the RFU in function of the time for arabinose and cAMP maximum:







We can observe that until 100 minutes the bacterium is growing (the absorbance increases), and we do not get any RFU signal. Then, the bacterium almost stops growing, and thus we begin to get a signal.

Sensitivity



The first question was: what is the sensitivity of our system?

Even if we do not know the exact value of all the parameters, we have enough information on them to be able to have a good evaluation of the sensitivity of our system.

To answer this question, we plot the evolution of the steady state of adenylate cyclase in function of [cAMPout], the initial concentration of cAMP, which is assumed to be constant. Here are the scripts that enable us to plot the graphs. We solve the differential equations to get the steady state, because if we wanted to solve a set of equations we would have had to give an initial point. If we had given 0, matlab would have stayed at this point, and we could not give another initial point without solving the equations.



Conclusion:

The sensitivity of our system is 10-6 mol.L-1 of initial cyclic AMP. When we introduced this quantity in the system, the bacteria will turn on. Then, the next question is to know when we are under this value in how much time we will be able to observe that one bacteria turned on.

Temporal evolution



To evaluate the time it will take to be able to detect a signal, we need to plot the evolution of the adenylate cyclase in the time for an initial concentration of cAMPout ≥ 10-6 mol.L-1. We first give the graph with cAMPout = 10-3 mol.L-1 :



We notice that after 200 minutes, we will be able to begin to detect a signal.

Then, we want to see the behavior of the system around the threshold. We give the evolution of the adenylate cyclase in the time in function with cAMPout = 10-6 mol.L-1 :



Here, to be able to begin to detect a signal, we should wait around 1300 minutes. So even if our system can detect cAMPout = 10-6 mol.L-1, we may not be able to say if it is a real detection or a false positive. We will be able to answer this question with the stochastic part.

Then we give one temporal evolution of the adenylate cyclase for cAMPout = 10-8 mol.L-1. It is bellow the threshold, but because of the basal values, we want to see exactly what happens.



We notice that Ca is expressed, and that half of the steady state is quickly reached (around 100 minutes). However the value of the concentration of the adenylate cyclase at steady state is too low to enable us to get a signal.

Conclusion: The more cAMPout is introduced in the system, the fastest the answer is. The quickest answer would be 200 minutes to reach half of the maximal expression of ca enabling us to get a signal. Because of the basal values, the adenylate cyclase is always expressed. Thus, we will make a steady state study of the system. This is useful seeing that with the sensitivity graph we could not see the low expression of adenylate cyclase and it is only in the temporal part that we could notice it, so we need a real study. Then, the stochastic part will be really important to be sure that there will not be too many false positives because of these basal values.

Steady state study



Here, we want to find the steady states of the system, to understand why our system works. At the steady states, the derivatives of our system are null:



If we solve this system of equations, we will get the steady states of the system. Because of the non linearity they are not obvious. However, we need to give an initial point to the matlab solver, and it will return only the closest steady state. Because we want to know all the steady states, we will use a different method.

Equation (1) gives us AraC in function of Ca. If we plot this function, we plot all the steady states of AraC, while Ca is varying. Equation (2) gives us Ca in function of AraC. If we plot this function, we plot all the steady states of Ca while AraC is varying. These functions are called the isoclines of the system. The intersections of these two graphs represent the steady states of the system.

Here you can find the scripts we worked with in this part. First, I give the isoclines with cAMPinit = 10-5 mol.L-1.


Isoclines of Ca and AraC with cAMPinit = 10-5 mol.L-1.


In black we have Ca isocline, in purple, AraC one, it will always be the same color code. We notice that there is only one steady state. It is high enough to enable us to get a signal.

Then I give the isoclines with cAMPinit = 10-6 mol.L-1.


Ca and AraC isoclines with cAMPinit = 10-6 mol.L-1.

It seems that there are two steady state. However, with a zoom around 0, we have:


Ca and Arac isoclines with cAMPinit = 10-6 mol.L-1, zoom around 0.

We notice that there is no low steady state, but only a high one.

Eventually, I give the isoclines with cAMPinit = 10-7 mol.L-1:


Ca and AraC isoclines with cAMPinit = 10-7 mol.L-1.

We have the same graph as previously. We zoom around 0 again:


Ca and AraC isoclines with cAMPinit = 10-7 mol.L-1, zoom around 0.

We can still see a high steady state, but this time there also is a low one!

To actually understand why our system works, we need to study the stability of our steady states. If they are stable when our system stops to one of these steady states, it will always stay there. It works like this:



To answer the question of the stability, we follow the classic dynamical system method. We first linearize the system around the steady states. We have:



The equilibrium points are
We compute the matrix



Thus, we have



Let us call λ1 and λ2 the eigenvalues of the matrix. Then, we get the solution of the system around the equilibrium points:



We give the computed values of the steady states and of the eigenvalues for:

- cAMPinit = 10-5 mol.L-1 :

AraC steady state = 0.167058129527727 10-4 mol.L-1

Ca seady states = 0.1837444563636 10-6 mol.L-1

λ1 = -0.006000000912526 λ2 = -0.005763188664176


- cAMPinit = 10-6 mol.L-1 :

AraC steady state = 0.166879570344986 10-4 mol.L-1

Ca seady states = 0.1832826298080 10-6 mol.L-1

λ1= -0.006000000910603   λ2= -0.005745344108236


- cAMPinit = 10-7 mol.L-1 :

AraC steady state = 0.182361098919416 10-6 mol.L-1

Ca seady states = 0.249177541683 10-9 mol.L-1

λ1= -0.006000006994365   λ2 = -0.002117175391388


The steady states are globally stable if and only if the real part of the proper values is < 0 .

Their real part are always < 0 ! The steady states are stable.

Conclusion: Now, we can be sure that our system will not always be turned on. In function of the quantity of initial cAMP, our system will stay at a low or a high steady state.

Conclusion

The sensitivity of our system is around 10-7 mol.L-1. To be able to know if it is satisfying, we need to link it with the signaling part. In addition, to know if our system is fast we need to link this part with the signaling too. That is what we are going to do in the next part.