Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
Type: Package
Package: sdmTMB
Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
Version: 0.6.0.9023
Version: 0.6.0.9024
Authors@R: c(
person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
role = c("aut", "cre"),
Expand Down
4 changes: 4 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -64,6 +64,10 @@

* Add `get_eao()` to calculate effective area occupied.

* Add option for `area` to be passed in as the name of a column in the
dataframe to be used for area weighting. Used in `get_index()`,
`get_cog()`, `get_eao()`, etc.

# sdmTMB 0.6.0

* Pass several arguments to `DHARMa::plotQQunif()`.
Expand Down
27 changes: 25 additions & 2 deletions R/index.R
Original file line number Diff line number Diff line change
Expand Up @@ -5,8 +5,9 @@
#' @param bias_correct Should bias correction be implemented [TMB::sdreport()]?
#' @param level The confidence level.
#' @param area Grid cell area. A vector of length `newdata` from
#' [predict.sdmTMB()] or a value of length 1, which will be repeated
#' internally to match.
#' [predict.sdmTMB()] *or* a value of length 1 which will be repeated
#' internally to match *or* a character value representing the column
#' used for area weighting.
#' @param silent Silent?
#' @param ... Passed to [TMB::sdreport()].
#'
Expand Down Expand Up @@ -107,6 +108,11 @@
#' }
#' @export
get_index <- function(obj, bias_correct = FALSE, level = 0.95, area = 1, silent = TRUE, ...) {
# if offset is a character vector, use the value in the dataframe
if (is.character(area)) {
area <- obj$data[[area]]
}

d <- get_generic(obj, value_name = "link_total",
bias_correct = bias_correct, level = level, trans = exp, area = area, ...)
names(d)[names(d) == "trans_est"] <- "log_est"
Expand All @@ -119,6 +125,12 @@ get_index <- function(obj, bias_correct = FALSE, level = 0.95, area = 1, silent
#' @export
get_cog <- function(obj, bias_correct = FALSE, level = 0.95, format = c("long", "wide"), area = 1, silent = TRUE, ...) {
if (bias_correct) cli_abort("Bias correction with get_cog() is currently disabled.")

# if offset is a character vector, use the value in the dataframe
if (is.character(area)) {
area <- obj$data[[area]]
}

d <- get_generic(obj, value_name = c("cog_x", "cog_y"),
bias_correct = bias_correct, level = level, trans = I, area = area, ...)
d <- d[, names(d) != "trans_est", drop = FALSE]
Expand All @@ -145,6 +157,12 @@ get_eao <- function(obj,
...
) {
if (bias_correct) cli_abort("Bias correction with get_eao() is currently disabled.")

# if offset is a character vector, use the value in the dataframe
if (is.character(area)) {
area <- obj$data[[area]]
}

d <- get_generic(obj, value_name = c("log_eao"),
bias_correct = bias_correct, level = level, trans = exp, area = area, ...)
names(d)[names(d) == "trans_est"] <- "log_est"
Expand All @@ -155,6 +173,11 @@ get_eao <- function(obj,
get_generic <- function(obj, value_name, bias_correct = FALSE, level = 0.95,
trans = I, area = 1, silent = TRUE, ...) {

# if offset is a character vector, use the value in the dataframe
if (is.character(area)) {
area <- obj$data[[area]]
}

reinitialize(obj$fit_obj)

if ((!isTRUE(obj$do_index) && value_name[1] == "link_total") || value_name[1] == "cog_x" || value_name[[1]] == "log_eao") {
Expand Down
5 changes: 3 additions & 2 deletions man/get_index.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

4 changes: 3 additions & 1 deletion man/sdmTMBcontrol.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

37 changes: 37 additions & 0 deletions tests/testthat/test-index.R
Original file line number Diff line number Diff line change
Expand Up @@ -114,6 +114,43 @@ test_that("Index integration with area vector works with extra time and possibly
expect_equal(ind$se - ind0$se[ind0$year %in% seq(2011, 2017, 2)], c(0, 0, 0, 0))
})

test_that("get_index works", {
skip_on_cran()

pcod_spde <- make_mesh(pcod, c("X", "Y"), n_knots = 50, type = "kmeans")
m <- sdmTMB(
data = pcod,
formula = density ~ 0 + as.factor(year),
spatiotemporal = "off", # speed
time = "year", mesh = pcod_spde,
family = tweedie(link = "log")
)

# add some jittered area data to qcs_grid for testing
qcs_grid$area <- runif(nrow(qcs_grid), 0.9, 1.1)
nd <- replicate_df(qcs_grid, "year", unique(pcod$year))

predictions <- predict(m, newdata = nd, return_tmb_object = TRUE)

# get predictions with area passed as vector
ind <- get_index(predictions, area = nd$area)
# get predictions with area as a named column
ind2 <- get_index(predictions, area = "area")
expect_equal(ind, ind2)

# get predictions with area passed as vector
eao <- get_eao(predictions, area = nd$area)
# get predictions with area as a named column
eao2 <- get_eao(predictions, area = "area")
expect_equal(eao, eao2)

# get predictions with area passed as vector
cog <- get_cog(predictions, area = nd$area)
# get predictions with area as a named column
cog2 <- get_cog(predictions, area = "area")
expect_equal(cog, cog2)
})

# test_that("get_index faster epsilon bias correction", {
# skip_on_cran()
#
Expand Down
Loading