forked from jonesor/Rage
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathREADME.Rmd
301 lines (236 loc) · 12.8 KB
/
README.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
---
output: github_document
---
<!-- README.md is generated from README.Rmd. Please edit that file -->
```{r, echo = FALSE}
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
fig.path = "man/figures/"
)
options(digits = 4)
```
# Rage
[](https://www.repostatus.org/#wip) [](https://travis-ci.org/jonesor/Rage) [](https://ci.appveyor.com/project/jonesor/Rage) [](https://codecov.io/github/jonesor/Rage?branch=devel)
Functions for calculating life history metrics from matrix population models (MPMs).
An R package for manipulating and analysing matrix population models (MPMs).
Includes functions for:
- deriving life history traits
- deriving life tables or life table components
- deriving vital rates
- perturbation analyses
- manipulating and transforming MPMs
## Installation
Install from GitHub with:
```{r, eval=FALSE}
# install.packages("remotes")
remotes::install_github("jonesor/Rage")
```
### Usage
```{r}
library(Rage)
```
### Loading an example MPM
The functions in Rage work on MPMs (or components of MPMs), so we'll start by
loading one of the example MPMs included in the Rage package (`mpm1`).
```{r}
library(Rage) # load Rage
data(mpm1) # load data object 'mpm1'
mpm1
```
The object `mpm1` is a list containing two elements: the growth/survival
component of the MPM (the __U__ matrix), and the sexual reproduction component
(the __F__ matrix). We can obtain the full MPM by adding the two components
together (__A__ = __U__ + __F__).
### Deriving life history traits from an MPM
One of the most common arguments among functions in Rage is `start`, which is
used to specify the stage class that represents the 'beginning of life' for the
purposes of calculation. Because the first stage class in `mpm1` is a 'seed'
stage, which we might consider functionally-distinct from the 'above-ground'
stages, we'll specify `start = 2` to set our starting stage class of interest to
the 'small' stage.
```{r}
life_expect(mpm1$matU, start = 2) # life expectancy
longevity(mpm1$matU, start = 2, lx_crit = 0.05) # longevity (age at lx = 0.05)
mature_age(mpm1$matU, mpm1$matF, start = 2) # mean age at first reproduction
mature_prob(mpm1$matU, mpm1$matF, start = 2) # prob survival to first repro
```
Some life history traits are independent of the starting stage class, in which
case we don't need to specify `start`.
```{r}
net_repro_rate(mpm1$matU, mpm1$matF) # net reproductive rate
gen_time(mpm1$matU, mpm1$matF) # generation time
```
Other life history traits are calculated from a life table rather than an MPM,
in which case we can first use the `mpm_to_` group of functions to derive the
necessary life table components.
```{r}
# first derive age-trajectories of survivorship (lx) and fecundity (mx)
lx <- mpm_to_lx(mpm1$matU, start = 2)
mx <- mpm_to_mx(mpm1$matU, mpm1$matF, start = 2)
# then calculate life history traits
entropy_k(lx) # Keyfitz' entropy
entropy_d(lx, mx) # Demetrius' entropy
shape_surv(lx) # shape of survival/mortality trajectory
shape_rep(lx) # shape of fecundity trajectory
```
### Life tables and the quasi-stationary distribution
Some MPMs are parameterized with a stasis loop at the maximum stage class, which
can lead to apparent plateaus in mortality or fertility trajectories derived
using age-from-stage methods. The function `qsd_converge()` can be used to
identify the time it takes for a cohort to reach the quasi-stationary
distribution (QSD). This quantity can then be used to subset age trajectories of
mortality or fertility to periods earlier than the QSD, so as to avoid
artefactual plateaus in mortality or fertility.
```{r, warning=FALSE, message=FALSE, fig.width=6.5, fig.height=4,cache=TRUE}
# derive life table from MPM
lt <- mpm_to_table(mpm1$matU, start = 2)
# calculate time to QSD
(q <- qsd_converge(mpm1$matU, start = 2))
# plot mortality trajectory w/ vertical line at time to QSD
par(mar = c(4.5, 4.5, 1, 1))
plot(qx ~ x, data = lt, type = "l", ylim = c(0, 0.65))
abline(v = q, lty = 2)
```
From the life table derived from `mpm1`, we can see a plateau in the mortality
rate (qx) beginning around age 5. However, this plateau corresponds to the QSD
and is therefore probably an artefact of the stasis loop rather than a
biological reality for the population represented by `mpm1`.
One approach to accounting for this artefactual plateau in subsequent life
history calculations is to limit our life table to the period prior to the QSD.
```{r}
# calculate the shape of the survival/mortality trajectory
shape_surv(lt$lx) # based on full lx trajectory
shape_surv(lt$lx[1:q]) # based on lx trajectory prior to the QSD
```
### Standardized vital rates
The transition rates that make up MPMs generally reflect products of two or more
vital rates (sometimes called ‘lower-level vital rates’). Assuming a
post-breeding census design, we can retroactively break apart each transition
rate into at least two vital rate components: survival, and ‘something’
conditional on survival. That ‘something’ might be growth, shrinkage, stasis,
dormancy, fecundity, or clonality.
##### Stage-specific vital rates (vector)
To summarize vital rates _within_ stage classes, we can use the `vr_vec_` group
of functions. We'll use the `exclude` argument here to exclude certain stage
classes ('seed' and 'dormant') from the calculation of certain vital rates (e.g.
we don't consider the large-to-dormant transition to actually represent
'growth').
```{r}
vr_vec_survival(mpm1$matU)
vr_vec_growth(mpm1$matU, exclude = c(1, 5))
vr_vec_shrinkage(mpm1$matU, exclude = 5)
vr_vec_stasis(mpm1$matU)
vr_vec_dorm_enter(mpm1$matU, dorm_stages = 5)
vr_vec_dorm_exit(mpm1$matU, dorm_stages = 5)
vr_vec_reproduction(mpm1$matU, mpm1$matF)
```
##### MPM-specific vital rates (scalar)
To summarize vital rates _across_ stage classes, we can use the `vr_` group of
functions. By default these functions take a simple average of the
stage-specific vital rates produced by the corresponding `vr_vec_` function.
However, here we'll demonstrate how to specify a _weighted_ average across
stages, based on the stable stage distribution at equilibrium (_w_).
```{r}
# derive full MPM (matA)
mpm1$matA <- mpm1$matU + mpm1$matF
# calculate stable stage distribution at equilibrium using popbio::stable.stage
library(popbio)
w <- popbio::stable.stage(mpm1$matA)
# calculate MPM-specific vital rates
vr_survival(mpm1$matU, exclude_col = c(1, 5), weights_col = w)
vr_growth(mpm1$matU, exclude = c(1, 5), weights_col = w)
vr_shrinkage(mpm1$matU, exclude = c(1, 5), weights_col = w)
vr_stasis(mpm1$matU, exclude = c(1, 5), weights_col = w)
vr_dorm_enter(mpm1$matU, dorm_stages = 5, weights_col = w)
vr_dorm_exit(mpm1$matU, dorm_stages = 5, weights_col = w)
vr_fecundity(mpm1$matU, mpm1$matF, weights_col = w)
```
Note how we've chosen to exclude the 'seed' and 'dormant' stage classes from our
vital rate summaries, because we consider these to be special classes (e.g.
'growth' from the 'seed' stage is really 'germination', which we may think of as
separate from somatic growth from 'small' to 'medium', or 'medium' to 'large').
### Perturbation analyses
The `perturb_matrix()` function measures the response of a demographic statistic
to perturbation of individual matrix elements (i.e. sensitivities and
elasticities). The `perturb_vr()` and `perturb_trans()` functions implement
perturbation analyses by vital rate type (survival, growth, etc.) and transition
type (stasis, retrogression, etc.), respectively.
```{r}
# matrix element perturbation
perturb_matrix(mpm1$matA, type = "sensitivity")
# vital rate perturbation
# (we use as.data.frame here for prettier printing)
as.data.frame(perturb_vr(mpm1$matU, mpm1$matF, type = "sensitivity"))
# transition type perturbation
as.data.frame(perturb_trans(mpm1$matU, mpm1$matF, type = "sensitivity"))
```
### Transforming MPMs
Rage includes a variety of functions that can be used to manipulate or transform
MPMs. For example, we can collapse an MPM to a smaller number of stage classes
using `mpm_collapse()`.
```{r}
# collapse 'small', 'medium', and 'large' stages into single stage class
col1 <- mpm_collapse(mpm1$matU, mpm1$matF, collapse = list(1, 2:4, 5))
col1$matA
```
The transition rates in the collapsed matrix are a weighted average of the
transition rates from the relevant stages of the original matrix, weighted by
the stable distribution at equilibrium. This process guarantees that the
collapsed MPM will retain the same population growth rate as the original.
However, other demographic and life history characteristics will not necessarily
be preserved.
```{r}
# compare population growth rate of original and collapsed MPM (preserved)
popbio::lambda(mpm1$matA)
popbio::lambda(col1$matA)
# compare net reproductive rate of original and collapsed MPM (not preserved)
net_repro_rate(mpm1$matU, mpm1$matF)
net_repro_rate(col1$matU, col1$matF)
```
## Complete list of functions
| Category | Function | Description |
|:---------------|:-------------|:-------------------------------------|
| Life history traits | `life_expect` | Life expectancy |
| | `longevity` | Longevity |
| | `net_repro_rate` | Net reproductive rate |
| | `gen_time` | Generation time |
| | `mature_age` | Age at reproductive maturity |
| | `mature_prob` | Probability of reaching reproductive maturity |
| | `mature_distrib` | Stage distribution of reproductive maturity |
| | `entropy_d` | Demetrius' entropy |
| | `entropy_k` | Keyfitz's entropy |
| | `shape_surv` | Shape of survival/mortality trajectory |
| | `shape_rep` | Shape of fecundity trajectory |
| Life table | `mpm_to_table` | MPM to life table |
| | `mpm_to_lx` | MPM to survivorship trajectory |
| | `mpm_to_px` | MPM to survival trajectory |
| | `mpm_to_hx` | MPM to mortality hazard trajectory |
| | `mpm_to_mx` | MPM to fecundity trajectory |
| | `lx_to_[px/hx]` | Convert from survivorship trajectory |
| | `px_to_[lx/hx]` | Convert from survival trajectory |
| | `hx_to_[lx/px]` | Convert from mortality hazard trajectory |
| | `qsd_converge` | Time to quasi-stationary distribution |
| Vital rates | `vr_[...]` | MPM-averaged vital rates
| | `vr_vec_[...]` | Stage-averaged vital rates
| | `vr_mat_[...]` | Survival-independent vital rates |
| Perturbation | `perturb_matrix` | Perturbation analysis of whole matrix |
| | `perturb_trans` | Perturbation analysis of transition types |
| | `perturb_vitals` | Perturbation analysis of vital rate types |
| | `perturb_stochastic` | Stochastic perturbation analysis |
| MPM transformation | `mpm_split` | Split MPM into survival and reproductive components |
| | `mpm_rearrange` | Rearrange MPM to segregate reproductive stages |
| | `mpm_collapse` | Collapse MPM to smaller number of stages |
| | `mpm_standardize` | Collapse MPM to standardized set of stages |
| | `standard_stages` | Group stages into standardized sets |
| | `repro_stages` | Identify reproductive stages |
| | `plot_life_cycle` | Plot a life cycle diagram |
=======
For a complete list of functions see the package [Reference](https://jonesor.github.io/Rage/reference/index.html) page.
## Contributions
All contributions are welcome. Please note that this project is released with a [Contributor Code of Conduct](CODE_OF_CONDUCT.md). By participating in this project you agree to abide by its terms.
There are numerous ways of contributing.
1. You can submit bug reports, suggestions etc. by [opening an issue](https://github.com/jonesor/Rcompadre/issues).
2. You can copy or fork the repository, make your own code edits and then send us a pull request. [Here's how to do that](https://jarv.is/notes/how-to-pull-request-fork-github/).
3. You can get to know us and join as a collaborator on the main repository.
4. You are also welcome to email us.