Team:SYSU-China/file/Project/Model/Chemicaldynamicalmodel.html

From 2014.igem.org

Revision as of 14:56, 11 January 2015 by JX (Talk | contribs)

Chemical Dynamics Model

Introduction

According to the design of experiment, how the system works from the chemical dynamics respective can be described using ordinary differential equations (ODEs) and statistics mechanics.

By solving the differential equations numerically, the replication of M13 phages can be simulated quantitatively. Moreover, how the interaction between antibody-antigen influence the number of progeny, which is a very important in evolution system, can be investigated quantitatively. Surprisingly, by introducing the selectable binding energy interval of two-hybrid system, such investigations not only unravel an unexpected limitation of two-hybrid, but also provide some practical suggestions for improving the function of our system.

Symbols and Definitions

<a class="fancybox" rel="group" href="SYSU_M_CKM_T1.png"><img src="SYSU_M_CKM_T1.png" alt="" /></a>

Model Analysis

To simulate the process of M13 replication and budding, many data about M13 is required. Limited by time and capacity, it’s impossible for us to investigate all the data needed in our experiment. As a result, data used in this model partly comes from our experiment, partly comes from related materials, and partly comes from crudely estimation.

Assumptions

To modeling with the help of chemical dynamics and statistics mechanics, several assumptions are required.
1)Since it takes some time for the synthesis of molecule, many process in biological system are time-discrete, such as the number of DNA replicated, the number of protein synthesized. In order to apply the ordinary equations, the effect of time-discrete is neglected.
2)Different biochemical process may happens at different part of the bacteria and the distribution of molecular in the bacteria is probably uneven. Also, it takes some time for a chemical reaction to reach the equilibrium state. For simplification, all of these are neglected in this model.

Derivation of Model and Estimation

This chemical dynamics model simulates the process from the replication of RF DNA to the budding of progeny phages

Replication of RF DNA

Double stranded RF DNA can generate ssDNAs by rolling circle replication according to pervious research. The newly generated ssDNA would be cyclized by P2 and replicate itself. A large amount of RF DNAs would stop their replication task and combined with P5 for further envelope when the copy of RF DNAs reaches 100-200. Thus the replication process of RF DNA can be divided into two steps.

According to the mechanism of rolling circle replication, if all the materials and DNA polymerase are abundant, the number of RF DNA will increase following Fibonacci Series. However, there are only 10-20 copies of DNAP III in E.coli. The lacking of DNAP III directly limits the RF DNA’s replication rate. Since the copy of DNAP III is generally a constant, we assume that the number of DNAP III that participated in rolling cycle replication (denote as numDNAPIII) is a constant, which is estimated to be 5. As for the replication rate, the polymerizing rate of DNAP III is around 1000bp/s and the M13 bacteriophage’s genome is around 7000bp when inserted by some functional genes. Taken other polymerizing steps into consideration, we set DNAP III’s polymerizing speed as 600bp/s for modification. Consequently, the differential equations describing the replication of RF DNA can be obtained:

[Equ.1]

<a class="fancybox" rel="group" href="SYSU_M_CKM_1.gif"><img src="SYSU_M_CKM_1.gif" style="width:200px; heigth:auto;margin-left:150px" alt="" /></a>


Where 2 stands for double stranded RF DNA.

Synthesis of Antibody & Antigen

In our system, λcI-Antigen(Ag in short) and Antibody-α(Ab in short) would be induced to express as RF DNA replicating. The time of induction are denoted as iTAb and iTAg respectively. To be more realistic, the degradation pathways are taken into consideration.

[Equ.2]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E2.gif"><img src="SYSU_M_CKM_E2.gif" style="width:300px; heigth:auto;margin-left:200px" alt="" /></a>

[Equ.3]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E3.gif"><img src="SYSU_M_CKM_E3.gif" style="width:300px; heigth:auto;margin-left:200px" alt="" /></a>


Where pRAb is the number of Ab synthesized in 1 second, dRAb is the ratio of antibody degraded in 1 second. The same for pRAg and dRAg.

The measurement in our experiment gives [Ag]/[Ab]=1.53±0.09. By qualitative theory of differential equation, the equilibrium number of [Ag] is proportional to pRAg. Thus, we set pRAb=1.5*pRAg. Since we fail to measure the concentration of Ag and Ab, the equilibrium number of Ag and Ab in each bacterium are estimated as 100 and 150 respectively. Ag and Ab are induced by one inductor, thus we have iTAg=iTAb

Two-hybrid System

Two-hybrid system is complicated from the view point of chemical dynamics. In our design, λcI gene and Ag gene are spliced and the fusion protein is expressed in E.coli. The subunit λcI can produce dimmers and each λcI has a site that can be bind by Or2, and the subunit Ag can bind with our target protein Ab-α. However, according to assumption 2, the equilibrium states of two-hybrid system can be calculated by introduction statistical mechanics. Specifically, we focus on the number of Ab that are linked with Or2 (whose copy number is 5)

First, consider a general ligand-receptor system. The probability of receptor that bound with ligand is given by Langmuir isotherm formula:

[Equ.4]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E4.gif"><img src="SYSU_M_CKM_E4.gif" style="width:180px; heigth:auto;margin-left:260px" alt="" /></a>


Where [c] is the concentration of ligand, [c0] is the reference concentration estimated to be here, εis the energy released in ligand-receptor combination, kB is Boltzmann constant and T is Kelvin’s temperature. kBT is about 0.026eV.

For dimerization of λcI, the energy released is about 11kcal/mol, or 0.48eV for each combination. Assuming each bacterium has about 100 λcI, the volume of E.coli is 10^-15L. Thus, the probability of λcI bound with another λcI is given by the following formula:

[Equ.5]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E5.gif"><img src="SYSU_M_CKM_E5.gif" style="width:300px; heigth:auto;margin-left:200px" alt="" /></a>


Thus, we can assume that all the λcI in the bacteria are dimerizated. So there are only two states for Or2: occupied by the dimer or not. Similarly, the binding energy between Or2 and dimer is about 12kcal/mol, so the binding probability is about 100%.

For the antigen-antibody interaction, denote their binding energy as ϵ_AA, their binding probability PAA is given by:

[Equ.6]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E6.gif"><img src="SYSU_M_CKM_E6.gif" style="width:180 px; heigth:auto;margin-left:260px" alt="" /></a>

With Equ(5) and Equ(6) we can calculate the number of α-subunit that linked with one Or2:

[Equ.7]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E7.gif"><img src="SYSU_M_CKM_E7.gif" style="width:100px; heigth:auto;margin-left:300px" alt="" /></a>

The background expression of GeneVIII is inevitable. Assume that the transcription rate of GeneVIII is proportional to the number of α-subunit around. Set the background number ofα-subunit (without Ag & Ab) as bgA, background transcription rate as bgPRRG8, and degradation rate of mRNA as dRRG8. We introduce an amplify coefficient that describe the function of two-hybrid system:

[Equ.8]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E8.gif"><img src="SYSU_M_CKM_E8.gif" style="width:120px; heigth:auto;margin-left:290px" alt="" /></a>


Thus, the number of mRNA of Gene8 can be calculated by:

[Equ.9]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E9.gif"><img src="SYSU_M_CKM_E9.gif" style="width:400px; heigth:auto;margin-left:150px" alt="" /></a>

Protein Synthesis and Budding

Gene VIII and gene III are continually expressed while RF DNA replicates itself. When the amount of RF DNAs’ copy reaches a certain level (150 copy in this model), newly generated ssDNA would be combined by P5, preventing it from cyclization. Thus ssDNA will accumulated in the bacterium. Later, when RF DNAs’ copies reaches 190 (assumed in this model), P8 starts to replace P5 on ssDNA, and the envelopment of progeny phages begins. When about 2700 P8 is assembled to a ssDNA (denote it as prototype I, PTI), P3s begin to be attached to the end of ssDNA. After some subsequent envelope process, the progenies are released.

According the life cycle of M13 mentioned above, we can derive the following differential equations:

[Equ.10]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E10.gif"><img src="SYSU_M_CKM_E10.gif" style="width:580px; heigth:auto;margin-left:60px" alt="" /></a>


where pRP8T describes the efficiency of prototype I construction.

[Equ.11]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E11.gif"><img src="SYSU_M_CKM_E11.gif" style="width:520px; heigth:auto;margin-left:90px" alt="" /></a>

Where 5 is the production rate of P8 estimated.

[Equ.12]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E12.gif"><img src="SYSU_M_CKM_E12.gif" style="width:560px; heigth:auto;margin-left:70px" alt="" /></a>


where aRP3T is a coefficient revealing the assembly rate of P3.

[Equ 13]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E13.gif"><img src="SYSU_M_CKM_E13.gif" style="width:420px; heigth:auto;margin-left:140px" alt="" /></a>


Where eTP3 is the time GeneIIII starts expression. Here, the assembly of P3 represents the subsequent process of envelopment, so the number of progeny released can be calculated by:

[Equ 14]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E14.gif"><img src="SYSU_M_CKM_E14.gif" style="width:240px; heigth:auto;margin-left:230px" alt="" /></a>


Result

Basic Properties

Antibody-antigen interaction is measured by their binding energy bEAA. We hope our design can select the mutated progeny that carry antibody with higher binding energy. As a result, we focus on how binding energy influences the replication, protein synthesis and budding process of M13 phages. Fig. 1 show how the number of each molecule changes with respect to time, with binding energy bEAA=15.5.

<a class="fancybox" rel="group" href="SYSU_M_CKM_BHB1.png"><img src="SYSU_M_CKM_BHB1.png" alt="" /></a> <a class="fancybox" rel="group" href="SYSU_M_CKM_BHB2.png"><img src="SYSU_M_CKM_BHB2.png" alt="" /></a> <a class="fancybox" rel="group" href="SYSU_M_CKM_BHB3.png"><img src="SYSU_M_CKM_BHB3.png" alt="" /></a>

<p1> Fig. 1 Antibody with Low binding energy </p1>

Fig. 1-B shows that the number of mRNA of GeneVIII rapidly increases till a new equilibrium concentration after induced by two-hybrid system, resulting a much large production rate of P8, as shown in Fig. 1-C. However, the production rate of P8 is not enough for the prototype I construction, suggested by the rapid descent of P8 and low new equilibrium concentration. Consequently, we can infer that the production rate of P8 is the dominant factor that constrains the number of progeny.

In contract, the main constraining factors changed when bEAA increases, as shown in Fig. 2. Here, the relatively low production rate of P3 impedes the later envelopment, as shown in Fig. 2-D.

<a class="fancybox" rel="group" href="SYSU_M_CKM_BLB1.png"><img src="SYSU_M_CKM_BLB1.png" alt="" /></a> <a class="fancybox" rel="group" href="SYSU_M_CKM_BLB2.png"><img src="SYSU_M_CKM_BLB2.png" alt="" /></a> <a class="fancybox" rel="group" href="SYSU_M_CKM_BLB3.png"><img src="SYSU_M_CKM_BLB3.png" alt="" /></a> <p1> Fig. 2 Antibody with High binding energy </p1>

The only difference between the parameters used in simulation for Fig. 1and Fig. 2 is the binding energy between antibody and antigen. Simulation successfully illustrates that phages, which carry DNA that encoding high binding energy antibody, can produce more progeny, which is right the assumption in Population Growth Model. Combined with the result of population growth model, we can draw a conclusion that antibody can be improved in each cycle of experiment.

Binding Energy & Progeny Number

To find the relation between binding energy bEAA and progeny number is the main reason for building this model. Here, progeny number of bEAA from 5 to 20kBT and t=1800s are shown as follows:

<a class="fancybox" rel="group" href="SYSU_M_CKM_Fig3.jpg"><img src="SYSU_M_CKM_Fig3.jpg" alt="" /></a> <p1> Fig.3 </p1>

For bEAA∈[10,18], the increase of binding energy will significantly increase the number of progeny. However, bEAA that out of [10,18] has no evidently influence on the number of progeny. Fig. 3 unravels a very important property of two-hybrid system: two-hybrid system can only discriminate protein-protein interaction whose binding energy lies in a certain interval. We denote this interval as selectable interval of two-hybrid system and the difference between its upper limit and lower limit as selection range.

Production Rate & Selectable Interval

Equ(6) describe the relation between binding probability, number of antibody and binding energy. Here, different production rate of antibody and antigen are applied in simulation while [Ab]=1.5[Ag] is hold.

<a class="fancybox" rel="group" href="SYSU_M_CKM_Fig4.png"><img src="SYSU_M_CKM_Fig4.png" alt="" /></a>

Fig.4 Parameter scan of [Ag]

As shown in Fig. 4, number of progeny is positively correlated to the number of antigen,[Ag]. As [Ag] increase, both the selectable interval and selection range change. To analysis how the interval and range change, we define the number of progeny (denote as numP) as a function of binding energy bEAA:

[Equ.15]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E15.gif"><img src="SYSU_M_CKM_E15.gif" style="width:140px; heigth:auto;margin-left:290px" alt="" /></a>


which can be obtained by interpolation. To quantitatively describe the selectable interval, we define is as the interval of binding energy bEAA that satisfying:

[Equ.16]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E16.gif"><img src="SYSU_M_CKM_E16.gif" style="width:80px; heigth:auto;margin-left:310px" alt="" /></a>


which can be calculated by differentiate the interpolation function. Selectable interval and selection range of different number of [Ag] are shown in Fig. 4 and Fig. 5:

<a class="fancybox" rel="group" href="SYSU_M_CKM_Fig5.png"><img src="SYSU_M_CKM_Fig5.png" alt="" /></a>

Fig. 5

<a class="fancybox" rel="group" href="SYSU_M_CKM_Fig6.png"><img src="SYSU_M_CKM_Fig6.png" alt="" /></a>

Fig. 6

Fig. 5 suggests that selectable interval shift to lower energy as [Ag] increases. So, when applying this system to evolve antibody with low bEAA, more inductors should be used increase [Ag] and [Ab] to make sure that the original bEAA is located in the selectable interval. In contract, to further improve an antibody, less inductor should be applied, so the selectable interval will shift to higher energy as required by high original bEAA. Besides, the selection range will shrink when [Ag] or [Ab] is either too large or too small, so more gradient experiments are required for tuning the parameters and thus optimizing our system.

Background Expression & Selectable Interval

Previously, background expression of GeneVIII is taken into consideration. Now we focus on how background expression influence the selection range and selectable interval. To do this, we should rewrite Equ (8) and Equ(9) first.

[Equ. 17]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E17.gif"><img src="SYSU_M_CKM_E17.gif" style="width:120px; heigth:auto;margin-left:290px" alt="" /></a>

[Equ.18]

<a class="fancybox" rel="group" href="SYSU_M_CKM_E18.gif"><img src="SYSU_M_CKM_E18.gif" style="width:380px; heigth:auto;margin-left:160px" alt="" /></a>


where uPRRG8 is the transcription rate when there are 1 α-subunit around GeneVIII. For consistency, uPRRG8 is tuned to have the same result with the previous result (bgA=0.1).

<a class="fancybox" rel="group" href="SYSU_M_CKM_Fig7.png"><img src="SYSU_M_CKM_Fig7.png" alt="" /></a> <p1> Fig.7 Parameter scan of background number of alpha-subunit of RNAP </p1>

Fig. 7 shows that the increase of background expression narrows the differences of phages with different ϵ¬_AA, which will strongly weaken the selection ability of our system. Moreover, the selectable interval is evidently shrunk by the increase of background expression, especially at high energy limit. Details are shown in Fig. 8 and Fig. 9.

<a class="fancybox" rel="group" href="SYSU_M_CKM_Fig8.png"><img src="SYSU_M_CKM_Fig8.png" alt="" /></a>

Fig.8

<a class="fancybox" rel="group" href="SYSU_M_CKM_Fig9.png"><img src="SYSU_M_CKM_Fig9.png" alt="" /></a>

Fig.9

Conclusion

By simulation the replication, protein synthesis and budding process of M13 phages, this chemical dynamical model investigate the properties of the system we design from a microscopic viewpoint. As illustrated by the model, phages with DNA that encoding high binding energy antibody can produce more progeny, which support the assumption of Population Growth Model.

Moreover, selectable interval is unraveled by our model unexpectedly, showing a very important limitation on the discriminating ability of two-hybrid system. Motivated by it, investigation on how [Ag], [Ab] and background expression influence selectable interval provide some suggestions for further experiments:

1) More inductor should be used when our system is applied to evolve antibody with low initial binding energy, while less inductor for the further improvement of antibody.
2) Background expression strongly impedes the discriminating ability of two-hybrid system and narrows the selectable interval. So it should be suppressed as much as possible.