Summary
AutoPop is an R-based joint-dynamic population simulation for diploids, triploids, and autotetraploids. For details on these methods, see Gaynor et al. 2023. To summarize, this package contains a stochastic stage-structured matrix population dynamics model for diploid, triploid, and autotetraploid perennial plants with overlapping generations. Specifically, this model includes three cytotypes (diploids, triploids, and autotetraploids), as well as two life stages, in regards to reproductive ability (immature and mature), for a total of six stages each represented by \(c_{i,j}\) where i = {2, 3, 4} and j = {1,2}.
Each time step includes three processes:
-
Reproduction
-
\(\color{green}{F_{i,k}(t)}\) is
the number of immature individuals of ploidal level \(i\) produced by mature individuals of
ploidal level \(k\).
-
\(\color{green}{F_{i,k}(t)}\) is
the number of immature individuals of ploidal level \(i\) produced by mature individuals of
ploidal level \(k\).
-
Survival
-
\(\color{red}{S_{i,j}(t)}\) is the
survival probability of individuals with ploidal level \(i\) and life stage \(j\).
-
\(\color{red}{S_{i,j}(t)}\) is the
survival probability of individuals with ploidal level \(i\) and life stage \(j\).
-
Maturation
- \(\color{blue}{M_{i}}\) is the maturation probability of immature individuals with ploidal level \(i\).
\[ = \left(\begin{array}{cc} c_{2,1}(t) \\ c_{3,1}(t) \\ c_{4,1}(t) \\ c_{2,2}(t) \\ c_{3,2}(t) \\ c_{4,2}(t) \\ \end{array}\right) \left(\begin{array}{cc} \color{red}{S_{2,1}(t)} & 0 & 0 & \color{green}{F_{2,2}(t)} & \color{green}{F_{2,3}(t)} & \color{green}{F_{2,4}(t)}\\ 0 & \color{red}{S_{3,1}(t)} & 0 & \color{green}{F_{3,2}(t)} & \color{green}{F_{3,3}(t)} & \color{green}{F_{3,4}(t)}\\ 0 & 0 & \color{red}{S_{4,1}(t)} & \color{green}{F_{4,2}(t)} & \color{green}{F_{4,3}(t)} & \color{green}{F_{4,4}(t)}\\ \color{blue}{M_{2}}& 0 & 0 & \color{red}{S_{2,2}(t)} & 0 & 0 \\ 0 & \color{blue}{M_{3}} & 0 & 0 & \color{red}{S_{3,2}(t)} & 0\\ 0 & 0 & \color{blue}{M_{4}} & 0 & 0 & \color{red}{S_{4,2}(t)} \\ \end{array}\right) \]
or in matrix notation: \[ \begin{equation} \mathbf{C}(t+1) = \mathbf{A}(t)\mathbf{C}(t), \end{equation} \]
Model Parameters
Table 1: Summary of model parameters and definitions. Color indicates processes: reproduction (green), survival (red), and maturation (blue).
| In-text | In-code | Definition |
|---|---|---|
| \(generations\) | \(generations\) | Number of generations to simulate. |
| \(init.pop\) | \(init.pop\) | Number of mature diploids in the initial founding population. |
| \(\color{green}{\lambda_{i}}\) | \(\color{green}{gnum.base}\) | Mean number of gametes per individual per cytotype. |
| \(\color{green}{d_{i}}\) | \(\color{green}{d}\) | Strength of density dependency on gamete production for each cytotype. |
| \(\color{green}{r_{j}}\) | \(\color{green}{rj}\) | Stage-specific density dependence. |
| \(\color{green}{\gamma_{i,k}}\) | \(\color{green}{gam.density.type}\) | Type-specific density dependence.Competitive effect of cytotype \(k\) on cytotype \(i\)’s reproduction. Options: “all” or “like-cytotype”. When equal to “all”, \(\gamma_{i,k} = 1\). Similarly, when equal to “like-cytotype”, \(\gamma_{i,k} = 0\) |
| \(\color{green}{b}\) | \(\color{green}{b}\) | Proportion of unreduced gamete formed by each diploid and tetraploid individual. |
| \(\color{green}{v}\) | \(\color{green}{cc}\) | Proportion of 3n gamete formation from each triploid individual. |
| \(\color{green}{s}\) | \(\color{green}{s}\) | Selfing rate. |
| \(\color{green}{mp}\) | \(\color{green}{mc}\) | Strength of mating proportion. |
| \(\color{green}{mate.lazy}\) | Prevents selfing from occurring during outcrossing. | |
| \(\color{red}{\rho}\) | \(\color{red}{rho}\) | Rank correlation between surivival of different cytotypes. |
| \(\color{red}{a_{i}}\) | \(\color{red}{aii.vec}\) | Cytotype-specific response to competition related to immature survival. |
| \(\color{red}{\alpha_{i,k}}\) | \(\color{red}{is.density.type}\) | Type-specific density dependence. Competitive effect of cytotype \(k\) on cytotype \(i\)’s survival. Options: “all” or “like-cytotype”. When equal to “all”, \(\alpha_{i,k} = 1\). Similarly, when equal to “like-cytotype”, \(\alpha_{i,k} = 0\) |
| \(\color{red}{y_{j}}\) | \(\color{red}{yj}\) | Stage-specific density dependence. |
| \(\color{red}{env.immature.survival}\) | Allows environmental stochasticity to influences immature survival. Defaults as “TRUE” | |
| \(\color{red}{\sigma}\) | \(\color{red}{env.sigma}\) | Proportion of environmental variance for immature individual. |
| \(\color{red}{\mu}\) | \(\color{red}{as.msurv}\) | Survival probability of an mature individual for each cytotype. |
| \(\color{red}{ci}\) | \(\color{red}{env.ci}\) | Proportion of environmental variance used to define mature survival rate per generation. |
| \(\color{blue}{as.matur}\) | \(\color{blue}{as.matur}\) | The probability of maturation from an immature stage to the mature stage for each cytotype. |
Simulation Example
To simulate multiple generations, we developed the function
gen.iter.f.choosy(). If you are exploring many parameter
sets, we recommend running model replicates in parallel using foreach and
doParallel.
Basic Example
# Load Packages
library(AutoPop)
# Example Simulation
example.sim <- gen.iter.f.choosy(generations = 10000, # Number of generations
init.pop = 100, # Number of mature diploids in initial pop
gnum.base = c(100, 100, 100), # Number of gametes per individual per cytotype
d = c(0.0001, 0.0009, 0.0001), # Strength of density dependency
gam.density.type = "all", # Type of density dependence for gamete production
b = 0.02, # Proportion of unreduced gamete formation
cc = 0.052, # Proportion of 3n gamete formation
s = 0.1, # Selfing rate
mc = 0.1, # Strength of mating choice
mate.lazy = FALSE, # Prevents selfing from occurring during outcrossing
rho = 0, # correlation coefficient
rj = c(1,1), # stage specific density dependence on gamete production
yj = c(1,1), # stage specific density dependence on immature surivial
env.ci = 0.1, # Proportion of environmental variance for mature individuals
env.immature.survival = TRUE,
env.sigma = 1, #Proportion of environmental variance for immature individuals
aii.vec = c(0.0005, 0.005, 0.0005), # Probability of survival of immature individuals
as.msurv = c(0.60, 0.06, 0.60), # Probability of survival of mature individuals
as.matur = c(0.4, 0.02, 0.4) # Maturation rate
)Parallel Example
Note, this example was designed to run on a SLURM-based cluster.
# Load Packages
library(foreach)
library(doParallel)
library(pryr)
# Detect CPIs available
## Eegister the cluster for using foreach
n.cpus <- as.numeric(Sys.getenv("SLURM_CPUS_PER_TASK"))
print(n.cpus)
cl <- makeCluster(n.cpus)
registerDoParallel(cl)
# Example Simulation
startparallel <- Sys.time()
temp <- foreach(iter = 1:500)%dopar%{
suppressWarnings(gen.iter.f.choosy(generations = 10000, # Number of generations
init.pop = 100, # Number of mature diploids in initial pop
gnum.base = c(100, 100, 100), # Number of gametes per individual per cytotype
d = c(0.0001, 0.0009, 0.0001), # Strength of density dependency
gam.density.type = "all", # Type of density dependence for gamete production
b = 0.02, # Proportion of unreduced gamete formation
cc = 0.052, # Proportion of 3n gamete formation
s = 0.1, # Selfing rate
mc = 0.1, # Strength of mating choice
mate.lazy = FALSE, # Prevents selfing from occurring during outcrossing
rho = 0, # correlation coefficient
rj = c(1,1), # stage specific density dependence on gamete production
yj = c(1,1), # stage specific density dependence on immature surivial
env.ci = 0.1, # Proportion of environmental variance for mature individuals
env.immature.survival = TRUE,
env.sigma = 1, #Proportion of environmental variance for immature individuals
aii.vec = c(0.0005, 0.005, 0.0005), # Probability of survival of immature individuals
as.msurv = c(0.60, 0.06, 0.60), # Probability of survival of mature individuals
as.matur = c(0.4, 0.02, 0.4) # Maturation rate
))
}
save(temp, file = "temp.RData")
## Benchmark
totalparallel <- Sys.time() - startparallel
print(totalparallel) # How much time did this simulation take?
print(pryr::mem_used()) # How much memory did this simulation use?