index.md

In levisc8/ipmr: Integral Projection Models

ipmr

ipmr is a package for implementing Integral Projection Models (IPMs) in R. It relies heavily on the mathematical syntax of the models, and does not try to abstract over the process of fitting vital rates. It is now relatively stable, though tweaks may be made. The News tab at the top of this page will contain more details on breaking changes as they are introduced to successive versions.

Below is a brief overview of how ipmr classifies different model types followed by examples of how to implement those types in this framework.

Installation

Install from CRAN:

install.packages("ipmr")

You can install the development version with the snippet below:

if(!require('remotes', quietly = TRUE)) {
  install.packages("remotes")
}

remotes::install_github("padrinoDB/ipmr", build_vignettes = TRUE)

Package scope

ipmr is intended to assist with IPM implementation and analysis. It is important to note that this package will not help with the process of fitting regression models for vital rates at all! That is a sufficiently different (and vast) topic that we decided it was not within the scope of this project. This will only help you turn those regression models into an IPM without shooting yourself in the foot. Furthermore, most of the documentation assumes a basic knowledge of IPM theory and construction. For those that are totally new to IPMs, it is strongly recommended to read a theoretical overview of the models first. Some favorites are Easterling et al. 2000, Merow et al. 2013, and Rees et al. 2014. The Introduction to ipmr vignette also contains a brief overview of IPM theory as well as a far more detailed introduction to this package. Thus, everything that follows here assumes you have some basic understanding of IPM theory, parameterized vital rate models, and are now ready to begin implementing your IPM.

Model classes

Once all parameters are estimated, the first step of defining a model in ipmr is to initialize the model using init_ipm(). This function has five arguments: sim_gen, di_dd, det_stoch, kern_param, and uses_age. We will ignore uses_age for now, because age-size models are less common and have their own vignette.

The combination of these arguments defines the type of IPM, and makes sure that the machinery for subsequent analyses works correctly. The possible entries for each argument are as follows:

• sim_gen: "simple"/"general"

  • A. simple: This describes an IPM with a single continuous state variable and no discrete stages.

  • B. general: This describes and IPM with either more than one continuous state variable, one or more discrete stages, or both of the above. Basically, anything other than an IPM with a single continuous state variable.

• di_dd: "di"/"dd"

  • A. di: This is used to denote a density-independent IPM.

  • B. dd: This is used to denote a density-dependent IPM.

• det_stoch: "det"/"stoch"

  • A. det: This is used to denote a deterministic IPM. If this is the third argument of init_ipm, kern_param must be left as NULL.

  • B. stoch: This is used to denote a stochastic IPM. If this is the third argument of init_ipm, kern_param must be specified. The two possibilities for the fourth are described next.

• kern_param: "kern"/"param" (Complete definitions found in Metcalf et al. 2015)

  • A. kern: This describes an IPM with discretely varying parameters such that their values are known before the model is specified. This is usually the case with models that estimate fixed and/or random year/site effects and for which defining a multivariate joint distribution to sample parameters from is not desirable/needed. These models can be a bit more computationally efficient than the param alternative because all kernels can be constructed before the iteration procedure begins, as opposed to requiring reconstruction for every single iteration.

  • B. param: This describes an IPM with parameters that are re-sampled from some distribution at each iteration of the model. For example, these could be a multivariate normal defined by covarying slopes and intercepts, or distributions of environmental variables that change from time to time. All that is required is that the parameters for the distribution are specified and that the function that generates the parameters at each iteration returns named lists that correspond to the parameter names in the model. Examples of this are available in the Introduction and General IPM vignettes.

The following possibilities are currently or will become available in ipmr (bold text denotes development progress):

• Simple, density independent models: Completed and ready

• "simple_di_det"

• "simple_di_stoch_kern"

• "simple_di_stoch_param"

• Simple, density dependent models: Completed, likely stable

• "simple_dd_det"

• "simple_dd_stoch_kern"

• "simple_dd_stoch_param"

• General, density independent models: Completed and ready

• "general_di_det"

• "general_di_stoch_kern"

• "general_di_stoch_param"

• General, density dependent models: Completed, likely stable

• "general_dd_det"

• "general_dd_stoch_kern"

• "general_dd_stoch_param"

Simple density-independent deterministic, simple kernel-resampled stochastic, and simple parameter resampled stochastic models (simple_di_det, simple_di_stoch_kern, simple_di_stoch_param) are described below as well as here. The general_* versions of these are described here. Density dependent versions are completed for simple and general models, and are probably stable, but have not been tested enough to be certain. A very brief (and incomplete) introduction is available here.

Below is an example of how to implement a simple IPM with ipmr.

Examples for simple IPM classes

Here is a simple model implemented with ipmr. This is a hypothetical plant species where plants can survive and grow (P(z′,z)), and reproduce sexually (F(z′,z)). There is no seed bank, and so there are no discrete traits in the model. We’ll use 4 regressions: survival (s(z)), growth (G(z′,z), fG), probability of reproducing (rr(z)), and number of seeds produced conditional on flowering (rs(z)). New recruits will be generated with a Gaussian distribution (denoted frd), which requires calculating the mean and standard deviation of new recruits from the data. For simplicity, we’ll assume there’s no maternal effect on recruit size. First, we’ll write out the functional forms for each component of the model:

1. n(z′,t+1) = ∫LUK(z′,z)n(z,t)dz

2. K(z′,z) = P(z′,z) + F(z′,z)

3. P(z′,z) = s(z) * G(z′,z)

4. Logit(s(z)) = αs + βs * z

5. G(z′,z) = fG(μG,σG)

6. μG = αG + βG * z

7. F(z′,z) = rr(z) * rs(z) * rd(z′)

8. Logit(rr(z)) = αrr + βrr * z

9. Log(rs(z)) = αrs + βrs * z

10. rd(z′) = frd(μrd,σrd)

Equation 1 describes how all the vital rates act on the initial trait distribution to produce a new one at t + 1. Equations 3-6 describe how existing individuals can survive, and if they survive, grow or shrink. Equations 7-10 describe how existing individuals create new individuals. In order to convert these functions to code, we usually fit regression models to our data. The following set of generalized linear models correspond to the functional forms described above:

1. Survival (s(z) / s): a generalized linear model w/ a logit link.

   • Example model formula: glm(surv ~ size_1, data = my_surv_data, family = binomial())

2. Growth (G(z′,z) / G): a linear model with a normal error distribution.

   • Example model formula: lm(size_2 ~ size_1, data = my_grow_data)

3. Pr(flowering) (rr(z) / r_r): a generalized linear model w/ a logit link.

   • Example model formula: glm(flower ~ size_1, data = my_repro_data, family = binomial())

4. Seed production (rs(z) / r_s): a generalized linear model w/ log link.

   • Example model formula: glm(seeds ~ size_1, data = my_flower_data, family = poisson())

5. Recruit size distribution (rd(z′) / r_d): a normal distribution w parameters mu_rd (mean) and sd_rd (standard deviation).

   • Example computations:

     • mu_rd = mean(seedling_data$size_2)

     • sd_rd = sd(seedling_data$size_2)

The example below assumes we’ve already fit our vital rate models from the raw data. In this example, the numbers are made up, but code that extracts the values you need from a real regression model is provided in comments.

# Load ipmr and get the parameter values. The data_list argument for define_kernel
# should hold every regression parameter and every constant used in the model.

library(ipmr)

my_data_list = list(s_int     = -2.2,  # coefficients(my_surv_mod)[1]
                    s_slope   = 0.25,  # coefficients(my_surv_mod)[2]
                    g_int     = 0.2,   # coefficients(my_grow_mod)[1]
                    g_slope   = 0.99,  # coefficients(my_grow_mod)[2]
                    sd_g      = 0.7,   # sd(resid(my_grow_mod))
                    r_r_int   = 0.003, # coefficients(my_pr_flower_mod)[1]
                    r_r_slope = 0.015, # coefficients(my_pr_flower_mod)[2]
                    r_s_int   = 0.45,  # coefficients(my_seed_mod)[1]
                    r_s_slope = 0.075, # coefficients(my_seed_mod)[2]
                    mu_fd     = 2,     # mean(recruit_data$size_next)
                    sd_fd     = 0.3)   # sd(recruit_data$size_next)

my_simple_ipm <- init_ipm(sim_gen = "simple",
                          di_dd   = "di",
                          det_stoch = "det")

my_simple_ipm <- define_kernel(

  proto_ipm = my_simple_ipm,

  # Name of the kernel

  name      = "P_simple",

  # The type of transition it describes (e.g. continuous - continuous, discrete - continuous).
  # These must be specified for all kernels!

  family    = "CC",

  # The formula for the kernel. We dont need to tack on the "z'/z"s here.  

  formula   = s * G,

  # A named set of expressions for the vital rates it includes. 
  # note the use of user-specified functions here. Additionally, each 
  # state variable has a stateVariable_1 and stateVariable_2, corresponding to
  # z and z' in the equations above. We don't need to define these variables ourselves,
  # just reference them correctly based on the way we've set up our model on paper.

  # Perform the inverse logit transformation to get survival probabilities
  # from your model. plogis from the "stats" package does this for us. 

  s         = plogis(s_int + s_slope * dbh_1), 

  # The growth model requires a function to compute the mean as a function of dbh.
  # The SD is a constant, so we don't need to define that in ... expression, 
  # just the data_list.

  G         = dnorm(dbh_2, mu_g, sd_g),
  mu_g      = g_int + g_slope * dbh_1,


  # Specify the constant parameters in the model in the data_list. 

  data_list = my_data_list,
  states    = list(c('dbh')),

  # If you want to correct for eviction, set evict_cor = TRUE and specify an
  # evict_fun. ipmr provides truncated_distributions() to help. This function
  # takes 2 arguments - the type of distribution, and the name of the parameter/
  # vital rate that it acts on.

  evict_cor = TRUE,
  evict_fun = truncated_distributions(fun    = 'norm',
                                      target = 'G')
  ) 

my_simple_ipm <- define_kernel(
  proto_ipm = my_simple_ipm,
  name      = 'F_simple',
  formula   = r_r * r_s * r_d,
  family    = 'CC',

  # Inverse logit transformation for flowering probability
  # (because we used a logistic regression)

  r_r       = plogis(r_r_int + r_r_slope * dbh_1),

  # Exponential function for seed progression 
  # (because we used a Poisson)

  r_s       = exp(r_s_int + r_s_slope * dbh_1),

  # The recruit size distribution has no maternal effect for size,
  # so mu_rd and sd_rd are constants. These get passed in the 
  # data_list

  r_d       = dnorm(dbh_2, mu_rd, sd_rd),
  data_list = my_data_list,
  states    = list(c('dbh')),

  # Again, we'll correct for eviction in new recruits by
  # truncating the normal distribution.

  evict_cor = TRUE,
  evict_fun = truncated_distributions(fun    = 'norm',
                                      target = 'r_d')
) 

# Next, we have to define the implementation details for the model. 
# We need to tell ipmr how each kernel is integrated, what state
# it starts on (i.e. z from above), and what state
# it ends on (i.e. z' above). In simple_* models, state_start and state_end will 
# always be the same, because we only have a single continuous state variable. 
# General_* models will be more complicated.

my_simple_ipm <- define_impl(
  proto_ipm = my_simple_ipm,
  make_impl_args_list(
    kernel_names = c("P_simple", "F_simple"),
    int_rule     = rep("midpoint", 2),
    state_start  = rep("dbh", 2),
    state_end    = rep("dbh", 2)
  )
) 

my_simple_ipm <- define_domains(
  proto_ipm = my_simple_ipm,
  dbh = c(0, # the first entry is the lower bound of the domain.
          50, # the second entry is the upper bound of the domain.
          100 # third entry is the number of meshpoints for the domain.
  ) 
) 

# Next, we define the initial state of the population. We must do this because
# ipmr computes everything through simulation, and simulations require a 
# population state.

my_simple_ipm <- define_pop_state(
  proto_ipm = my_simple_ipm,
  n_dbh = runif(100)
)

my_simple_ipm <- make_ipm(proto_ipm = my_simple_ipm)


lambda_ipmr <- lambda(my_simple_ipm)
w_ipmr      <- right_ev(my_simple_ipm)
v_ipmr      <- left_ev(my_simple_ipm)

# make_ipm_report works on either proto_ipms or made ipm objects

make_ipm_report(my_simple_ipm,
                title = "my_ipm_report",
                render_output = TRUE)

A more detailed introduction to ipmr is available in the Intro to ipmr article on this page. General IPM (ones with multiple continuous and/or discrete states) tutorials are also available there under the General IPMs article. Help with age × size models, density dependent models, and the parameter set index syntax are also available there. That’s all for now, happy modeling!

levisc8/ipmr documentation built on Feb. 22, 2023, 9:15 p.m.