For this first example, we use the hypothetical life history of species 1. First we have to set spip up to run with that life history.
spip has a large number of demographic parameters. Typically
spip is run as a command-line program in Unix. In CKMRpop, all that action goes on under the hood, but you still have to use the
spip parameters. This vignette is not about using
spip. For a short listing of all the
spip options, do this:
If you want a full, complete, long listing of all the
spip options, then you can do:
All of the “long-form” options to
spip are given on the Unix command line starting with two dashes, like
--fem-surv-probs. To set parameters within
CKMRpop to send to spip, you simply make a named list of input values. The names of the items in the list are the long-format option names without the leading two dashes. For an example, see the package data object
species_1_life_history, as described below.
These parameters are included in the package in the variable
species_1_life_history. It is named list of parameters to send to spip. The list names are the names of the spip options. It looks like this:
species_1_life_history #> $`max-age` #>  20 #> #> $`fem-surv-probs` #>  0.75 0.76 0.76 0.77 0.77 0.78 0.78 0.79 0.79 0.80 0.80 0.80 0.81 0.81 0.82 #>  0.82 0.82 0.82 0.82 0.82 #> #> $`male-surv-probs` #>  0.75 0.76 0.76 0.77 0.77 0.78 0.78 0.79 0.79 0.80 0.80 0.80 0.81 0.81 0.82 #>  0.82 0.82 0.82 0.82 0.82 #> #> $`fem-prob-repro` #>  0.00 0.00 0.00 0.02 0.09 0.36 0.75 0.94 0.99 1.00 1.00 1.00 1.00 1.00 1.00 #>  1.00 1.00 1.00 1.00 1.00 #> #> $`male-prob-repro` #>  0.00 0.00 0.00 0.00 0.01 0.05 0.22 0.64 0.91 0.98 1.00 1.00 1.00 1.00 1.00 #>  1.00 1.00 1.00 1.00 1.00 #> #> $`fem-asrf` #>  0.00 0.00 0.00 3.56 3.72 3.88 4.04 4.20 4.36 4.52 4.68 4.84 5.00 5.16 5.32 #>  5.48 5.64 5.80 5.96 6.12 #> #> $`male-asrp` #>  0.00 0.00 0.00 3.56 3.72 3.88 4.04 4.20 4.36 4.52 4.68 4.84 5.00 5.16 5.32 #>  5.48 5.64 5.80 5.96 6.12 #> #> $`offsp-dsn` #>  "negbin" #> #> $`fem-rep-disp-par` #>  0.7 #> #> $`male-rep-disp-par` #>  0.7 #> #> $`mate-fidelity` #>  0.3 #> #> $`sex-ratio` #>  0.5
We want to add instructions to those, telling spip how long to run the simulation, and what the initial census sizes should be.
So, first, we copy
species_1_life_history to a new variable,
Now, we can add things to SPD.
The number of new fish added each year is called the “cohort-size”. Once we know that, we can figure out what the stable age distribution would be given the survival rates, and we can use that as our starting point. There is a function in the package that helps with that:
# before we tell spip what the cohort sizes are, we need to # tell it how long we will be running the simulation SPD$`number-of-years` <- 100 # run the sim forward for 100 years # this is our cohort size cohort_size <- 300 # Do some matrix algebra to compute starting values from the # stable age distribution: L <- leslie_from_spip(SPD, cohort_size) # then we add those to the spip parameters SPD$`initial-males` <- floor(L$stable_age_distro_fem) SPD$`initial-females` <- floor(L$stable_age_distro_male) # tell spip to use the cohort size SPD$`cohort-size` <- paste("const", cohort_size, collapse = " ")
Spip let’s you specify what fraction of fish of different ages should be sampled in different years. Here we do something simple, and instruct spip to sample 1% of the fish of ages 1, 2, and 3 (after the episode of death, see the spip vignette…) every year from year 50 to 75.
samp_frac <- 0.03 samp_start_year <- 50 samp_stop_year <- 75 SPD$`discard-all` <- 0 SPD$`gtyp-ppn-fem-post` <- paste( samp_start_year, "-", samp_stop_year, " ", samp_frac, " ", samp_frac, " ", samp_frac, " ", paste(rep(0, SPD$`max-age` - 3), collapse = " "), sep = "" ) SPD$`gtyp-ppn-male-post` <- SPD$`gtyp-ppn-fem-post`
There are two function that do all this for you. The function
run_spip() runs spip in a temporary directory. After running spip, it also processes the output with a few shell scripts. The function returns the path to the temporary directory. You pass that temporary directory path into the function
slurp_spip() to read the output back into R. It looks like this:
Note that setting the seed allows you to get the same results from spip. If you don’t set the seed, that is fine. spip will be seeded by the next two integers in current random number sequence.
If you are doing multiple runs and you want them to be different, you should make sure that you don’t inadvertently set the seed to be the same each time.
Although during massive production simulations, you might not go back to every run and summarize it to see what it looks like, when you are parameterizing demographic simulations you will want to be able to quickly look at observed demographic rates and things. There are a few functions in CKMRpop that make this quick and easy to do.
This is just a convenience function to make a pretty plot so you can check to see what the population demographics look like:
This shows that the function
leslie_from_spip() does a good job of finding the initial population numbers that accord with the stable age distribution.
We can compute the survival rates like this:
That returns a list. One part of the list is a tibble with observed survival fractions. The first 40 rows look like this:
surv_rates$survival_tibble %>% slice(1:40) #> # A tibble: 40 x 7 #> year pop age sex n cohort surv_fract #> <int> <int> <int> <chr> <int> <int> <dbl> #> 1 20 0 20 female 1 0 0 #> 2 20 0 19 female 2 1 0.5 #> 3 21 0 20 female 1 1 0 #> 4 20 0 18 female 1 2 1 #> 5 21 0 19 female 1 2 1 #> 6 22 0 20 female 1 2 0 #> 7 20 0 17 female 2 3 1 #> 8 21 0 18 female 2 3 1 #> 9 22 0 19 female 2 3 1 #> 10 23 0 20 female 2 3 0 #> # … with 30 more rows
The second part of the list holds a plot with histograms of age-specific, observed survival rates across all years. The blue line is the mean over all years.
To compare these values to the parameter values for the simulation, you must pass those to the function:
Here, the red dashed line is the value chosen as the parameter for the simulations. The means are particularly different for the older age classes, which makes sense because there the total number of individuals in each of those year classes is smaller.
It makes sense to check that your simulation is delivering a reasonable distribution of offspring per year. This is the number of offspring that survive to just before the first prekill census. Keep in mind that, for super high-fecundity species, we won’t model every single larva, we just don’t start “keeping track of them” until they reach a stage that is recognizable in some way.
We make this summary from the pedigree information. In order to get the number of adults that were present, but did not produce any offspring, we also need to pass in the postkill census information. Also, to get lifetime reproductive output, we need to know how old individuals were when they died, so we also pass in the information about deaths.
To make all the summaries, we do:
Note that we are setting the lifetime reproductive output hexbin width to be suitable for this example.
The function above returns a list of plots, as follows:
Especially when dealing with viviparous species (like sharks and mammals) it is worth checking this to make sure that there aren’t some females having far too many offspring.
Especially with long-lived organisms, it can be instructive to see how lifetime reproductive output varies with age at death.
Yep, many individuals have no offspring, and you have more kids if you live longer.
Out of all the offspring born each year, we can tabulate the fraction that were born to males (or females) of each age. This summary shows a histogram of those values. The represent the distribution of the fractional contribution of each age group each year.
The blue vertical lines show the means over all years.
Some of the parameters in
spip affect the distribution of the number of mates that each individual will have. We can have a quick look at whether the distribution of number of mates (that produced at least one offspring) appears to be what we might hope it to be.
That gives us a tibble with a summary, like this:
head(mates$mate_counts) #> # A tibble: 6 x 6 #> sex year pop parent num_offs num_mates #> <chr> <int> <int> <chr> <int> <int> #> 1 female 20 0 F0_0_0 2 1 #> 2 female 20 0 F10_0_1 5 4 #> 3 female 20 0 F10_0_11 3 1 #> 4 female 20 0 F10_0_12 1 1 #> 5 female 20 0 F10_0_13 2 1 #> 6 female 20 0 F10_0_2 3 2
And also a plot: