如何在glmer模型及plot_model绘图中强制显示参考水平？

完全可以实现，核心思路是手动构造参考水平的统计数据，再与plot_model生成的绘图数据合并后重新绘图。以下分两种场景给出具体方案：

一、针对mtcars示例的实现代码 #

先加载所需包，运行模型后提取绘图数据，构造参考水平行再合并绘图：

library(sjPlot)
library(ggplot2)
library(dplyr)

# 运行原模型
carmodel <- glm(am ~ as.factor(cyl) + as.factor(gear) + as.factor(carb), 
                data = mtcars, family = binomial)

# 获取plot_model的原始绘图数据
plot_data <- plot_model(carmodel, type = "forest", return_data = TRUE)

# 构造参考水平的数据（OR=1，置信区间固定为[1,1]）
ref_data <- tibble(
  term = c("cyl: 4", "gear: 3", "carb: 1"),
  predicted = 1,
  conf.low = 1,
  conf.high = 1,
  group = c("as.factor(cyl)", "as.factor(gear)", "as.factor(carb)"),
  p.value = NA,
  std.error = 0
)

# 合并数据并设置因子顺序，确保参考水平在每个变量组最前面
combined_data <- bind_rows(ref_data, plot_data) %>%
  mutate(term = factor(term, levels = c(
    "cyl: 4", "cyl: 6", "cyl: 8",
    "gear: 3", "gear: 4", "gear: 5",
    "carb: 1", "carb: 2", "carb: 3", "carb: 4", "carb: 6", "carb: 8"
  )))

# 绘制带参考水平的森林图（风格对齐plot_model）
ggplot(combined_data, aes(x = predicted, y = term)) +
  geom_vline(xintercept = 1, linetype = "dashed", color = "gray50") +
  geom_point(aes(color = group), size = 3) +
  geom_errorbarh(aes(xmin = conf.low, xmax = conf.high, color = group), height = 0.2) +
  scale_x_log10(limits = c(0.1, 10)) +
  labs(title = "Odds Ratios for Automatic Transmission", 
       x = "Odds Ratio (OR)", y = "") +
  theme_minimal() +
  theme(legend.position = "bottom")

二、针对你的glmer分层模型的适配方案 #

逻辑和上面一致，只需根据你模型中变量的实际参考水平调整ref_data的内容：

library(sjPlot)
library(ggplot2)
library(dplyr)
library(lme4)

# 运行原分层模型
hmdl <- glmer(asww ~ bwo + dx_year + psadx + agedec + (1 | practiceid/npi),
              data = mdata, family = binomial(link = "logit"), 
              control = glmerControl(optimizer = "Nelder_Mead"))

# 获取plot_model的原始绘图数据
plot_data <- plot_model(hmdl, type = "forest", return_data = TRUE)

# 构造参考水平数据（替换为你模型中各分类变量的实际参考组）
ref_data <- tibble(
  term = c("bwo: 参考值", "dx_year: 参考值", "psadx: 参考值", "agedec: 参考值"),
  predicted = 1,
  conf.low = 1,
  conf.high = 1,
  group = c("bwo", "dx_year", "psadx", "agedec"),  # 对应模型中的变量名
  p.value = NA,
  std.error = 0
)

# 合并数据并设置因子顺序（确保参考水平在每组最前）
combined_data <- bind_rows(ref_data, plot_data) %>%
  mutate(term = factor(term, levels = c(
    "bwo: 参考值", "bwo: 其他值1", "bwo: 其他值2",  # 按实际变量水平补充
    "dx_year: 参考值", "dx_year: 其他值1", ...,
    psadx和agedec的对应水平
  )))

# 绘制带参考水平的森林图，保留你原plot_model的参数设置
ggplot(combined_data, aes(x = predicted, y = term)) +
  geom_vline(xintercept = 1, linetype = "dashed", color = "gray50") +
  geom_point(aes(color = group), size = 3) +
  geom_errorbarh(aes(xmin = conf.low, xmax = conf.high, color = group), height = 0.2) +
  scale_x_continuous(limits = c(0.7, 4)) +  # 对应你原axis.lim设置
  labs(title = "Predictors of active surveillance", 
       x = "Odds Ratio (OR)", y = "") +
  scale_color_brewer(palette = "Set1") +  # 对应你原colors设置
  theme_minimal() +
  theme(legend.position = "bottom")

关键说明 #

• plot_model默认不显示参考水平，因为参考组的回归系数为0（OR=exp(0)=1），无估计误差，所以不会出现在默认输出中。
• 手动构造的参考水平数据需与原数据的列名完全匹配，确保合并后绘图正常。
• 因子顺序的设置是为了让参考水平排在每个变量组的最前面，符合森林图的阅读习惯。

内容的提问来源于stack exchange，提问作者dr_coops