Elowitz.R

################################################################################
# A simulation of the Elowitz's repressilator.
#
# References
# * A Synthetic Oscillatory Network of Transcriptional Regulators; 
#   Michael Elowitz and Stanislas Leibler; Nature. 2000 Jan 20;403(6767):335-8.
#
# Author: Andrey Akinshin
################################################################################

library(deSolve)
library(ggplot2)
library(tidyr)
library(rgl)

Elowitz <- function(t, state, parameters) {
  with(as.list(c(state, parameters)), {
    dm1 <- alpha0 + alpha / (1 + p3^n) - m1
    dm2 <- alpha0 + alpha / (1 + p1^n) - m2
    dm3 <- alpha0 + alpha / (1 + p2^n) - m3
    dp1 <- -beta * (p1 - m1)
    dp2 <- -beta * (p2 - m2)
    dp3 <- -beta * (p3 - m3)
    list(c(dm1, dm2, dm3, dp1, dp2, dp3))
  })
}

parameters <- c(n = 2.0, alpha0 = 0.2, alpha = 200.0, beta = 3.0)
state      <- c(m1 = 100, m2 = 80, m3 = 50, p1 = 10, p2 = 10, p3 = 10)
times      <- seq(0, 50, by = 0.01)

out <- ode(y = state, times = times, func = Elowitz, parms = parameters)
df <- data.frame(
  time = out[,"time"],
  lacl = out[,"p1"],
  tetR = out[,"p2"],
  cl   = out[,"p3"]
)

plotElowitz <-
  df %>%
  gather("protein", "value", 2:4) %>%
  ggplot(aes(x=time, y=value, group=protein, colour=protein)) +
  geom_line() +
  theme_bw() +
  ggtitle("Repressilator") +
  xlab("Time (minutes)") +
  ylab("Proteins per cell")
plotElowitz

plot3d(df[,2:4], axes=T, type="l")

network.gv

digraph G {
  ratio="fill";
  size="8,4!";
  edge [
    arrowhead="tee";
    arrowsize=3
  ];
  node [
    width=1;
    height=1;
    fontsize=30;
    penwidth=4
  ];
  a[label="λ cl"; color=blue; fontcolor=blue]
  b[label="Lacl"; color=red; fontcolor=red]
  c[label="TetR"; color=orange; fontcolor=orange]
  d[label="GFP"; color=darkgreen; fontcolor=darkgreen]
  a -> b -> c -> {a, d}
  {rank=same; b; c; d}
}

network.png

Plot2d.png

Plot3d.png