"
Team:Cambridge/Modeling
From 2008.igem.org
(Prototype team page) |
|||
| Line 1: | Line 1: | ||
| - | + | = Modelling Reaction-Diffusion Systems = | |
| - | + | Mathematically, reaction-diffusion systems are coupled nonlinear differential equations that can be solved numerically. Our first step will be to implement/model a simple two-component system originally proposed by Turing. After validating our numerical method and computating our first patterns (img to the right), we will be thinking about a more realistic system describing the behaviour of the activator and inhibitor system that we intend to engineer with B.subtilis. This will include a detailed analysis of enzyme kinematics and we hope to deduce the parameter ranges, in which our B.subtilis construct will be able to form patterns, thus feeding back to our design decisions (promoter choice, rbs choice). | |
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | We will also take a significant investigative approach and ask whether Turing-like patterns can originate from systems that do not resemble reaction-diffusion. In particular, we would like to ask whether pattern formation can occur with a single signalling molecule only. | |
| + | == Numerical Routine == | ||
| + | I will soon upload my code. It's using a simple finite difference method. | ||
| + | == Introduction == | ||
| - | + | How does pattern formation occur? The main idea is to consider diffusion. Turing showed that pattern formation occurs if the system has a stable steady state without diffusion, but allows unstable states to develop when the diffusion term is added. This results in divergence (or, in a biological context, gene expression/cell differentiation) and subsequent pattern formation. In the following sections, we shall deduce conditions for pattern formation for a simple system. This analysis will apply analogously to subsequent systems as well. | |
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | We will work with a two-component system and we non-dimensionalise the variables. The resulting general form is: | |
| - | { | + | :<math>\frac{\partial A}{\partial t} = \gamma f(A,B) + \nabla^2 A</math> |
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | + | ||
| - | ( | + | |
| + | :<math>\frac{\partial B}{\partial t} = \gamma g(A,B) + d \nabla^2 B</math> | ||
| - | + | where <math>\gamma</math> (determining the scale) is a constant and d stands for the diffusion ratio. We note that non-dimensionalisation has the added advantage that we can now map a specific pattern onto a wide range of biological parameters, the easiest example being that two different pairs of diffusion rates will result in the same pattern if the respective diffusion ratios remain unchanged. | |
| - | + | == Turing System == | |
| + | (also known as Schnakenberg reaction) | ||
| + | |||
| + | === Case: Without diffusion, need stable steady state === | ||
| + | Consider the following system without diffusion terms: | ||
| + | |||
| + | :<math>\frac{\partial A}{\partial t} = \gamma (\alpha - A + A^2B) = \gamma f(A,B)</math> | ||
| + | :<math>\frac{\partial B}{\partial t} = \gamma (\beta - A^2B) = \gamma g(A,B)</math> | ||
| + | |||
| + | We determine the steady state solution and add a small pertubation <math>\tilde{A}, \tilde{B} </math> to linearise the system about the steady state in order to determine its stability. In Matrix notation, the linearised system can be written as: | ||
| + | |||
| + | :<math> | ||
| + | \begin{pmatrix} | ||
| + | \tilde{A_t} \\ | ||
| + | \tilde{B_t} \\ | ||
| + | \end{pmatrix} | ||
| + | = | ||
| + | \gamma | ||
| + | \begin{pmatrix} | ||
| + | f_A & f_B \\ | ||
| + | g_A & g_B \\ | ||
| + | \end{pmatrix} | ||
| + | . | ||
| + | \begin{pmatrix} | ||
| + | \tilde{A} \\ | ||
| + | \tilde{B} \\ | ||
| + | \end{pmatrix} | ||
| + | </math> | ||
| + | |||
| + | where the Jacobian is evaluated at the steady states of A and B. This is a set of coupled first-order ODEs and solutions are proportional to <math>exp(\lambda t)</math>. In order to have a stable steady state, we need the real part of <math>\lambda</math> to be negative. | ||
| + | |||
| + | This requires <math>f_A + g_B < 0</math> (need -ve trace) and <math>f_A g_B-f_B g_A > 0</math> (need +ve determinant), which we note as our first two conditions for pattern formation. | ||
| + | |||
| + | === Case: With diffusion, need unstable steady states === | ||
| + | Now we add diffusion to our system and require the resulting steady state to be unstable. After linearisation, the equation can be expressed as follows: | ||
| + | |||
| + | :<math> | ||
| + | \begin{align} | ||
| + | \begin{pmatrix} | ||
| + | \tilde{A_t} \\ | ||
| + | \tilde{B_t} \\ | ||
| + | \end{pmatrix} | ||
| + | &= | ||
| + | \gamma | ||
| + | \begin{pmatrix} | ||
| + | f_A & f_B \\ | ||
| + | g_A & g_B \\ | ||
| + | \end{pmatrix} | ||
| + | . | ||
| + | \begin{pmatrix} | ||
| + | \tilde{A} \\ | ||
| + | \tilde{B} \\ | ||
| + | \end{pmatrix} | ||
| + | |||
| + | + | ||
| + | |||
| + | \begin{pmatrix} | ||
| + | 1 & 0 \\ | ||
| + | 0 & d \\ | ||
| + | \end{pmatrix} | ||
| + | . | ||
| + | \nabla^2 | ||
| + | \begin{pmatrix} | ||
| + | \tilde{A} \\ | ||
| + | \tilde{B} \\ | ||
| + | \end{pmatrix} | ||
| + | \\ | ||
| + | & \equiv \gamma \mathbf{M}.\mathbf{w} + \mathbf{D}.\nabla^2\mathbf{w} | ||
| + | \end{align} | ||
| + | </math> | ||
| + | |||
| + | In the following steps, we consider solutions of this particular form <math>\mathbf{w}(\mathbf{x},t) = \sum_{k} c_{k} e^{\lambda t} \mathbf{v(x)}_k </math> | ||
| + | |||
| + | with <math> \mathbf{v}_k </math> satisfying <math> \nabla^2 \mathbf{v}_k + k^2\mathbf{v}_k = 0</math>. | ||
| + | |||
| + | |||
| + | After substituting into the linearised system, we can see that for each index <b>k</b>, we require: | ||
| + | |||
| + | :<math>det(\lambda \mathbf{I} - \gamma \mathbf{M} + k^2\mathbf{D}.\mathbf{v}_k) = 0</math>. | ||
| + | |||
| + | |||
| + | We now need to determine the conditions for unstable steady states, i.e. conditions on the parameters for which <math>Re(\lambda)>0</math> hold. The calculations will be standard and similarly to the previous case. In this case, we need the trace and the determinant of the matrix <math>\gamma \mathbf{M}-k^2\mathbf{D.v}_k</math> to infer stability properties as previously. A detailed exposition can be found in <i>Mathematical Biology</i> by Murray (and on which our approach is based). Thus, the conditions are: | ||
| + | |||
| + | |||
| + | :<math> | ||
| + | \begin{align} | ||
| + | dk^4-\gamma k^2(df_u+g_v)+\gamma^2(f_A g_B-g_A f_B) &< 0 \\ | ||
| + | \gamma (f_A + g_B) - k^2(1+d) &> 0 \\ | ||
| + | f_A+g_B &< 0 \\ | ||
| + | f_A g_B - f_A g_B &> 0\\ | ||
| + | \end{align} | ||
| + | </math> | ||
| + | |||
| + | |||
| + | with the last two lines being required by the previous case. Now, after calculating the steady state solutions for f and g, we just need to substitute in the functions into these inequalities and we will arrive at a further set of inequalities involving the parameters <math>\alpha, \beta, \gamma</math> only: | ||
| + | |||
| + | |||
| + | :<math> | ||
| + | \begin{align} | ||
| + | 0 < \beta - \alpha &< (\alpha+\beta)^3 \\ | ||
| + | d(\beta-\alpha) &> (\alpha+\beta)^3 \\ | ||
| + | (d(\beta-\alpha)-(\alpha+\beta)^3)^2 &> 4d(\alpha+\beta)^4\\ | ||
| + | \end{align} | ||
| + | </math> | ||
| + | |||
| + | |||
| + | These are a set of conditions that our parameters must satisfy in order for this particular system to develop patterns. Given a different system (i.e. different f and g), we now have a quick and easy way of findings conditions for our parameters to satisfy. | ||
| + | |||
| + | This will be of particular use in our project. Given a (two-component) system, we can now infer which parameters to choose (i.e. which promoters or binding sites to use) in order for the system to develop patterns. The bacillus model in the following section will be slightly more complicated and will involve more components, but we will try to use the same approach in order to find a set of parameters for which the bacillus construct can, in theory, develop patterns. | ||
| + | |||
| + | =Bacillus Model= | ||
| + | |||
| + | =Beyond Turing= | ||
| + | <math>Insert formula here</math> | ||
Revision as of 13:59, 4 August 2008
Contents |
Modelling Reaction-Diffusion Systems
Mathematically, reaction-diffusion systems are coupled nonlinear differential equations that can be solved numerically. Our first step will be to implement/model a simple two-component system originally proposed by Turing. After validating our numerical method and computating our first patterns (img to the right), we will be thinking about a more realistic system describing the behaviour of the activator and inhibitor system that we intend to engineer with B.subtilis. This will include a detailed analysis of enzyme kinematics and we hope to deduce the parameter ranges, in which our B.subtilis construct will be able to form patterns, thus feeding back to our design decisions (promoter choice, rbs choice).
We will also take a significant investigative approach and ask whether Turing-like patterns can originate from systems that do not resemble reaction-diffusion. In particular, we would like to ask whether pattern formation can occur with a single signalling molecule only.
Numerical Routine
I will soon upload my code. It's using a simple finite difference method.
Introduction
How does pattern formation occur? The main idea is to consider diffusion. Turing showed that pattern formation occurs if the system has a stable steady state without diffusion, but allows unstable states to develop when the diffusion term is added. This results in divergence (or, in a biological context, gene expression/cell differentiation) and subsequent pattern formation. In the following sections, we shall deduce conditions for pattern formation for a simple system. This analysis will apply analogously to subsequent systems as well.
We will work with a two-component system and we non-dimensionalise the variables. The resulting general form is:
- <math>\frac{\partial A}{\partial t} = \gamma f(A,B) + \nabla^2 A</math>
- <math>\frac{\partial B}{\partial t} = \gamma g(A,B) + d \nabla^2 B</math>
where <math>\gamma</math> (determining the scale) is a constant and d stands for the diffusion ratio. We note that non-dimensionalisation has the added advantage that we can now map a specific pattern onto a wide range of biological parameters, the easiest example being that two different pairs of diffusion rates will result in the same pattern if the respective diffusion ratios remain unchanged.
Turing System
(also known as Schnakenberg reaction)
Case: Without diffusion, need stable steady state
Consider the following system without diffusion terms:
- <math>\frac{\partial A}{\partial t} = \gamma (\alpha - A + A^2B) = \gamma f(A,B)</math>
- <math>\frac{\partial B}{\partial t} = \gamma (\beta - A^2B) = \gamma g(A,B)</math>
We determine the steady state solution and add a small pertubation <math>\tilde{A}, \tilde{B} </math> to linearise the system about the steady state in order to determine its stability. In Matrix notation, the linearised system can be written as:
- <math>
\begin{pmatrix}
\tilde{A_t} \\
\tilde{B_t} \\
\end{pmatrix} = \gamma \begin{pmatrix}
f_A & f_B \\ g_A & g_B \\
\end{pmatrix} . \begin{pmatrix}
\tilde{A} \\
\tilde{B} \\
\end{pmatrix} </math>
where the Jacobian is evaluated at the steady states of A and B. This is a set of coupled first-order ODEs and solutions are proportional to <math>exp(\lambda t)</math>. In order to have a stable steady state, we need the real part of <math>\lambda</math> to be negative.
This requires <math>f_A + g_B < 0</math> (need -ve trace) and <math>f_A g_B-f_B g_A > 0</math> (need +ve determinant), which we note as our first two conditions for pattern formation.
Case: With diffusion, need unstable steady states
Now we add diffusion to our system and require the resulting steady state to be unstable. After linearisation, the equation can be expressed as follows:
- <math>
\begin{align} \begin{pmatrix}
\tilde{A_t} \\
\tilde{B_t} \\
\end{pmatrix} &= \gamma \begin{pmatrix}
f_A & f_B \\ g_A & g_B \\
\end{pmatrix} . \begin{pmatrix}
\tilde{A} \\
\tilde{B} \\
\end{pmatrix}
+
\begin{pmatrix}
1 & 0 \\ 0 & d \\
\end{pmatrix} . \nabla^2 \begin{pmatrix}
\tilde{A} \\
\tilde{B} \\
\end{pmatrix} \\ & \equiv \gamma \mathbf{M}.\mathbf{w} + \mathbf{D}.\nabla^2\mathbf{w} \end{align} </math>
In the following steps, we consider solutions of this particular form <math>\mathbf{w}(\mathbf{x},t) = \sum_{k} c_{k} e^{\lambda t} \mathbf{v(x)}_k </math>
with <math> \mathbf{v}_k </math> satisfying <math> \nabla^2 \mathbf{v}_k + k^2\mathbf{v}_k = 0</math>.
After substituting into the linearised system, we can see that for each index k, we require:
- <math>det(\lambda \mathbf{I} - \gamma \mathbf{M} + k^2\mathbf{D}.\mathbf{v}_k) = 0</math>.
We now need to determine the conditions for unstable steady states, i.e. conditions on the parameters for which <math>Re(\lambda)>0</math> hold. The calculations will be standard and similarly to the previous case. In this case, we need the trace and the determinant of the matrix <math>\gamma \mathbf{M}-k^2\mathbf{D.v}_k</math> to infer stability properties as previously. A detailed exposition can be found in Mathematical Biology by Murray (and on which our approach is based). Thus, the conditions are:
- <math>
\begin{align}
dk^4-\gamma k^2(df_u+g_v)+\gamma^2(f_A g_B-g_A f_B) &< 0 \\ \gamma (f_A + g_B) - k^2(1+d) &> 0 \\ f_A+g_B &< 0 \\ f_A g_B - f_A g_B &> 0\\
\end{align} </math>
with the last two lines being required by the previous case. Now, after calculating the steady state solutions for f and g, we just need to substitute in the functions into these inequalities and we will arrive at a further set of inequalities involving the parameters <math>\alpha, \beta, \gamma</math> only:
- <math>
\begin{align}
0 < \beta - \alpha &< (\alpha+\beta)^3 \\ d(\beta-\alpha) &> (\alpha+\beta)^3 \\ (d(\beta-\alpha)-(\alpha+\beta)^3)^2 &> 4d(\alpha+\beta)^4\\
\end{align} </math>
These are a set of conditions that our parameters must satisfy in order for this particular system to develop patterns. Given a different system (i.e. different f and g), we now have a quick and easy way of findings conditions for our parameters to satisfy.
This will be of particular use in our project. Given a (two-component) system, we can now infer which parameters to choose (i.e. which promoters or binding sites to use) in order for the system to develop patterns. The bacillus model in the following section will be slightly more complicated and will involve more components, but we will try to use the same approach in order to find a set of parameters for which the bacillus construct can, in theory, develop patterns.
Bacillus Model
Beyond Turing
<math>Insert formula here</math>
