Team:Freiburg/Project/modeling
From 2013.igem.org
Lisaschmunk (Talk | contribs) |
|||
Line 225: | Line 225: | ||
- | <p id="h3"> Results </p> | + | <p id="h3"> Fitting Procedure and Results </p> |
+ | |||
+ | <p>Assuming the given ODE and using the fminsearch-function implemented in matlab with various initial parameter vectors the fitting process results in one optimal parameter composition (Table.1). During the fitting process not only the kinetic parameters were assumed as variable, but also the initial concentrations, because our first measured point is at time point 6 h. However the model is written in a way that the change between flexible and variable parameters is easily to perform, because of an additional vector (qfit). This vector contains boolean values depending on wether the parameter is fixed or flexible during the fitting process. Moreover another parameter is required for adjusting the absoulte values of the duplicates.</p> | ||
+ | <p> | ||
+ | The model reflects the general construction of the network (Fig.2). As assumed Cas seems to converge assymptotically to a stable state and there is an exponential increase in the SEAP concentration. | ||
+ | The model also shows a potential behaviour of the not measured components, the free tracr/crRNA-complex and the gene recognition complex. The tracr/crRNA-complex seems to rise linearly, however the small curvature underlines the assumption of a steady state. There is no free gene recognition complex during the measured time periode, which means that after tracr/crRNA binding the complex immediatly binds to DNA.</p> | ||
Line 231: | Line 236: | ||
<table class="gelpic"> | <table class="gelpic"> | ||
<tbody><tr> | <tbody><tr> | ||
- | <td> <img id="bild" src=https://static.igem.org/mediawiki/2013/1/1b/Freiburg2013_Akt500_fehler1_1.png> </td> | + | <td> <img id="bild" src="https://static.igem.org/mediawiki/2013/f/fc/Freiburg2013-Act500_fertig.png"> </td> |
+ | </tr> | ||
+ | <tr> | ||
+ | <td> <b>Fig. 4: Modeling Result:</b><br> Shown are the experimental results (purple square) in comparison to the model prediction values (cyan cross) for SEAP and Cas, as well as the model prediction for not measured components </td> | ||
+ | </tr> | ||
+ | </tbody></table> | ||
+ | </div> | ||
+ | |||
+ | <div> | ||
+ | <table class="gelpic"> | ||
+ | <tbody><tr> | ||
+ | <td> <img id="bild" src="https://static.igem.org/mediawiki/2013/1/1b/Freiburg2013_Akt500_fehler1_1.png"> </td> | ||
</tr> | </tr> | ||
<tr> | <tr> | ||
Line 238: | Line 254: | ||
</tbody></table> | </tbody></table> | ||
</div> | </div> | ||
+ | |||
+ | |||
+ | |||
+ | |||
+ | |||
<p id="h3"> The Code Files </p> | <p id="h3"> The Code Files </p> |
Revision as of 13:34, 2 October 2013
Modeling uniCAS
Introduction
We used a thermodynamic approach to model and characterize our system. It is based on various ordinary differential equations that describe the behaviour of our network. Due to the limited measurment possibilities and the unwritten law, that you should at least measure half of the number of components of your network we started by using a small network with a limited amount of different components.
The Components are: the protein of Cas9, the target protein (SEAP), and two RNAs (tracrRNA and crRNA).
The two RNAs (tracrRNA and crRNA) are modeled as one complex, therefore there is no possibility to distinguish between them. There might be some differences in their expressions, especially because of the different promoters (crRNA expressed under U6-promoter; tracrRNA expressed under h1-promoter), however the model won't show them.
The Network
uniCAS-VP16
Our network includes four different components uniCas-VP16, a RNA complex (cr/tracr RNA), a RNA-uniCas-VP16 complex and the Secreted alkaline phosphatase (SEAP). UniCas-VP16 binds the RNA-complex and the whole complex binds the DNA, which leads to the production of SEAP.
Fig. 1: Transcriptional Activation via uniCAS-VP16: The uniCAS-VP16 fusion protein is guided to the desired DNA sequence by a co-expressed crRNA, which binds a trRNA. The binding leads to an expression of SEAP. |
Setting up the ODE
According to the graphical reaction network the ODE can be set up.
Cas9 is constitutively expressed by the CBh promoter and degraded proportional to the current concentration. It is used to build the DNA recognition complex and produced during complex decay.
d[Cas]/dt = k1 - k2 ∗ [Cas] - k5 ∗ [Cas] ∗ [tracr,crRNA] + [tracr,crRNA,Cas]
The RNA-complex is build linearly. The production constant can be seen as production constant of the lower expressed RNA, because this expression regulates the complex building. It is assumed that the RNA is degraded after DNA recognition complex decay and therefore the complex decay does not lead to more RNA.
d[tracr,crRNA]/dt = k3 - k4 ∗ [tracr,crRNA] - k2 ∗ [tracr,crRNA] ∗ [Cas]
The complex is build, when Cas9 and RNA meets and degraded proportional to the current DNA recognition complex-concentration.
d[tracr,crRNA,Cas]/dt = k5 ∗ [Cas] ∗ [tracr,crRNA] - k6 ∗ [tracr,crRNA,Cas]
There is a leaky SEAP production and one that depends on the current concentration of the Cas9/RNA Complex. This dependency is assumed to follow the Michaelis-Menten-Rules. Because of the long half time (T2 > 500 h) of SEAP we can neglect the SEAP decay.(Müller et al., 2013 )
d[SEAP]/dt = k7 + k8∗ ([tracr,crRNA,Cas]/(k9+[tracr,crRNA,Cas])))
The parameters are: | |
---|---|
k1: linear production rate of Cas9 k2: Cas9 degradation rate k3: cr/trRNA production rate k4: cr/trRNA degradation rate k5: cr/trRNA /Cas9 complex building rate |
k6: cr/trRNA /Cas9 degradation rate k7: SEAPs leaky production rate k8: Complex dependent SEAP production rate k9: Michaelis-Menten constant |
Finding the parameters
Up to this point, it is only a theoretical description, which fits to every four component system connected in the given way. The reason is, that the parameter values are not yet defined and therefore there is a n-dimensional hypothesis space which contains the hypothesis which fits to the Cas9-system (n is the number of parameters). At least you assume this. Finding this hypothesis means finding the right parameters.
To find these parameters we use a method called maximum likelihood. The maximum likelihood hypothesis is the hypothesis which has the highest probability to generate the measured data. It is shown (Müller et al., 2013 ), that using the maximum likelihood approach and assuming gaussion noise in the data leads to a least-square error minimization problem.
A minimization problem is an optimization problem. You search for parameters (p0) for which holds, that the value of the function (f) at the point of the parameters is smaller than all other values. (f(p0)<=f(p)). In three dimensions the function can be thought as a landscape and minimization is finding the deepest valley. Depending on the method you use different problems arise. The most common problem is finding only a local minimum and not the global one.
Fig. 2: Example of a minimization problem. Shown is a 3D landscape. Depending on the start position (the initial parameters), the found minimum is either a local or the global one. |
To avoid this and to be sure to have found a global minimum we started our minimization procedure using different start values for our parameters. To sample these parameters we used the latin hypercube sampling on a logarithmic scale.
Fig. 3: Illustration of the latin hypercube sampling in a two dimensional parameter space. The number of initial parameter vectors is 5. Therefore the parameter space is divided in 25 subspaces. Shown is one possible parameter combination. |
N is set as the number of different initial parameter settings and the parameter range is divided by N. For the initial parameter the values are chosen so that there is only one parameter in each row and column.
The resulting errors we plotted in an increasing order to be sure to have found a global minimum.
Data generation
Cas is quantified by using Western blot and we used SEAP as target protein that can be quantified by a SEAP assay. For more detailed information refer our modeling notebook.
Fitting Procedure and Results
Assuming the given ODE and using the fminsearch-function implemented in matlab with various initial parameter vectors the fitting process results in one optimal parameter composition (Table.1). During the fitting process not only the kinetic parameters were assumed as variable, but also the initial concentrations, because our first measured point is at time point 6 h. However the model is written in a way that the change between flexible and variable parameters is easily to perform, because of an additional vector (qfit). This vector contains boolean values depending on wether the parameter is fixed or flexible during the fitting process. Moreover another parameter is required for adjusting the absoulte values of the duplicates.
The model reflects the general construction of the network (Fig.2). As assumed Cas seems to converge assymptotically to a stable state and there is an exponential increase in the SEAP concentration. The model also shows a potential behaviour of the not measured components, the free tracr/crRNA-complex and the gene recognition complex. The tracr/crRNA-complex seems to rise linearly, however the small curvature underlines the assumption of a steady state. There is no free gene recognition complex during the measured time periode, which means that after tracr/crRNA binding the complex immediatly binds to DNA.
Fig. 4: Modeling Result: Shown are the experimental results (purple square) in comparison to the model prediction values (cyan cross) for SEAP and Cas, as well as the model prediction for not measured components |
Fig. 4: Different error values plotted in increasing order. |
The Code Files