Model averaged models with polynomial terms have no confidence intervals

[MarcRieraDominguez]

Hi! I have notice that ggpredict will not plot confidence intervals for model-averaged objects (otabined with package MuMIn), when the model includes a polinomial term declared with poly(x, 2, raw = TRUE). Could the functionality be extended to including polynomials in averaged models? Thank you for your time!

See below a reproducible example.

options(na.action = "na.fail")

library(ggeffects)
#> Warning: package 'ggeffects' was built under R version 4.2.3
library(MuMIn)
#> Warning: package 'MuMIn' was built under R version 4.2.3
library(tidyverse)
#> Warning: package 'ggplot2' was built under R version 4.2.2
#> Warning: package 'tibble' was built under R version 4.2.3
#> Warning: package 'tidyr' was built under R version 4.2.3
#> Warning: package 'purrr' was built under R version 4.2.3
#> Warning: package 'dplyr' was built under R version 4.2.3
#> Warning: package 'stringr' was built under R version 4.2.3

mtcars$am <- factor(mtcars$am)
mtcars$disp.sq <- mtcars$disp * mtcars$disp

ggplot(aes(x = mpg, y = disp, col = am), data = mtcars) +
  geom_jitter() +
  geom_smooth()
#> `geom_smooth()` using method = 'loess' and formula = 'y ~ x'

lm(disp ~ mpg*am + gear, mtcars) %>% 
  dredge() %>% 
  model.avg(fit = TRUE) %>% 
  ggpredict(terms = c("mpg[all]", "am")) %>% 
  plot()
#> Fixed term is "(Intercept)"

lm(disp ~ poly(mpg, 2, raw = TRUE)*am + gear, mtcars) %>% 
  dredge() %>% 
  model.avg(fit = TRUE) %>% 
  ggpredict(terms = c("mpg[all]", "am")) %>% 
  plot()
#> Fixed term is "(Intercept)"
#> Warning: Some model terms could not be found in model data.
#>   You probably need to load the data into the environment.
#> Error: Confidence intervals could not be computed.
#> Warning: Some model terms could not be found in model data.
#>   You probably need to load the data into the environment.

sessionInfo()
#> R version 4.2.0 (2022-04-22 ucrt)
#> Platform: x86_64-w64-mingw32/x64 (64-bit)
#> Running under: Windows 10 x64 (build 22000)
#> 
#> Matrix products: default
#> 
#> locale:
#> [1] LC_COLLATE=Spanish_Spain.utf8  LC_CTYPE=Spanish_Spain.utf8   
#> [3] LC_MONETARY=Spanish_Spain.utf8 LC_NUMERIC=C                  
#> [5] LC_TIME=Spanish_Spain.utf8    
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#>  [1] forcats_0.5.1     stringr_1.5.0     dplyr_1.1.2       purrr_1.0.1      
#>  [5] readr_2.1.2       tidyr_1.3.0       tibble_3.2.1      ggplot2_3.4.0    
#>  [9] tidyverse_1.3.1   MuMIn_1.47.5      ggeffects_1.3.1.5
#> 
#> loaded via a namespace (and not attached):
#>  [1] lubridate_1.8.0  lattice_0.20-45  assertthat_0.2.1 digest_0.6.29   
#>  [5] utf8_1.2.2       R6_2.5.1         cellranger_1.1.0 backports_1.4.1 
#>  [9] reprex_2.0.2     stats4_4.2.0     evaluate_0.15    httr_1.4.3      
#> [13] highr_0.9        pillar_1.9.0     rlang_1.1.1      readxl_1.4.0    
#> [17] rstudioapi_0.13  Matrix_1.4-1     rmarkdown_2.14   labeling_0.4.2  
#> [21] splines_4.2.0    munsell_0.5.0    broom_0.8.0      compiler_4.2.0  
#> [25] modelr_0.1.8     xfun_0.40        pkgconfig_2.0.3  mgcv_1.8-40     
#> [29] htmltools_0.5.6  insight_0.19.5   tidyselect_1.2.0 fansi_1.0.3     
#> [33] crayon_1.5.1     tzdb_0.3.0       dbplyr_2.1.1     withr_2.5.0     
#> [37] grid_4.2.0       nlme_3.1-157     jsonlite_1.8.0   gtable_0.3.0    
#> [41] lifecycle_1.0.3  DBI_1.1.2        magrittr_2.0.3   scales_1.2.1    
#> [45] datawizard_0.7.0 cli_3.6.1        stringi_1.7.6    farver_2.1.0    
#> [49] fs_1.5.2         snakecase_0.11.0 xml2_1.3.3       ellipsis_0.3.2  
#> [53] generics_0.1.2   vctrs_0.6.3      sjlabelled_1.2.0 tools_4.2.0     
#> [57] glue_1.6.2       hms_1.1.1        fastmap_1.1.0    yaml_2.3.5      
#> [61] colorspace_2.0-3 rvest_1.0.2      knitr_1.39       haven_2.5.0

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

[MarcRieraDominguez]

Hi! I forgot to mention in the original issue that ggpredict will return confidence bands if model-averaged polynomial terms are declared with I(x^2), as shown in the reprex below. This is however more cumbersome than the more compact poly(x, 2, raw = TRUE). Note that I used dependency chain (dc) and subset to mimic the behaviour of poly(x, 2, raw = TRUE) while using I(x^2) (basically, that the linear and quadratic terms are always together).

library(ggeffects)
library(MuMIn)
#> Warning: package 'MuMIn' was built under R version 4.1.3
library(tidyverse)
#> Warning: package 'tidyverse' was built under R version 4.1.3
#> Warning: package 'ggplot2' was built under R version 4.1.3
#> Warning: package 'tibble' was built under R version 4.1.3
#> Warning: package 'tidyr' was built under R version 4.1.3
#> Warning: package 'readr' was built under R version 4.1.3
#> Warning: package 'purrr' was built under R version 4.1.3
#> Warning: package 'dplyr' was built under R version 4.1.3
#> Warning: package 'stringr' was built under R version 4.1.3
#> Warning: package 'forcats' was built under R version 4.1.3

options(na.action = "na.fail")

mtcars$am <- factor(mtcars$am)

mod.avg.i <-
  lm(disp ~ mpg + I(mpg^2) + am + gear, mtcars) %>% 
  dredge(subset = dc(mpg, I(mpg^2))) %>% 
  subset(!( has(mpg) & !has(I(mpg^2)))) %>% 
  model.avg(fit = TRUE)
#> Fixed term is "(Intercept)"

mod.avg.poly <-
  lm(disp ~ poly(mpg, 2, raw = TRUE) + am + gear, mtcars) %>% 
  dredge() %>% 
  model.avg(fit = TRUE)
#> Fixed term is "(Intercept)"

# ggeffects uses full-averaging if I recall correctly
lapply(list(mod.avg.poly, mod.avg.i), function(x) coef(x, full = TRUE))
#> [[1]]
#>               (Intercept)                      gear poly(mpg, 2, raw = TRUE)1 
#>               954.5036243               -10.5936478               -51.2140535 
#> poly(mpg, 2, raw = TRUE)2                       am1 
#>                 0.7936609               -10.0217995 
#> 
#> [[2]]
#> (Intercept)        gear         mpg    I(mpg^2)         am1 
#> 954.5036243 -10.5936478 -51.2140535   0.7936609 -10.0217995

plot(ggpredict(mod.avg.i, terms = c("mpg[all]", "am")))

plot(ggpredict(mod.avg.poly, terms = c("mpg[all]", "am")))
#> Warning: Some model terms could not be found in model data.
#>   You probably need to load the data into the environment.
#> Error: Confidence intervals could not be computed.
#> Warning: Some model terms could not be found in model data.
#>   You probably need to load the data into the environment.

sessionInfo()
#> R version 4.1.0 (2021-05-18)
#> Platform: x86_64-w64-mingw32/x64 (64-bit)
#> Running under: Windows 10 x64 (build 19045)
#> 
#> Matrix products: default
#> 
#> locale:
#> [1] LC_COLLATE=Spanish_Spain.1252  LC_CTYPE=Spanish_Spain.1252   
#> [3] LC_MONETARY=Spanish_Spain.1252 LC_NUMERIC=C                  
#> [5] LC_TIME=Spanish_Spain.1252    
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#>  [1] forcats_0.5.1   stringr_1.5.0   dplyr_1.1.1     purrr_1.0.1    
#>  [5] readr_2.1.2     tidyr_1.3.0     tibble_3.2.1    ggplot2_3.4.2  
#>  [9] tidyverse_1.3.1 MuMIn_1.46.0    ggeffects_1.3.1
#> 
#> loaded via a namespace (and not attached):
#>  [1] lubridate_1.8.0  lattice_0.20-44  assertthat_0.2.1 digest_0.6.29   
#>  [5] utf8_1.2.2       R6_2.5.1         cellranger_1.1.0 backports_1.4.1 
#>  [9] reprex_2.0.1     stats4_4.1.0     evaluate_0.15    httr_1.4.2      
#> [13] highr_0.9        pillar_1.9.0     rlang_1.1.0      readxl_1.4.0    
#> [17] rstudioapi_0.13  Matrix_1.3-3     rmarkdown_2.13   labeling_0.4.2  
#> [21] munsell_0.5.0    broom_0.8.0      compiler_4.1.0   modelr_0.1.8    
#> [25] xfun_0.30        pkgconfig_2.0.3  htmltools_0.5.2  insight_0.19.1  
#> [29] tidyselect_1.2.0 fansi_1.0.3      crayon_1.5.1     tzdb_0.3.0      
#> [33] dbplyr_2.1.1     withr_2.5.0      grid_4.1.0       nlme_3.1-152    
#> [37] jsonlite_1.8.0   gtable_0.3.0     lifecycle_1.0.3  DBI_1.1.2       
#> [41] magrittr_2.0.3   scales_1.2.0     datawizard_0.6.3 cli_3.6.0       
#> [45] stringi_1.7.6    farver_2.1.0     fs_1.5.2         xml2_1.3.3      
#> [49] ellipsis_0.3.2   generics_0.1.2   vctrs_0.6.1      tools_4.1.0     
#> [53] glue_1.6.2       hms_1.1.1        fastmap_1.1.0    yaml_2.3.5      
#> [57] colorspace_2.0-3 rvest_1.0.2      knitr_1.38       haven_2.5.0

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