m-clark / mixedup

An R package for extracting results from mixed models that are easy to use and viable for presentation.
http://m-clark.github.io/mixedup
MIT License
64 stars 2 forks source link

extract_ranef() fails for random slopes on factor predictors with gam #12

Closed m-clark closed 4 years ago

m-clark commented 4 years ago

This was found via gammit for doing random categorical effects. It can be overcome by simply creating another grouping variable that is the interaction of the categorical predictor and grouping variable for the random effect (specifically model_w_wm at that link). This would also make clear what coefficients go with which effects. However it'd be nice to not to have to worry with it if specified as follows, but the level names are not retained in the smooth elements/attributes if you do it this way.

library(mgcv)
library(mixedup)
library(lme4)
library(dplyr)

data(sleepstudy)

ga_model <- gam(Reaction ~ Days + s(Subject, bs = "re") + s(Days, Subject, bs = "re"),
                data = sleepstudy %>% 
                    mutate(Days = factor(case_when(
                        Days < 2 ~ "x",
                        Days < 5 ~ "y",
                        TRUE ~ "z"
                    )))
                 ##
)

extract_random_coefs(ga_model) # error
m-clark commented 4 years ago

This works now but is perhaps less than ideal, as mgcv doesn't save the factor levels with the coefficients. They are ordered correctly (as checked by lme4). I still would recommend the 'interaction term' approach linked above instead of this random slope approach. I also updated the extract_random_coefs (which attempts to combine the fixed and random effects) documentation to reflect that it (at present) will typically not work on some models. This would be one where it will only produce the random intercepts, but not the random slopes, as there isn't a way to easily attach the fixed effects (I may revisit that issue in the future).

library(mgcv)
#> Loading required package: nlme
#> This is mgcv 1.8-31. For overview type 'help("mgcv-package")'.
library(mixedup)
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following object is masked from 'package:nlme':
#> 
#>     collapse
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union

data(sleepstudy)
#> Warning in data(sleepstudy): data set 'sleepstudy' not found

ga_model <-
  gam(
    Reaction ~ Days + s(Subject, bs = "re") + s(Days, Subject, bs = "re"),
    data = lme4::sleepstudy %>%
      mutate(Days = factor(
        case_when(Days < 2 ~ "x",
                  Days < 5 ~ "y",
                  TRUE ~ "z")
      ))
  )

extract_random_effects(ga_model) %>% tail()
#> # A tibble: 6 x 7
#>   group_var effect group value    se lower_2.5 upper_97.5
#>   <chr>     <chr>  <chr> <dbl> <dbl>     <dbl>      <dbl>
#> 1 Subject   Days   371    5.46  15.5     -24.9       35.8
#> 2 Subject   Days   371    1.04  14.9     -28.1       30.2
#> 3 Subject   Days   371   -7.12  14.3     -35.1       20.8
#> 4 Subject   Days   372   -1.53  15.5     -31.9       28.8
#> 5 Subject   Days   372    2.98  14.9     -26.2       32.1
#> 6 Subject   Days   372    5.83  14.3     -22.1       33.8

Created on 2020-06-24 by the reprex package (v0.3.0)