1
votes

I'm fitting a y ~ v + m + s + m:s + (1|subunit) model with lmer(). s is a continuous variable interacting with m, a categorical factor with 3 levels: A, B, and C.

Fitting the model uses A as the reference level for factor m:

fit_ref_A <- lmer(y ~ 1 + v + m*s + (1|subunit), data=df)

The parameter estimates for B and C can just be calculated using the estimate for A and the differences for B and C. I'm interested in extracting the confidence intervals.

Running confint() gives the confidence interval for the slope of variable s at A. I'm also interested in the confidence intervals of the slopes of s at B and C, not just the confidence intervals for the differences from the slope at A. Is there a way to extract this from fit_ref_A? So far the only thing I've been able to figure out is to relevel with B as the reference, fit a new fit, then relevel with C as the reference, and fit the third fit.

Question: Is there a way to extract everything (especially the confidence intervals) from fit_ref_A?

Code:

library(lme4)
# create the dataset, unbalanced at the lowest stratum ( 2 repeats for m==A instead of 3)
set.seed(2)
s_levels <- 1:5
m_levels <- c("A", "B", "C")
v_levels <- c("L2", "L3", "L4")
reps <- 1:3
df <- expand.grid(rep=reps, s=s_levels, m=m_levels, v=v_levels)
df$subunit <- as.factor(paste(df$v,"-",df$m,"-",df$s, sep=""))
df$y <- rnorm(nrow(df), 0, 1)
df <- subset(df, !(rep==3 & m=="A"))  # drop the 3rd repeat for m=="A"
table(df$m)  # shows 30 for A, 45 for B, 45 for C as expected

# fit 3 different models, with three different reference levels for 'm'
fit_ref_A <- lmer(y ~ 1 + v + m*s + (1|subunit), data=df)

df$m <- relevel(df$m, ref = "B")
fit_ref_B <- lmer(y ~ 1 + v + m*s + (1|subunit), data=df)

df$m <- relevel(df$m, ref = "C")
fit_ref_C <- lmer(y ~ 1 + v + m*s + (1|subunit), data=df)

# Calculate the confidence intervals for the continuous variable s at the three
# different levels for categorical factor m.  Must use 3 separate fits.

cis_at_m_reference_A <- confint(fit_ref_A)
cis_at_m_reference_B <- confint(fit_ref_B)
cis_at_m_reference_C <- confint(fit_ref_C)

cis_at_m_reference_A["s",]
cis_at_m_reference_B["s",]
cis_at_m_reference_C["s",]

# Any way to just extract all three from fit_ref_A?
1
The slope at B is the s coef + the (m == B):s coef, and similarly for the other values. You don't need to relevel and re-fit the model, you can add the coefficients you have. - Gregor Thomas
@GregorThomas Right. The point estimates are easy. I'm interested in the confidence limits, which cannot be gotten in the same easy way. I've edited my question to make it more clear that I want to extract the confidence limits. - MichiganWater

1 Answers

2
votes

You can get approximate CIs with Gaussian error propagation:

sum(fixef(fit_ref_A)[c("s", "mB:s")]) + 
  c(-1.96, 1.96) * sqrt(sum(vcov(fit_ref_A)[c("s", "mB:s"), c("s", "mB:s")]))
#[1] -0.3346310  0.1863014

Or you could bootstrap:

myboot <- bootMer(fit_ref_A, function(x) {
  cf <- fixef(x)
  c(sA = cf[["s"]], sB = cf[["s"]] + cf[["mB:s"]], sC = cf[["s"]] + cf[["mC:s"]])
}, nsim = 1e4, seed = 42)

apply(myboot$t, 2, quantile, probs = c(0.025, 0.975))
#              sA         sB         sC
#2.5%  -0.4022927 -0.3415690 -0.3969831
#97.5%  0.2041610  0.1858731  0.1266355