Possible bug with estimate_grouplevel() for brms model

[strengejacke]

Fixes #229

[strengejacke]

Here's the current implementation:

library(modelbased)
# lme4 model
data(mtcars)
model <- lme4::lmer(mpg ~ hp + (1 | carb), data = mtcars)

estimate_grouplevel(model)
#> Group | Level | Parameter   | Coefficient |   SE |        95% CI
#> ----------------------------------------------------------------
#> carb  | 1     | (Intercept) |        0.41 | 0.84 | [-1.24, 2.05]
#> carb  | 2     | (Intercept) |        0.11 | 0.78 | [-1.42, 1.65]
#> carb  | 3     | (Intercept) |       -0.32 | 0.94 | [-2.16, 1.51]
#> carb  | 4     | (Intercept) |       -0.78 | 0.78 | [-2.31, 0.75]
#> carb  | 6     | (Intercept) |        0.09 | 1.00 | [-1.87, 2.05]
#> carb  | 8     | (Intercept) |        0.50 | 1.00 | [-1.47, 2.46]

estimate_grouplevel(model, type = "total")
#> Group | Level | Parameter   | Coefficient
#> -----------------------------------------
#> carb  | 1     | (Intercept) |       30.18
#> carb  | 2     | (Intercept) |       29.88
#> carb  | 3     | (Intercept) |       29.45
#> carb  | 4     | (Intercept) |       28.99
#> carb  | 6     | (Intercept) |       29.87
#> carb  | 8     | (Intercept) |       30.27

model <- insight::download_model("brms_sigma_3")

estimate_grouplevel(model)
#> Component   | Group | Level | Parameter   |    Median |        95% CI
#> ---------------------------------------------------------------------
#> conditional | Group | G1    | Intercept   |      1.12 | [-0.65, 5.70]
#> conditional | Group | G1    | Petal.Width | -9.46e-03 | [-1.01, 0.94]
#> conditional | Group | G2    | Intercept   |     -0.57 | [-2.46, 3.64]
#> conditional | Group | G2    | Petal.Width |      0.07 | [-0.80, 0.94]
#> conditional | Group | G3    | Intercept   |     -0.45 | [-2.31, 4.13]
#> conditional | Group | G3    | Petal.Width |     -0.14 | [-1.21, 0.52]
#> sigma       | Group | G1    | Intercept   |      0.15 | [-1.05, 1.64]
#> sigma       | Group | G1    | Petal.Width |      0.01 | [-1.06, 1.62]
#> sigma       | Group | G2    | Intercept   |     -0.18 | [-1.84, 0.92]
#> sigma       | Group | G2    | Petal.Width |     -0.01 | [-1.19, 1.21]
#> sigma       | Group | G3    | Intercept   |     -0.06 | [-1.57, 1.60]
#> sigma       | Group | G3    | Petal.Width |     -0.04 | [-1.29, 1.02]

estimate_grouplevel(model, type = "total")
#> Component   | Group | Level | Parameter   | Median |         95% CI
#> -------------------------------------------------------------------
#> conditional | Group | G1    | Intercept   |   3.22 | [ 2.99,  3.44]
#> conditional | Group | G1    | Petal.Width |   0.83 | [ 0.06,  1.69]
#> conditional | Group | G2    | Intercept   |   1.49 | [ 0.97,  1.94]
#> conditional | Group | G2    | Petal.Width |   0.96 | [ 0.62,  1.35]
#> conditional | Group | G3    | Intercept   |   1.62 | [ 1.04,  2.22]
#> conditional | Group | G3    | Petal.Width |   0.67 | [ 0.38,  0.95]
#> sigma       | Group | G1    | Intercept   |  -1.04 | [-1.42, -0.71]
#> sigma       | Group | G1    | Petal.Width |   0.16 | [-0.81,  1.63]
#> sigma       | Group | G2    | Intercept   |  -1.51 | [-2.71, -0.62]
#> sigma       | Group | G2    | Petal.Width |   0.08 | [-0.57,  0.97]
#> sigma       | Group | G3    | Intercept   |  -1.32 | [-2.54, -0.29]
#> sigma       | Group | G3    | Petal.Width |   0.02 | [-0.48,  0.64]

model <- insight::download_model("brms_zi_4")

estimate_grouplevel(model)
#> Component     | Group   | Level | Parameter | Median |          95% CI
#> ----------------------------------------------------------------------
#> conditional   | ID      | 1     | Intercept |  -0.23 | [ -1.16,  0.79]
#> conditional   | ID      | 1     | zg        |   0.16 | [ -0.01,  0.32]
#> conditional   | ID      | 2     | Intercept |  -0.15 | [ -0.98,  0.86]
#> conditional   | ID      | 2     | zg        |   0.25 | [  0.09,  0.41]
#> conditional   | ID      | 3     | Intercept |   0.14 | [ -0.74,  1.15]
#> conditional   | ID      | 3     | zg        |   0.08 | [ -0.05,  0.23]
#> conditional   | ID      | 4     | Intercept |  -0.65 | [ -1.55,  0.36]
#> conditional   | ID      | 4     | zg        |   0.37 | [  0.28,  0.45]
#> conditional   | persons | 1     | Intercept |  -0.76 | [ -3.22,  1.04]
#> conditional   | persons | 1     | xb        |   1.57 | [  1.17,  1.99]
#> conditional   | persons | 2     | Intercept |  -0.40 | [ -2.78,  1.37]
#> conditional   | persons | 2     | xb        |   1.30 | [  1.10,  1.48]
#> conditional   | persons | 3     | Intercept |  -0.57 | [ -3.06,  1.23]
#> conditional   | persons | 3     | xb        |   1.32 | [  1.14,  1.51]
#> conditional   | persons | 4     | Intercept |   0.51 | [ -1.82,  2.35]
#> conditional   | persons | 4     | xb        |   0.89 | [  0.79,  0.98]
#> zero_inflated | ID      | 1     | Intercept |  -2.04 | [-16.78,  2.35]
#> zero_inflated | ID      | 1     | nofish1   |  -0.17 | [ -6.88,  3.39]
#> zero_inflated | ID      | 1     | zg        |  -7.23 | [-25.28, -1.96]
#> zero_inflated | ID      | 2     | Intercept |  -3.46 | [-25.28,  0.86]
#> zero_inflated | ID      | 2     | nofish1   |  -0.23 | [ -5.31,  2.32]
#> zero_inflated | ID      | 2     | zg        |  -5.33 | [-16.05, -1.84]
#> zero_inflated | ID      | 3     | Intercept |  -3.46 | [-22.62,  1.15]
#> zero_inflated | ID      | 3     | nofish1   |  -0.14 | [ -6.21,  3.03]
#> zero_inflated | ID      | 3     | zg        |  -7.13 | [-21.74, -2.65]
#> zero_inflated | ID      | 4     | Intercept |  -3.40 | [-21.72,  0.87]
#> zero_inflated | ID      | 4     | nofish1   |  -0.08 | [ -4.47,  4.18]
#> zero_inflated | ID      | 4     | zg        |  -3.74 | [-13.44, -0.83]

estimate_grouplevel(model, type = "total")
#> Component   | Group   | Level | Parameter | Median |          95% CI
#> --------------------------------------------------------------------
#> conditional | ID      | 1     | Intercept |  -0.89 | [ -2.70,  1.48]
#> conditional | ID      | 1     | zg        |   0.16 | [ -0.01,  0.32]
#> conditional | ID      | 2     | Intercept |  -0.81 | [ -2.59,  1.61]
#> conditional | ID      | 2     | zg        |   0.25 | [  0.09,  0.41]
#> conditional | ID      | 3     | Intercept |  -0.53 | [ -2.28,  1.86]
#> conditional | ID      | 3     | zg        |   0.08 | [ -0.05,  0.23]
#> conditional | ID      | 4     | Intercept |  -1.32 | [ -3.15,  1.09]
#> conditional | ID      | 4     | zg        |   0.37 | [  0.28,  0.45]
#> conditional | persons | 1     | Intercept |  -1.45 | [ -2.56, -0.41]
#> conditional | persons | 1     | xb        |   1.57 | [  1.17,  1.99]
#> conditional | persons | 2     | Intercept |  -1.07 | [ -2.15, -0.11]
#> conditional | persons | 2     | xb        |   1.30 | [  1.10,  1.48]
#> conditional | persons | 3     | Intercept |  -1.27 | [ -2.33, -0.26]
#> conditional | persons | 3     | xb        |   1.32 | [  1.14,  1.51]
#> conditional | persons | 4     | Intercept |  -0.13 | [ -1.20,  0.74]
#> conditional | persons | 4     | xb        |   0.89 | [  0.79,  0.98]
#> zi          | ID      | 1     | Intercept |  -5.98 | [-22.84,  0.57]
#> zi          | ID      | 1     | zg        |  -7.23 | [-25.28, -1.96]
#> zi          | ID      | 2     | Intercept |  -7.43 | [-29.46, -0.72]
#> zi          | ID      | 2     | nofish1   |  -0.23 | [ -5.31,  2.32]
#> zi          | ID      | 3     | Intercept |  -7.23 | [-28.82,  0.27]
#> zi          | ID      | 3     | zg        |  -7.13 | [-21.74, -2.65]
#> zi          | ID      | 4     | Intercept |  -7.13 | [-30.06,  0.16]
#> zi          | ID      | 4     | nofish1   |  -0.08 | [ -4.47,  4.18]
#> zi          | persons | 1     | Intercept |  -3.56 | [-14.59,  4.05]
#> zi          | persons | 2     | Intercept |  -3.56 | [-14.59,  4.05]
#> zi          | persons | 3     | Intercept |  -3.56 | [-14.59,  4.05]
#> zi          | persons | 4     | Intercept |  -3.56 | [-14.59,  4.05]

^{Created on 2025-03-04 with reprex v2.1.1}

[DominiqueMakowski]

Rather than making a special case in visualization_recipe_internal to accommodate no CIs, how about in visualization_recipe.estimate_grouplevel we add 2 cols CI_low and CI_high equal to the estimate?

[strengejacke]

Rather than making a special case in visualization_recipe_internal to accommodate no CIs, how about in visualization_recipe.estimate_grouplevel we add 2 cols CI_low and CI_high equal to the estimate?

Were just minor changes. You can still refactor, if you like.

[DominiqueMakowski]

It's fine :) thanks

[strengejacke]

Plots are not fully working, we need to split by Component.

[strengejacke]

library(modelbased)

data(mtcars)
model <- lme4::lmer(mpg ~ hp + (1 | carb), data = mtcars)

estimate_grouplevel(model) |> plot()

estimate_grouplevel(model, type = "total") |> plot()
#> Warning: Removed 6 rows containing missing values or values outside the scale range
#> (`geom_segment()`).

model <- insight::download_model("brms_sigma_3")

estimate_grouplevel(model) |> plot()

estimate_grouplevel(model, type = "total") |> plot()

^{Created on 2025-03-04 with reprex v2.1.1}

[strengejacke]

library(modelbased)
model <- insight::download_model("brms_zi_4")

estimate_grouplevel(model) |> plot()

estimate_grouplevel(model, type = "total") |> plot()

^{Created on 2025-03-04 with reprex v2.1.1}

[strengejacke]

Ok, plotting should work.
Now let's find the 1e+8 edge-cases for brms ;-)

Still to do: the final decision on the argument option for type (here) resp. effects (in parameters).

_{I still prefer "total" because I see no need to add a "random_" prefix... 😬}

[DominiqueMakowski]

I still prefer "total" because I see no need to add a "random_" prefix... 😬

Yeah in general I'm not a fan of underscores, but in this case it's warranted:

• "total" is misleading: people might assume that it means "everything". It is prone of being interpreted outside the context of random effects
• "total" doesn't imply that only random effects will be returned
• In parameters, this feature is really an "extension" / modification of what "random" does, hence in this very specific case the underscores makes sense, it's like I want "random" using a non-default modified "_" to get the total random

That said, in modelbased, because the function is already only about random effects, it can be "total" as it was before, so that there are no breaking changes.