Getting delta_st values from a random walk random field model

[seananderson]

Paraphrased question via email:

When predicting on a model with random walk random fields, the epsilon_st column gives the realized spatiotemporal values in space over time. How would I get the 'delta_st' values as described in the model description vignette that represent the MVN draws each year?

[seananderson]

You'd have to difference each year's values from the previous ones to get the MVN draws. Using predictions on a grid:

library(ggplot2)
library(dplyr)
library(sdmTMB)
theme_set(theme_light())

mesh <- make_mesh(pcod, c("X", "Y"), cutoff = 20)
fit <- sdmTMB(
  density ~ s(depth),
  data = pcod, mesh = mesh,
  spatiotemporal = "rw", spatial = "off",
  time = "year",
  family = tweedie(link = "log")
)
#> Detected irregular time spacing with an AR(1) or random walk process.
#> Consider filling in the missing time slices with `extra_time`.
#> `extra_time = c(2006, 2008, 2010, 2012, 2014, 2016)`

grid <- replicate_df(qcs_grid, "year", unique(pcod$year))
predictions <- predict(fit, newdata = grid)

ggplot(predictions, aes(X, Y, fill = epsilon_st)) +
  facet_wrap(~year) +
  geom_raster() +
  scale_fill_gradient2() +
  coord_fixed()

# get first difference of epsilon_st:
yrs <- sort(unique(pcod$year))
out <- list()
for (i in seq(2, length(yrs))) {
  last_yr <- filter(predictions, year == yrs[i - 1])
  this_yr <- filter(predictions, year == yrs[i])
  delta_st <- this_yr$epsilon_st - last_yr$epsilon_st
  out[[i]] <- data.frame(year = yrs[i], X = this_yr$X, Y = this_yr$Y, delta_st = delta_st)
}
d <- bind_rows(out)

ggplot(d, aes(X, Y, fill = delta_st)) +
  facet_wrap(~year) +
  geom_raster() +
  scale_fill_gradient2() +
  coord_fixed()

^{Created on 2023-02-17 with reprex v2.0.2}