Help me to understand, how to build a simple hierarchical model. Suppose we have a categorical variable, and we want to model Intercept and Slope for each level. Since level “b” has almost no data available, we would want the model to pull information for this level from the global means.
library(brms)
dd = data.frame(y = numeric(), x = numeric(), c = character())
dd = rbind(dd,
{thisX = rnorm(1000); data.frame(y = -4 + 1 * thisX, x = thisX, c = "a")},
{thisX = rnorm(10); data.frame(y = 2 + 1 * thisX, x = thisX, c = "b")},
{thisX = rnorm(1000); data.frame(y = 4 + 3 * thisX, x = thisX, c = "c")}
)
model <- brm(
y~ 1 + x + (1 + x | c),
data = dd,
family = gaussian(),
prior= c(
prior(normal(0, 0.1), class = "b"),
prior(gamma(1, 10), class = "sd")
)
)
fixef(model)
ranef(model)
Expectation:
Global Intercept ~0 Slope ~2
Global+a: Intercept =-4 Slope =1
Global+b: Intercept ~0.1 Slope ~1.9 (because not enough data is available)
Global+c: Intercept =4 Slope =3
What am I missing? I tried setting strong priors on the random sd, and tried replacing random component with the cross-effect, both to no avail. Thank you in advance for any hints.
Recognized by R Language Collective
4