Dynamical Systems in Biology

Tutorial 1

Plotting with ggplot2

First, get ggplot2 from CRAN

install.packages("ggplot2")

and load it

library(ggplot2)
theme_set(theme_bw())

For illustration purposes, we load a data set from another package. In this exercise sheet we are going to produce our own data.frames.

data(Oxboys, package = "nlme")
mydata <- subset(Oxboys, Subject %in% as.character(1:9))

The data set contains columns Subject, age, height and Occasion. The standard procedure to produce a plot with ggplot2 is the following

1. Generate the structure of the plot: ggplot(data, aes(x = ..., y = ..., ...)) from a data set data and a mapping returned by the aes() function.
2. Add graphical elements to the plot: ... + geom_line() + geom_point() + ...

Graphical elements can be lines, points, bars, etc. but also instructions how to distribute plots, e.g., facet_wrap() or facet_grid to arrange the plot in subplots.

We use the Oxboys dataset to generate some example plots

Grouping of different occasions in one plot

ggplot(mydata, aes(x = age, y = height, group = Occasion, color = Occasion)) + geom_point()

ggplot(mydata, aes(x = age, y = height, group = Subject, color = Subject)) + geom_line()

Wrapping occasions up in a panels (facets)

ggplot(mydata, aes(x = age, y = height)) + facet_wrap(~Subject) + geom_point() + geom_line()

ggplot(Oxboys, aes(x = height)) + facet_wrap(~Occasion) + geom_histogram()

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Please visit http://docs.ggplot2.org for more examples and the complete documentation of ggplot2.

Exercise 1: Plotting Verhulst dynamics

A solution of the Verhulst ODE \[ \dot x = rx\left( 1-\frac{x}{K} \right), \quad x(0)=x_0 \qquad (1) \] is given by \[ x(t) = \frac{Kx_0}{x_0 + (K-x_0)e^{-rt}} \qquad (2) \]

• Write a function f <- function(t, p) for Eq. (2) returning a data.frame with columns t, x, x0, K, r.
• Plot \(x(t)\) on the interval \([0,1]\) for \(K=0.5, r=4\) and a set of values for \(x_0\) in [0,1] (use lapply() and do.call()). Group the plot by values of \(x_0\).
• In addition, vary \(K\) within [0.2,1] (see expand.grid()) and plot the curves in panels separated by their \(K\) values.

• Make a grid plot by varying both, parameters \(K\) and \(r\).

• Write a function VH <- function(x, p) for Eq. (1) returning a data.frame with columns x, dx.
• Plot \(\dot{x}(t)\) on the interval \(x\in [0,1]\) for \(K=0.75, r=4, x_0 = 1\).

• Add arrows to the plot along the \(x\)-axis pointing from \(x\) to \(x+\dot x {\rm d}t\) with \({\rm d}t=0.1\). Hint: take a look at the examples of geom_segment().

Exercise 2: Data simulation and error bars

Consider Eq. (2) on the interval \(t\in [0,1]\) given \(x_0 = 0.1, K = 0.5, r = 4\). A data point is to be simulated by adding a normally distributed random number \(r\) with mean 0 and standard deviation 0.05, \(r\sim\mathcal{N}(0,0.05^2)\) (see rnorm), to the function value \(x\).

• Simulate data for 10 timepoints in [0,1] generating a data.frame with columns t, x, sigma.
• Generate a data.frame for 100 timepoints in [0,1] for \(x(t)\) without noise.
• Plot the datasets, including data points (geom_point()), error bars (geom_errorbar()), curves \(x(t)\) (geom_line()) and error bands around the curves (geom_ribbon()).

Cathedral exercise

What distinguishes our cathedral from all other German Gothic cathedrals?