bumbl {bumbl}

Estimate colony growth, switch point, and decay parameters

Description

Fits generalized linear models that assume bumblebee colonies will switch from growth to gyne production at some point, \tau. This allows for a different switchpoint (\tau) for each colony, chosen by maximum likelihood methods.

Usage

bumbl(
  data,
  t,
  formula,
  family = gaussian(link = "log"),
  colonyID = NULL,
  augment = FALSE,
  keep.model = FALSE,
  tau_optim_maxit = 100,
  ...
)

Arguments

data

a dataframe or tibble with a column for colony ID (as a character or factor), a column for some measure of time (numeric), and a column for some measure of colony growth (numeric), at minimum.

t

the unquoted column name of the time variable.

formula

a formula with the form response ~ time + covariates where response is your measure of colony growth, time is whatever measure of time you have (date, number of weeks, etc.) and covariates are any optional co-variates you want to fit at the colony level.

family

a description of the error distribution and link function. This is passed to glm() except in the case of family = "negbin", which causes MASS::glm.nb() to be used to fit a negative binomial GLM.

colonyID

the unquoted column name of the colony ID variable. This is required, so to run bumbl() on a single colony, you must create a dummy column with a colony ID.

augment

when FALSE, bumbl returns a summary dataframe with one row for each colonyID. When TRUE, it returns the original data with additional columns containing model coefficients. Cannot be used in conjunction with keep.model = TRUE.

keep.model

If TRUE, then the output will contain a list-column with the models for each colony. This may be useful for extracting statistics and performing model diagnostics not provided by bumbl(). Learn more about working with list columns with vignette("nest", package = "tidyr").

tau_optim_maxit

passed to optim() which is used to find the optimal change point, tau. Mostly used for testing purposes (to force convergence errors), but could be increased if optimal switchpoint doesn't converge.

...

additional arguments passed to glm() or MASS::glm.nb().

Details

Colony growth is modeled as increasing exponentially until the colony switches from producing workers to producing reproductive individuals (drones and gynes), at which time the workers die and gynes leave the colony, causing the colony to decline. The switch point, \tau, may vary among colonies. bumbl() finds the value of \tau that maximizes likelihood and this "winning" model is used to calculate statistics returned in the output. This function works by fitting generalized linear models (GLMs) to modified colony growth data. Because of this, the assumptions for GLMs apply, namely independence and homogeneity of variance. See vignette("bumbl", package = "bumbl") for more details on the underlying math of the model.

Value

A data.frame with the additional class bumbldf containing a summary of the data with a row for every colony and the following columns:

• converged indicates whether the winning model converged.

• tau is the switchpoint, in the same units as t, for each colonyID. The colony grows for \tau weeks, then begins to decline in week \tau + 1.

• logN0 is the intercept of the growth function. It reflects actual initial colony size, if the colony initially grows exponentially. It would also be lower if there were a few weeks lag before growth started in the field.

• logLam is the average (log-scale) colony growth rate (i.e., rate of weight gain per unit t) during the growth period.

• decay reflects the rate of decline during the decline period. Equivalent to ln(\delta) - ln(\lambda) (see vignette for more in-depth explanation).

• logNmax is the maximum weight reached by each colony. It is a function of tau, logN0 and logLam

• Additional columns are coefficients for any covariates supplied in the formula

When augment = TRUE, the original data are returned with these columns as well as fitted values (.fitted) residuals (.resid) and standard error (.se.fit). When keep.model = TRUE a list-column with the glm models for each colony is returned as well.

Note

This function assumes there is a switchpoint and does not test whether the switchpoint model is significantly better than a log-linear model. As a result, it may estimate a switchpoint even if the data do not represent a true switchpoint. See the vignette for an example of how to extract the GLMs—you could compare them to a simpler log-linear model without the switchpoint by AIC or a likelihood ratio test to test the significance of the switchpoint.

References

Crone EE, Williams NM (2016) Bumble bee colony dynamics: quantifying the importance of land use and floral resources for colony growth and queen production. Ecology Letters 19:460–468. https://doi.org/10.1111/ele.12581

See Also

plot.bumbldf()

Examples

bumbl(bombus, colonyID = colony, t = week, formula = d.mass ~ week)