Skip to main content
edited tags
Link
kjetil b halvorsen
  • 82.8k
  • 32
  • 201
  • 663
Source Link

Is it possible to extract confidence intervals for non-reference levels in a lm() model in R? Is there a multiple comparison issue?

I've got a dataset with two categorical predictors (v and m) and a continuous predictor (s). I'm modeling the response with an interaction between the continuous predictor and one of the categorical predictors. y ~ v + m + s + m:s I'd like to get the confidence intervals for the slope of s for each level of m. summary(fit) gives the slope at the reference level of m (A), but for levels B and C it gives the differences from the slope at A. Same thing for confint(fit). I know I can just add the estimate of the difference between B and A to get B, but I don't know what to do for the confidence intervals other than relevel and re-fit the data.

Is there a (relatively easy?) way to get the parameter estimates and confidence intervals from a fit, rather than re-fitting using a new reference for each level of a parameter?

Given that I do extract the confidence intervals, is there any issue with multiple-comparisons and having to correct? Honestly, I'm more or less considering these confidence intervals as approximations to credible intervals (flat priors on the coefficients and all that), but I do want to know if the widths of the confidence intervals would be affected by what I intend to do.

Example code below.

# How to get confidence intervals for the parameters, not just contrasts

set.seed(123)

s_avg <- 1:5
m_avg <- c("A", "B", "C")  # three years
v_avg <- c("L2", "L3", "L4") # three fields per year
df <- expand.grid(s=s_avg, m=m_avg, v=v_avg)
df$y <- as.numeric(as.factor(df$v)) + rnorm(nrow(df), mean=0, sd=1)

fit <- lm(
  formula = y ~ v + m + s + s:m ,
  data=df)

# The "s" estimate and CIs apply to the slope under modification A
summary(fit)
confint(fit)

# I can manually calculate get the slope of modification B by adding the mB:s
# estimate. How do I get the confidence intervals?  I can do it by changing the
# reference level and re-fitting the model, but I'm guessing there's a way to
# just extract the estimates and confidence intervals for the parameters (not
# the contrasts) without re-fitting. Since s is a continuous predictor, I can't
# just use level means coding, as far as I know.

# fit in order to get the confidence interval for the slpoe of s at modification B
df$m <- relevel(df$m, ref = "B")

fit2 <- lm(
  formula = y ~ v + m + s + s:m ,
  data=df)

summary(fit2)
confint(fit2)


# fit in order to get the confidence interval for the slpoe of s at modification C
df$m <- relevel(df$m, ref = "C")

fit3 <- lm(
  formula = y ~ v + m + s + s:m ,
  data=df)

summary(fit3)
confint(fit3)


# Plot the data and model -------------------------------------------------

library(ggplot2)
library(sjPlot)

ggplot(data = df ,
       aes(x = s, y = y)) +
  geom_point () + 
  facet_wrap( ~ v + m)

plot_model(
  model = fit ,
  type="pred" ,
  terms=c("s", "m", "v") ,
  ci.lvl = 0.95
)