class: center, middle, title-slide .title[ # Population Projection Models (PPMs) ] .author[ ### Sarah Cubaynes for the Team ] .date[ ### last updated: 2023-06-01 ] --- # What we've learned so far ### Obj. 1: Assess current and past trends in population abundance - see class 1 by Aurélien -- ### Obj. 2: Estimate demographic parameters and identify causes of variation - see class 2 by Olivier -- ### Obj. 3: Analyze population viability to inform management decisions - some hints about this now! --- ## Temporal dynamics of a virtual population <img src="matrixmodels2_files/figure-html/unnamed-chunk-1-1.png" style="display: block; margin: auto;" /> --- ## What are the risks of population decline ? <img src="matrixmodels2_files/figure-html/unnamed-chunk-2-1.png" style="display: block; margin: auto;" /> --- # Population Viability Analysis (PVA) + Evaluation of **data and models** to determine the **likelihood** that a population will **persist for some arbitrarily chosen time** into the future -- + Use of **quantitative methods** to predict the **likely future status** of a population or collection of populations of conservation concern -- + **Tentative assessments** based upon what we now know (**rather than as iron-clad predictions **of population fate) -- + Main references Caswell, H. 2003. Matrix Population Models: Construction, Analysis, and Interpretation . Sinauer. Morris et al. (1999) A practical handbook for population viability analysis. The Nature Conservancy. --- class: center, middle # Important things to keep in mind --- ## Amount of data available .pull-left[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-3-1.png" style="display: block; margin: auto;" /> ] .pull-right[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-4-1.png" style="display: block; margin: auto;" /> ] --- ## Temporal variance .pull-left[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-5-1.png" style="display: block; margin: auto;" /> ] .pull-right[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-6-1.png" style="display: block; margin: auto;" /> ] --- ## Precision and robustness of estimates .pull-left[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-7-1.png" style="display: block; margin: auto;" /> ] .pull-right[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-8-1.png" style="display: block; margin: auto;" /> ] --- ## How to assess population viability from census data ? <img src="img/pop_size.png" width="60%" style="display: block; margin: auto;" /> --- class: center, middle # 1- Count-based PVA --- class: center, middle ## 1-1 Deterministic exponential growth model --- ### Start simple, let's assume exponential growth at a constant rate of change .pull-left[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-10-1.png" style="display: block; margin: auto;" /> ] -- .pull-right[ - `\(\lambda\)` gives the **proportional change** in population size `$$N_{2}=\lambda \cdot N_{1}$$` `$$N_{3}=\lambda \cdot N_{2}$$` `$$\cdot \cdot \cdot$$` `$$N_{t+1}=\lambda \cdot N_{t}$$` `$$\lambda = \displaystyle{\frac{N_{t+1}}{N_{t}}}$$` - `\(\lambda\)` is the **population growth rate** ] --- ## Population growth rate `\(\lambda\)` <img src="matrixmodels2_files/figure-html/unnamed-chunk-11-1.png" style="display: block; margin: auto;" /> --- ## Growth rate and rate of increase - `\(\lambda = \displaystyle{\frac{N_{t+1}}{N_{t}}}\)` is log-normally distributed -- - **Population growth rate** `\(\lambda \sim \mbox{Lognormal}(\mu,\sigma^{2})\)` -- - Easier to work with the **rate of increase** `$$r = \log(\lambda) = \log\left(\frac{N_{t+1}}{N_{t}}\right) = \log(N_{t+1}) - \log(N_{t})$$` -- - Which is normally distributed (back to it later) `$$r \sim N(\mu,\sigma^{2})$$` -- - if `\(r> 0\)` the population is growing --- ## To project the population in this simple case - Calculate the rate of increase `$$r = \log(\lambda) = \log(N_{t+1}) - \log(N_{t})$$` -- - Iterate the population starting from initial population `\(N_{1}\)` -- - After `\(t\)` time steps, population size is simply obtained as `$$log (N_{t} ) = log(N_{1}) + r \cdot t$$` --- ## Example of the Yellowstone grizzly bear population From Morris & Doak (2002). *Quantitative conservation biology: Theory and practice of population viability analysis*. Massachusetts, USA. .pull-left[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-12-1.png" width="85%" /> ] .pull-right[ <img src="img/bear.jpg" width="60%" style="display: block; margin: auto;" /> ] --- ## Step 1: Calculate mean rate of increase ```r N <- grizzly$N # pop size each year Years <- grizzly$year # years of sampling # average rate of increase over the years rt <- log( N[-1]/N[-length(N)] ) # log( Nt+1 / Nt ) mu <- mean(rt) mu ``` ``` [1] 0.02134027 ``` - remember we said that `\(r_{t}\)` is normally distributed `$$r_{t} \sim N(\mu,\sigma^{2}_{t})$$` -- - mean rate of increase `\(\mu > 0\)`, so `\(\lambda > 1\)`, on average the grizzly population is growing --- ## Step 2: Project the population You only need to specify * Initial population `\(N_{1}=44\)` * Time steps to project over: `\(t = 38\)` years (nb. years - 1) * Mean rate of increase `\(\mu = 0.021\)` -- * Project the population by iteration using `$$log(N_{t}) = log(N_{1}) + \mu \cdot t$$` --- ## Step 3: Evaluate the quality of your projection .pull-left[ <img src="matrixmodels2_files/figure-html/unnamed-chunk-15-1.png" style="display: block; margin: auto;" /> ] -- .pull-right[ - average population trajectory is interesting but not sufficient - only using the **mean rate of increase** is not recommended - except for exponentially growing populations - when environmental variance is low ] --- # Environmental variance - Remember that the rate of increase `\(r_{t}\)` is normally distributed `$$r_{t} \sim N(\mu,\sigma^{2}_{t})$$` - `\(\mu\)` is the **mean rate of increase** - `\(\sigma^{2}\)` is the **environmental variance** (we ignored it until now) -- - If the mean rate of increase `\(\mu < 0\)`, extinction will certainly occur. -- - A population can still decline or go extinct even if the mean rate of increase `\(\mu>0\)`, because of environmental variance `\(\sigma^{2}\)` --- # Environmental stochasticity - Unpredictable fluctuations in environmental conditions -- - Stochastic fluctuations in survival and reproductive rates over the years reduce long-run population growth rate -- - The chance of occurrence of consecutive periods of unfavourable environmental conditions **can drive small populations to extinctions** -- - Better to not ignore environmental variance when assessing extinction risks --- class: center, middle ## 1-2 Stochastic count-based PVA --- ## Environmental stochasticity - **Stochastic PVAs** project the population using a series of annual stochastic rates of increase to account for environmental variance -- - Each projection is unique, you need to run stochastic projection several times (at least 100 times, 500 is good in most cases) -- - Using all population trajectories, you can calculate interesting quantities such as **extinction probability**, **distribution of population sizes** at various times --- # Step 1: Calculate `\(\mu\)` and `\(\sigma^{2}\)` from the data ```r # rate of increase over years logN <- log(N[-1]/N[-length(N)]) # log(Nt+1) - log(Nt) #mean rate of increase *mu <- mean(logN) #environmental variance *sigma2 <- var(logN) ``` ``` mu sigma2 1 0.02134027 0.01305092 ``` -- - `\(\mu >0\)` so on average the population is growing -- - `\(\sigma^{2} = 0.013\)` reflects low inter-annual variance in the rate of increase --- ## Get confidence intervals for `\(\mu\)` and `\(\sigma_{2}\)` ```r ## 95% confidence interval for mu c(mu - qnorm(0.975) * sd(logN) /sqrt(39), mu + qnorm(0.975) * sd(logN) /sqrt(39)) ``` ``` [1] -0.01451363 0.05719416 ``` - Confidence interval of mean rate of increase encompasses 0, therefore we cannot rule out a potential risk of decline ! -- ```r ## Confidence interval for sigma 2 df1 <- length(logN) - 1 df1 * sigma2 / qchisq(c(.975, .025), df = df1) ``` ``` [1] 0.008674359 0.021844393 ``` --- ## Step 2: Set initial population size, number of time steps, and number of runs ```r n1 <- grizzly$N[1] # initial pop. T <- 50 # time iterations to project over, projection depth runs <- 500 # number of repetitions (pop. trajectories) ``` --- ## Step 3: Set a quasi-extinction threshold - Some simulated populations might go extinct or drop below a viability threshold ```r Ne <- 30 # threshold for minimum viable pop. ``` -- **Quasi-extinction threshold** - 1 female or a minimum viable population (genetic drift, demographic stochasticity) - Can also be the lowest level of abundance at which it remains feasible to attempt intervention to prevent further decline. --- ## Step 4: Project the population while accounting for observed variation in growth rate - **Principle** Generate time series of stochastic rates of increase by sampling from a Normal distribution `$$r_{1,...,T-1} \sim N(\mu,\sigma^{2})$$` -- - You can now iterate the population using the generated stochastic rates of increase `$$log (N_{t} ) = log(N_{1}) + r_{t} \cdot t$$` -- - Repeat until last time step, for each run ! --- ## Step 4: Run the stochastic projections ```r stoch.pop <- matrix(NA,T,runs) # to store results stoch.pop[1,] <- n1 # initil population # let's project the population *for (i in 1:runs){ * for (t in 2:T){ # Draw r from normal using estimates of mu and sigma2 r <- rnorm(n = 1, mean = mu, sd = sqrt(sigma2)) # back-transform to get lambda and get pop. size lambda <- exp(r) #project one time step from the current pop size stoch.pop[t,i] <- stoch.pop[(t-1),i] * lambda # leave the loop if pop <= threshold if(stoch.pop[t,i] <= Ne){ stoch.pop[t,i] <- 0 i < i+1} } } ``` --- ## Step 4: Run the stochastic projections ```r # let's project the population for (i in 1:runs){ # loop over repetitions for (t in 2:T){ # loop over years # Draw r from normal using estimates of mu and sigma2 * r <- rnorm(n = 1, mean = mu, sd = sqrt(sigma2)) # back-transform to get lambda and get pop. size lambda <- exp(r) #project one time step from the current pop size stoch.pop[t,i] <- stoch.pop[(t-1),i] * lambda # leave the loop if pop <= threshold if(stoch.pop[t,i] <= Ne){ stoch.pop[t,i] <- 0 i < i+1} } } ``` --- ## Step 4: Run the stochastic projections ```r # let's project the population for (i in 1:runs){ # loop over repetitions for (t in 2:T){ # loop over years # Draw r from normal using estimates of mu and sigma2 r <- rnorm(n = 1, mean = mu, sd = sqrt(sigma2)) # back-transform to get lambda and get pop. size * lambda <- exp(r) #project one time step from the current pop size * stoch.pop[t,i] <- stoch.pop[(t-1),i] * lambda # leave the loop if pop <= threshold if(stoch.pop[t,i] <= Ne){ stoch.pop[t,i] <- 0 i < i+1} } } ``` --- ## Step 4: Run the stochastic projections ```r # let's project the population for (i in 1:runs){ # loop over repetitions for (t in 2:T){ # loop over years # Draw r from normal using estimates of mu and sigma2 r <- rnorm(n = 1, mean = mu, sd = sqrt(sigma2)) # back-transform to get lambda and get pop. size lambda <- exp(r) #project one time step from the current pop size stoch.pop[t,i] <- stoch.pop[(t-1),i] * lambda # leave the loop if pop <= threshold * if(stoch.pop[t,i] <= Ne){ * stoch.pop[t,i] <- 0 * i < i+1} } } ``` --- ## Step 5: Examine the results * Each simulated population has different trajectory <img src="matrixmodels2_files/figure-html/unnamed-chunk-27-1.png" style="display: block; margin: auto;" /> --- ## Step 5: Examine the results Distribution of population sizes at the last time step <img src="matrixmodels2_files/figure-html/unnamed-chunk-28-1.png" style="display: block; margin: auto;" /> --- ## Step 5: Examine the results Compare predictions with data <img src="matrixmodels2_files/figure-html/unnamed-chunk-29-1.png" style="display: block; margin: auto;" /> --- ## Step 6: Quantify extinction risks - The average population growth rate doesn’t do a good job at predicting what most population realizations will do -- - What are the chances that the population will go extinct at various times? + **Extinction risk** + **Time to extinction** --- ## Step 6: Quantify extinction risks - **Ultimate extinction probability** = percentage of trajectories (over the 500 runs) reaching the extinction threshold at some point (over T years) ```r Pr.ext <- sum(lastN <= Ne) / runs # prob. to reach the quasi-extinction threshold Pr.ext ``` ``` [1] 0.196 ``` --- ## Step 6: Quantify extinction risks - **Cumulative extinction risk** over the years <img src="matrixmodels2_files/figure-html/unnamed-chunk-31-1.png" style="display: block; margin: auto;" /> --- ## Step 7: Quantify time to extinction - Mean time to extinction is an overestimation (because of few pop. growing fast) - **Median time to extinction** is a better measure ```r # Time to reach extinction for extinct pop. maxt <- NULL # empty vector to store results for (i in 1:runs){ # loop over repetitions N <- stoch.pop[,i] # max time N > threshold * maxt[i] <- max(which(N>0)) } # time at extinction for pseudo-extinct populations *time.ext <- maxt[maxt < T] median(time.ext) ``` ``` [1] 12.5 ``` --- ## Step 7: Quantify time to extinction <img src="matrixmodels2_files/figure-html/unnamed-chunk-33-1.png" style="display: block; margin: auto;" /> --- ## Step 8: Perturb the model ### Interesting to evaluate the sensitivity of the results to changes in: - Initial population size - Extinction threshold - Amount of environmental variance - Projection depth (number of years) --- # Live demo on grizzly bears From Morris & Doak (2002). *Quantitative conservation biology: Theory and practice of population viability analysis*. Massachusetts, USA. <img src="img/bear.jpg" width="40%" style="display: block; margin: auto;" /> --- ## Count-based extinction analyses rely on strong hypotheses - Exhaustive counts (no sampling error) -- - No density-dependence (exponential growth) -- - Only source of variation is environmental stochasticity + no demographic stochasticity, + no trends in mean or variance over time + uncorrelated environment among successive years -- - Moderate environmental variability (no catastrophe, no bonanzas) --- # Advantages - Simplicity (data requirement at least 10 censuses and calculation) - Work relatively well when assumptions are met - Assess model quality by hindcasting - Always keep a critical eye on your results, assess the quality of your data and sampling protocols --- # Limitations - `\(\lambda\)` is only a summary of the population dynamics -- - No info about the mechanisms governing the dynamics -- - No hints about which management action might be most efficient? - Is it better to act on survival? Fecundity? Which age classes govern population dynamics? --- # How do we assess population viability from demographic parameters ? <img src="img/PPM.png" width="60%" style="display: block; margin: auto;" /> --- class: center, middle # 2- PVA using Matrix Projection Models (MPMs) --- # Let's assume no migrations for now - Survival and fecundity rates are enough to fully describe the population dynamics `$$N_{t+1} = N_{t} * F + N_{t} * S$$` + F = fecundity + S = survival (= 1- mortality) --- # Demographic parameters (see class 2) - Demographic parameters are often **heterogeneous across the life cycle**: survival and fecundity vary with age and/or stage in many species. <img src="img/oak-tree.jpg" width="60%" style="display: block; margin: auto;" /> --- # Matrix Population Models (MPMs) - Incorporate vital rates that are heterogeneous across the life cycle - Project the population based on the matrix summarizing the age- or stage- dependent demographic parameters --- ## The Barn swallow (*Hirundo rustica*) example <img src="img/barnswallowlifecycle.png" width="60%" style="display: block; margin: auto;" /> --- # Demographic parameters Survival rates: - `\(S_{0}\)` first-year survival - `\(S_{1}\)` juvenile survival (1 yo) - `\(S_{2}\)` adult survival (2+ yo) Fecundity: - `\(F_{1}\)` number of females produced by a juvenile female - `\(F_{2}\)` number of females produced by an adult female --- # Timing of data collection ? <img src="img/PreBreeding.png" width="100%" style="display: block; margin: auto;" /> + Pre-breeding census --- # Timing of data collection ? <img src="img/PostBreeding.png" width="100%" style="display: block; margin: auto;" /> + Post-breeding census --- # Step 1: write the agenda of events + Let's first consider a pre-breeding census <img src="img/AgendaPreB1.png" width="100%" style="display: block; margin: auto;" /> --- # Step 1: write the agenda of events + Let's first consider a pre-breeding census <img src="img/AgendaPreB2.png" width="100%" style="display: block; margin: auto;" /> --- # Step 1: write the agenda of events + Let's first consider a pre-breeding census <img src="img/AgendaPreB3.png" width="100%" style="display: block; margin: auto;" /> -- + Note that newborn are not observed directly! --- # Step 2: Translate into life cycle graph or 'Caswell representation' - A trick is to go 'up the arrows' <img src="img/preBLifeCycle.png" width="35%" style="display: block; margin: auto;" /> --- # Step 3: Translate into equations - Link `\(N_{(t+1)}\)` to `\(N_{(t)}\)` via survival and fertility rates - A trick is to read the parameters going up the arrows `\(N_{(1,t+1)} = F_{1} \cdot S_{0} \cdot N_{(1,t)}+ F_{2} \cdot S_{0} \cdot N_{(2,t)}\)` `\(N_{(2,t+1)} = S_{1} . N_{(1,t)} + S_{2} . N_{(2,t)}\)` --- # Step 4: Arrange in a matrix format + Called the **transition matrix**, or the **projection matrix** `$$N_{(t+1)} = \begin{bmatrix}F_{1} \cdot S_{0} & F_{2} \cdot S_{0}\\ S_{1} & S_{2} \end{bmatrix} \cdot N_{(t)}$$` --- # Step 4: Arrange in a matrix format + With `\(S_{0} = 0.2\)`, `\(S_{1} = 0.5\)` and `\(S_{2} = 0.65\)` + `\(F_{1} = 3/2\)` and `\(F_{2} = 6/2\)` `$$N_{(t+1)} = \begin{bmatrix}0.3 & 0.6\\ 0.5 & 0.65 \end{bmatrix} \cdot N_{(t)}$$` --- ## What is the difference with a post-breeeding census ? -- <img src="img/PostBreeding.png" width="100%" style="display: block; margin: auto;" /> + Post-breeding census --- # Step 1: write the agenda of events <img src="img/AgendaPostB1.png" width="100%" style="display: block; margin: auto;" /> --- # Step 1: write the agenda of events <img src="img/AgendaPostB2.png" width="100%" style="display: block; margin: auto;" /> --- # Step 1: write the agenda of events <img src="img/AgendaPostB3.png" width="100%" style="display: block; margin: auto;" /> -- + Note that newborn are now observed! --- # Step 2: Translate into life cycle graph or 'Caswell representation' <img src="img/postBLifeCycle.png" width="70%" style="display: block; margin: auto;" /> --- # Step 3: Translate into equations `\(N_{(0,t+1)} = S_{0} \cdot F_{1} \cdot N_{(0,t)} + S_{1} \cdot F_{2} \cdot N_{(1,t)}+ S_{2} \cdot F_{2} \cdot N_{(2,t)}\)` `\(N_{(1,t+1)} = S_{0} \cdot N_{(0,t)}\)` `\(N_{(2,t+1)} = S_{1} \cdot N_{(1,t)} + S_{2} \cdot N_{(2,t)}\)` -- # Step 4: Arrange in matrix format `$$N_{(t+1)} = \begin{bmatrix} S_{0} \cdot F_{1} & S_{1} \cdot F_{2} & S_{2} \cdot F _{2} \\ S_{0} & 0 & 0 \\ 0 & S_{1} & S_{2}\\ \end{bmatrix} \cdot N_{(t)}$$` --- # Why is the matrix format interesting ? - Easier to read than multiple equations - Intrinsic numeric features (back to it later) - Work the same way for complex life cycles --- # Examples <img src="img/Compadre.png" width="80%" style="display: block; margin: auto;" /> --- ## Several sites: Barn swallow example <img src="img/twosites.png" width="80%" style="display: block; margin: auto;" /> --- ## Variable age at first reproduction: Slender-billed gull example <img src="img/recruitment.png" width="80%" style="display: block; margin: auto;" /> --- ## Transition among stages: Peony example <img src="img/stages.png" width="80%" style="display: block; margin: auto;" /> --- # The projection matrix - How many age classes / stages in the life cycle? - Age at first reproduction? - Maximum age at death fixed or not? - Pre- or post-breeding census? <img src="img/Leslie.png" width="60%" style="display: block; margin: auto;" /> --- class: center, middle # 1- Deterministic MPMs --- ## Step 5: Project the population + Write the transition matrix in R: ```r A.swallow <- matrix(c(0.3, 0.6, 0.5, 0.65), # pre-breeding leslie matrix nrow = 2, byrow = TRUE) A.swallow ``` ``` [,1] [,2] [1,] 0.3 0.60 [2,] 0.5 0.65 ``` --- ## Step 5: Project the population + Start from an initial population `\(n_{0}\)` ```r n0 <- c(50,30) # vector with initial population n0 ``` ``` [1] 50 30 ``` --- ## Step 5: Project the population + Project to the next time step by iteration `$$N_{t+1} = A \cdot N_{t}$$` + equivalent to in R: ```r n1 <- A.swallow %*% n0 # matrix product n1 ``` ``` [,1] [1,] 33.0 [2,] 44.5 ``` --- ## Step 5: Project the population + Project the population over 10 years: `$$N_{t+1} = A^{t} \cdot N_{0}$$` + equivalent to in R: ```r require(matrixcalc) t <- 10 # matrix product *n10 <- matrix.power(A.swallow,t) %*% n0 n10 ``` ``` [,1] [1,] 53.82434 [2,] 67.28043 ``` --- ## Step 5: Project the population + Or using built-in functions from package `popbio`: ```r library(popbio) # load package t <- 11 # project the population *results <- pop.projection(A.swallow,n0,iterations = t) results ``` --- ## Step 6: Examine the results Let's plot the projected population dynamics <img src="matrixmodels2_files/figure-html/unnamed-chunk-59-1.png" style="display: block; margin: auto;" /> --- ## Convergence to a stable distribution Transient dynamics: + Depends on the initial population + Damping ratio measures how fast the population converges toward equilibrium -- Stationary phase: + Independent of initial conditions -- + Depends on the transition matrix only -- + Constant growth rate = stable or **asymptotic population growth rate** ( `\(\lambda\)` ) -- + Constant proportion of individuals per age/stage = **stable age/stage structure** --- ## Step 6: Examine the results ```r results <- pop.projection(A.swallow, n0, iterations = t) names(results) ``` ``` [1] "lambda" "stable.stage" "stage.vectors" "pop.sizes" [5] "pop.changes" ``` --- ## Step 6: Examine the results - First element contains the average growth rate in stationary phase = **asymptotic growth rate** ```r results$lambda ``` ``` [1] 1.05 ``` --- ## Step 6: Examine the results - Second element contains the **stable age/stage structure** = % of each age/stage in the population in stationary phase ```r results$stable.stage ``` ``` [1] 0.4444444 0.5555556 ``` --- ## Step 6: Examine the results - Third element contains **pop. sizes per age/stage** class at each time step ```r results$stage.vectors ``` ``` 0 1 2 3 4 5 6 7 8 [1,] 50 33.0 36.600 38.23500 40.16625 42.17261 44.28144 46.49549 48.82027 [2,] 30 44.5 45.425 47.82625 50.20456 52.71609 55.35177 58.11937 61.02533 9 10 [1,] 51.26128 53.82434 [2,] 64.07660 67.28043 ``` --- ## Step 6: Examine the results - Fourth element contains **total pop. size** at each time step ```r results$pop.sizes ``` ``` [1] 80.00000 77.50000 82.02500 86.06125 90.37081 94.88870 99.63320 [8] 104.61486 109.84560 115.33788 121.10477 ``` --- ## Step 6: Examine the results - Fifth element is `\(\lambda_{t}\)` the **rate of change** at each time step ```r results$pop.changes ``` ``` [1] 0.968750 1.058387 1.049208 1.050076 1.049993 1.050001 1.050000 1.050000 [9] 1.050000 1.050000 ``` --- ## Quantities in stationary phase Calculated directly from the transition matrix: ```r lambda(A.swallow) #stable population growth rate ``` ``` [1] 1.05 ``` --- ## Quantities in stationary phase Calculated directly from the transition matrix: ```r stable.stage(A.swallow) # stable age/stage structure ``` ``` [1] 0.4444444 0.5555556 ``` --- ## Quantities in stationary phase Calculated directly from the transition matrix: ```r # relative contribution of each age/stage to the next generation reproductive.value(A.swallow) ``` ``` [1] 1.0 1.5 ``` --- ## Quantities in stationary phase Calculated directly from the transition matrix: ```r generation.time(A.swallow) # average time between generations ``` ``` [1] 4.150924 ``` ```r # average age of mothers at birth of their daughters ``` --- # Step 7: Sensitivity analysis Let's perturb the model: what happens if female adult survival is reduced by 50%? ```r A.swallow.modified <- matrix(c(0.3, 0.6, 0.5, 0.65/2), # pre-breeding leslie matrix nrow = 2, byrow = TRUE) A.swallow.modified ``` ``` [,1] [,2] [1,] 0.3 0.600 [2,] 0.5 0.325 ``` --- # Step 7: Sensitivity analysis ```r n0 <- c(100,100) t <- 10 results.modified <- pop.projection(A.swallow.modified, n0, iterations = t) ``` --- # Step 7: Sensitivity analysis <img src="matrixmodels2_files/figure-html/unnamed-chunk-72-1.png" style="display: block; margin: auto;" /> --- # Step 7: Sensitivity analysis ```r lambda(A.swallow.modified) ``` ``` [1] 0.8603652 ``` + The population is now declining + Adult survival is important to population growth --- # Step 7: Sensitivity analysis + What happens if juvenile survival is reduced ? -- + What happens if fecundity or chick survival is reduced ? -- + Which demographic parameter(s) drive the population dynamics ? --- # Objectives of sensitivity analyses + Measure the **impact of a change in a specific demographic parameter on population dynamics** -- + **Sensitivity** measures **absolute change** (e.g. -0.1 in parameter) `$$\frac{\delta \lambda}{\delta \theta}$$` -- + **Elasticity** measures **relative change** (e.g. -0.1% change in parameter) `$$\frac{\delta \lambda}{\delta \theta} \cdot \frac{\theta}{\lambda}$$` -- + Both are useful, elasticity is better to compare parameters that are on different scales (change in survival versus fertility) --- # Step 7: Sensitivity analysis ```r swallow.param <- list(s0 = 0.20, s1 = 0.5, s2 = 0.65, f1 = 3/2, f2 = 6/2) swallow.equation <- expression( s0 * f1, s0 * f2, s1, s2) VS <- vitalsens(swallow.equation, swallow.param) ``` --- # Step 7: Sensitivity analysis ```r VS <- vitalsens(swallow.equation, swallow.param) VS ``` ``` estimate sensitivity elasticity s0 0.20 1.82608696 0.34782609 s1 0.50 0.52173913 0.24844720 s2 0.65 0.65217391 0.40372671 f1 1.50 0.06956522 0.09937888 f2 3.00 0.08695652 0.24844720 ``` --- # Step 7: Sensitivity analysis <img src="matrixmodels2_files/figure-html/unnamed-chunk-76-1.png" style="display: block; margin: auto;" /> --- # Implications for management + Identify key parameters for population management - best strategy here is to reduce adult mortality -- + Evaluate the impact of relative management actions - actions focused on fecundity (e.g. nest protection) will have a reduced impact -- + see package 'Rramas' for implementation of specific management actions in the simulations --- # Live demo on the barn swallow <img src="img/swallow.jpg" width="328" style="display: block; margin: auto;" /> --- ## Assumptions and limitations of deterministic MPMs + One sex model: Females drive the demography, males are not limiting -- + Synchronous breeding: discrete time -- + No density-dependence (exponential growth or decay) -- + Demographic parameters are constant in time -- + No environmental stochasticity -- + No demographic stochasticity --- # Interest of deterministic models + Project the average population trajectory -- + Species living in stable environments (such as protected areas) -- + Exponential growth (recolonization, recovering from over-exploitation) -- + Sensitivity analyses provide useful information to identify key parameters for population management --- class: center, middle # 2- Stochastic MPMs --- ## Environmental stochasticity - Unpredictable fluctuations in environmental conditions -- - Demographic rates vary over the years -- - Stochastic fluctuations in survival and reproductive rates reduce long-run population growth rate -- - The chance of occurrence of consecutive periods of unfavourable environmental conditions can drive small populations to extinctions -- - **MPMs in stochastic environments** project the population using a series of annual stochastic transition matrices (instead of one transition matrix) --- ## MPMs with environmental stochasticity **Principle** - generate a random sequence of demographic parameters -- - pile up the stochastic transition matrices -- <img src="img/stochmat1.png" width="1632" style="display: block; margin: auto;" /> --- ## Several options exist to include annual variability on demographic parameters ### Random annual variation around mean values <img src="img/stochdemo_mat.png" width="2307" style="display: block; margin: auto;" /> --- ## Several options exist to include annual variability on demographic parameters ### Catastrophic events <img src="img/stochmat3.png" width="1632" style="display: block; margin: auto;" /> --- ## MPMs with environmental stochasticity **Principle** - generate a random sequence of demographic parameters - pile up the stochastic transition matrices - iterate the population `$$N_{t+1} = A_{t} \cdot N_{t}$$` --- ## MPMs with environmental stochasticity **Principle** - generate a random sequence of demographic parameters - pile up the stochastic transition matrices - iterate the population `$$N_{t+1} = A_{t} \cdot N_{t}$$` - repeat several times -- - derive several quantities of interest (back to it in a bit) --- ## Analyze results of the simulations * population sizes at the last iteration -- * stochastic growth rate -- * probabilities of extinction -- * etc ... --- ## Example on the crested newt (*Triturus cristatus*) <img src="img/triton.jpg" width="50%" style="display: block; margin: auto;" /> --- ## Crested newt life cycle and transition matrix <img src="img/tritionlifecycle.png" width="2560" style="display: block; margin: auto;" /> --- ## The crested newt example High inter-annual variation in demographic parameters -- - Mean adult survival = 0.52 varied between 0.22 et 0.74 in 8 years of monitoring -- - Mean adult fecundity = 3.07 juvenile females / adult female varied form 0.31 to 5.40 in 10 years -- - Bad years with dry pond in spring (about 1 out of 3 years) induces quasi-complete failure of reproduction --- ## Step 1: Start from the deterministic transition matrix ``` [,1] [,2] [,3] [1,] 0.00 1.03152 1.2894 [2,] 0.48 0.00000 0.0000 [3,] 0.00 0.52000 0.5200 ``` --- ## Step 2: Draw annual values for each demographic parameters by sampling around the mean ```r library(truncnorm) ns <- 1500 s0stoch <- rtruncnorm(n=ns, a=0, b=1, mean=s0, sd=0.15) s1stoch <- rtruncnorm(n=ns, a=0, b=1, mean=s1,sd=0.15) s2stoch <- rtruncnorm(n=ns, a=0, b=1, mean=s2,sd=0.15) s3stoch <- rtruncnorm(n=ns, a=0, b=1, mean=s3,sd=0.15) fstoch <- rnorm(n=ns, mean=f,sd=0.5) ``` alternatively you can use the logit function --- ## Step 2: Draw annual values for each demographic parameters from a normal distribution <img src="matrixmodels2_files/figure-html/unnamed-chunk-85-1.png" width="40%" style="display: block; margin: auto;" /> --- ## Step 3: Pile up the stochastic transition matrices ```r # Create a list of stochastic transition matrices A.newtSE <- list() # fill in by sampling from distrib. of demo. param. for(i in 1:ns){ A.newtSE[[i]] <- matrix( c( 0, * sample(alpha2,1) * sample(fstoch,1) * sample(s0stoch,1), * sample(alpha3,1) * sample(fstoch,1) * sample(s0stoch,1), * sample(s1stoch,1), 0, 0, * 0, sample(s2stoch,1), sample(s3stoch,1) ) , nrow = 3, ncol = 3) } ``` --- ## Step 4: Project the population using the stochastic transition matrices ```r T <- 30 *runSE <- stoch.projection(A.newtSE, n0 = c(50, 50, 50), tmax = T, nreps = 1000, verbose = FALSE) runSE ``` --- ## Step 4: Project the population using the stochastic transition matrices ```r runSE <- stoch.projection(A.newtSE, # initial population * n0 = c(50, 50, 50), tmax = T, nreps = 1000, verbose = FALSE) runSE ``` --- ## Step 4: Project the population using the stochastic transition matrices ```r runSE <- stoch.projection(A.newtSE, n0 = c(50, 50, 50), # number of time steps to project over * tmax = T, nreps = 1000, verbose = FALSE) runSE ``` --- ## Step 4: Project the population using the stochastic transition matrices ```r runSE <- stoch.projection(A.newtSE, n0 = c(50, 50, 50), tmax = T, # number of repetitions * nreps = 1000, verbose = FALSE) head(runSE) ``` ``` [,1] [,2] [,3] [1,] 58.52637 152.56630 225.88360 [2,] 1770.05211 2815.25599 5284.29485 [3,] 98.21278 307.17791 258.54641 [4,] 12.53178 79.21635 83.97932 [5,] 161.18859 237.06524 308.14630 [6,] 405.44206 2167.06492 1570.76964 ``` --- ## Step 5: Examine the results Distribution of population sizes after `\(T=30\)` years <img src="matrixmodels2_files/figure-html/unnamed-chunk-91-1.png" width="40%" style="display: block; margin: auto;" /> --- ## Step 5: Examine the results Long-run stochastic growth rate `\(\lambda_{s}\)`: ```r lambdastoch <- stoch.growth.rate(A.newtSE, maxt = 5000, verbose = FALSE) names(lambdastoch) ``` ``` [1] "approx" "sim" "sim.CI" ``` --- ## Step 5: Examine the results Long-run stochastic growth rate `\(\lambda_{s}\)` on a log-scale: ```r lambdastoch$approx # by Tuljapukar's approximation ``` ``` [1] 0.03045988 ``` ```r lambdastoch$sim # by simulation ``` ``` [1] 0.03188456 ``` ```r lambdastoch$sim.CI # with confidence interval ``` ``` [1] 0.02638590 0.03738323 ``` --- ## Step 5: Examine the results Long-run stochastic growth rate `\(\lambda_{s}\)`: ```r exp(lambdastoch$approx) # exponentiate to get stochastic growth rate ``` ``` [1] 1.030929 ``` --- ## Step 5: Examine the results Probability of extinction: ```r proba.ext <- stoch.quasi.ext(A.newtSE, # initial population size * n0 = c(50,50, 50), Nx = 30, nreps = 1000, tmax=50, maxruns = 10, verbose = FALSE) ``` --- ## Step 5: Examine the results Probability of extinction: ```r proba.ext <- stoch.quasi.ext(A.newtSE, # initial population size * n0 = c(50, 50, 50), Nx = 30, nreps = 1000, tmax = 50, maxruns = 10, verbose = FALSE) ``` --- ## Step 5: Examine the results Probability of extinction: ```r proba.ext <- stoch.quasi.ext(A.newtSE, n = c(50, 50, 50), # quasi-extinction threshold * Nx = 30, nreps = 1000, tmax = 50, maxruns = 10, verbose = FALSE) ``` --- ## Step 5: Examine the results Probability of extinction: ```r proba.ext <- stoch.quasi.ext(A.newtSE, n = c(50, 50, 50), Nx = 30, # number of runs * nreps = 1000, tmax = 50, maxruns = 10, verbose = FALSE) ``` --- ## Step 5: Examine the results Probability of extinction: ```r proba.ext <- stoch.quasi.ext(A.newtSE, n = c(50, 50, 50), Nx = 30, nreps = 1000, # number of time steps * tmax = 50, maxruns = 10, verbose = FALSE) ``` --- ## Step 5: Examine the results Probability of extinction: ```r proba.ext <- stoch.quasi.ext(A.newtSE, n = c(50, 50, 50), Nx = 30, nreps = 1000, tmax = 50, # repeat to get robust estimates * maxruns = 10, verbose = FALSE) ``` --- ## Step 5: Examine the results Probability of extinction: <img src="matrixmodels2_files/figure-html/unnamed-chunk-101-1.png" width="40%" style="display: block; margin: auto;" /> --- ## Step 5: Examine the results Probability of extinction: ```r proba.ext.mean <- apply(proba.ext,1,mean) proba.ext.mean[20] # in 20 years ``` ``` [1] 0.0027 ``` --- ## Several options to generate stochastic transition matrices Include annual variability on demographic parameters: - Random annual variation around mean values - **Catastrophic events** --- ## Step 1: Define transition matrix in good vs bad years ```r fgood <- 3.07 # fecundity in normal years fbad <- 0 # fecundity in years with Spring dryness of the pond ``` --- ## Step 1: Define transition matrix in good vs bad years ```r # Transition matrix in normal years A.newtGood <- matrix( c(0, s0 * alpha2 * fgood, s0 * alpha3 * fgood, s1, 0, 0, 0, s2, s3) , nrow = 3, ncol = 3, byrow = TRUE) # Transition matrix in bad years A.newtBad <- matrix( c(0, s0 * alpha2 * fbad, s0 * alpha3 * fbad, s1, 0, 0, 0, s2, s3) , nrow = 3, ncol = 3, byrow = TRUE) ``` --- ## Step 1: Define transition matrix in good vs bad years ``` $good [,1] [,2] [,3] [1,] 0.00 1.03152 1.2894 [2,] 0.48 0.00000 0.0000 [3,] 0.00 0.52000 0.5200 $bad [,1] [,2] [,3] [1,] 0.00 0.00 0.00 [2,] 0.48 0.00 0.00 [3,] 0.00 0.52 0.52 ``` --- ## Step 2: Define frequency of catastrophic events ```r # Spring dyness of the pond occurs every 3 years in average freqbad <- 1/3 ``` --- ## Step 3: Project the population Project using the 'good' or 'bad' transition matrix with probabillity defined above: ```r *stochCATA <- stoch.projection(A.newtCATA, prob = c( (1-freqbad), (freqbad)), n0 = c(50, 50, 50), tmax = 100, nreps = 1000, verbose = FALSE) head(stochCATA) ``` ``` [,1] [,2] [,3] [1,] 1.662253e-03 0.000000e+00 6.703674e-04 [2,] 9.818829e-03 3.201417e-03 4.600103e-03 [3,] 0.000000e+00 0.000000e+00 2.646197e-03 [4,] 0.000000e+00 0.000000e+00 7.350913e-04 [5,] 6.182643e-01 0.000000e+00 2.880107e-01 [6,] 3.033455e-06 2.556872e-06 1.223357e-06 ``` --- # Step 3: Project the population Project using the 'good' or 'bad' transition matrix with probabillity defined above: ```r stochCATA <- stoch.projection(A.newtCATA, # frequency of castastrophic events * prob = c( (1-freqbad), (freqbad)), n0 = c(50, 50, 50), tmax = 100, nreps = 1000, verbose = FALSE) stochCATA ``` --- # Step 3: Project the population Project using the 'good' or 'bad' transition matrix with probabillity defined above: ```r stochCATA <- stoch.projection(A.newtCATA, prob = c( (1-freqbad), (freqbad)), # initial population * n0 = c(50, 50, 50), tmax = 100, nreps = 1000, verbose = FALSE) stochCATA ``` --- # Step 3: Project the population Project using the 'good' or 'bad' transition matrix with probabillity defined above: ```r stochCATA <- stoch.projection(A.newtCATA, prob = c( (1-freqbad), (freqbad)), n0 = c(50, 50, 50), # number of time steps * tmax = 100, nreps = 1000, verbose = FALSE) stochCATA ``` --- # Step 3: Project the population Project using the 'good' or 'bad' transition matrix with probabillity defined above: ```r stochCATA <- stoch.projection(A.newtCATA, prob = c( (1-freqbad), (freqbad)), n0 = c(50, 50, 50), tmax = 100, # number of replicates * nreps = 1000, verbose = FALSE) stochCATA ``` --- # Step 4: Examine the results Long-run stochastic growth rate: ```r lambdaCATA <- stoch.growth.rate(A.newtCATA, prob = c( (1-freqbad), (freqbad)), maxt = 5000, verbose = FALSE) exp(lambdaCATA$approx) ``` ``` [1] 0.8871908 ``` --- # Step 4: Examine the results Probability of extinction: ```r proba.extCATA <- stoch.quasi.ext(A.newtCATA, prob = c( (1-freqbad), (freqbad)), n = c(50, 15, 50), # initial population Nx = 30, # Quasi-extinction threshold nreps = 1000, # nb. of replicates tmax = 50, # nb. of time steps maxruns = 10, # nb of repetitions verbose = TRUE) ``` ``` Calculating extinction probability for run 1 ``` ``` Calculating extinction probability for run 2 ``` ``` Calculating extinction probability for run 3 ``` ``` Calculating extinction probability for run 4 ``` ``` Calculating extinction probability for run 5 ``` ``` Calculating extinction probability for run 6 ``` ``` Calculating extinction probability for run 7 ``` ``` Calculating extinction probability for run 8 ``` ``` Calculating extinction probability for run 9 ``` ``` Calculating extinction probability for run 10 ``` ```r proba.extCATA ``` ``` [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 [2,] 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 [3,] 0.037 0.037 0.030 0.042 0.026 0.042 0.032 0.030 0.032 0.034 [4,] 0.100 0.127 0.104 0.122 0.109 0.115 0.098 0.108 0.107 0.106 [5,] 0.155 0.195 0.166 0.190 0.185 0.187 0.164 0.166 0.167 0.169 [6,] 0.223 0.268 0.233 0.255 0.249 0.243 0.234 0.231 0.229 0.250 [7,] 0.285 0.334 0.306 0.308 0.301 0.320 0.290 0.299 0.295 0.323 [8,] 0.336 0.384 0.358 0.372 0.360 0.372 0.323 0.369 0.341 0.377 [9,] 0.406 0.423 0.416 0.414 0.419 0.427 0.372 0.420 0.403 0.427 [10,] 0.459 0.475 0.475 0.469 0.478 0.479 0.435 0.470 0.454 0.483 [11,] 0.515 0.523 0.533 0.525 0.519 0.540 0.485 0.520 0.514 0.539 [12,] 0.560 0.578 0.585 0.563 0.573 0.584 0.549 0.562 0.565 0.593 [13,] 0.601 0.628 0.625 0.601 0.610 0.639 0.592 0.613 0.619 0.644 [14,] 0.643 0.668 0.675 0.647 0.661 0.673 0.644 0.668 0.660 0.680 [15,] 0.680 0.697 0.713 0.688 0.703 0.704 0.684 0.701 0.691 0.705 [16,] 0.717 0.734 0.742 0.712 0.737 0.735 0.723 0.729 0.719 0.741 [17,] 0.744 0.766 0.771 0.740 0.760 0.759 0.753 0.766 0.751 0.760 [18,] 0.776 0.787 0.805 0.764 0.787 0.790 0.778 0.798 0.780 0.785 [19,] 0.802 0.812 0.830 0.797 0.810 0.813 0.806 0.820 0.818 0.800 [20,] 0.828 0.839 0.844 0.819 0.836 0.833 0.830 0.850 0.836 0.824 [21,] 0.851 0.856 0.867 0.836 0.850 0.848 0.858 0.872 0.853 0.845 [22,] 0.868 0.867 0.885 0.851 0.867 0.866 0.868 0.893 0.870 0.868 [23,] 0.884 0.885 0.895 0.862 0.877 0.883 0.883 0.909 0.880 0.883 [24,] 0.893 0.893 0.913 0.871 0.895 0.893 0.897 0.918 0.895 0.896 [25,] 0.899 0.911 0.923 0.882 0.903 0.906 0.909 0.926 0.911 0.902 [26,] 0.907 0.920 0.929 0.900 0.915 0.922 0.917 0.936 0.922 0.915 [27,] 0.923 0.929 0.933 0.912 0.930 0.929 0.924 0.942 0.931 0.926 [28,] 0.943 0.938 0.941 0.922 0.941 0.946 0.936 0.948 0.936 0.935 [29,] 0.949 0.948 0.949 0.932 0.953 0.950 0.943 0.955 0.944 0.940 [30,] 0.953 0.952 0.955 0.942 0.962 0.954 0.951 0.960 0.950 0.946 [31,] 0.958 0.961 0.962 0.951 0.971 0.961 0.954 0.967 0.955 0.951 [32,] 0.965 0.964 0.963 0.962 0.976 0.965 0.959 0.971 0.964 0.956 [33,] 0.966 0.967 0.966 0.970 0.977 0.967 0.962 0.976 0.966 0.963 [34,] 0.973 0.971 0.970 0.975 0.980 0.969 0.970 0.979 0.969 0.971 [35,] 0.974 0.977 0.974 0.979 0.981 0.971 0.974 0.981 0.973 0.977 [36,] 0.978 0.980 0.976 0.980 0.983 0.974 0.975 0.985 0.977 0.981 [37,] 0.979 0.985 0.980 0.982 0.984 0.977 0.977 0.985 0.980 0.982 [38,] 0.980 0.988 0.982 0.985 0.986 0.978 0.980 0.987 0.985 0.985 [39,] 0.985 0.990 0.986 0.988 0.987 0.981 0.982 0.989 0.989 0.988 [40,] 0.986 0.992 0.987 0.990 0.991 0.982 0.983 0.990 0.990 0.990 [41,] 0.991 0.992 0.988 0.992 0.992 0.987 0.985 0.991 0.990 0.993 [42,] 0.991 0.994 0.990 0.993 0.995 0.991 0.986 0.992 0.990 0.993 [43,] 0.992 0.994 0.991 0.993 0.997 0.994 0.987 0.992 0.992 0.994 [44,] 0.993 0.995 0.993 0.994 0.997 0.995 0.988 0.992 0.993 0.996 [45,] 0.994 0.995 0.997 0.995 0.997 0.996 0.988 0.992 0.993 0.996 [46,] 0.996 0.996 0.999 0.995 0.998 0.996 0.989 0.994 0.994 0.997 [47,] 0.998 0.996 0.999 0.996 0.999 0.998 0.992 0.994 0.995 0.997 [48,] 0.998 0.996 0.999 0.999 0.999 0.999 0.992 0.994 0.995 0.997 [49,] 0.998 0.996 0.999 0.999 0.999 0.999 0.992 0.994 0.995 0.997 [50,] 0.998 0.998 0.999 0.999 0.999 0.999 0.992 0.995 0.995 0.997 ``` --- # Step 4: Examine the results Probability of extinction: <img src="matrixmodels2_files/figure-html/unnamed-chunk-114-1.png" width="40%" style="display: block; margin: auto;" /> --- # Step 4: Examine the results Probability of extinction: ```r proba.ext.mean <- apply(proba.extCATA,1,mean) proba.ext.mean[20] # in 20 years ``` ``` [1] 0.8339 ``` --- ## Step 5: Interesting to evaluate the sensitivity of the results to changes in: - Amount of environmental variance -- - Frequency of catastrophic events -- - Initial population size -- - Extinction threshold -- - Number of time steps --- # Live demo on crested newts <img src="img/triton.jpg" width="50%" style="display: block; margin: auto;" /> --- # How is PVA useful ? + Quantify rate of population change over time -- + Estimate extinction risks (used by IUCN) -- + Identify key parameters for population management -- + Evaluate and compare the relative impact of population management actions --- ## Useful packages * Package 'popbio' (by Stubben, Milligan & Nantel) * Package ‘Rramas’ (by de la Cruz) * Package 'popdemo' (by Stott, Hodgson & Townley) * Package 'lefko3' (by Shefferson & Ehrlen) --- ## Assumptions and limitations - **One sex model: Females drive the demography, males are not limiting** (see e.g. Jenouvrier et al (2010) Mating behavior, population growth, and the operational sex ratio: a periodic two-sex model approach. The American Naturalist) -- - **Synchronous breeding: discrete time** (see e.g. Bacaer (2009). Periodic matrix population models: growth rate, basic reproduction number, and entropy. Bulletin of Mathematical Biology) -- - **No density-dependence** (see e.g. Caswell, Takada & Hunter (2004) Sensitivity analysis of equilibrium in density‐dependent matrix population models. Ecology Letters) --- ## Assumptions and limitations - **No demographic stochasticity** (see e.g. Vindenes, Engen & Sæther (BE.)2008) Individual heterogeneity in vital parameters and demographic stochasticity. The American Naturalist) -- - **No trend or change of mean in environmental variance** (possible to include environmental covariates) --- ## Always remember - we are making assumptions about demographic parameters and therefore environmental conditions in the future -- - short-term studies = underestimation of variance in demographic rates while extinction risks increase with increased temporal variance in pop size -- - keep a critical eyes on results: efficiency of PVAs are debated in the literature --- ## Pay attention to model construction **Further efforts may be required to educate biologists on the construction of MPMs** -- 3 errors commonly encountered in published MPMs (1) 34% of published studies failing to include survival in the fertility coefficient (2) 62% introducing a one-year delay in age at first reproduction (3) 53% incorrectly for estimating the asymptotic population growth rate or its sensitivity" *from Kendall et al (2019) Persistent problems in the construction of matrix population models in Ecological Modelling* --- ## In which cases constant MPMs are useful despite its limitations? - Good knowledge of the species biology, life cycle and estimate of demographic parameters are critical -- - Sensitivity analyses remain a powerful tool to identify key demographic parameters and evalute management actions. -- - Especially for big populations (little impact of demographic stochasticity) -- - When survival and fertility are best structured by age or stage --- # Other (more complicated) models exist - MPMs including demographic stochasticity for small populations - Density-dependent MPMs - Two-sex models for species with skewed reproductive success - Integral projection models for size-structured demography - Continuous-time models - Multi-species models... -- - the best model is not always the most complex one, it depends on the **species**, and the **data at hand** --- # Must-read references Morris, William, et al. (1999). *A practical handbook for population viability analysis.* The Nature Conservancy. Caswell, H. (2000). *Matrix population models (Vol. 1)*. Sunderland, MA: Sinauer. Brook, Barry W., et al. (2000). Predictive accuracy of population viability analysis in conservation biology. *Nature* 404.6776: 385-387. Beissinger, Steven R., and Dale R. McCullough, eds. (2002). *Population viability analysis*. University of Chicago Press. Reed, J. Michael, et al. (2002). Emerging issues in population viability analysis. *Conservation Biology* 16.1: 7-19.