# Filename: EcoSuccessionPlot.R
# R script to
# - Print out a table showing landscape structure over time
# A function to perform matrix exponentiation
"%^%"<-function(A,n){
if(n==1) A else {B<-A; for(i in (2:n)){A<-A%*%B}}; A
}
# Enter the transfer matrix
T = matrix(c( 0.94, 0.05, 0.01,
0.02, 0.86, 0.12,
0.01, 0.06, 0.93 ), ncol = 3)
# Enter the initial state
x0 = matrix(c( 1, 0, 0 ), ncol = 1)
# We will create a matrix x that has three columns.
# Each column will contain time series data for one class.
# Each row will correspond to a time step.
x = matrix(rep(0, 201*3), ncol = 3)
x[1, ] = x0 # Data for time step t = 0
# Use for loop to generate times series data
for (t in 1:200){
x[t+1, ] = T%^%t %*% x0 # Data for time step t
}
# Time series information for proportion underwater is
# in the first column
u = x[ ,1]
# Time series information for proportion saturated but
# not underwater is in the second column
s = x[ ,2]
# Time series information for proportion dry is in the
# thid column
d = x[ , 3]
# Generate plot
time = seq(1, 200, by = 1)
plot(u,
col = "red",
type = "l",
xlab = "Time step t",
ylab = "Proportion of Wetlands",
xlim = c(1, 200),
ylim = c(0, 1))
par(new = TRUE)
plot(s,
col = "green",
type = "l",
xlab = "",
ylab = "",
xlim = c(1, 200),
ylim = c(0, 1),
xaxt = "n",
yaxt = "n")
par(new = TRUE)
plot(d,
col = "blue",
type = "l",
xlab = "",
ylab = "",
xlim = c(1, 200),
ylim = c(0, 1),
xaxt = "n",
yaxt = "n")
legend("topright",
c("u", "s", "d"),
col = c("red", "green", "blue"),
lty = c(1, 1, 1))