在glmmTMB中建模不规则时间自相关时如何绘制随机效应

问题说明 #

当使用glmmTMB的ou(times + 0 | YearF)结构为不规则时间序列建模自相关后，调用ranef()或plot_model()会返回每个times水平对应的随机效应，而非我们需要的按年份分组的随机截距。这是因为OU自相关结构本质是为每个组内的时间点构建了相关的随机过程，而非单一的组级截距。

解决方案 #

方案1：调整模型结构，同时包含随机截距与自相关过程 #

直接在模型中同时定义组水平的随机截距和组内的OU自相关过程，这样就能单独提取并绘制组级随机截距：

# 重新构建模型，同时包含组随机截距和OU自相关
mod3 <- glmmTMB(N ~ DayYear + YearF + 
                  (1 | YearF) +  # 年份组的随机截距
                  ou(times + 0 | YearF),  # 组内时间自相关结构
                family = nbinom2,
                data = df3)

# 提取组水平随机截距
glmmTMB::ranef(mod3, condVar = T)$cond$YearF
broom.mixed::tidy(mod3, effects = "ran_vals")

# 绘制随机截距图
library(sjPlot)
plot_model(mod3, type = "re", terms = "(1|YearF)")

这种方法的优势是模型结构清晰，明确分离了组间截距差异和组内时间自相关，结果解释性更强。

方案2：从OU模型的随机效应推导组水平截距 #

如果无法重新建模，可以通过计算每个组内所有时间点随机效应的均值，近似得到组级截距，再手动绘图：

# 提取OU模型的所有随机效应
re_ou <- ranef(mod2, condVar = T)$cond$YearF

# 转换为数据框并拆分组与时间信息
re_df <- data.frame(
  group_time = rownames(re_ou),
  re_value = re_ou[, 1],
  se = sqrt(diag(attr(re_ou, "postVar")))
) %>%
  tidyr::separate(group_time, into = c("YearF", "times"), sep = ":")

# 按年份分组计算平均随机效应与标准误
group_re <- re_df %>%
  dplyr::group_by(YearF) %>%
  dplyr::summarise(
    avg_re = mean(re_value),
    avg_se = sqrt(mean(se^2))  # 近似组水平标准误
  )

# 手动绘制随机截距图
library(ggplot2)
ggplot(group_re, aes(x = YearF, y = avg_re)) +
  geom_point(size = 3) +
  geom_errorbar(aes(ymin = avg_re - 1.96*avg_se, ymax = avg_re + 1.96*avg_se), width = 0.2) +
  labs(title = "OU模型组水平平均随机截距", y = "平均随机截距", x = "年份") +
  theme_minimal()

关键提示 #

• ou(times + 0 | group)结构的核心是为组内每个时间点建模相关随机过程，每个时间点对应一个随机效应，这和单纯的组随机截距是完全不同的模型组件。
• 优先选择方案1，它的模型设定更符合“组截距+自相关”的需求，结果的统计意义更明确；方案2属于近似方法，仅适用于无法调整模型的场景。

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