Team:SJTU-BioX-Shanghai/Modelling

Home

The team

• Introduction
• Members
• The album

The project

• Idea & Background
• Project design
• Network
• Modelling
• Protocols

Results

• Overview
• Expanding study
• Parts & Devices
• Characterization
• Human practice

Calendar

• Events
• Lablog

Others

• Tools for Wiki
• Biological safety
• Team collaboration
• Acknowledgements

Contents 1 Motivation 2 Introduction 3 Assumptions and hypothesis 4 Model development 4.1 Single gene expression 4.2 relE cleavage process 4.3 relB-Lon hydrolysis 4.4 relE and relB polymerization 4.5 Modeling the full system 5 Simulation 5.1 oscillate with periodicity 5.2 damped oscillations 5.3 equilibrium states 6 Conclusion 7 References & Attachment

Motivation

From the background of our system, one can see that our model is designed to regulate the bio-clock of E.coli. However, this system is more complicated than one can think in that it consists of a variety of reactions, which is quite hard for us to analysis without mathematics. So it is indispensible to justify that the concentration of relE will oscillate with a specific periodicity. Moreover, understanding the important parameters of the system allows one to choose appropriate promoters and rbses to make our system work. Therefore, we have created the mathematical models to prove our project is feasible, to guide our experiment, and to predict the length of E.coli’s sleeping time.

Introduction

First of all, let's pick out the three key genes, which names relE, relB and Lon, from the original network. Then the interaction between them can be expressed as a graph as follows.

Actually, there are four types of reactions in these systems, which are highlighted in different colors.

• The arrow marked in green is a single gene expression
• The arrow marked in yellow is the relE cleavage process
• The arrow marked in blue is the relB-Lon hydrolysis
• the arrow marked in purple is the relE and relB polymerization

From the diagram mentioned above, obviously almost each two genes or its mRNAs will have a bio-chemical reaction.

Assumptions and hypothesis

1. The number of ribosomes in E.coli is sufficiently large at all times
2. The cleavage reaction caused by relE is an enzyme process
3. Any protein and its related mRNAs generate at a constant rate
4. Either protein or mRNA degrades in direct ratio to its concentration

In order to make the quantitative analysis, we should simplify our model by making the above-mentioned assumptions. The first one is the vital assumption in our system. Since a cell may have many thousands of ribosomes, the number being positively correlated with growth rate^[1], we consider this assumption reasonable. In fact we can overlook the opposite situation, which will lead us to consider the piecewise ODEs (ordinary differential equations).Unfortunately, such kind of equations cannot be solved by Matlab. The second assumption is supported by The Bacterial Toxin RelE Displays Codon-Specific Cleavage of mRNAs in the Ribosomal A Site^[2] while the last two assumptions are illustrated by Stochastic simulations Application to molecular networks^[3].

Model development

Single gene expression

According to assumption 3 and 4, we have made, we will have the following six equations:

Where α₁, α₃and α₅ represent the strength of promoters. α₂, α₄ and α₆represent the strength of rbs. β_i(i=1,2...6) represents the degradation rate of the associated proteins or mRNAs.

relE cleavage process

Based on assumption 1 and 2, we choose the classical Michaelis-Menten equations to describe the relE cleavage process:

Where k_m is the Michaelis constant, which is only related with the affinity of the enzyme (relE) for the substrate. k_cati(i=1,2,3) is related with the cleavage rate of relE.

relB-Lon hydrolysis

For relB-Lon hydrolysis, we have the following equations:

Where k_m is the Michaelis constant, which is related with the affinity of the Lon. k_cat0 is related with both relB and Lon.

relE and relB polymerization

For relE and relB polymerization, we considered it as the second-order reaction:

where k represents the reaction coefficient.

Modeling the full system

Based on the former analysis, we can generalize our model into the following six ordinary differential equations:

The following table shows the description of the parameters.

However, since we proposed these equations, we began to realize some problems for us before we can get satisfactory results from them:

1. From the table above, we can see that there are heaps of parameters for us to find out. And some of them are even impossible for us to get.
2. We are only interested in the amounts of some proteins, mainly the relB and relE. However, there are many other things, such as the various mRNA, for us to cope within the equations.
3. The complexity of the relations among different elements adds great difficulties for us to calculate and analyze.

So, we decided to only focus on the relations between two of the three proteins, which mostly interests us.

Having the above notion in our minds, we proposed our simplified new model. We keep the majority of the terms in our new model, except the relE cleavage term. Since we are not pretty sure about the principle of relE cleavage process, we choose the classical hill equations. Now, we can write down our new model in the following three ordinary differential equations:

Now, we can see the new model is much simpler for Matlab to cope with. However, we still should find the values of some parameters from a variety of databases. They are listed as below:

Since most parameters can be found in several papers, what we should do is to choose the appropriate parameters, such as the promoter strength and the degradation rate, to make the system oscillate. Fortunately, Matlab can handle all these works, which will be presented in the next section.

Simulation

At the beginning of the simulation part, let's replace the parameters with simple characters.

• a: relE promoter strength
• b: relE degradation rate
• c: Reaction coefficient(relE and relB)
• d: relB promoter strength
• e: relB degradation rate
• f: reaction coefficient(Lon and relB)
• g: Lon promoter strength
• h: Lon degradation rate

Now, we are ready to calculate the right value of these eight parameters by Matlab. Remember our target is to select the proper parameters to make the system oscillate. If so, the bacteria will fall asleep when relE rises to a high level and the bacteria will be active when relE falls down to a low level.

We set the following value for the parameters, to see if they satisfy our requirement.

a = {11，16，21，26，31};b = {1，6};c = {1，6，11};d = {21，26，31};

e = {1，6};f = 46;g = {21，26，31};h = {1，6}

At first we find that the majority of the systems with specific parameters are asymptotic stable. In order to figure out the reason why they cannot show us oscillations, we consider the oscillation criterion.

How to judge a system that cannot oscillate?

With the help of the criterion, we can remove many situations which will lead to asymptotic stable. This will reduce the burden of matlab to some degree. Matlab tells us the following values parameter value can help the system to reach an oscillation state.

And these parameter value listed below will lead to a damping oscillation state.

After several minutes, one can see three situations with specific parameter values.

oscillate with periodicity

damped oscillations

equilibrium states

Conclusion

Since we have three types of promoters (j23110,j23118,j23116) with their relative strength 1：0.6：0.07, we are now able to choose the proper one to make our system oscillate. All the well-performing situation are listed as below:

From these pictures, it is obvious that a larger value of “a” means a higher probability of oscillation. In addition, a less value of “d” means a higher probability of oscillation. Moreover, the value of “g” plays an important role in oscillation system. We find that when we set the value of “g” lower than 13, then this system will fall into chaos. When we set the value of “g” from 18 to 22, we observed the oscillations. When we set the value of “g” higher than 25, the system will reduce to damping oscillations and reach the equilibrium state immediately. In conclusion, we can set g=22 to have the best effect. According to our former discussion, we are now able to choose the appropriate promoters to make our system oscillate.

References & Attachment

[1] Michael T. Madigan, John M. Martinko, Brock biology of microorganism 11th edition, Pearson Education International 2006, ISBN:0-13-196893-9

[2] Kim Pedersen, Andrey V. Zavialov,Michael Yu. Pavlov, Johan Elf,Kenn Gerdes, and Mans Ehrenberg. The Bacterial Toxin RelE Displays Codon-Specific Cleavage of mRNAs in the Ribosomal A Site. Cell, Vol. 112, 131–140, January 10, 2003

[3] Didier Gonze, Stochastic simulations Application to molecular networks March 27, 2007

---

Since this site doesn't support Tex format for formulas, we provide our
source document File:SJTU09 Modeling.zip for you to download.