General tools (beyond Boreal Caribou)

bbousims contains more flexible, general functionality todo the following for any set of time periods and stages: generate survival and fecundity rates; generate survival, birth and age process matrices; project population; allocate individuals from population into groups; sample a proportion of groups in a given period each year.

Survival and fecundity rates

bbs_survival() and bbs_fecundity() can be used to stochastically generate survival and fecundity rates over time. Survival varies by period (e.g., month, season) with options to set the intercept, year trend, annual random effect, period random effect and period within annual random effect for each stage, on the log-odds scale.

The rates are provided as a vector with each element corresponding to each stage.

set.seed(1)
# 2 stages
survival <- bbs_survival(
  intercept = logit(c(0.94, 0.98)),
  trend = c(0, 0.3),
  annual_sd = rep(0.05, 2),
  period_sd = rep(0.1, 2),
  nyear = 3,
  nperiod_within_year = 4
)
survival$eSurvival
#> , , 1
#> 
#>           [,1]      [,2]      [,3]
#> [1,] 0.9409750 0.9431849 0.9403914
#> [2,] 0.9423532 0.9445145 0.9417824
#> [3,] 0.9414638 0.9436565 0.9408848
#> [4,] 0.9364146 0.9387841 0.9357890
#> 
#> , , 2
#> 
#>           [,1]      [,2]      [,3]
#> [1,] 0.9840585 0.9873763 0.9900676
#> [2,] 0.9821995 0.9858986 0.9889015
#> [3,] 0.9803428 0.9844216 0.9877351
#> [4,] 0.9770252 0.9817796 0.9856468

Fecundity varies by year, with options to set the intercept, year trend and annual random effect for each stage, on the log-odds scale. Setting intercept to NA will force all rates in that stage to be 0.

set.seed(1)
fecundity <- bbs_fecundity(
  intercept = c(NA, logit(0.4)),
  trend = c(0, -0.2),
  annual_sd = c(0, 0.1),
  nyear = 3
)
fecundity$eFecundity
#>      [,1]      [,2]
#> [1,]    0 0.3850636
#> [2,]    0 0.3573003
#> [3,]    0 0.2913105

Process matrices

Survival and fecundity rate arrays can be converted into process matrices for use in matrix model population projection. In this example there are two stages representing recruit and adult.
male_recruit_stage is set to NULL to indicate that there is only one recruit stage.

survival_mat <- bbs_matrix_survival_period(survival$eSurvival)
birth_mat <- bbs_matrix_birth_year(fecundity$eFecundity,
  female_recruit_stage = 1,
  male_recruit_stage = NULL
)
# first period, first year
survival_mat[, , 1, 1]
#>          [,1]      [,2]
#> [1,] 0.940975 0.0000000
#> [2,] 0.000000 0.9840585

# first year
birth_mat[, , 1]
#>      [,1]      [,2]
#> [1,]    1 0.1925318
#> [2,]    0 1.0000000

bbs_matrix_age() is used to generate an age process matrix, which does not vary in time. The input vector denotes which stage that stage will age into (i.e., stage 1 ages to stage 2, stage 2 collects).

age_mat <- bbs_matrix_age(c(2, 2))
age_mat
#>      [,1] [,2]
#> [1,]    0    0
#> [2,]    1    1

Population projection

bbousims contains a custom infix operator %*b% to perform matrix multiplication with stochasticity. Instead of multiplying the population at each stage by the corresponding rate, %*b% draws from a binomial distribution, where size is the population and probability is the rate of interest.

For example, with an initial population vector and survival process matrix for the first period of the first year

pop0 <- c(25, 52)
survival_mat1 <- survival_mat[, , 1, 1]

regular (deterministic) matrix multiplication yields

survival_mat1 %*% pop0
#>          [,1]
#> [1,] 23.52438
#> [2,] 51.17104

and stochastic matrix multiplication yields

set.seed(1)
survival_mat1 %*b% pop0
#>      [,1]
#> [1,]   24
#> [2,]   52

bbs_population() is used to project population forward in time from survival, age and birth process matrices. Survival occur at the end of each period, and survival, ageing and birth occur at the end of each year, in that order.

population <- bbs_population(pop0,
  birth = birth_mat,
  age = age_mat,
  survival = survival_mat
)
population
#>      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13]
#> [1,]   25   23   23   20   10   10   10   10   15    13    12    12    15
#> [2,]   52   50   48   47   66   65   63   62   68    68    68    68    79

Group allocation

bbs_population_groups() is used to randomly allocate individuals at each time step into groups. Group sizes are randomly drawn from a gamma-poisson distribution with user-specified lambda and theta (dispersion parameter) values. Minimum and maximum (i.e., as a maximum proportion of the total population) group sizes can be set. Group sizes are drawn until the cumulative size exceeds the total individuals. The remaining individuals comprise the final group. If the remaining number of individuals is < minimum group size specified by the user, these will be added to the previous group.

The output is a list of lists containing the individuals (identified by stage) in each group at each time period.

groups <- bbs_population_groups(population,
  group_size_lambda = 20,
  group_size_theta = 0
)

# first time period
groups[[1]]
#> [[1]]
#>  [1] 2 2 2 2 2 2 2 1 2 2 2 2 2 1 2 2
#> 
#> [[2]]
#>  [1] 1 2 2 2 2 2 1 2 1 1 2 1 1 2
#> 
#> [[3]]
#>  [1] 2 1 2 2 1 2 2 2 1 1 2 1 2 2 2 2 2 1 2
#> 
#> [[4]]
#>  [1] 2 2 1 1 2 2 1 2 2 1 2 2 2 2 1 2 2
#> 
#> [[5]]
#>  [1] 1 2 2 2 1 1 2 1 1 2 1

bbs_population_groups_pairs() behaves similarly but keeps calf-cow pairs (or equivalent for different animals) together in groups whenever possible.

bbs_population_groups_survey() allocates individuals into groups once a year (month of composition survey relative to biological year start) and allows the user to set the sample coverage (i.e., the proportion of groups observed). Groups are sampled randomly.

groups <- bbs_population_groups_survey(population,
  group_size_lambda = 20,
  month_composition = 3L,
  group_coverage = 0.2
)

# sampled groups in first time period
groups[[1]]
#> [[1]]
#>  [1] 2 1 2 1 2 2 2 2 2 1 2
#> 
#> [[2]]
#> [1] 1 2 1 2 1