how to obtain p-value (or CIs) for correlation of random effects in a GLMM (lme4)

Question

I want to test for a correlation between the random effects of a GLMM model calculated in lme4. I have already been suggested to conduct a likelihood ratio comparison of a model with and without the random correlation. That is indeed significant but I wanted to ask if there is any way to get the confidence intervals or p-values of this correlation from the model.

(Specifically, I have compared a model with the random effects structure (1 + X1 + X2 || group) against (1 + X1 + X2 | group) but the problem is that in the second model also the correlation with the intercept is included and I want to specifically test for the significance of the correlation between X1 and X2. Unfortunately, a model with (1 + X1 | group) + (1 + X2 | group) does not converge)

Any help would be appreciated

Ben Bolker · Answer 1 · 2023-02-20T00:59:44.790

You can use confint() to get likelihood profile confidence intervals. P-values would be harder; you could do parametric bootstrapping but it would be slow.

set.seed(101)
dd <- data.frame(x = rnorm(1000), y = rnorm(1000),
                 g = factor(sample(1:20, size = 1000, replace = TRUE)))
library(lme4)
dd$z <- simulate(~ x + y + (1 + x + y | g),
                 newdata = dd,
                 newparams = list(beta = rep(1, 3),
                                  sigma = 1,
                                  theta = rep(1, 6)))[[1]]

m <- lmer(z ~ x + y + (1 + x + y | g),
          data = dd)

In the confint() call below, parm = "theta_" means "all covariance parameters". You could use parm = c(2, 3, 5) to select only the correlation parameters, but you'd have to read ?profile.merMod and think carefully to figure out the correct indices ...


cc <- confint(m, parm = "theta_", oldNames = FALSE)

Results give you 95% (by default) CIs for all of the covariance parameters. In this example, the x/y slope correlation is significant but the correlations between (intercept and x) and (intercept and y) aren't. (Note that the correlations aren't necessarily invariant to reparameterizing the model, in particular centering or otherwise shifting the predictors will change the answers ...)

cc
                          2.5 %    97.5 %
sd_(Intercept)|g     0.38142602 0.7451417
cor_x.(Intercept)|g -0.15990774 0.6492967
cor_y.(Intercept)|g -0.01148283 0.7294138
sd_x|g               0.67205037 1.2800681
cor_y.x|g            0.53404483 0.9116571
sd_y|g               0.83378353 1.5742580
sigma                0.94201110 1.0311559

how to obtain p-value (or CIs) for correlation of random effects in a GLMM (lme4)

1 Answers1