# multiple.matrices.r
# work with matrices that vary among years
# and work out stochastic growth rate by 2 methods:
#  1) Tuljapurkar's approximation method and 
#  2) simulation-based projection over many time steps

# Clear any variables in work space
rm(list=ls(all=TRUE))

# load package that has many useful functions
require(popbio)

# load data for exercise, which are:
#  projection matrices for a tropical understory herb
#  (Calathea ovandensis) for plots 1-4 and years 1982-1985 
#  and the pooled matrix = 17 matrices.
# Each matrix constructed using a post-breeding census with 
#  8 size classes: seed, seedling, juvenile, pre-reproductive, 
#  and 4 reproductive classes divided by leaf area.
data(calathea)
calathea
calathea$pooled = NULL # drop pooled matrix 

# check that this worked
names(calathea)

## Calculate and examine growth rates for each plot by year

# calculate lambda's for each of the 16 matrices still stored in 'calathea'
x<-sapply(calathea, lambda) 

# store and then view lambdas in plot by year matrix
x<-matrix(x, nrow=4, byrow=TRUE, dimnames= list(paste("plot", 1:4,sep="."), 1982:1985))
round(x,3)

# create matrix plot of lambda values stored in 'x'
matplot2(x, type='b', ylab='Annual Population Growth Rate', 
         main='Calathea growth rates', lsort=FALSE)

## Work with data for plot 2 from 4 years

# examine mean matrix and SD & CV of values for matrix elements
p3=calathea[9:12]
round(mean(p3),3)
std.dev=sqrt(var2(p3))
round(std.dev,2)
coef.var=std.dev/mean(p3)
round(coef.var,2)

# examine lambda values by year for plot 3
sapply(p3,lambda) # annual lambdas
mean(sapply(p3,lambda)) # mean of annual lambdas

# conduct stochastic projection
sgr<-stoch.growth.rate(p3)
exp(sgr$approx) # exp(log stochastic growth rate by Tuljapukar's approximation)
exp(sgr$sim)    # exp(log stochastic growth rate from simulation)
exp(sgr$sim.CI) # exponentiate the confidence bounds of simulated log stoch. rate