mgcv拟合GAM随机效应：估计值提取与效应图解读方法

开发者阿泽

2026-8-28

你拟合模型的参考代码如下：

gam_fit <- gam(y ~ s(age) + s(region, bs='re'), 
               data = my_data, 
               method = 'REML')

mgcv本身不需要依赖第三方包就可以提取随机效应估计值，只是没有直接提供和lme4同名的ranef()快捷泛型，常用两种原生实现方式：

直接从模型系数提取：bs='re'定义的随机效应项本质是带岭惩罚的特殊光滑项，对应位置的模型系数就是各分组水平的随机效应BLUP（最佳线性无偏预测）值。先定位随机效应光滑项的位置，再按参数索引切片即可，返回值顺序和分组变量的因子水平顺序完全一致，参考代码：
```
# 定位s(region)光滑项的索引
re_loc <- which(vapply(gam_fit$smooth, `[[`, character(1), "label") == "s(region)")
# 提取对应位置的随机效应值
re_est <- gam_fit$coefficients[
  gam_fit$smooth[[re_loc]]$first.para : gam_fit$smooth[[re_loc]]$last.para
]
# 匹配分组水平名称
names(re_est) <- levels(my_data$region)
```
按项预测提取：如果需要拿到原始数据中每一行观测对应的随机效应预测值，直接调用预测函数指定对应项即可，加se.fit=TRUE还能同步返回随机效应的标准误：
```
re_pred <- predict(gam_fit, type = "terms", terms = "s(region)", se.fit = TRUE)
```

对上述gam_fit调用默认plot()方法时，连续光滑项s(age)会展示自变量与效应值的拟合曲线，而随机效应项s(region)不会输出效应变化曲线，会自动生成随机效应估计值的正态Q-Q诊断图：

横轴是理论正态分布（高斯分布）的分位数，纵轴是各分组水平的随机效应估计值。
这是纯粹的模型诊断图，核心作用是验证随机效应的正态性假设：如果图中所有点大致沿对角参考线分布，说明随机效应符合模型预设的正态分布假设；如果点出现明显的S形弯曲、尾部大量偏离参考线的离群点，说明正态性假设不成立，需要考虑调整随机效应结构、检查异常分组、或对响应/分组变量做预处理。
注意不要用这个图解读随机效应的实际大小或方向，它不具备效应解释功能，仅做诊断用。如果需要展示各分组随机效应的具体取值，用问题1中提取到的re_est单独绘制点图、柱状图即可。

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