gammit/mgcv随机效应预测异常：结果仅匹配组随机截距

问题概述 #

长期使用gammit::predict_gamm()（基于mgcv::predict.gam()）执行GAMM/GLMM模型预测时，遇到以下异常：

• 指定re_form = c("s(Country_fact)")纳入目标随机效应时，预测值完全等于对应Country_fact组的随机截距系数，完全不随测试集中exp_A、exp_B等其他变量变化；
• 指定re_form = NA不纳入任何随机效应时，预测值能正常随自变量变化；
• 直接使用mgcv::predict.gam()时，无论是否通过re.form指定随机效应，预测结果均一致，且与gammit::predict_gamm()不带随机效应的结果匹配。

复现代码 #

# 创建训练数据集
set.seed(123)
dat <- tibble(
  Reference_fact = as.factor(rep(seq(1,10,1),5)),
  Country_fact = as.factor(rep(stri_rand_strings(10,2),5)),
  exp_A = runif(50,0,500),
  exp_B = runif(50,0,20),
  resp = exp_A^2 + 0.1*exp_B * 0.3*exp_A*exp_B + rnorm(50, 0, 10)
)

# 定义模型公式
model_spec <- c("resp ~ 
                 s(exp_A) +
                 s(exp_B) +
                 s(Reference_fact, bs = 're') + 
                 s(Country_fact, bs = 're')"
)

# 拟合GAMM模型
fit <- mgcv::gam(formula(str_replace_all(model_spec, "[\r\n]", "")),
                 method = 'REML',
                 family = 'gaussian',
                 data=dat)

# 创建测试数据集
set.seed(123)
dat_pred <- tibble(
  Reference_fact = as.factor(rep(seq(1,10,1),5)),
  Country_fact = as.factor(rep(stri_rand_strings(10,2),5)), 
  exp_A = runif(50,0,500),
  exp_B = runif(50,0,20)
)

# 带随机效应的预测（异常结果）
p <- gammit::predict_gamm(
  fit, 
  dat_pred, 
  re_form = c("s(Country_fact)"),
  keep_prediction_data = TRUE, 
  newdata.guaranteed = FALSE, 
  se.fit = FALSE) 

原因分析 #

1. gammit::predict_gamm()的参数逻辑缺陷
   该函数的re_form参数设计与mgcv::predict.gam()的re.form存在差异。mgcv中re.form用于指定要保留的随机效应（支持公式或字符向量格式），但gammit在处理re_form = c("s(Country_fact)")时，仅提取了该随机效应的截距值，完全忽略了固定效应部分（s(exp_A)、s(exp_B)等平滑项）的贡献，导致预测值仅等于随机截距。

2. mgcv::predict.gam()的参数使用误区
   你测试中mgcv结果一致，大概率是因为re.form参数传递方式错误。mgcv中指定保留特定随机效应需使用公式形式（如re.form = ~s(Country_fact)），而非字符向量；若未正确指定，函数会默认包含所有随机效应或按默认逻辑处理，导致结果无差异。

3. gammit包的潜在bug
   从复现结果看，当指定具体随机效应项时，函数未将固定效应预测值与随机效应值相加，这属于gammit::predict_gamm()内部计算逻辑的错误，导致输出仅为随机效应部分。

解决方案 #

• 优先使用mgcv原生预测函数
  若无需gammit的额外功能，直接调用mgcv::predict.gam()，注意正确设置re.form：

  • 包含所有随机效应：默认不指定re.form
  • 仅包含Country_fact的随机效应：re.form = ~s(Country_fact)
  • 排除所有随机效应：re.form = NA

• 调整gammit的参数格式
  若必须使用gammit，尝试将re_form改为公式形式~s(Country_fact)而非字符向量，验证是否能正常结合固定效应与随机效应。若问题仍存在，建议向gammit包的开发者提交issue报告该bug。

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