Add and document ggeffects ggpredict support

[seananderson]

Moving discussion here to its own thread.

• ggeffects::ggeffect() (marginal effects based on the fixed effect covariance matrix, doesn't work for splines or breakpoints) was already working in the 'dev' branch
• ggeffects::ggpredict() (conditional predictions, a bit slower because it calls TMB, works for splines and breakpoints) should work if you install this fork and branch of ggeffects:

remotes::install_github("seananderson/ggeffects", ref = "sdmTMB")

Still need to get delta models working.

Examples:

library(magrittr)
library(sdmTMB)
library(ggeffects)
pcod_2011$fyear <- as.factor(pcod_2011$year)
mesh <- make_mesh(pcod_2011, c("X", "Y"), cutoff = 15)
fit <- sdmTMB(
  present ~ poly(depth_scaled, 2) + fyear,
  data = pcod_2011,
  spatial = "off",
  family = binomial(link = "logit")
)

ggpredict(fit, "depth_scaled [all]") %>% plot()
#> Loading required namespace: ggplot2

ggeffect(fit, "depth_scaled [all]") %>% plot()

ggeffect(fit, c("depth_scaled [all]", "fyear")) %>% plot()

fit <- sdmTMB(
  present ~ s(depth_scaled) + fyear,
  data = pcod_2011,
  spatial = "off",
  family = binomial(link = "logit")
)
ggpredict(fit, c("depth_scaled [all]", "fyear")) %>% plot()

^{Created on 2022-06-01 by the reprex package (v2.0.1)}

[jindivero]

Awesome. If you want me to add anything to the ggeffects vignette just let me know.

[ericward-noaa]

Going back through this today, I realized that most of the examples aren't working for me. Probably should flag as needing refreshing these examples. Specifically I get

ggpredict(fit, "depth_scaled [all]") %>% plot()
Error in `ggplot2::geom_line()`:
! Problem while computing aesthetics.
ℹ Error occurred in the 1st layer.
Caused by error in `.data[["predicted"]]`:
! Column `predicted` not found in `.data`.
Run `rlang::last_trace()` to see where the error occurred.
Warning messages:
1: Could not compute predictions for model of class `sdmTMB`. 
2: In `[<-.factor`(`*tmp*`, seq_along(fac[[n]]), value = 1:4) :
  invalid factor level, NA generated
3: Could not compute standard errors or confidence intervals because the model and
  variance-covariance matrices are non-conformable. This can sometimes happen when the
  `data` used to make predictions fails to include all the levels of a factor variable
  or all the interaction components. 

[seananderson]

Did you install this fork and start a fresh session?

remotes::install_github("seananderson/ggeffects", ref = "sdmTMB")

It still seems to be working for me. We should get this into a pull request I guess. Ideally it would work on delta models first? Or maybe just go ahead.

[ericward-noaa]

I realized I had some TMB conflict, but after fixing that everything works (I was using that branch). I think probably fine to put this into a pull request though -- and sort out the delta issues later

[seananderson]

Pull request: https://github.com/strengejacke/ggeffects/pull/393

[seananderson]

The pull request has been merged in ggeffects and ggpredict() should be working for non-delta models. Here is a proposed interface for delta models that doesn't require any changes to ggeffects (other than maybe making the delta model error in ggpredict() more helpful if called incorrectly).

I'm trying to decide if the name set_delta_model() is right. I think this might be the slickest way rather than a ton of bespoke code within ggeffects/insight. It basically just adds an attribute that various methods then pick up on to know which model to select for formulas, families, etc.

https://github.com/pbs-assess/sdmTMB/blob/b10ce3d884538440a1b5922a2f9a1d7f9e21cba9/R/utils.R#L524-L554

e.g. https://github.com/pbs-assess/sdmTMB/blob/b10ce3d884538440a1b5922a2f9a1d7f9e21cba9/R/methods.R#L147-L158

Longterm this should let us phase out visreg, which has an inconsistent interface, requires a ton of complex code within sdmTMB to work, and does not seem to be actively maintained.

Examples:

library(sdmTMB)

# non-delta
fit <- sdmTMB(density ~ poly(depth_scaled, 2), data = pcod_2011,
  spatial = "off", family = tweedie())

ggeffects::ggpredict(fit, "depth_scaled [all]") |> 
  plot()

# delta:
fit <- sdmTMB(density ~ poly(depth_scaled, 2), data = pcod_2011,
  spatial = "off", family = delta_gamma())

# binomial part:
set_delta_model(fit, model = 1) |>
  ggeffects::ggpredict("depth_scaled [all]") |> 
  plot()

# gamma part:
set_delta_model(fit, model = 2) |>
  ggeffects::ggpredict("depth_scaled [all]") |> 
  plot()

# combined:
set_delta_model(fit, model = NA) |>
  ggeffects::ggpredict("depth_scaled [all]") |> 
  plot()

^{Created on 2023-10-24 with reprex v2.0.2}