Edinburgh/DivisionPopper/Modelling

From 2007.igem.org

< Edinburgh | DivisionPopper(Difference between revisions)

Latest revision as of 17:27, 17 December 2007

This page contains the mathematical and computational models we constructed in order to analyze the behaviour of the Division PoPper system. You can find the models of three constructs (the proof of concept, the Division PoPper, the composition of the Division PoPper with a counter device), using different approaches to modelling (with Ordinary Differential Equations and with Gillespie simulation on a CTMC). Look at the table of contents below for the section you are interested in.

The proof of concept system

As explained in the design section, the first part of our project consists in the realization of a proof of concept system. This is a smaller and less complex construct that we use to test our main assumption: the DNA region in between the two Dif sites is "flipped" (reversed) at each cell division. The model should help us in answering the following questions:

Is the system working well with standard E. coli division time?
Is the system working well with normal YFP, GFP Expression and degradation rate?
How can we tune these parameters in order to have meaningful tests from the construct in the lab?

Since the mechanisms and the values that characterize the system are sometimes not known or time consuming/difficult to discovery in the lab, we use the models as "in silico" tests. Thus our models should be dependent on (functions of) the following quantities:

Cell division frequency (units: minutes): is the time elapsed from one division to the next one.
Flipping duration (units: minutes): is the time elapsed from starting the flipping to the complete execution.
YFP/GFP Expression rate (units: protein per second): the rate at which fluorescent proteins are expressed.
YFP/GFP Degradation rate (units: protein per second): the rate at which fluorescent proteins degrade.

In order to construct the proper model, we need to identify the dynamics of the system. Since the real biological mechanisms are too complex, several approximations and abstractions are introduced at this step. In this sense, we can simplify the dynamics of the system by observing that our construct can be in one of three states during the time evolution: Forward, Flipping and Backward phase. These three phases happen in an order controlled by cell division. In detail:

Configuration	Processes
	The system starts in this phase. It is the phase in which the promoter is oriented in the 5' to 3' direction, thus controlling the expression of the YFP coding region. In this phase in the cell there is expression of YFP and no expression of GFP.
	The system is in this phase when the recombinase activity has been initialized and not yet concluded. In this phase in the cell we assume there is no expression of YFP and no expression of GFP because the two coding regions are not controlled by the promoter regulatory effect.
	This is the phase in which the promoter is oriented in the 3' to 5' direction, thus controlling the expression of the GFP coding region. In this phase in the cell there is expression of GFP and no expression of YFP.

These are the three states (phases) in which the proof of concept system can be during its time course evolution. The system passes through these phases according to cell division timing and following the pattern: forward phase, flipping phase, backward phase, flipping phase and so on. A picture clarifies better the time course behaviour:

The arrows mark the points in time in which a cell division starts.

Modelling a single cell

As first step, we are interested in investigating the behaviour of the system at the level of a single cell. In particular, we are interested in following the time course of behaviour of one single cell line. This means that in our model we start observing a individual cell behaviour and at each cell division we follow one of the two daughter cells and so on. We make the assumption of only one plasmid present in each cell.

For that we construct two different approaches: using Ordinary Differential Equations and Stochastic Process Algebra. The first is deterministic (the behaviour of the system is defined only by the values of the variables) and treats the population of molecules as a continuous value. It gives us a general understanding of the "average" behaviour of the time evolution of a cell. The second is stochastic (the behaviour of the system is defined by the value of normal and random variables) and treats the population of molecules as a discrete value. It gives us (theoretically) a better understanding of the different trajectories in the behaviour of individual cells and can help in the future for simulating the behaviour of a population.

ODEs (deterministic and continuous)

We define a set of ordinary differential equations that represents the behaviour of a individual cell during time. First, we define a set of equations for each one of the three phases (Forward, Flipping, Backward) and then we construct a unified set of equations by using two support functions. The concentrations in which we are interested are obviously GFP and YFP concentrations. The biological mechanisms (or processes) we need to represent are: gene expression of proteins and degradation of proteins. We assume gene expression of a protein to be in on or off state and, when in the on state, to be expressed at a constant rate. We assume degradation of a protein to be rate dependent on the amount of protein present. Thus the variables are:

YFP expression rate (mol/sec) called YFPexp
GFP expression rate (mol/sec) called GFPexp
YFP degradation rate (1/sec) called YFPdeg
GFP degradation rate (1/sec) called GFPdeg

For each phase the set of equations are:

Phase	Equations
Forward Phase:	d[YFP]/dt = YFPexp - YFPdeg[YFP] d[GFP]/dt = - GFPdeg**[GFP]
Flipping Phase:	d[YFP]/dt = - YFPdeg[YFP] d[GFP]/dt = - GFPdeg**[GFP]
Backward Phase:	d[YFP]/dt = - YFPdeg[YFP] d[GFP]/dt = GFPexp* - GFPdeg*[GFP]

In order to construct a model that represents the sequence of phases we introduce two support functions called respectively ForwardPhase(t) and BackwardPhase(t). They are time depending: ForwardPhase(t) returns value 1 when the system is in the Forward phase and 0 elsewhere, BackwardPhase(t) returns value 1 when the system is in the Backward Phase and 0 elsewhere.

Then the whole system is defined by equations:

d[YFP]/dt = YFPexp*ForwardPhase(t) - YFPdeg*[YFP]

d[GFP]/dt = GFPexp*BackwardPhase(t) - GFPdeg*[GFP]

The data for the expression and degradation rates and the have been taken from peer-reviewed publications as reference values and summarized in the table below:

Data	Values
YFP expression rate	0.0001
GFP expression rate	0.0001
YFP degradation	0.05 (fast because of LVA tag)
GFP degradation	0.05 (fast because of LVA tag)
Flipping duration	5 minutes (estimated)
Division duration	80 minutes (slowed division rate)

We implemented this model using the [http://www.sbtoolbox.org/ Systems Biology Toolbox for MATLAB] and the related code can be downloaded here. With the parameter stated above we simulated the system obtaining this graph:

After having run a sensitivity analysis on the most influencing parameters we realized that our system is strongly dependent on the fast degradation of GFP and YFP. By putting a slower values of YFP and GFP degradation (0.005) we obtain this graph:

CTMC (stochastic and discrete)

Obviously the behaviour simulated by the ODEs system accounts as an "average" behaviour: we do not expect all cells to be exactly in the same state during time. Moreover we know that physically populations are not continues values but change in a discrete way. In order to deal with the implicit stochasticity and diversity of biological system, in the last years stochastic approaches have been investigated. In our case we use an approach based on [http://en.wikipedia.org/wiki/Process_calculi Process Algebras] (as formalism to define the system) and [http://en.wikipedia.org/wiki/CTMC Continuous Time Markov Chain (CTMC)] (as mathematical model). The Process Algebra we use to define the system is [http://www.dcs.ed.ac.uk/pepa/ PEPA], developed in Edinburgh University. For each species, we define a set of PEPA processes, each one representing a possible population level. Then we add some stochastic actions representing the event of expressing a protein or degrading a protein. The PEPA definitions are:

YFP(I)= (YFPexp,rYFPexp).YFP(I+1) + (YFPdeg,rYFPdeg).YFP(I-1) 
for 0<I<maxYFP, where maxYFP is the maximal number of molecules of YFP

GFP(I)= (GFPexp,rGFPexp).GFP(I+1) + (GFPdeg,rGFPdeg).GFP(I-1)
for 0<J<maxGFP, where maxGFP is the maximal number of molecules of GFP

FORWARDPHASE= (YFPexp,rYFPexp).FORWARDPHASE + (ChangePhaseToFlipping,rChangePhaseToFlipping).FLIPPINGPHASE1
FLIPPINGPHASE1= (ChangePhaseToBacward,rChangePhaseToBacward).BACKWARDPHASE 
FLIPPINGPHASE2= (ChangePhaseToForward,rChangePhaseToForward).FORDWARDPHASE
BACKWARDPHASE= (GFPexp,rGFPexp).BACKWARDPHASE + (ChangePhaseToFlipping,rChangePhaseToFlipping).FLIPPINGPHASE2

The exponential rates are have been derived from the ODEs related parameter or reasoned.

The system initial configuration is:

YFP(initYFP)<YFPexp>FORWARDPHASE<GFPexp>GFP(initGFP)

where initYFP is the initial number of YFP molecules and initiGFP is the initial number of GFP molecules

Using CTMC we are forced to use a stochastic event also for division time. We assume the division time and the flipping time to be exponentially distributed. From the initial configuration is possible then to derive the underlying CTMC of the model. Because of the simplicity of the system, we have implemented it directly in [http://www.prismmodelchecker.org/ PRISM] (the code can be downloaded here) and then simulated the CTMC (using the [http://en.wikipedia.org/wiki/Gillespie_algorithm Gillespie Algorithm]). A single simulation looks like:

in detail when switching from one phase to the other.

on a longer run that consider more than one division.

NB: This approach can be used to model more complex systems (pathways and gene networks for examples). A graphical notation used as human interface to process algebra, the theoretical background, guidelines and case studies of this approach are the Master Dissertation topic of one of the team member (Luca Gerosa). The Dissertation PDF can be viewed or downloaded from here.

Results

The ODEs and the stochastic models, not surprisingly, show the same kind of behaviour regarding to the expression and degradation processes. By running sensitivity analysis and comparing simulation runs with different parameter values, it is possible to observe that the most important factor for the functioning of the Proof of Concept construct is the degradation rate of the two fluoresent proteins. This has to be made as fast as possible, for example by addition of rapid-degradation protein tails. The result of a slow protein degradation rate is the continuously high presence of fluoresence for both colours, leading to an inability to observe flipping behaviours. Regarding the cell division duration, the simulation with a duration between 40-80 minutes and a flipping phase of 0-20 minutes has showed to not be influencing adversely the device behaviour. The stochastic simulation permitted the analysis of the behaviour in the case of highly random division frequency in the colony. Although the exponential random distribution seems quite unrealistic for this kind of parameter, it is still interesting to consider that the functions of the test construct don't seem to be altered except in the case of very rapid division.