Add offset, similar to implementation in sdmTMB

R/fitting.R

@@ -5,7 +5,10 @@
 #' through time. It is also possible to fit a spatial random fields model
 #' without a time component.
 #'
-#' @param formula The model formula.
+#' @param formula The model formula. This is generally similar to formulas used
+#' in other packages, such as [lm()] or [glm()] with the exception of how the offset
+#' is handled. To include the offset, a named column of the dataframe must be
+#' `offset`, and the formula needs to include `+ offset`.
 #' @param data A data frame.
 #' @param time A character object giving the name of the time column. Leave
 #'   as `NULL` to fit a spatial GLMM without a time element.
@@ -219,6 +222,9 @@ glmmfields <- function(formula, data, lon, lat,
   y <- model.response(mf, "numeric")
   fixed_intercept <- ncol(X) == 0
 
+  offset <- as.vector(model.offset(mf))
+  if (is.null(offset)) offset <- rep(0, length(y))
+
   if (is.null(time)) {
     data$null_time_ <- 1
     time <- "null_time_"
@@ -264,6 +270,7 @@ glmmfields <- function(formula, data, lon, lat,
       prior_rw_sigma = parse_t_prior(prior_rw_sigma),
       prior_beta = parse_t_prior(prior_beta),
       prior_phi = parse_t_prior(prior_phi),
+      offset = offset,
       cov_func = switch(covariance,
         exponential = 0L,
         `squared-exponential` = 1L,

inst/stan/glmmfields.stan

@@ -8,6 +8,7 @@ data {
   int<lower=0> binomialN[N];
   real y[N]; // y for normal and gamma obs. model
   int y_int[N]; // y for NB2 or poisson or binomial obs. model
+  real offset[N]; // optional offset, is 0 if not included
   real prior_gp_theta[3];
   real prior_gp_sigma[3];
   real prior_sigma[3];
@@ -126,6 +127,7 @@ transformed parameters {
         y_hat[i] = X[i] * B + spatialEffects[yearID[i], stationID[i]] + yearEffects[yearID[i]];
       }
     }
+    y_hat[i] = y_hat[i] + offset[i]; // optional offset, additive in link space
   }
 
   if (obs_model==0) {