## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(echo = TRUE)
library(aae.pop)

## ----eval = FALSE-------------------------------------------------------------
#  install.packages("aae.pop")

## ----eval = FALSE-------------------------------------------------------------
#  install.packages("remotes")
#  remotes::install_github("aae-stats/aae.pop")

## -----------------------------------------------------------------------------
popmat <- rbind(
  c(0,    0,    2,    4,    7),  # reproduction from 3-5 year olds
  c(0.25, 0,    0,    0,    0),  # survival from age 1 to 2
  c(0,    0.45, 0,    0,    0),  # survival from age 2 to 3
  c(0,    0,    0.70, 0,    0),  # survival from age 3 to 4
  c(0,    0,    0,    0.85, 0)   # survival from age 4 to 5
)

## -----------------------------------------------------------------------------
new_offspring <- c(2, 4, 7)
transition_probabilities <- c(0.25, 0.45, 0.70, 0.85)
popmat <- matrix(0, nrow = 5, ncol = 5)
popmat[reproduction(popmat, dims = 3:5)] <- new_offspring
popmat[transition(popmat)] <- transition_probabilities

## -----------------------------------------------------------------------------
popdyn <- dynamics(popmat)

## -----------------------------------------------------------------------------
plot(popdyn)

## -----------------------------------------------------------------------------
sims <- simulate(popdyn)

## ----dev = "png", fig.alt = "Line plot showing a single simulated trajectory initialised with a random initial condition. A single blue line shows deterministic changes in population abundance through time."----
plot(sims, col = "#2171B5")

## ----dev = "png"--------------------------------------------------------------
initials <- c(100, 50, 20, 10, 5)  # some initial conditions
sims <- simulate(
  popdyn,
  nsim = 100,
  init = initials,
  options = list(ntime = 20)
)

## -----------------------------------------------------------------------------
sims <- simulate(popdyn, nsim = 1000)
summary(sims)

## -----------------------------------------------------------------------------
pr_extinct(sims)
emps(sims)
risk_curve(sims, n = 10)

## -----------------------------------------------------------------------------
pr_extinct(sims, subset = 3:5, times = 40:50)
emps(sims, subset = 3:5, times = 40:50)
risk_curve(sims, subset = 3:5, times = 40:50, n = 10)

## -----------------------------------------------------------------------------
pr_extinct(sims, threshold = 100, subset = 3:5, times = 40:50)

## -----------------------------------------------------------------------------
risk_curve(sims, threshold = c(0, 10, 50, 100, 1000))

## -----------------------------------------------------------------------------
exps(sims, fun_within = quantile, fun_among = median, probs = 0.95)

## -----------------------------------------------------------------------------
# subset the population to adults 
sims <- subset(sims, subset = 3:5)

# drop the first 10 generations (the drop = FALSE
#   argument is a safeguard that keeps the third array
#   dimension when filtering to a single time step)
sims <- sims[, , 11:51, drop = FALSE]

# sum abundances over all classes, which
#   gives a matrix (2D array) with replicates
#   in rows and time steps in columns
abundance <- apply(sims, c(1, 3), sum)

# and calculate median over all trajectories, which
#   requires keeping the second dimension (time steps)
#   while iterating over the first
apply(abundance, 2, median)