Include MAD_SD of random parameters in `mod$ses`?

[emstruong]

Summary:

I'm trying to easily get the MAD_SD for the random parameters of my model, but they seem to be missing from the ses list of the model for at least stan_glmer. Could they be added? It seems like the fixed effects and random values are included... This is related to https://github.com/bbolker/broom.mixed/issues/156#issue-2662483581

Reproducible Steps:

library(rstanarm)
#> Loading required package: Rcpp
#> This is rstanarm version 2.32.1
#> - See https://mc-stan.org/rstanarm/articles/priors for changes to default priors!
#> - Default priors may change, so it's safest to specify priors, even if equivalent to the defaults.
#> - For execution on a local, multicore CPU with excess RAM we recommend calling
#>   options(mc.cores = parallel::detectCores())
fit <- stan_glmer(mpg ~ wt + (1|cyl) + (1+wt|gear), data = mtcars,
                  iter = 500, chains = 2)
#> 
#> SAMPLING FOR MODEL 'continuous' NOW (CHAIN 1).
#> Chain 1: 
#> Chain 1: Gradient evaluation took 4e-05 seconds
#> Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.4 seconds.
#> Chain 1: Adjust your expectations accordingly!
#> Chain 1: 
#> Chain 1: 
#> Chain 1: Iteration:   1 / 500 [  0%]  (Warmup)
#> Chain 1: Iteration:  50 / 500 [ 10%]  (Warmup)
#> Chain 1: Iteration: 100 / 500 [ 20%]  (Warmup)
#> Chain 1: Iteration: 150 / 500 [ 30%]  (Warmup)
#> Chain 1: Iteration: 200 / 500 [ 40%]  (Warmup)
#> Chain 1: Iteration: 250 / 500 [ 50%]  (Warmup)
#> Chain 1: Iteration: 251 / 500 [ 50%]  (Sampling)
#> Chain 1: Iteration: 300 / 500 [ 60%]  (Sampling)
#> Chain 1: Iteration: 350 / 500 [ 70%]  (Sampling)
#> Chain 1: Iteration: 400 / 500 [ 80%]  (Sampling)
#> Chain 1: Iteration: 450 / 500 [ 90%]  (Sampling)
#> Chain 1: Iteration: 500 / 500 [100%]  (Sampling)
#> Chain 1: 
#> Chain 1:  Elapsed Time: 0.402 seconds (Warm-up)
#> Chain 1:                0.212 seconds (Sampling)
#> Chain 1:                0.614 seconds (Total)
#> Chain 1: 
#> 
#> SAMPLING FOR MODEL 'continuous' NOW (CHAIN 2).
#> Chain 2: 
#> Chain 2: Gradient evaluation took 1.7e-05 seconds
#> Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.17 seconds.
#> Chain 2: Adjust your expectations accordingly!
#> Chain 2: 
#> Chain 2: 
#> Chain 2: Iteration:   1 / 500 [  0%]  (Warmup)
#> Chain 2: Iteration:  50 / 500 [ 10%]  (Warmup)
#> Chain 2: Iteration: 100 / 500 [ 20%]  (Warmup)
#> Chain 2: Iteration: 150 / 500 [ 30%]  (Warmup)
#> Chain 2: Iteration: 200 / 500 [ 40%]  (Warmup)
#> Chain 2: Iteration: 250 / 500 [ 50%]  (Warmup)
#> Chain 2: Iteration: 251 / 500 [ 50%]  (Sampling)
#> Chain 2: Iteration: 300 / 500 [ 60%]  (Sampling)
#> Chain 2: Iteration: 350 / 500 [ 70%]  (Sampling)
#> Chain 2: Iteration: 400 / 500 [ 80%]  (Sampling)
#> Chain 2: Iteration: 450 / 500 [ 90%]  (Sampling)
#> Chain 2: Iteration: 500 / 500 [100%]  (Sampling)
#> Chain 2: 
#> Chain 2:  Elapsed Time: 0.308 seconds (Warm-up)
#> Chain 2:                0.135 seconds (Sampling)
#> Chain 2:                0.443 seconds (Total)
#> Chain 2:
#> Warning: There were 1 divergent transitions after warmup. See
#> https://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
#> to find out why this is a problem and how to eliminate them.
#> Warning: Examine the pairs() plot to diagnose sampling problems
#> Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
#> Running the chains for more iterations may help. See
#> https://mc-stan.org/misc/warnings.html#bulk-ess
#> Warning: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
#> Running the chains for more iterations may help. See
#> https://mc-stan.org/misc/warnings.html#tail-ess
fit$ses
#>           (Intercept)                    wt  b[(Intercept) cyl:4] 
#>             3.1357937             0.9469535             2.0728322 
#>  b[(Intercept) cyl:6]  b[(Intercept) cyl:8] b[(Intercept) gear:3] 
#>             1.4681377             1.9467202             0.3974553 
#>          b[wt gear:3] b[(Intercept) gear:4]          b[wt gear:4] 
#>             0.2688166             0.3394349             0.2915114 
#> b[(Intercept) gear:5]          b[wt gear:5] 
#>             0.3644062             0.2846547

^{Created on 2024-11-15 with reprex v2.1.1}

RStanARM Version:

The version of the rstanarm package you are running (e.g., from packageVersion("rstanarm")): ‘2.32.1’

R Version:

The version of R you are running (e.g., from getRversion()): ‘4.4.1’

Operating System:

Your operating system (e.g., OS X 10.11.3): Ubuntu 22.04

[jgabry]

Thanks for pointing that out. I think for now the easiest way to get the MAD for all parameters is just by computing it from the posterior draws. You could do:

apply(as.matrix(fit), 2, mad)

Does that work for you?

[emstruong]

So that definitely works for my current purposes--though a fix would still be nice for broom.mixed functionality. Was there any particular technical reason why the MAD for the random parameters was omitted or was it just a lapse?

[jgabry]

Good question, I'm not sure, it's been so long! I'll have to take a look at the source code and refresh my memory when I have a chance. I don't have time at the moment, so I'll leave this open for now until I can get to it.