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
67 stars 2 forks source link

issue with new mixedup or new mgcv #44

Closed BenoitLondon closed 1 year ago

BenoitLondon commented 1 year ago

Hi, first thank you for that great package! Unfortunately, it does not work for my use case anymore :( Hope you can help :D Here's a reprex:

packageVersion("mixedup")
# [1] ‘0.4.0’
packageVersion("mgcv")
# [1] ‘1.8.41’
mtcarsf <- mtcars
mtcarsf$cyl <- factor(mtcarsf$cyl)
mod <- mgcv::bam(mpg ~ s(cyl, bs = "re", by = I(log(wt))), data = mtcarsf)
coefs <- mixedup::extract_random_effects(mod)
# Error in `[.data.frame`(model$model, , re_names[i]) : 
# undefined columns selected

digging into the code vn does not look to contain a valid column anymore when grabbed here

$vn
[1] "eraI(log(expected_dog_goals))"

if no by term in the spline, it works fine

> mod <- mgcv::bam(mpg ~ s(cyl, bs = "re"), data = mtcarsf)
> coefs <- mixedup::extract_random_effects(mod)
> coefs
# A tibble: 3 × 7
  group_var effect    group value    se lower_2.5 upper_97.5
  <chr>     <chr>     <chr> <dbl> <dbl>     <dbl>      <dbl>
1 cyl       Intercept 4      6.00  3.42     -0.71      12.7 
2 cyl       Intercept 6     -0.72  3.45     -7.48       6.04
3 cyl       Intercept 8     -5.28  3.41    -12.0        1.41
m-clark commented 1 year ago

I don't think it has to do with the latest mixedup, it's just that I don't have functionality for 'by' stuff (except maybe by accident). But I think it's also that you can do this more clearly as a random slope of log_wt across cyl approach

From your primary model:

> coef(mod)
        (Intercept) s(cyl):I(log(wt)).1 s(cyl):I(log(wt)).2 s(cyl):I(log(wt)).3 
          33.933241           -9.382861          -12.449715          -13.648564 

Differently expressed:

> mod <- mgcv::bam(mpg ~ s(log_wt, cyl, bs = "re"), data = mtcarsf)
> coefs <- mixedup::extract_random_effects(mod)
> coefs
# A tibble: 3 × 7
  group_var effect group  value    se lower_2.5 upper_97.5
  <chr>     <chr>  <chr>  <dbl> <dbl>     <dbl>      <dbl>
1 cyl       log_wt 4      -9.38  2.66     -14.6      -4.17
2 cyl       log_wt 6     -12.4   2.08     -16.5      -8.37
3 cyl       log_wt 8     -13.6   1.64     -16.9     -10.4 

Same result. 😄 I also think that model would be more clear to you in the future as well.

Other notes:

Hope that helps, but let me know if you have other issues.

m-clark commented 1 year ago

Going ahead and closing but feel free to reopen if needed.

BenoitLondon commented 1 year ago

thanks that all makes sense.