User:James Chappell /R RK4 Modeling
<html> <style type="text/css"> .firstHeading {display: none;} </style> </html>
James Chappell
Imperial College London Synthetic Biology Lab
Modelling With ODEs
RK4
WORK IN PROGRESS!!!!!
The RK4 tool is a tool that can be used to model ordinary differential equations (ODE). Like other ODE solvers you can output the results of a simulation into graph or tables.
Example 1 : Simple Constitutive Expression
Ordinary Differential Equations of constitutive gene expression
[math]\displaystyle{
\begin{alignat}{1}
\frac{d[mRNA]}{dt} & = k_{1} - d_{1}[mRNA] \\
\frac{d[Protein]}{dt} & = k_{2}[mRNA] - d_{2}[Protein] \\
\end{alignat}
}[/math]
- mRNA = mRNA species
- Protein = Protein species
- k1 = Rate of transcription
- k2 = Rate of translation
- d1 = Rate of mRNA degradation
- d2 = Rate of Protein degradation
(Images taken fromImperial College Synthetic Biology Course)
R CODE
library(odesolve)
#recalls the rk4 allgorithm - need to have install packages to be able to recall and use this function
rm(list=objects())
#Function removes any previous objects from previous simulations using this or other functions, means always need to define any objects for each simulation
(constitutiveexpression <- function(t, x, parms) {
#Define the name of the model, in this case "constitutiveexpression"
mRNA <- x[1]
protein <- x[2]
#Name the response variables, which are the outputs of the model. So in this case the two outputs of this model are the mRNA and protein species, which will both be solved over time given the parameters defined below
with(as.list(parms),{
dmRNA <- k1 - d1*mRNA
dprotein <- k2*mRNA - d2*protein
#Define the differential equations in the format given above
res<-c(dmRNA, dprotein)
list(res)
})}
#Define the species which will be solved, these are the species which have been defined previous as x[1][2][...] and combines these vectors into a list called res</nowiki>
times <-seq(0, 1000, length=1001)
parms <-c(k1=0.69 , k2=f ,d1=0.23 ,d2=0.012)
#Generate a time series over which to solve the equations (1000 in steps of 1)
#Set the parameter values in parms
y <- xstart <-c (mRNA=1, protein=1)
#Set the starting value for response variables
Plotting a Graph
output <- as.data.frame(rk4(xstart, times, constitutiveexpression, parms))
plot(output$time, output$protein, ylim=c(0,2000), ylab="protein (molecules)",xlab="time (min)",type="n") }