Skip to content

Commit cf46f2b

Browse files
seanandersonseananderson
committed
Let emmeans work with smoothers included
1 parent 08faf82 commit cf46f2b

File tree

4 files changed

+206
-4
lines changed

4 files changed

+206
-4
lines changed

‎DESCRIPTION

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,7 @@
11
Type: Package
22
Package: sdmTMB
33
Title: Spatial and Spatiotemporal SPDE-Based GLMMs with 'TMB'
4-
Version: 0.7.4.9000
4+
Version: 0.7.4.9001
55
Authors@R: c(
66
person(c("Sean", "C."), "Anderson", , "sean@seananderson.ca",
77
role = c("aut", "cre"),

‎NEWS.md

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,9 @@
11
# sdmTMB (development version)
22

3+
* Fix `emmeans` support for models with smoothers (`s()` terms). Previously,
4+
`emmeans` would fail with "Non-conformable elements in reference grid" when
5+
smoothers were included in the model formula.
6+
37
# sdmTMB 0.7.4
48

59
## Minor improvements and fixes

‎R/emmeans.R

Lines changed: 13 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -70,11 +70,21 @@ emm_basis.sdmTMB <- function(object, trms, xlev, grid, ...) {
7070
contrasts <- contrasts[names(contrasts) %in% all.vars(terms(object))]
7171
m <- model.frame(trms, grid, na.action = stats::na.pass, xlev = xlev)
7272
X <- model.matrix(trms, m, contrasts.arg = contrasts)
73-
bhat <- fixef(object)
73+
74+
# Get coefficients using existing fixef method and filter to parametric terms only
75+
# This fixes the issue with smoothers by only including parametric terms
76+
all_bhat <- fixef(object)
77+
78+
# Only keep coefficients that correspond to columns in the design matrix X
79+
bhat <- all_bhat[names(all_bhat) %in% colnames(X)]
80+
# Ensure coefficient order matches design matrix column order
81+
bhat <- bhat[colnames(X)]
82+
7483
if (length(bhat) < ncol(X)) {
7584
kept <- match(names(bhat), dimnames(X)[[2]])
76-
bhat <- NA * X[1, ]
77-
bhat[kept] <- fixef(object)
85+
full_bhat <- NA * X[1, ]
86+
full_bhat[kept] <- bhat
87+
bhat <- full_bhat
7888
modmat <- model.matrix(
7989
trms, model.frame(object),
8090
contrasts.arg = contrasts

‎tests/testthat/test-emmeans.R

Lines changed: 188 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,188 @@
1+
test_that("emmeans works with parametric models", {
2+
skip_if_not_installed("emmeans")
3+
library(emmeans)
4+
5+
# Basic parametric model without spatial component
6+
fit <- sdmTMB(
7+
present ~ as.factor(year),
8+
data = pcod_2011,
9+
spatial = "off",
10+
family = binomial()
11+
)
12+
13+
# Test emmeans creation
14+
emm <- emmeans::emmeans(fit, ~ year)
15+
expect_s4_class(emm, "emmGrid")
16+
17+
# Test pairwise comparisons
18+
pairs_result <- pairs(emm) # pairs is a method, not namespaced
19+
expect_s4_class(pairs_result, "emmGrid")
20+
21+
# Test response scale
22+
emm_resp <- emmeans::emmeans(fit, ~ year, type = "response")
23+
expect_s4_class(emm_resp, "emmGrid")
24+
25+
# Check that results have correct dimensions
26+
emm_summary <- summary(emm)
27+
expect_equal(nrow(emm_summary), length(unique(pcod_2011$year)))
28+
expect_true(all(c("year", "emmean", "SE", "df", "lower.CL", "upper.CL") %in% names(emm_summary)))
29+
})
30+
31+
test_that("emmeans works with continuous covariates", {
32+
skip_if_not_installed("emmeans")
33+
library(emmeans)
34+
35+
fit <- sdmTMB(
36+
present ~ as.factor(year) + depth_scaled,
37+
data = pcod_2011,
38+
spatial = "off",
39+
family = binomial()
40+
)
41+
42+
# Test emmeans with continuous covariate
43+
emm <- emmeans::emmeans(fit, ~ year)
44+
expect_s4_class(emm, "emmGrid")
45+
46+
# Test emtrends for slopes
47+
emtrends_result <- emmeans::emtrends(fit, ~ year, var = "depth_scaled")
48+
expect_s4_class(emtrends_result, "emmGrid")
49+
50+
# Test interaction
51+
fit_int <- sdmTMB(
52+
present ~ as.factor(year) * depth_scaled,
53+
data = pcod_2011,
54+
spatial = "off",
55+
family = binomial()
56+
)
57+
58+
emm_int <- emmeans::emmeans(fit_int, ~ year)
59+
expect_s4_class(emm_int, "emmGrid")
60+
})
61+
62+
test_that("emmeans works with smoothers", {
63+
skip_if_not_installed("emmeans")
64+
library(emmeans)
65+
66+
# Model with smoother - this was previously broken
67+
fit_smooth <- sdmTMB(
68+
present ~ as.factor(year) + s(depth_scaled),
69+
data = pcod_2011,
70+
spatial = "off",
71+
family = binomial()
72+
)
73+
74+
# Test that emmeans works with smoothers
75+
emm_smooth <- emmeans::emmeans(fit_smooth, ~ year)
76+
expect_s4_class(emm_smooth, "emmGrid")
77+
78+
# Test pairwise comparisons with smoothers
79+
pairs_smooth <- pairs(emm_smooth)
80+
expect_s4_class(pairs_smooth, "emmGrid")
81+
82+
# Test response scale with smoothers
83+
emm_smooth_resp <- emmeans::emmeans(fit_smooth, ~ year, type = "response")
84+
expect_s4_class(emm_smooth_resp, "emmGrid")
85+
86+
# Check that results are reasonable
87+
emm_summary <- summary(emm_smooth)
88+
expect_equal(nrow(emm_summary), length(unique(pcod_2011$year)))
89+
expect_true(all(is.finite(emm_summary$emmean)))
90+
expect_true(all(emm_summary$SE > 0))
91+
})
92+
93+
test_that("emmeans works with multiple smoothers", {
94+
skip_if_not_installed("emmeans")
95+
library(emmeans)
96+
97+
# Model with multiple smoothers
98+
fit_multi_smooth <- sdmTMB(
99+
present ~ as.factor(year) + s(depth_scaled) + s(X),
100+
data = pcod_2011,
101+
spatial = "off",
102+
family = binomial()
103+
)
104+
105+
# Test that emmeans still works with multiple smoothers
106+
emm_multi <- emmeans::emmeans(fit_multi_smooth, ~ year)
107+
expect_s4_class(emm_multi, "emmGrid")
108+
109+
# Test pairwise comparisons
110+
pairs_multi <- pairs(emm_multi)
111+
expect_s4_class(pairs_multi, "emmGrid")
112+
113+
# Results should be sensible
114+
emm_summary <- summary(emm_multi)
115+
expect_true(all(is.finite(emm_summary$emmean)))
116+
})
117+
118+
test_that("emmeans results are consistent between parametric and smooth models", {
119+
skip_if_not_installed("emmeans")
120+
library(emmeans)
121+
122+
# Parametric model
123+
fit_param <- sdmTMB(
124+
present ~ as.factor(year) + depth_scaled,
125+
data = pcod_2011,
126+
spatial = "off",
127+
family = binomial()
128+
)
129+
130+
# Smooth model
131+
fit_smooth <- sdmTMB(
132+
present ~ as.factor(year) + s(depth_scaled),
133+
data = pcod_2011,
134+
spatial = "off",
135+
family = binomial()
136+
)
137+
138+
# Get emmeans for both
139+
emm_param <- emmeans::emmeans(fit_param, ~ year)
140+
emm_smooth <- emmeans::emmeans(fit_smooth, ~ year)
141+
142+
# Results should be reasonably similar (smooth vs linear depth effect)
143+
param_est <- summary(emm_param)$emmean
144+
smooth_est <- summary(emm_smooth)$emmean
145+
146+
# They shouldn't be identical but should be correlated
147+
expect_true(cor(param_est, smooth_est) > 0.8)
148+
})
149+
150+
test_that("emmeans error handling", {
151+
skip_if_not_installed("emmeans")
152+
library(emmeans)
153+
154+
# Delta model should produce informative error
155+
fit_delta <- sdmTMB(
156+
density ~ as.factor(year),
157+
data = pcod_2011,
158+
spatial = "off",
159+
family = delta_gamma()
160+
)
161+
162+
expect_error(
163+
emmeans::emmeans(fit_delta, ~ year),
164+
"Delta models not yet supported"
165+
)
166+
})
167+
168+
test_that("emmeans works with binomial family and smoothers", {
169+
skip_if_not_installed("emmeans")
170+
library(emmeans)
171+
172+
# Test binomial with smoother - this is the main use case that was broken
173+
fit <- sdmTMB(
174+
present ~ as.factor(year) + s(depth_scaled),
175+
data = pcod_2011,
176+
spatial = "off",
177+
family = binomial()
178+
)
179+
180+
# Test emmeans works
181+
emm <- emmeans::emmeans(fit, ~ year)
182+
expect_s4_class(emm, "emmGrid")
183+
184+
# Test that we can extract meaningful results
185+
emm_summary <- summary(emm)
186+
expect_true(nrow(emm_summary) == 4) # 4 years
187+
expect_true(all(c("year", "emmean") %in% names(emm_summary)))
188+
})

0 commit comments

Comments
 (0)