User:James Chappell /R RK4 Modeling

From OpenWetWare
Jump to: navigation, search

<html> <style type="text/css"> .firstHeading {display: none;} </style> </html>

James Chappell

Imperial College London Synthetic Biology Lab

Modelling With ODEs



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

CellDesigner Constitutive Gene Expression Network.png

Ordinary Differential Equations of constitutive gene expression

  • 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)



#recalls the rk4 allgorithm - need to have install packages to be able to recall and use this function


#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


dmRNA <- k1 - d1*mRNA

dprotein <- k2*mRNA - d2*protein

#Define the differential equations in the format given above

res<-c(dmRNA, dprotein)



#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 <-, times, constitutiveexpression, parms))

plot(output$time, output$protein, ylim=c(0,2000), ylab="protein (molecules)",xlab="time (min)",type="n") }