knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  eval = TRUE
)
library(popcycles)

Introduction

Code to implement Caswell & Trevisan periodic matrix models. Code may yet fully updated; function calls are being to be updated and checked (e.g., functions to calc sensitivity)

Caswell and Trevisan 1994 Ecology Sensitivity analysis of periodic matrix models

This paper develope the theory and mechanics for extending sensitivity analysis to periodic matrix models.

Notation:

Sensivitiy of periodic matrix product A

"to calc the sensitivity of lambda to change in each element of Bh, first cyclically permute the individual matrics until B()h appears first in the series (ie at the right hand end of the product matrics). Scond calcualte the sensitivity matrix S(ha) corresponding to this product. Third premultiply it by the TRANSPOSE of the product of the indvidual matrices, EXCLUDING B(h)" Caswell and Treviasan 1994, pg xx (check this quote)

Example 1: Toy Example

Two made up 2 x 2 matrices, Bs and Bw.

Bs <- matrix(data = c(2.0, 0.5,
                      0.5, 0.5),
             byrow = TRUE,
             nrow = 2)

Bw <- matrix(data = c(0.2, 5.0,
                      0.5, 0.5),
             byrow = TRUE,
             nrow = 2)

Periodic matrix product of Bw and Bs, As.

As <- Bw%*%Bs
As

Lambda

calc_lam(As)

Note that Bw * Bs and Bs * Bw result in different matrices, but the same lambda

Bs%*%Bw

Check that lambda is the same

calc_lam(Bs%*%Bw) == calc_lam(Bw%*%Bs)

Sensitivities of Bs

Ssa <- calc_S(As)$S
B.term <- t(Bw)

B.term%*%Ssa

Sensitivities of Bw

Aw <- Bs%*%Bw
Swa <- calc_S(Aw)$S
B.term <- t(Bs)

B.term%*%Swa

One liners

Ss <- round(t(Bw)%*%calc_S(Bw%*%Bs)$S, digits = 2)
Sw <- round(t(Bs)%*%calc_S(Bs%*%Bw)$S, digits = 2)

Real example: McFadden 1991 Alcyonium data

McFadden 1991. A COMPARATIVE DEMOGRAPHIC ANALYSIS OF CLONAL REPRODUCTION IN A TEMPERATE SOFT CORAL. Ecology.

McFadden Table 1 report number of colonies dieing etc.

63/177

Caswell looks at population T1 Data in McFadden table 2

For the "winter" period there are 2 matrices, which she averaged

"I first averaged those transition probabilities that were measured for thesame period in different years at a site to determine a mean transition matrix for that period (e.g., for May-August at sites B1 and B2 and September-May at sites T1 and T2). I then computed the product of the transition matrices for the two different periods at each site. For this calculation, 4- and 8-mo periods were treated as equal (e.g., 6-mo) lengths"

Matrices

# Winter 1:September 1985-May 1986
Bw1 <- matrix(
data = 
c(0,      0,      0,      0,    0,
  0.007,    0.13,   0.12,   0.04,   0.1,
  0,        0.26,   0.37,   0.2,    0.1,
  0,        0.17,   0.35,   0.47,   0.24,
  0,       0.17,    0.09,   0.29,   0.76),
byrow = TRUE,nrow =5
)

#Summer
Bs <- matrix(
  data = 
    c(0,   0,       0.22,     0.45, 4.09,
      0,     0.58,  0.18,     0.13, 0.05,
      0,     0.17,  0.76,     0.34, 0.21,
      0,     0,     0,      0.47,   0.36,
      0,     0,     0,      0.03,   0.43),
      byrow = TRUE,nrow =5
    )


#Winter 2
Bw2 <- matrix(
  data = 
    c(0,      0,       0,     0,    0,
      0.007,    0.3,    0.14,   0.16,   0,
      0,        0.13,   0.21,   0.28,   0.06,
      0,        0.17,   0.31,   0.23,   0.17,
      0,        0.09,   0.18,   0.52,   0.87),
byrow = TRUE,nrow =5
)

Average winter matrix

Bw <- (Bw1+Bw2)/2

Bw

Multiply matrices

Aw <- Bs%*%Bw
Aw

Sensitivites

McFadden reports her S in Table 6. Caswell reproduces (this is equation 10 in Caswell)

round(calc_S(Aw)$E*100,0)

Elast by hand

round((Aw/calc_lam(Aw))*calc_S(Aw)$S*100,0)

Sensitivities as per caswel and Trevisan 1994 Ecology

Winter

Aw <- Bs%*%Bw            
Sha <- calc_S(Aw)$S
B.term <- t(Bs)
Swa <-B.term%*%Sha

Lw <- calc_lam(Aw)

#Calculate elasticity
e.w <- round((Bw/Lw)*Swa*100,0)  #NOTE: the last is calcualted using the component matrix b, not final matrix A!!!

# e = matrix of interest / overall lambda * Special periodic sensitivity
# where overall lambda is NOT lambda fo the matrix of interest but of all of them mult together

#check: should sum to ~100
sum(round((Bw/Lw)*Swa*100,0) )

1+3+3+3+14+7+15
0+2+9+10+18

Summer

As <- Bw%*%Bs
Sha <- calc_S(As)$S
b.term <- t(Bw)
Ssa <-b.term%*%Sha

Ls <- calc_lam(As)

#Calculate elasticity
e.s <- round((Bs/Ls)*Ssa*100,0)  #NOTE: the elast is calculated using the 
                                 #      component matrix b, 
                                 #      not final matrix A!!!

The elasticity

e.s
#check: should sum to ~100
sum(e.s,0)

apply(e.s,2,sum)


brouwern/popcycles documentation built on June 10, 2021, 3:23 p.m.