get_index: Extract a relative biomass/abundance index, center of...

In sdmTMB: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'

Extract a relative biomass/abundance index, center of gravity, effective area occupied, or weighted average

Description

Extract a relative biomass/abundance index, center of gravity, effective area occupied, or weighted average

Usage

get_index(obj, bias_correct = TRUE, level = 0.95, area = 1, silent = TRUE, ...)

get_index_split(
  obj,
  newdata,
  bias_correct = FALSE,
  nsplit = 1,
  level = 0.95,
  area = 1,
  silent = FALSE,
  predict_args = list(),
  ...
)

get_cog(
  obj,
  bias_correct = FALSE,
  level = 0.95,
  format = c("long", "wide"),
  area = 1,
  silent = TRUE,
  ...
)

get_weighted_average(
  obj,
  vector,
  bias_correct = FALSE,
  level = 0.95,
  area = 1,
  silent = TRUE,
  ...
)

get_eao(obj, bias_correct = FALSE, level = 0.95, area = 1, silent = TRUE, ...)

Arguments

obj	Output from predict.sdmTMB() with return_tmb_object = TRUE (the usual case). Alternatively, if sdmTMB() was called with do_index = TRUE, or if using get_index_split(), an object from sdmTMB().
bias_correct	Should bias correction be implemented via TMB::sdreport()? Bias correction accounts for the non-linear transformation of random effects when calculating the index. Recommended to be TRUE for final analyses, but can be set to FALSE for faster calculation while experimenting with models. See Thorson and Kristensen (2016) in the References.
level	The confidence level.
area	Grid cell area for area weighting the index. Can be: (1) a numeric vector of length nrow(newdata) with area for each grid cell, (2) a single numeric value to apply to all grid cells, or (3) a character value giving the column name in newdata containing areas.
silent	Silent?
...	Passed to TMB::sdreport().
newdata	New data (e.g., a prediction grid by year) to pass to predict.sdmTMB() in the case of get_index_split().
nsplit	The number of splits to do the calculation in. For memory intensive operations (large grids and/or models), it can be helpful to do the prediction, area integration, and bias correction on subsets of time slices (e.g., years) instead of all at once. If nsplit > 1, this will usually be slower but with reduced memory use.
predict_args	A list of arguments to pass to predict.sdmTMB() in the case of get_index_split().
format	Long or wide.
vector	A numeric vector of the same length as the prediction data, containing the values to be averaged (e.g., depth, temperature).

Value

For get_index(): A data frame with columns for time, estimate (area-weighted total abundance or biomass), lower and upper confidence intervals, log estimate, and standard error of the log estimate.

For get_cog(): A data frame with columns for time, estimate (center of gravity: the abundance-weighted mean x and y coordinates), lower and upper confidence intervals, and standard error of center of gravity coordinates.

For get_eao(): A data frame with columns for time, estimate (effective area occupied: the area required if the population was spread evenly at the arithmetic mean density), lower and upper confidence intervals, log EAO, and standard error of the log EAO estimates.

For get_weighted_average(): A data frame with columns for time, estimate (weighted average of the provided vector, weighted by predicted density), lower and upper confidence intervals, and standard error of the estimates.

References

Geostatistical model-based indices of abundance (along with many newer papers):

Shelton, A.O., Thorson, J.T., Ward, E.J., and Feist, B.E. 2014. Spatial semiparametric models improve estimates of species abundance and distribution. Canadian Journal of Fisheries and Aquatic Sciences 71(11): 1655–1666. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1139/cjfas-2013-0508")}

Thorson, J.T., Shelton, A.O., Ward, E.J., and Skaug, H.J. 2015. Geostatistical delta-generalized linear mixed models improve precision for estimated abundance indices for West Coast groundfishes. ICES J. Mar. Sci. 72(5): 1297–1310. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1093/icesjms/fsu243")}

Geostatistical model-based centre of gravity:

Thorson, J.T., Pinsky, M.L., and Ward, E.J. 2016. Model-based inference for estimating shifts in species distribution, area occupied and centre of gravity. Methods Ecol Evol 7(8): 990–1002. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1111/2041-210X.12567")}

Geostatistical model-based effective area occupied:

Thorson, J.T., Rindorf, A., Gao, J., Hanselman, D.H., and Winker, H. 2016. Density-dependent changes in effective area occupied for sea-bottom-associated marine fishes. Proceedings of the Royal Society B: Biological Sciences 283(1840): 20161853. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1098/rspb.2016.1853")}

Bias correction:

Thorson, J.T., and Kristensen, K. 2016. Implementing a generic method for bias correction in statistical models using random effects, with spatial and population dynamics examples. Fisheries Research 175: 66–74. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1016/j.fishres.2015.11.016")}

See Also

get_index_sims()

Examples

library(ggplot2)

# use a small number of knots for this example to make it fast:
mesh <- make_mesh(pcod, c("X", "Y"), n_knots = 60)

# fit a spatiotemporal model:
m <- sdmTMB(
 data = pcod,
 formula = density ~ 0 + as.factor(year),
 time = "year", mesh = mesh, family = tweedie(link = "log")
)

# prepare a prediction grid:
nd <- replicate_df(qcs_grid, "year", unique(pcod$year))

# Note `return_tmb_object = TRUE` and the prediction grid:
predictions <- predict(m, newdata = nd, return_tmb_object = TRUE)

# biomass index:
ind <- get_index(predictions, bias_correct = TRUE)
ind
ggplot(ind, aes(year, est)) + geom_line() +
  geom_ribbon(aes(ymin = lwr, ymax = upr), alpha = 0.4) +
  ylim(0, NA)

# do that in 2 chunks
# only necessary for very large grids to save memory
# will be slower but save memory
# note the first argument is the model fit object:
ind <- get_index_split(m, newdata = nd, nsplit = 2, bias_correct = TRUE)

# center of gravity:
cog <- get_cog(predictions, format = "wide")
cog
ggplot(cog, aes(est_x, est_y, colour = year)) +
  geom_point() +
  geom_linerange(aes(xmin = lwr_x, xmax = upr_x)) +
  geom_linerange(aes(ymin = lwr_y, ymax = upr_y)) +
  scale_colour_viridis_c()

# effective area occupied:
eao <- get_eao(predictions)
eao
ggplot(eao, aes(year, est)) + geom_line() +
  geom_ribbon(aes(ymin = lwr, ymax = upr), alpha = 0.4) +
  ylim(0, NA)

# weighted average (e.g., depth-weighted by biomass):
wa <- get_weighted_average(predictions, vector = nd$depth)
wa
ggplot(wa, aes(year, est)) + geom_line() +
  geom_ribbon(aes(ymin = lwr, ymax = upr), alpha = 0.4)

sdmTMB documentation built on Nov. 26, 2025, 9:06 a.m.