From f5b4bcd72a255021219d5acb2f10b856cf7dd74c Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Tue, 18 Jun 2024 15:07:47 -0700
Subject: [PATCH 01/21] Fix truncated families #350

---
 DESCRIPTION                           |   2 +-
 NEWS.md                               |   6 ++
 R/families.R                          |  24 ++++-
 R/fit.R                               |  10 ++
 R/tmb-sim.R                           |   6 +-
 src/sdmTMB.cpp                        |  59 ++++++++++--
 tests/testthat/test-truncated-dists.R | 127 ++++++++++++++++++++++++++
 7 files changed, 219 insertions(+), 15 deletions(-)
 create mode 100644 tests/testthat/test-truncated-dists.R

diff --git a/DESCRIPTION b/DESCRIPTION
index acc06556..ec16ae5c 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9001
+Version: 0.6.0.9002
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/NEWS.md b/NEWS.md
index f7c4743e..6d89efcb 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,11 @@
 # sdmTMB (development version)
 
+* Fix bug in prediction from `delta_truncated_nbinom1()` and 
+  `delta_truncated_nbinom2()` families. The positive component
+  needs to be transformed to represent the mean of the *un*truncated
+  distribution first before multiplying by the probability of a non-zero.
+  Thanks to @tom-peatman #350
+
 * Add `get_eao()` to calculate effective area occupied.
 
 # sdmTMB 0.6.0
diff --git a/R/families.R b/R/families.R
index 7b2b8ae4..60e5b205 100644
--- a/R/families.R
+++ b/R/families.R
@@ -235,11 +235,21 @@ truncated_nbinom2 <- function(link = "log") {
     stats <- stats::make.link(linktemp)
   else if (is.character(link))
     stats <- stats::make.link(link)
-
+  linkinv <- function(eta, phi = NULL) {
+    s1 <- eta
+    if (is.null(phi)) phi <- .phi
+    s2 <- logspace_add(0, s1 - log(phi)) # log(1 + mu/phi)
+    log_nzprob <- logspace_sub(0, -phi * s2)
+    exp(eta) / exp(log_nzprob)
+  }
   structure(list(family = "truncated_nbinom2", link = linktemp, linkfun = stats$linkfun,
-    linkinv = stats$linkinv), class = "family")
+    linkinv = linkinv), class = "family")
 }
 
+logspace_sub <- function (lx, ly) lx + log1mexp(lx - ly)
+logspace_add <- function (lx, ly) pmax(lx, ly) + log1p(exp(-abs(lx - ly)))
+log1mexp <- function(x) ifelse(x <= log(2), log(-expm1(-x)), log1p(-exp(-x)))
+
 #' @export
 #' @examples
 #' truncated_nbinom1(link = "log")
@@ -253,9 +263,15 @@ truncated_nbinom1 <- function(link = "log") {
     stats <- stats::make.link(linktemp)
   else if (is.character(link))
     stats <- stats::make.link(link)
-
+  linkinv <- function(eta, phi = NULL) {
+    mu <- exp(eta)
+    if (is.null(phi)) phi <- .phi
+    s2 <- logspace_add(0, log(phi)) # log(1 + phi)
+    log_nzprob <- logspace_sub(0, -mu / phi * s2) # 1 - prob(0)
+    mu / exp(log_nzprob)
+  }
   structure(list(family = "truncated_nbinom1", link = linktemp, linkfun = stats$linkfun,
-    linkinv = stats$linkinv), class = "family")
+    linkinv = linkinv), class = "family")
 }
 
 #' @param df Student-t degrees of freedom fixed value parameter.
diff --git a/R/fit.R b/R/fit.R
index 538486fd..c53eeb68 100644
--- a/R/fit.R
+++ b/R/fit.R
@@ -1501,6 +1501,16 @@ sdmTMB <- function(
   sd_report <- TMB::sdreport(tmb_obj, getJointPrecision = get_joint_precision)
   conv <- get_convergence_diagnostics(sd_report)
 
+  ## save params that families need to grab from environments:
+  if (any(family$family %in% c("truncated_nbinom1", "truncated_nbinom2"))) {
+    phi <- exp(tmb_obj$par[["ln_phi"]])
+    if (delta) {
+      assign(".phi", phi, environment(out_structure[["family"]][[2]][["linkinv"]]))
+    } else {
+      assign(".phi", phi, environment(out_structure[["family"]][["linkinv"]]))
+    }
+  }
+
   out_structure$tmb_obj <- tmb_obj
   out <- c(out_structure, list(
     model      = tmb_opt,
diff --git a/R/tmb-sim.R b/R/tmb-sim.R
index 7cac7dc8..08cc8e7a 100644
--- a/R/tmb-sim.R
+++ b/R/tmb-sim.R
@@ -271,7 +271,11 @@ sdmTMB_simulate <- function(formula,
   d[["omega_s"]] <- if (all(s$omega_s_A != 0)) s$omega_s_A
   d[["epsilon_st"]] <- if (all(s$epsilon_st_A_vec != 0)) s$epsilon_st_A_vec
   d[["zeta_s"]] <- if (all(s$zeta_s_A != 0)) s$zeta_s_A
-  d[["mu"]] <- family$linkinv(s$eta_i)
+  if (any(family$family %in% c("truncated_nbinom1", "truncated_nbinom2"))) {
+    d[["mu"]] <- family$linkinv(s$eta_i, phi = phi)
+  } else {
+    d[["mu"]] <- family$linkinv(s$eta_i)
+  }
   d[["eta"]] <- s$eta_i
   d[["observed"]] <- s$y_i
   d <- do.call("data.frame", d)
diff --git a/src/sdmTMB.cpp b/src/sdmTMB.cpp
index 7c3891f0..c76fab17 100644
--- a/src/sdmTMB.cpp
+++ b/src/sdmTMB.cpp
@@ -98,6 +98,30 @@ Type Link(Type eta, int link)
   return out;
 }
 
+// Modified from glmmTMB:
+/* log-prob of non-zero value in conditional distribution  */
+template<class Type>
+Type calc_log_nzprob(Type mu, Type phi, int family) {
+  Type ans, s1, s2;
+  switch (family) {
+  case truncated_nbinom1_family:
+    s2 = logspace_add(Type(0), log(phi));      // log(1. + phi(i)
+    ans = logspace_sub(Type(0), -mu / phi * s2); // 1-prob(0)
+    break;
+  case truncated_nbinom2_family:
+    s1 = log(mu);
+    // s2 := log( 1. + mu(i) / phi(i) )
+    s2 = logspace_add(Type(0), s1 - log(phi));
+    ans = logspace_sub(Type(0), -phi * s2);
+    break;
+  // case truncated_poisson_family:
+  //   ans = logspace_sub(Type(0), -mu);  // log(1-exp(-mu(i))) = P(x>0)
+  //   break;
+  default: ans = Type(0);
+  }
+  return ans;
+}
+
 // ------------------ Main TMB template ----------------------------------------
 
 template <class Type>
@@ -755,21 +779,21 @@ Type objective_function<Type>::operator()()
   Type s1, s2, s3, lognzprob, tmp_ll, ll_1, ll_2, p_mix, mix_ratio, tweedie_p, s1_large, s2_large;
 
   // calcs for mix distr. first:
-  int mix_model;
+  int pos_model;
   if (n_m > 1) {
-    mix_model = 1;
+    pos_model = 1;
   } else {
-    mix_model = 0;
+    pos_model = 0;
   }
   vector<Type> mu_i_large(n_i);
-  switch (family(mix_model)) {
+  switch (family(pos_model)) {
   case gamma_mix_family:
   case lognormal_mix_family:
   case nbinom2_mix_family: {
     p_mix = invlogit(logit_p_mix); // probability of larger event
     mix_ratio = exp(log_ratio_mix) + Type(1.); // ratio of large:small values, constrained > 1.0
     for (int i = 0; i < n_i; i++) {
-      mu_i_large(i) = exp(log(mu_i(i, mix_model)) + log(mix_ratio));  // mean of large component = mean of smaller * ratio
+      mu_i_large(i) = exp(log(mu_i(i, pos_model)) + log(mix_ratio));  // mean of large component = mean of smaller * ratio
     }
     ADREPORT(logit_p_mix);
     ADREPORT(log_ratio_mix);
@@ -1186,12 +1210,12 @@ Type objective_function<Type>::operator()()
     // for families that implement mixture models, adjust proj_eta by
     // proportion and ratio of means
     // (1 - p_mix) * mu_i(i,m) + p_mix * (mu(i,m) * mix_ratio);
-    switch (family(mix_model)) {
+    switch (family(pos_model)) {
       case gamma_mix_family:
       case lognormal_mix_family:
       case nbinom2_mix_family:
-        proj_eta.col(mix_model) = log((1. - p_mix) * exp(proj_eta.col(mix_model)) + // regular part
-               p_mix * exp(proj_eta.col(mix_model)) * mix_ratio); //large part
+        proj_eta.col(pos_model) = log((1. - p_mix) * exp(proj_eta.col(pos_model)) + // regular part
+               p_mix * exp(proj_eta.col(pos_model)) * mix_ratio); //large part
         break;
       default:
         break;
@@ -1237,6 +1261,17 @@ Type objective_function<Type>::operator()()
     vector<Type> mu_combined(n_p);
     mu_combined.setZero();
 
+    int truncated_dist;
+    switch (family(pos_model)) {
+      case truncated_nbinom1_family:
+      case truncated_nbinom2_family:
+        truncated_dist = 1;
+        break;
+      default:
+        truncated_dist = 0;
+        break;
+    }
+
     if (calc_index_totals || calc_cog || calc_eao) {
       // ------------------ Derived quantities ---------------------------------
       Type t1;
@@ -1249,7 +1284,13 @@ Type objective_function<Type>::operator()()
             // Type R1 = Type(1.) - exp(-exp(proj_eta(i,0)));
             // Type R2 = exp(proj_eta(i,0)) / R1 * exp(proj_eta(i,1))
             mu_combined(i) = exp(proj_eta(i,0) + proj_eta(i,1)); // prevent numerical issues
-          } else {
+          } else if (truncated_dist) {
+            t1 = InverseLink(proj_eta(i,0), link(0));
+            // convert from mean of *un-truncated* to mean of *truncated* distribution
+            Type log_nzprob = calc_log_nzprob(exp(proj_eta(i,1)), phi(1), family(1));
+            t2 = exp(proj_eta(i,1)) / exp(log_nzprob);
+            mu_combined(i) = t1 * t2;
+          } else  {
             t1 = InverseLink(proj_eta(i,0), link(0));
             t2 = InverseLink(proj_eta(i,1), link(1));
             mu_combined(i) = t1 * t2;
diff --git a/tests/testthat/test-truncated-dists.R b/tests/testthat/test-truncated-dists.R
new file mode 100644
index 00000000..978e51c1
--- /dev/null
+++ b/tests/testthat/test-truncated-dists.R
@@ -0,0 +1,127 @@
+test_that("Delta-truncated NB2 works with various types of predictions", {
+  skip_on_cran()
+  skip_if_not_installed("glmmTMB")
+
+  set.seed(1)
+  d0 <- data.frame(y = rnbinom(n = 1e3, mu = 3, size = 5))
+
+  fit <- sdmTMB(
+    y ~ 1,
+    data = d0,
+    spatial = FALSE,
+    family = delta_truncated_nbinom2()
+  )
+  fit$sd_report
+
+  fit2 <- glmmTMB::glmmTMB(
+    y ~ 1,
+    ziformula = ~ 1,
+    data = d0,
+    family = glmmTMB::truncated_nbinom2()
+  )
+
+  p <- predict(fit)
+  p2 <- predict(fit2)
+  expect_equal(p2, p$est2, tolerance = 1e-4)
+
+  p <- predict(fit, type = "response")
+  p2 <- predict(fit2, type = "response")
+  expect_equal(p2, p$est, tolerance = 1e-4)
+
+  p <- predict(fit)
+  p2 <- predict(fit2, type = "zlink")
+  expect_equal(-p2, p$est1, tolerance = 1e-4)
+
+  p <- predict(fit, type = "response")
+  p2 <- predict(fit2, type = "conditional")
+  expect_equal(p2, p$est2, tolerance = 1e-4)
+})
+
+test_that("Delta-truncated NB1 works with various types of predictions", {
+  skip_on_cran()
+  skip_if_not_installed("glmmTMB")
+
+  set.seed(1)
+  d0 <- data.frame(y = rnbinom(n = 1e3, mu = 3, size = 5))
+
+  fit <- sdmTMB(
+    y ~ 1,
+    data = d0,
+    spatial = FALSE,
+    family = delta_truncated_nbinom1()
+  )
+  fit$sd_report
+
+  fit2 <- glmmTMB::glmmTMB(
+    y ~ 1,
+    ziformula = ~ 1,
+    data = d0,
+    family = glmmTMB::truncated_nbinom1()
+  )
+
+  p <- predict(fit)
+  p2 <- predict(fit2)
+  expect_equal(p2, p$est2, tolerance = 1e-4)
+
+  p <- predict(fit, type = "response")
+  p2 <- predict(fit2, type = "response")
+  expect_equal(p2, p$est, tolerance = 1e-4)
+
+  p <- predict(fit)
+  p2 <- predict(fit2, type = "zlink")
+  expect_equal(-p2, p$est1, tolerance = 1e-4)
+
+  p <- predict(fit, type = "response")
+  p2 <- predict(fit2, type = "conditional")
+  expect_equal(p2, p$est2, tolerance = 1e-4)
+})
+
+test_that("Truncated NB2 works with various types of predictions", {
+  skip_on_cran()
+  skip_if_not_installed("glmmTMB")
+  set.seed(1)
+  d0 <- data.frame(y = rnbinom(n = 1e3, mu = 3, size = 5))
+  d0 <- subset(d0, y > 0)
+  fit <- sdmTMB(
+    y ~ 1,
+    data = d0,
+    spatial = FALSE,
+    family = truncated_nbinom2()
+  )
+  fit2 <- glmmTMB::glmmTMB(
+    y ~ 1,
+    data = d0,
+    family = glmmTMB::truncated_nbinom2()
+  )
+  p <- predict(fit)
+  p2 <- predict(fit2)
+  expect_equal(p2, p$est, tolerance = 1e-4)
+  p <- predict(fit, type = "response")
+  p2 <- predict(fit2, type = "response")
+  expect_equal(p2, p$est, tolerance = 1e-4)
+})
+
+test_that("Truncated NB1 works with various types of predictions", {
+  skip_on_cran()
+  skip_if_not_installed("glmmTMB")
+  set.seed(1)
+  d0 <- data.frame(y = rnbinom(n = 1e3, mu = 3, size = 5))
+  d0 <- subset(d0, y > 0)
+  fit <- sdmTMB(
+    y ~ 1,
+    data = d0,
+    spatial = FALSE,
+    family = truncated_nbinom1()
+  )
+  fit2 <- glmmTMB::glmmTMB(
+    y ~ 1,
+    data = d0,
+    family = glmmTMB::truncated_nbinom1()
+  )
+  p <- predict(fit)
+  p2 <- predict(fit2)
+  expect_equal(p2, p$est, tolerance = 1e-4)
+  p <- predict(fit, type = "response")
+  p2 <- predict(fit2, type = "response")
+  expect_equal(p2, p$est, tolerance = 1e-4)
+})

From 75b4ecc0ff1d93f5e6c3b9fb44fb00ca17dc0172 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Tue, 18 Jun 2024 15:25:54 -0700
Subject: [PATCH 02/21] Fix indexes for truncated NB1/NB2 #350

---
 src/sdmTMB.cpp                        |  8 +++-
 tests/testthat/test-truncated-dists.R | 64 +++++++++++++++++++++++++++
 2 files changed, 71 insertions(+), 1 deletion(-)

diff --git a/src/sdmTMB.cpp b/src/sdmTMB.cpp
index c76fab17..3adee6db 100644
--- a/src/sdmTMB.cpp
+++ b/src/sdmTMB.cpp
@@ -1297,7 +1297,13 @@ Type objective_function<Type>::operator()()
           }
           total(proj_year(i)) += mu_combined(i) * area_i(i);
         } else { // non-delta model
-          mu_combined(i) = InverseLink(proj_eta(i,0), link(0));
+          if (truncated_dist) {
+            // convert from mean of *un-truncated* to mean of *truncated* distribution
+            Type log_nzprob = calc_log_nzprob(exp(proj_eta(i,0)), phi(0), family(0));
+            mu_combined(i) = exp(proj_eta(i,0)) / exp(log_nzprob);
+          } else  {
+            mu_combined(i) = InverseLink(proj_eta(i,0), link(0));
+          }
           total(proj_year(i)) += mu_combined(i) * area_i(i);
         }
       }
diff --git a/tests/testthat/test-truncated-dists.R b/tests/testthat/test-truncated-dists.R
index 978e51c1..2039e41a 100644
--- a/tests/testthat/test-truncated-dists.R
+++ b/tests/testthat/test-truncated-dists.R
@@ -125,3 +125,67 @@ test_that("Truncated NB1 works with various types of predictions", {
   p2 <- predict(fit2, type = "response")
   expect_equal(p2, p$est, tolerance = 1e-4)
 })
+
+test_that("Truncated NB1/2 indexes are right", {
+  skip_on_cran()
+  set.seed(1)
+  d0 <- data.frame(y = rnbinom(n = 1e3, mu = 3, size = 5))
+  d0$year <- 1
+
+  # D-TNB2
+  fit <- sdmTMB(
+    y ~ 1,
+    data = d0,
+    time = "year",
+    spatiotemporal = "off",
+    spatial = FALSE,
+    family = delta_truncated_nbinom2()
+  )
+  p <- predict(fit, newdata = d0, return_tmb_object = TRUE)
+  x <- get_index(p, area = 1/nrow(d0))
+  pp <- predict(fit, type = "response")
+  expect_equal(x$est, pp$est[1], tolerance = 1e-3)
+
+  # D-TNB1
+  fit <- sdmTMB(
+    y ~ 1,
+    data = d0,
+    time = "year",
+    spatiotemporal = "off",
+    spatial = FALSE,
+    family = delta_truncated_nbinom1()
+  )
+  p <- predict(fit, newdata = d0, return_tmb_object = TRUE)
+  x <- get_index(p, area = 1/nrow(d0))
+  pp <- predict(fit, type = "response")
+  expect_equal(x$est, pp$est[1], tolerance = 1e-3)
+
+  # TNB2
+  dd <- subset(d0, y > 0)
+  fit <- sdmTMB(
+    y ~ 1,
+    data = dd,
+    time = "year",
+    spatiotemporal = "off",
+    spatial = FALSE,
+    family = truncated_nbinom2()
+  )
+  p <- predict(fit, newdata = d0, return_tmb_object = TRUE)
+  x <- get_index(p, area = 1/nrow(d0))
+  pp <- predict(fit, type = "response")
+  expect_equal(x$est, pp$est[1], tolerance = 1e-3)
+
+  # TNB1
+  fit <- sdmTMB(
+    y ~ 1,
+    data = dd,
+    time = "year",
+    spatiotemporal = "off",
+    spatial = FALSE,
+    family = truncated_nbinom2()
+  )
+  p <- predict(fit, newdata = d0, return_tmb_object = TRUE)
+  x <- get_index(p, area = 1/nrow(d0))
+  pp <- predict(fit, type = "response")
+  expect_equal(x$est, pp$est[1], tolerance = 1e-3)
+})

From 8ff58fe522609a7812a69cc7367ba01f3406bdc4 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Tue, 18 Jun 2024 15:33:42 -0700
Subject: [PATCH 03/21] utils::globalVariables(".phi")

---
 R/families.R | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/R/families.R b/R/families.R
index 60e5b205..c6994652 100644
--- a/R/families.R
+++ b/R/families.R
@@ -222,6 +222,8 @@ nbinom1 <- function(link = "log") {
   add_to_family(x)
 }
 
+utils::globalVariables(".phi") ## avoid R CMD check NOTE
+
 #' @export
 #' @examples
 #' truncated_nbinom2(link = "log")

From 4464143e509b7eedfb1511ee0aeed56df7e1d843 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Wed, 19 Jun 2024 10:10:56 -0700
Subject: [PATCH 04/21] Fix exponentiate in tidy() #353

And make conf.int = TRUE default
---
 R/tidy.R           | 5 +++--
 man/tidy.sdmTMB.Rd | 2 +-
 2 files changed, 4 insertions(+), 3 deletions(-)

diff --git a/R/tidy.R b/R/tidy.R
index e1dff011..0d42765a 100644
--- a/R/tidy.R
+++ b/R/tidy.R
@@ -45,7 +45,7 @@
 #' tidy(fit, "ran_vals")
 
 tidy.sdmTMB <- function(x, effects = c("fixed", "ran_pars", "ran_vals"), model = 1,
-                 conf.int = FALSE, conf.level = 0.95, exponentiate = FALSE,
+                 conf.int = TRUE, conf.level = 0.95, exponentiate = FALSE,
                  silent = FALSE, ...) {
   effects <- match.arg(effects)
   assert_that(is.logical(exponentiate))
@@ -131,7 +131,6 @@ tidy.sdmTMB <- function(x, effects = c("fixed", "ran_pars", "ran_vals"), model =
   b_j_se <- se$b_j[!fe_names == "offset", drop = TRUE]
   fe_names <- fe_names[!fe_names == "offset"]
   out <- data.frame(term = fe_names, estimate = b_j, std.error = b_j_se, stringsAsFactors = FALSE)
-  if (exponentiate) out$estimate <- trans(out$estimate)
 
   if (x$tmb_data$threshold_func > 0) {
     if (x$threshold_function == 1L) {
@@ -152,6 +151,8 @@ tidy.sdmTMB <- function(x, effects = c("fixed", "ran_pars", "ran_vals"), model =
     out$conf.low <- as.numeric(trans(out$estimate - crit * out$std.error))
     out$conf.high <- as.numeric(trans(out$estimate + crit * out$std.error))
   }
+  out$estimate <- trans(out$estimate)
+  if (exponentiate) out$std.error <- NULL
 
   out_re <- list()
   log_name <- c("log_range")
diff --git a/man/tidy.sdmTMB.Rd b/man/tidy.sdmTMB.Rd
index a9ae62a6..d7d03ce3 100644
--- a/man/tidy.sdmTMB.Rd
+++ b/man/tidy.sdmTMB.Rd
@@ -8,7 +8,7 @@
   x,
   effects = c("fixed", "ran_pars", "ran_vals"),
   model = 1,
-  conf.int = FALSE,
+  conf.int = TRUE,
   conf.level = 0.95,
   exponentiate = FALSE,
   silent = FALSE,

From 51d7d455b135271f235b837a2fd3e7dadcef54e4 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Wed, 19 Jun 2024 10:12:01 -0700
Subject: [PATCH 05/21] Bump version [skip ci]

---
 DESCRIPTION | 2 +-
 NEWS.md     | 3 +++
 2 files changed, 4 insertions(+), 1 deletion(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index ec16ae5c..24979c49 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9002
+Version: 0.6.0.9003
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/NEWS.md b/NEWS.md
index 6d89efcb..501aa2a3 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,8 @@
 # sdmTMB (development version)
 
+* Fig bug in `exponentiate` argument for `tidy()`. Set `conf.int = TRUE` as
+  default. #353
+
 * Fix bug in prediction from `delta_truncated_nbinom1()` and 
   `delta_truncated_nbinom2()` families. The positive component
   needs to be transformed to represent the mean of the *un*truncated

From d2d54f67cb9757ae052e5263dc9be2efbadc0b61 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Wed, 19 Jun 2024 11:21:58 -0700
Subject: [PATCH 06/21] Fix emmeans bug

---
 R/tidy.R | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/R/tidy.R b/R/tidy.R
index 0d42765a..71fe69e0 100644
--- a/R/tidy.R
+++ b/R/tidy.R
@@ -151,7 +151,8 @@ tidy.sdmTMB <- function(x, effects = c("fixed", "ran_pars", "ran_vals"), model =
     out$conf.low <- as.numeric(trans(out$estimate - crit * out$std.error))
     out$conf.high <- as.numeric(trans(out$estimate + crit * out$std.error))
   }
-  out$estimate <- trans(out$estimate)
+  # must wrap in as.numeric() otherwise I() leaves 'AsIs' class that affects emmeans package
+  out$estimate <- as.numeric(trans(out$estimate))
   if (exponentiate) out$std.error <- NULL
 
   out_re <- list()

From f63bd84d7cf89b98eb9da1ee02748248ab5b9a43 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Wed, 3 Jul 2024 11:35:07 -0700
Subject: [PATCH 07/21] Fix plot.sdmTMB docs [352]

---
 DESCRIPTION      |  2 +-
 R/mesh.R         | 23 ++---------------------
 man/make_mesh.Rd |  6 ++----
 3 files changed, 5 insertions(+), 26 deletions(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index 24979c49..17d71bb1 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9003
+Version: 0.6.0.9004
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/R/mesh.R b/R/mesh.R
index 8b97c69d..9af2a3d0 100644
--- a/R/mesh.R
+++ b/R/mesh.R
@@ -216,30 +216,11 @@ binary_search_knots <- function(loc_xy,
 #' @param ... Passed to [graphics::plot()].
 #'
 #' @importFrom graphics points
-#' @return `plot.sdmTMBmesh()`: A plot of the mesh and data points. If
-#'   \pkg{ggplot2} is installed, a \pkg{ggplot2} object is
-#'   returned, otherwise a base graphics R plot is returned. To make your own,
-#'   pass `your_mesh$mesh` to `inlabru::gg()`.
+#' @return `plot.sdmTMBmesh()`: A plot of the mesh and data points. To make your
+#'   own \pkg{ggplot2} version, pass `your_mesh$mesh` to `inlabru::gg()`.
 #' @rdname make_mesh
 #' @export
 plot.sdmTMBmesh <- function(x, ...) {
-  # r1 <- requireNamespace("inlabru", quietly = TRUE)
-  # r2 <- requireNamespace("ggplot2", quietly = TRUE)
-  # if (r1 && r2) {
-  #   dat <- data.frame(
-  #     x = x$loc_xy[,1,drop=TRUE],
-  #     y = x$loc_xy[,2,drop=TRUE]
-  #   )
-  #   ggplot2::ggplot() +
-  #     # inlabru::gg(x$mesh, ext.color = "grey20", ext.linewidth = 0.5, edge.color = "grey50") +
-  #     ggplot2::coord_sf() +
-  #     fmesher::geom_fm(data = x$mesh) +
-  #     ggplot2::geom_point(
-  #       data = dat,
-  #       mapping = ggplot2::aes(x = .data$x, y = .data$y), alpha = 0.4, pch = 20, colour = "#3182BD") +
-  #     # ggplot2::coord_fixed() +
-  #     ggplot2::labs(x = x$xy_cols[[1]], y = x$xy_cols[[2]])
-  # } else {
     plot(x$mesh, main = NA, edge.color = "grey60", asp = 1, ...)
     points(x$loc_xy, pch = 21, cex = 0.3, col = "#00000080")
     points(x$loc_centers, pch = 20, col = "red")
diff --git a/man/make_mesh.Rd b/man/make_mesh.Rd
index 9ea51aea..d993b535 100644
--- a/man/make_mesh.Rd
+++ b/man/make_mesh.Rd
@@ -61,10 +61,8 @@ Distance to extend non-convex hull from data.}
 from \code{fmesher_func} (default is \code{\link[fmesher:fm_mesh_2d]{fmesher::fm_mesh_2d_inla()}}). See
 \code{mesh$mesh$n} for the number of vertices.
 
-\code{plot.sdmTMBmesh()}: A plot of the mesh and data points. If
-\pkg{ggplot2} is installed, a \pkg{ggplot2} object is
-returned, otherwise a base graphics R plot is returned. To make your own,
-pass \code{your_mesh$mesh} to \code{inlabru::gg()}.
+\code{plot.sdmTMBmesh()}: A plot of the mesh and data points. To make your
+own \pkg{ggplot2} version, pass \code{your_mesh$mesh} to \code{inlabru::gg()}.
 }
 \description{
 Construct an SPDE mesh for use with sdmTMB.

From 5a600beaa58346eba35e985da1da942ca2244695 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 26 Aug 2024 11:14:52 -0700
Subject: [PATCH 08/21] Move articles from web_only to articles #363

---
 .Rbuildignore                                              | 1 +
 vignettes/{web_only => articles}/basic-intro.Rmd           | 0
 vignettes/{web_only => articles}/bayesian.Rmd              | 0
 vignettes/{web_only => articles}/cross-validation.Rmd      | 0
 vignettes/{web_only => articles}/delta-models.Rmd          | 0
 vignettes/{web_only => articles}/ggeffects.Rmd             | 0
 vignettes/{web_only => articles}/index-standardization.Rmd | 0
 vignettes/{web_only => articles}/poisson-link.Rmd          | 0
 vignettes/{web_only => articles}/pretty-plots.Rmd          | 0
 vignettes/{web_only => articles}/residual-checking.Rmd     | 0
 vignettes/{web_only => articles}/spatial-trend-models.Rmd  | 0
 vignettes/{web_only => articles}/threshold-models.Rmd      | 0
 vignettes/{web_only => articles}/visreg.Rmd                | 0
 13 files changed, 1 insertion(+)
 rename vignettes/{web_only => articles}/basic-intro.Rmd (100%)
 rename vignettes/{web_only => articles}/bayesian.Rmd (100%)
 rename vignettes/{web_only => articles}/cross-validation.Rmd (100%)
 rename vignettes/{web_only => articles}/delta-models.Rmd (100%)
 rename vignettes/{web_only => articles}/ggeffects.Rmd (100%)
 rename vignettes/{web_only => articles}/index-standardization.Rmd (100%)
 rename vignettes/{web_only => articles}/poisson-link.Rmd (100%)
 rename vignettes/{web_only => articles}/pretty-plots.Rmd (100%)
 rename vignettes/{web_only => articles}/residual-checking.Rmd (100%)
 rename vignettes/{web_only => articles}/spatial-trend-models.Rmd (100%)
 rename vignettes/{web_only => articles}/threshold-models.Rmd (100%)
 rename vignettes/{web_only => articles}/visreg.Rmd (100%)

diff --git a/.Rbuildignore b/.Rbuildignore
index 41f8b9ec..cac90758 100644
--- a/.Rbuildignore
+++ b/.Rbuildignore
@@ -21,3 +21,4 @@
 ^CRAN-SUBMISSION$
 ^vignettes/web_only$
 revdep/
+^vignettes/articles$
diff --git a/vignettes/web_only/basic-intro.Rmd b/vignettes/articles/basic-intro.Rmd
similarity index 100%
rename from vignettes/web_only/basic-intro.Rmd
rename to vignettes/articles/basic-intro.Rmd
diff --git a/vignettes/web_only/bayesian.Rmd b/vignettes/articles/bayesian.Rmd
similarity index 100%
rename from vignettes/web_only/bayesian.Rmd
rename to vignettes/articles/bayesian.Rmd
diff --git a/vignettes/web_only/cross-validation.Rmd b/vignettes/articles/cross-validation.Rmd
similarity index 100%
rename from vignettes/web_only/cross-validation.Rmd
rename to vignettes/articles/cross-validation.Rmd
diff --git a/vignettes/web_only/delta-models.Rmd b/vignettes/articles/delta-models.Rmd
similarity index 100%
rename from vignettes/web_only/delta-models.Rmd
rename to vignettes/articles/delta-models.Rmd
diff --git a/vignettes/web_only/ggeffects.Rmd b/vignettes/articles/ggeffects.Rmd
similarity index 100%
rename from vignettes/web_only/ggeffects.Rmd
rename to vignettes/articles/ggeffects.Rmd
diff --git a/vignettes/web_only/index-standardization.Rmd b/vignettes/articles/index-standardization.Rmd
similarity index 100%
rename from vignettes/web_only/index-standardization.Rmd
rename to vignettes/articles/index-standardization.Rmd
diff --git a/vignettes/web_only/poisson-link.Rmd b/vignettes/articles/poisson-link.Rmd
similarity index 100%
rename from vignettes/web_only/poisson-link.Rmd
rename to vignettes/articles/poisson-link.Rmd
diff --git a/vignettes/web_only/pretty-plots.Rmd b/vignettes/articles/pretty-plots.Rmd
similarity index 100%
rename from vignettes/web_only/pretty-plots.Rmd
rename to vignettes/articles/pretty-plots.Rmd
diff --git a/vignettes/web_only/residual-checking.Rmd b/vignettes/articles/residual-checking.Rmd
similarity index 100%
rename from vignettes/web_only/residual-checking.Rmd
rename to vignettes/articles/residual-checking.Rmd
diff --git a/vignettes/web_only/spatial-trend-models.Rmd b/vignettes/articles/spatial-trend-models.Rmd
similarity index 100%
rename from vignettes/web_only/spatial-trend-models.Rmd
rename to vignettes/articles/spatial-trend-models.Rmd
diff --git a/vignettes/web_only/threshold-models.Rmd b/vignettes/articles/threshold-models.Rmd
similarity index 100%
rename from vignettes/web_only/threshold-models.Rmd
rename to vignettes/articles/threshold-models.Rmd
diff --git a/vignettes/web_only/visreg.Rmd b/vignettes/articles/visreg.Rmd
similarity index 100%
rename from vignettes/web_only/visreg.Rmd
rename to vignettes/articles/visreg.Rmd

From 075440b2d31bf82aadac67770d428f43e4821dea Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 26 Aug 2024 12:18:34 -0700
Subject: [PATCH 09/21] Remove 'web_only' article links

---
 DESCRIPTION                         | 2 +-
 NEWS.md                             | 2 +-
 R/dharma.R                          | 2 +-
 R/predict.R                         | 2 +-
 README.Rmd                          | 2 +-
 README.md                           | 2 +-
 man/dharma_residuals.Rd             | 2 +-
 man/predict.sdmTMB.Rd               | 2 +-
 vignettes/articles/delta-models.Rmd | 2 +-
 9 files changed, 9 insertions(+), 9 deletions(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index 17d71bb1..7cc976e5 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -109,5 +109,5 @@ Config/testthat/parallel: true
 Encoding: UTF-8
 LazyData: true
 Roxygen: list(markdown = TRUE)
-RoxygenNote: 7.3.1
+RoxygenNote: 7.3.2
 SystemRequirements: GNU make
diff --git a/NEWS.md b/NEWS.md
index 501aa2a3..b3db0f70 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -61,7 +61,7 @@
 # sdmTMB 0.5.0
 
 * Overhaul residuals vignette ('article') 
-  <https://pbs-assess.github.io/sdmTMB/articles/web_only/residual-checking.html>
+  <https://pbs-assess.github.io/sdmTMB/articles/residual-checking.html>
   including brief intros to randomized quantile residuals, simulation-based
   residuals, 'one-sample' residuals, and uniform vs. Gaussian residuals.
 
diff --git a/R/dharma.R b/R/dharma.R
index a3f29323..ca66b210 100644
--- a/R/dharma.R
+++ b/R/dharma.R
@@ -27,7 +27,7 @@
 #'
 #' @details
 #'
-#' See the [residuals vignette](https://pbs-assess.github.io/sdmTMB/articles/web_only/residual-checking.html).
+#' See the [residuals vignette](https://pbs-assess.github.io/sdmTMB/articles/residual-checking.html).
 #'
 #' Advantages to these residuals over the ones from the [residuals.sdmTMB()]
 #' method are (1) they work with delta/hurdle models for the combined
diff --git a/R/predict.R b/R/predict.R
index c28b96a5..b8172471 100644
--- a/R/predict.R
+++ b/R/predict.R
@@ -45,7 +45,7 @@
 #' @param mcmc_samples See `extract_mcmc()` in the
 #'   \href{https://github.com/pbs-assess/sdmTMBextra}{sdmTMBextra} package for
 #'   more details and the
-#'   \href{https://pbs-assess.github.io/sdmTMB/articles/web_only/bayesian.html}{Bayesian vignette}.
+#'   \href{https://pbs-assess.github.io/sdmTMB/articles/bayesian.html}{Bayesian vignette}.
 #'   If specified, the predict function will return a matrix of a similar form
 #'   as if `nsim > 0` but representing Bayesian posterior samples from the Stan
 #'   model.
diff --git a/README.Rmd b/README.Rmd
index ca07128c..95180c7c 100644
--- a/README.Rmd
+++ b/README.Rmd
@@ -509,7 +509,7 @@ See [`?sdmTMBpriors`](https://pbs-assess.github.io/sdmTMB/reference/priors.html)
 
 ### Bayesian MCMC sampling with Stan
 
-The fitted model can be passed to the tmbstan package to sample from the posterior with Stan. See the [Bayesian vignette](https://pbs-assess.github.io/sdmTMB/articles/web_only/bayesian.html).
+The fitted model can be passed to the tmbstan package to sample from the posterior with Stan. See the [Bayesian vignette](https://pbs-assess.github.io/sdmTMB/articles/bayesian.html).
 
 ### Turning off random fields
 
diff --git a/README.md b/README.md
index 341bcc1d..b55a7061 100644
--- a/README.md
+++ b/README.md
@@ -729,7 +729,7 @@ for more details.
 
 The fitted model can be passed to the tmbstan package to sample from the
 posterior with Stan. See the [Bayesian
-vignette](https://pbs-assess.github.io/sdmTMB/articles/web_only/bayesian.html).
+vignette](https://pbs-assess.github.io/sdmTMB/articles/bayesian.html).
 
 ### Turning off random fields
 
diff --git a/man/dharma_residuals.Rd b/man/dharma_residuals.Rd
index d531b460..91523c80 100644
--- a/man/dharma_residuals.Rd
+++ b/man/dharma_residuals.Rd
@@ -50,7 +50,7 @@ around \code{\link[DHARMa:createDHARMa]{DHARMa::createDHARMa()}} to facilitate i
 expected distribution. This is \emph{not} the default.
 }
 \details{
-See the \href{https://pbs-assess.github.io/sdmTMB/articles/web_only/residual-checking.html}{residuals vignette}.
+See the \href{https://pbs-assess.github.io/sdmTMB/articles/residual-checking.html}{residuals vignette}.
 
 Advantages to these residuals over the ones from the \code{\link[=residuals.sdmTMB]{residuals.sdmTMB()}}
 method are (1) they work with delta/hurdle models for the combined
diff --git a/man/predict.sdmTMB.Rd b/man/predict.sdmTMB.Rd
index c33cb8dd..1745cd03 100644
--- a/man/predict.sdmTMB.Rd
+++ b/man/predict.sdmTMB.Rd
@@ -78,7 +78,7 @@ both sets of predictions are returned.}
 \item{mcmc_samples}{See \code{extract_mcmc()} in the
 \href{https://github.com/pbs-assess/sdmTMBextra}{sdmTMBextra} package for
 more details and the
-\href{https://pbs-assess.github.io/sdmTMB/articles/web_only/bayesian.html}{Bayesian vignette}.
+\href{https://pbs-assess.github.io/sdmTMB/articles/bayesian.html}{Bayesian vignette}.
 If specified, the predict function will return a matrix of a similar form
 as if \code{nsim > 0} but representing Bayesian posterior samples from the Stan
 model.}
diff --git a/vignettes/articles/delta-models.Rmd b/vignettes/articles/delta-models.Rmd
index bfb4006a..f15b21dc 100644
--- a/vignettes/articles/delta-models.Rmd
+++ b/vignettes/articles/delta-models.Rmd
@@ -47,7 +47,7 @@ Here and with other delta models, the `link1` and `link2` can be omitted and lef
 2. Delta-lognormal: `family = delta_lognormal()`.
 This fits a binomial presence-absence model (i.e., `binomial(link = "logit")`) and then a model for the positive catches only with a lognormal observation distribution and a log link (i.e., `lognormal(link = "log")`
 
-3. Poisson-link delta-Gamma or delta-lognormal. See the [Poisson-link delta model vignette](https://pbs-assess.github.io/sdmTMB/articles/web_only/poisson-link.html).
+3. Poisson-link delta-Gamma or delta-lognormal. See the [Poisson-link delta model vignette](https://pbs-assess.github.io/sdmTMB/articles/poisson-link.html).
 
 4. Delta-truncated-negative-binomial: `family = delta_truncated_nbinom1()` or `family = delta_truncated_nbinom2()`.
 This fits a binomial presence-absence model (`binomial(link = "logit")`) and a `truncated_nbinom1(link = "log")` or `truncated_nbinom1(link = "log")` distribution for positive catches.

From 1418f42a2b905f1fcd9febf235103133161a5b2b Mon Sep 17 00:00:00 2001
From: Eric Ward <5046884+ericward-noaa@users.noreply.github.com>
Date: Fri, 13 Sep 2024 07:58:26 -0700
Subject: [PATCH 10/21] Add warning messages when components are collapsing to
 0

---
 R/tmb-sim.R | 15 +++++++++++++++
 1 file changed, 15 insertions(+)

diff --git a/R/tmb-sim.R b/R/tmb-sim.R
index 08cc8e7a..a95d0a07 100644
--- a/R/tmb-sim.R
+++ b/R/tmb-sim.R
@@ -271,6 +271,21 @@ sdmTMB_simulate <- function(formula,
   d[["omega_s"]] <- if (all(s$omega_s_A != 0)) s$omega_s_A
   d[["epsilon_st"]] <- if (all(s$epsilon_st_A_vec != 0)) s$epsilon_st_A_vec
   d[["zeta_s"]] <- if (all(s$zeta_s_A != 0)) s$zeta_s_A
+  
+    # Warnings for pieces collapsing to 0. The sum() > 0 piece is needed because parameters that are input are turned to 0x0 matrices
+  if (sum(sigma_O) > 0) {
+    if(is.null(d[["omega_s"]])) cli::cli_alert_info(paste("Warning: spatial field is collapsing to 0 but sigma_O was specified.",
+                                                        "Please check the spatial range and cutoff distance used to construct the mesh."))
+  }
+  if (sum(sigma_E) > 0) {
+    if(is.null(d[["epsilon_st"]])) cli::cli_alert_info(paste("Warning: spatiotemporal field is collapsing to 0 but sigma_E was specified.",
+                                                        "Please check the spatial range and cutoff distance used to construct the mesh."))
+  }
+  if (sum(sigma_Z) > 0) {
+    if(is.null(d[["zeta_s"]])) cli::cli_alert_info(paste("Warning: spatially varying coefficient field is collapsing to 0 but sigma_Z was specified.",
+                                                               "Please check the spatial range and cutoff distance used to construct the mesh."))
+  }
+  
   if (any(family$family %in% c("truncated_nbinom1", "truncated_nbinom2"))) {
     d[["mu"]] <- family$linkinv(s$eta_i, phi = phi)
   } else {

From cb83a62207a2d62a4e1f7a8b2fa02d1535a2e703 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Fri, 13 Sep 2024 13:46:51 -0600
Subject: [PATCH 11/21] Fix #370

---
 R/tmb-sim.R                          | 37 ++++++++++++++--------------
 tests/testthat/test-tmb-simulation.R | 25 +++++++++++++++++++
 2 files changed, 44 insertions(+), 18 deletions(-)

diff --git a/R/tmb-sim.R b/R/tmb-sim.R
index a95d0a07..724e884d 100644
--- a/R/tmb-sim.R
+++ b/R/tmb-sim.R
@@ -268,24 +268,25 @@ sdmTMB_simulate <- function(formula,
   if (!is.null(fit$time)) d[[fit$time]] <- data[[fit$time]]
   d[[mesh$xy_cols[1]]] <- data[[mesh$xy_cols[1]]]
   d[[mesh$xy_cols[2]]] <- data[[mesh$xy_cols[2]]]
-  d[["omega_s"]] <- if (all(s$omega_s_A != 0)) s$omega_s_A
-  d[["epsilon_st"]] <- if (all(s$epsilon_st_A_vec != 0)) s$epsilon_st_A_vec
-  d[["zeta_s"]] <- if (all(s$zeta_s_A != 0)) s$zeta_s_A
-  
-    # Warnings for pieces collapsing to 0. The sum() > 0 piece is needed because parameters that are input are turned to 0x0 matrices
-  if (sum(sigma_O) > 0) {
-    if(is.null(d[["omega_s"]])) cli::cli_alert_info(paste("Warning: spatial field is collapsing to 0 but sigma_O was specified.",
-                                                        "Please check the spatial range and cutoff distance used to construct the mesh."))
-  }
-  if (sum(sigma_E) > 0) {
-    if(is.null(d[["epsilon_st"]])) cli::cli_alert_info(paste("Warning: spatiotemporal field is collapsing to 0 but sigma_E was specified.",
-                                                        "Please check the spatial range and cutoff distance used to construct the mesh."))
-  }
-  if (sum(sigma_Z) > 0) {
-    if(is.null(d[["zeta_s"]])) cli::cli_alert_info(paste("Warning: spatially varying coefficient field is collapsing to 0 but sigma_Z was specified.",
-                                                               "Please check the spatial range and cutoff distance used to construct the mesh."))
-  }
-  
+
+  d[["omega_s"]] <- if (sum(sigma_O) > 0) s$omega_s_A
+  d[["epsilon_st"]] <- if (sum(sigma_E) > 0) s$epsilon_st_A_vec
+  d[["zeta_s"]] <- if (sum(sigma_Z) > 0) s$zeta_s_A
+
+  # # Warnings for fields collapsing to 0
+  # info_collapse <- function(sig, vec, .par, .name) {
+  #   if (sum(sig) > 0 && all (vec == 0)) {
+  #     msg <- paste0("The ", .name, " has been returned as all zeros although ", .par,
+  #       " was specified as > 0. Try making your mesh finer, e.g., with a lower ",
+  #       "`cutoff` or a higher number of knots. Triangle edge length needs to be ",
+  #       "lower than the range size (distance correlation is effectively independent.")
+  #     cli::cli_alert_info(msg)
+  #   }
+  # }
+  # info_collapse(sigma_O, s$omega_s_A, "sigma_O", "spatial field")
+  # info_collapse(sigma_E, s$epsilon_st_A_vec, "sigma_E", "spatiotemporal field")
+  # info_collapse(sigma_Z, s$zeta_s_A, "sigma_Z", "spatially varying coefficient field")
+
   if (any(family$family %in% c("truncated_nbinom1", "truncated_nbinom2"))) {
     d[["mu"]] <- family$linkinv(s$eta_i, phi = phi)
   } else {
diff --git a/tests/testthat/test-tmb-simulation.R b/tests/testthat/test-tmb-simulation.R
index fe8851f5..d6c41874 100644
--- a/tests/testthat/test-tmb-simulation.R
+++ b/tests/testthat/test-tmb-simulation.R
@@ -183,3 +183,28 @@ test_that("simulate.sdmTMB returns the right length", {
   s <- simulate(m, nsim = 2)
   expect_equal(nrow(s), nrow(pcod))
 })
+
+test_that("coarse meshes with zeros in simulation still return fields #370", {
+  set.seed(123)
+  predictor_dat <- data.frame(
+    X = runif(100), Y = runif(100),
+    a1 = rnorm(100), year = rep(1:2, each = 50))
+  mesh <- sdmTMB::make_mesh(predictor_dat, xy_cols = c("X", "Y"), n_knots = 30)
+  sim_dat <- sdmTMB::sdmTMB_simulate(
+    formula = ~ 1 + a1,
+    data = predictor_dat,
+    time = "year",
+    mesh = mesh,
+    family = gaussian(),
+    range = 0.5,
+    sigma_E = 0.1,
+    phi = 0.1,
+    sigma_O = 0.2,
+    seed = 42,
+    B = c(0.2, -0.4) # B0 = intercept, B1 = a1 slope
+  )
+  nm <- names(sim_dat)
+  expect_true("omega_s" %in% nm)
+  expect_true("epsilon_st" %in% nm)
+  expect_false("zeta_s" %in% nm)
+})

From 4b785e06e845e48f6390bf87e64414c061323e31 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Fri, 20 Sep 2024 13:34:34 -0700
Subject: [PATCH 12/21] Fix offset prediction in cross validation #372

---
 DESCRIPTION                  |  2 +-
 NEWS.md                      |  3 +++
 R/cross-val.R                | 10 ++++++----
 tests/testthat/test-offset.R | 22 ++++++++++++++++++++++
 4 files changed, 32 insertions(+), 5 deletions(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index 7cc976e5..c69e2f1d 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9004
+Version: 0.6.0.9005
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/NEWS.md b/NEWS.md
index b3db0f70..18375ff4 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,8 @@
 # sdmTMB (development version)
 
+* Fix passing of `offset` argument through in `sdmTMB_cv()`. Before it was being
+  omitted in the prediction (i.e., set to 0). #372
+
 * Fig bug in `exponentiate` argument for `tidy()`. Set `conf.int = TRUE` as
   default. #353
 
diff --git a/R/cross-val.R b/R/cross-val.R
index 06abe9af..6e871f13 100644
--- a/R/cross-val.R
+++ b/R/cross-val.R
@@ -287,10 +287,10 @@ sdmTMB_cv <- function(
     cli_abort("`weights` cannot be specified within sdmTMB_cv().")
   }
   if ("offset" %in% names(dot_args)) {
-    .offset <- eval(dot_args$offset)
-    if (parallel && !is.character(.offset) && !is.null(.offset)) {
-      cli_abort("We recommend using a character value for 'offset' (indicating the column name) when applying parallel cross validation.")
+    if (!is.character(dot_args$offset)) {
+      cli_abort("Please use a character value for 'offset' (indicating the column name) for cross validation.")
     }
+    .offset <- eval(dot_args$offset)
   } else {
     .offset <- NULL
   }
@@ -369,7 +369,9 @@ sdmTMB_cv <- function(
 
     # FIXME: only use TMB report() below to be faster!
     # predict for withheld data:
-    predicted <- predict(object, newdata = cv_data, type = "response")
+    predicted <- predict(object, newdata = cv_data, type = "response",
+      offset = if (!is.null(.offset)) cv_data[[.offset]] else rep(0, nrow(cv_data)))
+
     cv_data$cv_predicted <- predicted$est
     response <- get_response(object$formula[[1]])
     withheld_y <- predicted[[response]]
diff --git a/tests/testthat/test-offset.R b/tests/testthat/test-offset.R
index aa6a636f..d19420d7 100644
--- a/tests/testthat/test-offset.R
+++ b/tests/testthat/test-offset.R
@@ -138,6 +138,28 @@ test_that("Offset prediction matches glm()", {
   # p_glmmTMB <- predict(fit_glmmTMB, newdata = dat)
   # expect_equal(p$est, unname(p_glmmTMB))
 })
+
+test_that("offset gets passed through cross validation as expected #372", {
+  dat <- subset(dogfish, catch_weight > 0)
+  expect_error(
+    x <- sdmTMB_cv(catch_weight ~ 1,
+      data = dat,
+      family = Gamma("log"), offset = log(dat$area_swept), spatial = "off"
+    ),
+    "offset"
+  )
+  set.seed(1)
+  x <- sdmTMB_cv(catch_weight ~ 1,
+    data = dat, family = Gamma("log"),
+    offset = "area_swept", spatial = "off",
+    mesh = make_mesh(dat, c("X", "Y"), cutoff = 10), k_folds = 2
+  )
+  y <- x$data[, c("catch_weight", "cv_predicted")]
+  # plot(y$catch_weight, y$cv_predicted)
+  # if offset is applied, will have unique values because an intercept-only model:
+  expect_true(length(unique(y$cv_predicted)) == 684L)
+})
+
 # #
 # # offset/prediction setting checks:
 #

From 1be189e3a81006851ec012799223399242bbe8bf Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 23 Sep 2024 13:14:22 -0700
Subject: [PATCH 13/21] Add message on offset newdata prediction #372

---
 DESCRIPTION                  |  2 +-
 NEWS.md                      |  3 +++
 R/predict.R                  |  9 +++++++++
 tests/testthat/test-offset.R | 16 ++++++++++++++++
 4 files changed, 29 insertions(+), 1 deletion(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index c69e2f1d..a3e767fd 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9005
+Version: 0.6.0.9006
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/NEWS.md b/NEWS.md
index 18375ff4..411966ac 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,8 @@
 # sdmTMB (development version)
 
+* Add informative message if fitting with an offset but predicting with offset
+  argument left at NULL on newdata. #372
+
 * Fix passing of `offset` argument through in `sdmTMB_cv()`. Before it was being
   omitted in the prediction (i.e., set to 0). #372
 
diff --git a/R/predict.R b/R/predict.R
index b8172471..1981545b 100644
--- a/R/predict.R
+++ b/R/predict.R
@@ -495,6 +495,15 @@ predict.sdmTMB <- function(object, newdata = NULL,
       cli_abort("`area` should be of the same length as `nrow(newdata)` or of length 1.")
     }
 
+    # newdata, null offset in predict, and non-null in fit #372
+    if (isFALSE(nd_arg_was_null) && is.null(offset) && !all(object$offset == 0)) {
+      msg <- c(
+        "Fitted object contains an offset but the offset is `NULL` in `predict.sdmTMB()`.",
+        "Prediction will proceed assuming the offset vector is 0 in the prediction.",
+        "Specify an offset vector in `predict.sdmTMB()` to override this.")
+      cli_inform(msg)
+    }
+
     if (!is.null(offset)) {
       if (nrow(proj_X_ij[[1]]) != length(offset))
         cli_abort("Prediction offset vector does not equal number of rows in prediction dataset.")
diff --git a/tests/testthat/test-offset.R b/tests/testthat/test-offset.R
index d19420d7..b8448429 100644
--- a/tests/testthat/test-offset.R
+++ b/tests/testthat/test-offset.R
@@ -160,6 +160,22 @@ test_that("offset gets passed through cross validation as expected #372", {
   expect_true(length(unique(y$cv_predicted)) == 684L)
 })
 
+test_that("predicting on newdata with a non-null offset in fit but a null offset in predict informs the user appropriately", {
+  dat <- subset(dogfish, catch_weight > 0)
+  fit <- sdmTMB(
+    catch_weight ~ 1,
+    data = dat,
+    family = Gamma("log"),
+    offset = "area_swept",
+    spatial = "off"
+  )
+  pred <- predict(fit)
+  pred <- predict(fit, offset = rep(0, nrow(dat)))
+  pred <- predict(fit, newdata = qcs_grid, offset = rep(0, nrow(qcs_grid)))
+  pred <- predict(fit, newdata = qcs_grid)
+  expect_message({pred <- predict(fit, newdata = qcs_grid)}, regexp = "offset")
+})
+
 # #
 # # offset/prediction setting checks:
 #

From 870ad8ce9c0c794d38b6d42b66b9edce3338064f Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 23 Sep 2024 13:27:11 -0700
Subject: [PATCH 14/21] Add error msg if coord NAs in make_mesh #365

---
 DESCRIPTION                |  2 +-
 NEWS.md                    |  2 ++
 R/mesh.R                   |  8 ++++++++
 tests/testthat/test-mesh.R | 13 +++++++++++++
 4 files changed, 24 insertions(+), 1 deletion(-)
 create mode 100644 tests/testthat/test-mesh.R

diff --git a/DESCRIPTION b/DESCRIPTION
index a3e767fd..667d1951 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9006
+Version: 0.6.0.9007
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/NEWS.md b/NEWS.md
index 411966ac..27124a55 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,7 @@
 # sdmTMB (development version)
 
+* Add helpful error message if some coordinates in make_mesh() are NA. #365
+
 * Add informative message if fitting with an offset but predicting with offset
   argument left at NULL on newdata. #372
 
diff --git a/R/mesh.R b/R/mesh.R
index 9af2a3d0..0a071ed1 100644
--- a/R/mesh.R
+++ b/R/mesh.R
@@ -88,6 +88,14 @@ make_mesh <- function(data, xy_cols,
     cli_abort(msg)
   }
 
+  all_x_non_na <- sum(is.na(data[[xy_cols[[1]]]])) == 0L
+  all_y_non_na <- sum(is.na(data[[xy_cols[[2]]]])) == 0L
+  if (!all_x_non_na || !all_y_non_na) {
+    msg <- c("Some coordinates in `xy_cols` were NA.", "
+      Remove or fix these rows before proceeding.")
+    cli_abort(msg)
+  }
+
   if (max(data[[xy_cols[1]]]) > 1e4 || max(data[[xy_cols[2]]] > 1e4)) {
     msg <- paste0(
       "The x or y column values are fairly large. ",
diff --git a/tests/testthat/test-mesh.R b/tests/testthat/test-mesh.R
new file mode 100644
index 00000000..b3d069ad
--- /dev/null
+++ b/tests/testthat/test-mesh.R
@@ -0,0 +1,13 @@
+test_that("make_mesh returns informative error if some coordinates are NA", {
+  d <- data.frame(x = c(NA, runif(10)), y = c(NA, runif(10)))
+  expect_error({make_mesh(d, xy_cols = c("x", "y"), cutoff = 1)}, regexp = "NA")
+
+  d <- data.frame(x = c(1, runif(10)), y = c(NA, runif(10)))
+  expect_error({make_mesh(d, xy_cols = c("x", "y"), cutoff = 1)}, regexp = "NA")
+
+  d <- data.frame(x = c(NA, runif(10)), y = c(1, runif(10)))
+  expect_error({make_mesh(d, xy_cols = c("x", "y"), cutoff = 1)}, regexp = "NA")
+
+  d <- data.frame(x = c(1, runif(10)), y = c(1, runif(10)))
+  mesh <- make_mesh(d, xy_cols = c("x", "y"), cutoff = 1)
+})

From e868be28ac98a1592ae33de0443a1ef527aeeaba Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 23 Sep 2024 13:52:28 -0700
Subject: [PATCH 15/21] Allow coef to work for delta models #351

---
 DESCRIPTION                   |  2 +-
 NEWS.md                       |  3 +++
 R/methods.R                   | 14 +++++++++++---
 man/coef.sdmTMB.Rd            | 21 +++++++++++++++++++++
 tests/testthat/test-methods.R | 12 ++++++++++++
 5 files changed, 48 insertions(+), 4 deletions(-)
 create mode 100644 man/coef.sdmTMB.Rd

diff --git a/DESCRIPTION b/DESCRIPTION
index 667d1951..30269120 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9007
+Version: 0.6.0.9008
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/NEWS.md b/NEWS.md
index 27124a55..3256deda 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,8 @@
 # sdmTMB (development version)
 
+* Add `model` (linear predictor number) argument to coef() method. Also,
+  write documentation for `?coef.sdmTMB`. #351
+
 * Add helpful error message if some coordinates in make_mesh() are NA. #365
 
 * Add informative message if fitting with an offset but predicting with offset
diff --git a/R/methods.R b/R/methods.R
index 58ec2163..314577ab 100644
--- a/R/methods.R
+++ b/R/methods.R
@@ -45,12 +45,20 @@ fitted.sdmTMB <- function(object, ...) {
 #'
 #' @param object The fitted sdmTMB model object
 #' @param complete Currently ignored
+#' @param model Linear predictor for delta models. Defaults to the first
+#'   linear predictor.
 #' @param ... Currently ignored
 #' @importFrom stats coef
 #' @export
-#' @noRd
-coef.sdmTMB <- function(object, complete = FALSE, ...) {
-  x <- tidy(object)
+coef.sdmTMB <- function(object, complete = FALSE, model = 1, ...) {
+  if (is_delta(object)) {
+    assert_that(length(model) == 1L)
+    model <- as.integer(model)
+    assert_that(model %in% c(1L, 2L))
+    msg <- paste0("Returning coefficients from linear predictor ", model, " based on the `model` argument.")
+    cli_inform(msg)
+  }
+  x <- tidy(object, model = model)
   out <- x$estimate
   names(out) <- x$term
   out
diff --git a/man/coef.sdmTMB.Rd b/man/coef.sdmTMB.Rd
new file mode 100644
index 00000000..d0060e66
--- /dev/null
+++ b/man/coef.sdmTMB.Rd
@@ -0,0 +1,21 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/methods.R
+\name{coef.sdmTMB}
+\alias{coef.sdmTMB}
+\title{Get fixed-effect coefficients}
+\usage{
+\method{coef}{sdmTMB}(object, complete = FALSE, model = 1, ...)
+}
+\arguments{
+\item{object}{The fitted sdmTMB model object}
+
+\item{complete}{Currently ignored}
+
+\item{model}{Linear predictor for delta models. Defaults to the first
+linear predictor.}
+
+\item{...}{Currently ignored}
+}
+\description{
+Get fixed-effect coefficients
+}
diff --git a/tests/testthat/test-methods.R b/tests/testthat/test-methods.R
index 6b31ba3b..99a86ffe 100644
--- a/tests/testthat/test-methods.R
+++ b/tests/testthat/test-methods.R
@@ -27,6 +27,18 @@ test_that("coef and vcov and confint work", {
   expect_true(grepl("Estimate", colnames(x))[3])
 })
 
+test_that("coef works with delta models and informs as needed", {
+  skip_on_cran()
+  fit <- sdmTMB(
+    density ~ depth,
+    data = pcod_2011, spatial = "off",
+    family = delta_gamma()
+  )
+  expect_message(x <- coef(fit), regexp = "model")
+  expect_message(x <- coef(fit, model = 1), regexp = "model")
+  expect_message(x <- coef(fit, model = 2), regexp = "model")
+})
+
 test_that("various methods work", {
   skip_on_cran()
   fit <- sdmTMB(

From 242079d30d2e0761d326503b92e99dd2682b1699 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 23 Sep 2024 13:52:41 -0700
Subject: [PATCH 16/21] Skip 2 new tests on CRAN for speed

---
 tests/testthat/test-offset.R | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/tests/testthat/test-offset.R b/tests/testthat/test-offset.R
index b8448429..bb3e9120 100644
--- a/tests/testthat/test-offset.R
+++ b/tests/testthat/test-offset.R
@@ -140,6 +140,7 @@ test_that("Offset prediction matches glm()", {
 })
 
 test_that("offset gets passed through cross validation as expected #372", {
+  skip_on_cran()
   dat <- subset(dogfish, catch_weight > 0)
   expect_error(
     x <- sdmTMB_cv(catch_weight ~ 1,
@@ -161,6 +162,7 @@ test_that("offset gets passed through cross validation as expected #372", {
 })
 
 test_that("predicting on newdata with a non-null offset in fit but a null offset in predict informs the user appropriately", {
+  skip_on_cran()
   dat <- subset(dogfish, catch_weight > 0)
   fit <- sdmTMB(
     catch_weight ~ 1,

From f175db8218b242454a4a43dd33e474a850c856d7 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 23 Sep 2024 14:04:54 -0700
Subject: [PATCH 17/21] Add coef.sdmTMB in pkgdown index

---
 _pkgdown.yml | 9 +++++----
 1 file changed, 5 insertions(+), 4 deletions(-)

diff --git a/_pkgdown.yml b/_pkgdown.yml
index beae68bc..6adb1f35 100644
--- a/_pkgdown.yml
+++ b/_pkgdown.yml
@@ -32,15 +32,16 @@ reference:
     Core tools for model fitting, prediction, and inspection.
   contents:
   - sdmTMB
-  - predict.sdmTMB
-  - tidy.sdmTMB
   - sanity
-  - sdmTMBcontrol
-  - run_extra_optimization
+  - tidy.sdmTMB
+  - predict.sdmTMB
   - residuals.sdmTMB
   - dharma_residuals
+  - sdmTMBcontrol
+  - run_extra_optimization
   - replicate_df
   - set_delta_model
+  - coef.sdmTMB
 
 - title: 'Families'
   desc: |

From 11506d92f9185ab48289919db3fd1e24b92f9529 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 23 Sep 2024 14:47:51 -0700
Subject: [PATCH 18/21] Add AUC and TSS to crossval vignette #268

---
 .github/workflows/pkgdown.yaml          |  1 +
 NEWS.md                                 |  2 +
 vignettes/articles/cross-validation.Rmd | 78 +++++++++++++++++++++++++
 3 files changed, 81 insertions(+)

diff --git a/.github/workflows/pkgdown.yaml b/.github/workflows/pkgdown.yaml
index 3d1d9f05..a4721902 100644
--- a/.github/workflows/pkgdown.yaml
+++ b/.github/workflows/pkgdown.yaml
@@ -40,6 +40,7 @@ jobs:
           remotes::install_cran("cowplot")
           remotes::install_cran("rnaturalearth")
           remotes::install_cran("rnaturalearthdata")
+          remotes::install_cran("pROC")
           install.packages("Matrix", type = "source")
           install.packages("TMB", type = "source")
         shell: Rscript {0}
diff --git a/NEWS.md b/NEWS.md
index 3256deda..bca9cc51 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,7 @@
 # sdmTMB (development version)
 
+* Add AUC and TSS examples to cross validation vignette. #268
+
 * Add `model` (linear predictor number) argument to coef() method. Also,
   write documentation for `?coef.sdmTMB`. #351
 
diff --git a/vignettes/articles/cross-validation.Rmd b/vignettes/articles/cross-validation.Rmd
index 2eb4669e..c453bd27 100644
--- a/vignettes/articles/cross-validation.Rmd
+++ b/vignettes/articles/cross-validation.Rmd
@@ -221,4 +221,82 @@ If we had just wanted to use the predictions from the first fold onto the 10% te
 weights <- sdmTMB_stacking(model_list, include_folds = 1)
 ```
 
+# Calculating measures of predictive skill for binary data
+
+For delta models, or models of presence-absence data, several measures of predictive ability are available.
+These are applicable to cross validation, although we demonstrate them here first in a non-cross validation context for simplicity.
+
+A first commonly used diagnostic is the AUC (Area Under the Curve), which quantifies the ability of a model to discriminate between the two classes; this is done from the Receiver Operating Characteristic (ROC) curve, which plots the true positive rate vs. false positive rate.
+There are several packages to calculate AUC in R, but this can be done with the `pROC` package, where inputs are a vector of 0s and 1s (or factor equivalents) in the raw data, and a vector of estimated probabilities (generated from a call to `predict()`, as shown below).
+The `plogis()` function is needed to convert estimated values in logit space to probabilities in natural (zero to one) space. 
+
+```{r roc}
+mesh <- make_mesh(pcod, c("X", "Y"), cutoff = 10)
+fit <- sdmTMB(present ~ s(depth), data = pcod, mesh = mesh)
+pred <- predict(fit) # presence-absence model
+roc <- pROC::roc(pcod$present, plogis(pred$est))
+auc <- pROC::auc(roc)
+auc
+```
+
+With a delta model, two estimated values are returned, so only the first would be used. E.g.,
+
+```{r}
+fit <- sdmTMB(density ~ 1, data = pcod, 
+  mesh = mesh, family = delta_gamma())
+pred <- predict(fit)
+
+# the first linear predictor is the binomial component (est1):
+roc <- pROC::roc(pcod$present, plogis(pred$est1))
+auc <- pROC::auc(roc)
+auc
+```
+
+If we wanted to apply this in the context of cross validation, we could do it like this:
+
+```{r, eval=FALSE}
+x <- sdmTMB_cv(
+  present ~ s(depth), data = pcod, spatial = "off",
+  mesh = mesh, family = binomial(), k_folds = 2
+)
+roc <- pROC::roc(x$data$present, plogis(x$data$cv_predicted))
+auc <- pROC::auc(roc)
+auc
+```
+
+AUC may be sensitive to imbalances in the data, however, and alternative metrics may better approximate skill.
+Here we highlight an example of using true skill score (implemented in packages such as SDMtune):
+
+```{r}
+mesh <- make_mesh(pcod, c("X", "Y"), cutoff = 10)
+fit <- sdmTMB(present ~ 1, data = pcod, 
+  mesh = mesh, family = binomial())
+```
+
+Next, we can generate predicted probabilities and classes using a threshold of 0.5 as an example:
+
+```{r}
+pred <- predict(fit)
+pred$p <- plogis(pred$est)
+pred$pred_01 <- ifelse(pred$p < 0.5, 0, 1)
+```
+
+Next we create a confusion matrix and calculate the true skill score:
+
+```{r}
+conmat <- table(pred$pred_01, pred$present)
+true_neg <- conmat[1, 1]
+false_neg <- conmat[1, 2]
+false_pos <- conmat[2, 1]
+true_pos <- conmat[2, 2]
+
+# Calculate TSS:
+true_pos_rate <- true_pos / (true_pos + false_neg)
+true_neg_rate <- true_neg / (true_neg + false_pos)
+TSS <- true_pos_rate + true_neg_rate - 1
+TSS
+```
+
+In some cases, reporting the true negative or true positive rate might be of interest in addition to TSS.
+
 # References

From 796eb08f9dbf9fb980ab19c76a2282011055d43c Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 23 Sep 2024 14:48:27 -0700
Subject: [PATCH 19/21] Add progress bar in simulate.sdmTMB() #346

Also silent = FALSE by default
---
 DESCRIPTION            |  2 +-
 NEWS.md                |  2 ++
 R/tmb-sim.R            | 23 +++++++++++++----------
 man/simulate.sdmTMB.Rd |  2 +-
 4 files changed, 17 insertions(+), 12 deletions(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index 30269120..e287cdf3 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9008
+Version: 0.6.0.9009
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/NEWS.md b/NEWS.md
index bca9cc51..5815ca76 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,7 @@
 # sdmTMB (development version)
 
+* Add progress bar to `simulate.sdmTMB()`. #346
+
 * Add AUC and TSS examples to cross validation vignette. #268
 
 * Add `model` (linear predictor number) argument to coef() method. Also,
diff --git a/R/tmb-sim.R b/R/tmb-sim.R
index 724e884d..c54ff2a9 100644
--- a/R/tmb-sim.R
+++ b/R/tmb-sim.R
@@ -381,7 +381,7 @@ sdmTMB_simulate <- function(formula,
 simulate.sdmTMB <- function(object, nsim = 1L, seed = sample.int(1e6, 1L),
                             type = c("mle-eb", "mle-mvn"),
                             model = c(NA, 1, 2),
-                            re_form = NULL, mcmc_samples = NULL, silent = TRUE, ...) {
+                            re_form = NULL, mcmc_samples = NULL, silent = FALSE, ...) {
   set.seed(seed)
   type <- tolower(type)
   type <- match.arg(type)
@@ -421,18 +421,21 @@ simulate.sdmTMB <- function(object, nsim = 1L, seed = sample.int(1e6, 1L),
     new_par <- mcmc_samples
   }
 
-  # do the sim
+  # do the simulation
+  if (!silent) cli::cli_progress_bar("Simulating...", total = nsim)
+  ret <- list()
   if (!is.null(mcmc_samples)) { # we have a matrix
-    ret <- lapply(seq_len(nsim), function(i) {
-      if (!silent) cat("-")
-      newobj$simulate(par = new_par[, i, drop = TRUE], complete = FALSE)$y_i
-    })
+    for (i in seq_len(nsim)) {
+      if (!silent) cli::cli_progress_update()
+      ret[[i]] <- newobj$simulate(par = new_par[, i, drop = TRUE], complete = FALSE)$y_i
+    }
   } else {
-    ret <- lapply(seq_len(nsim), function(i) {
-      if (!silent) cat("-")
-      newobj$simulate(par = new_par, complete = FALSE)$y_i
-      })
+    for (i in seq_len(nsim)) {
+      if (!silent) cli::cli_progress_update()
+      ret[[i]] <- newobj$simulate(par = new_par, complete = FALSE)$y_i
+    }
   }
+  cli::cli_progress_done()
 
   if (isTRUE(object$family$delta)) {
     if (is.na(model[[1]])) {
diff --git a/man/simulate.sdmTMB.Rd b/man/simulate.sdmTMB.Rd
index e3cda1e8..2d8b25ca 100644
--- a/man/simulate.sdmTMB.Rd
+++ b/man/simulate.sdmTMB.Rd
@@ -12,7 +12,7 @@
   model = c(NA, 1, 2),
   re_form = NULL,
   mcmc_samples = NULL,
-  silent = TRUE,
+  silent = FALSE,
   ...
 )
 }

From ed977c76726082b6819374afb123d79d2bd830de Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Mon, 23 Sep 2024 15:15:26 -0700
Subject: [PATCH 20/21] Add print method for sdmTMB_cv() #319

---
 DESCRIPTION                            |  2 +-
 NAMESPACE                              |  1 +
 NEWS.md                                |  2 ++
 R/cross-val.R                          | 28 +++++++++++++++++++++++++-
 tests/testthat/test-cross-validation.R |  1 +
 5 files changed, 32 insertions(+), 2 deletions(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index e287cdf3..032b519c 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Type: Package
 Package: sdmTMB
 Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
-Version: 0.6.0.9009
+Version: 0.6.0.9010
 Authors@R: c(
     person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
            role = c("aut", "cre"),
diff --git a/NAMESPACE b/NAMESPACE
index 8c094a0e..ec805523 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -14,6 +14,7 @@ S3method(nobs,sdmTMB)
 S3method(plot,sdmTMBmesh)
 S3method(predict,sdmTMB)
 S3method(print,sdmTMB)
+S3method(print,sdmTMB_cv)
 S3method(ranef,sdmTMB)
 S3method(residuals,sdmTMB)
 S3method(simulate,sdmTMB)
diff --git a/NEWS.md b/NEWS.md
index 5815ca76..7cc10b26 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,7 @@
 # sdmTMB (development version)
 
+* Add print method for `sdmTMB_cv()` output. #319
+
 * Add progress bar to `simulate.sdmTMB()`. #346
 
 * Add AUC and TSS examples to cross validation vignette. #268
diff --git a/R/cross-val.R b/R/cross-val.R
index 6e871f13..0264d5c6 100644
--- a/R/cross-val.R
+++ b/R/cross-val.R
@@ -453,7 +453,7 @@ sdmTMB_cv <- function(
   pdHess <- vapply(out, `[[`, "pdHess", FUN.VALUE = logical(1L))
   max_grad <- vapply(out, `[[`, "max_gradient", FUN.VALUE = numeric(1L))
   converged <- all(pdHess)
-  list(
+  out <- list(
     data = data,
     models = models,
     fold_loglik = fold_cv_ll,
@@ -462,9 +462,35 @@ sdmTMB_cv <- function(
     pdHess = pdHess,
     max_gradients = max_grad
   )
+  `class<-`(out, "sdmTMB_cv")
 }
 
 log_sum_exp <- function(x) {
   max_x <- max(x)
   max_x + log(sum(exp(x - max_x)))
 }
+
+#' @export
+#' @import methods
+print.sdmTMB_cv <- function(x, ...) {
+  nmods <- length(x$models)
+  nconverged <- sum(x$converged)
+  cat(paste0("Cross validation of sdmTMB models with ", nmods, " folds.\n"))
+  cat("\n")
+  cat("Summary of the first fold model fit:\n")
+  cat("\n")
+  print(x$models[[1]])
+  cat("\n")
+  cat("Access the rest of the models in a list element named `models`.\n")
+  cat("E.g. `object$models[[2]]` for the 2nd fold model fit.\n")
+  cat("\n")
+  cat(paste0(nconverged, " out of ", nmods, " models are consistent with convergence.\n"))
+  cat("Figure out which folds these are in the `converged` list element.\n")
+  cat("\n")
+  cat(paste0("Out-of-sample log likelihood for each fold: ", paste(round(x$fold_loglik, 2), collapse = ", "), ".\n"))
+  cat("Access these values in the `fold_loglik` list element.\n")
+  cat("\n")
+  cat("Sum of out-of-sample log likelihoods:", round(x$sum_loglik, 2), "\n")
+  cat("More positive values imply better out-of-sample prediction.\n")
+  cat("Access this value in the `sum_loglik` list element.\n")
+}
diff --git a/tests/testthat/test-cross-validation.R b/tests/testthat/test-cross-validation.R
index e7f5abfe..5b1be663 100644
--- a/tests/testthat/test-cross-validation.R
+++ b/tests/testthat/test-cross-validation.R
@@ -12,6 +12,7 @@ test_that("Basic cross validation works", {
     data = d, mesh = spde,
     family = tweedie(link = "log"), time = "year", k_folds = 2
   )
+  print(x)
   expect_equal(class(x$sum_loglik), "numeric")
   expect_equal(x$sum_loglik, sum(x$data$cv_loglik))
   expect_equal(x$sum_loglik, sum(x$fold_loglik))

From 3b7b1bcf103a097ebb9b4e2243b63ac06505e733 Mon Sep 17 00:00:00 2001
From: Sean Anderson <sean@seananderson.ca>
Date: Tue, 24 Sep 2024 21:01:41 -0700
Subject: [PATCH 21/21] https cranlogs [skip ci]

---
 README.md | 2 +-
 header.md | 2 +-
 2 files changed, 2 insertions(+), 2 deletions(-)

diff --git a/README.md b/README.md
index b55a7061..7147a42e 100644
--- a/README.md
+++ b/README.md
@@ -10,7 +10,7 @@
 [![Documentation](https://img.shields.io/badge/documentation-sdmTMB-orange.svg?colorB=E91E63)](https://pbs-assess.github.io/sdmTMB/)
 [![R-CMD-check](https://github.com/pbs-assess/sdmTMB/workflows/R-CMD-check/badge.svg)](https://github.com/pbs-assess/sdmTMB/actions)
 [![Codecov test coverage](https://codecov.io/gh/pbs-assess/sdmTMB/branch/main/graph/badge.svg)](https://app.codecov.io/gh/pbs-assess/sdmTMB?branch=main)
-[![downloads](http://cranlogs.r-pkg.org/badges/sdmTMB)](https://cranlogs.r-pkg.org/)
+[![downloads](https://cranlogs.r-pkg.org/badges/sdmTMB)](https://cranlogs.r-pkg.org/)
 <!-- badges: end -->
 
 sdmTMB is an R package that fits spatial and spatiotemporal GLMMs (Generalized Linear Mixed Effects Models) using Template Model Builder ([TMB](https://github.com/kaskr/adcomp)), [R-INLA](https://www.r-inla.org/), and Gaussian Markov random fields. One common application is for species distribution models (SDMs). See the [documentation site](https://pbs-assess.github.io/sdmTMB/) and a preprint:
diff --git a/header.md b/header.md
index 0d9224a5..a2499511 100644
--- a/header.md
+++ b/header.md
@@ -9,7 +9,7 @@
 [![Documentation](https://img.shields.io/badge/documentation-sdmTMB-orange.svg?colorB=E91E63)](https://pbs-assess.github.io/sdmTMB/)
 [![R-CMD-check](https://github.com/pbs-assess/sdmTMB/workflows/R-CMD-check/badge.svg)](https://github.com/pbs-assess/sdmTMB/actions)
 [![Codecov test coverage](https://codecov.io/gh/pbs-assess/sdmTMB/branch/main/graph/badge.svg)](https://app.codecov.io/gh/pbs-assess/sdmTMB?branch=main)
-[![downloads](http://cranlogs.r-pkg.org/badges/sdmTMB)](https://cranlogs.r-pkg.org/)
+[![downloads](https://cranlogs.r-pkg.org/badges/sdmTMB)](https://cranlogs.r-pkg.org/)
 <!-- badges: end -->
 
 sdmTMB is an R package that fits spatial and spatiotemporal GLMMs (Generalized Linear Mixed Effects Models) using Template Model Builder ([TMB](https://github.com/kaskr/adcomp)), [R-INLA](https://www.r-inla.org/), and Gaussian Markov random fields. One common application is for species distribution models (SDMs). See the [documentation site](https://pbs-assess.github.io/sdmTMB/) and a preprint: