r - 访问 r 包 `gamm4` 生成的广义加法混合模型 (GAMM) 的 AIC

Question

我正在使用 r 软件包 gamm4对具有八个预测变量的二元响应变量运行 gamm，其中两个使用薄板样条进行平滑。该模型运行良好，我可以验证它已从 summery() 生成结果

这是我的数据的子集: structure(list(DATE = c("1/1/1990", "1/2/1990", "1/3/1990", "1/4/1990", "1/5/1990", "1/6/1990"), RESPONSE = c(0L, 0L, 0L, 0L, 0L, 0L), PREDICTOR1 = c(100L, 80L, 60L, 100L, 100L, 100L), PREDICTOR2 = c(5.3, 3.8, 2.5, 2.2, 2.8, 4), PREDICTOR3 = c(0.001016002, 0, 0, 0.001778004, 0.021590043, 0), PREDICTOR4 = c(315, 90, 90, 315, 360, 225), PREDICTOR5 = c(1019.326087, 1028.770833, 1027.920833, 1024.4625, 1022.25, 1018.891667), PREDICTOR6 = c(17.36086957, 12.1, 13.02083333, 16.5375, 16.48333333, 19.6375), PREDICTOR7 = c(13.47916667, 12.775, 12.525, 12.3875, 12.725, 13.09166667), PREDICTOR8 = c(0.28335247, 0.212789468, 0.239630838, 0.386550722, 0.221811468, 0.379079404), random_effect1 = c(1L, 1L, 1L, 2L, 2L, 2L),

这是我的模型公式:

model_1 <- gamm(RESPONSE ~ PREDICTOR1 + PREDICTOR2 + PREDICTOR3 + s(PREDICTOR4) + 
           PREDICTOR5 + s(PREDICTOR6) + PREDICTOR7 + PREDICTOR8, data = data, 
           random = random_effect1, family = binomial())

但是，当我尝试在 gamm 结果对象上运行 AIC 时，它失败了

aic_model1 <- AIC(model_1)
Error in UseMethod("logLik") :
no applicable method for 'logLik' applied to an object of class
 "c('gamm', 'list')"

我尝试分解 AIC 部分并使用 predict.gam 计算拟合值和

deviance_model <- 
  -2 * sum(log(dbinom(model_1$RESPONSE, size = 1, prob = fitted_values)))

这两个都会生成 NULL结果。 我也尝试过AIC.gam和AIC.gamm这些失败了。 问题(尽管我可能弄错了)似乎是对象的类 gamm, list 。有谁知道如何解决这个问题？提前感谢您的帮助!

Answer 1

首先，您没有在代码中使用gamm4()。包 gamm4 提供了一个名为 gamm4() 的函数，而不是 gamm()。

如果您愿意使用混合模型公式中的 AIC(请阅读 ?mgcv::logLik.gam，我在下面显示的内容将产生边际可能性 -基于 AIC，它整合了随机效应，因此 [波动位] 变得平滑)，那么您可以将 AIC() 函数应用于由 gamm4()

返回的对象

# use example from ?gamm4::gamm4
set.seed(0) 
dat <- gamSim(1,n=400,scale=2) ## simulate 4 term additive truth

## Now add 20 level random effect `fac'...
dat$fac <- fac <- as.factor(sample(1:20,400,replace=TRUE))
dat$y <- dat$y + model.matrix(~fac-1)%*%rnorm(20)*.5
br <- gamm4(y~s(x0)+x1+s(x2),data=dat,random=~(1|fac))

现在将 AIC() 应用于 br$mer，以便混合模型类的 logLik() 方法 ("merMod") 用于:

> AIC(br$mer)                                                                
[1] 1811.787

您无法使用具有非高斯族的 gamm() 模型来执行此操作，因为这些模型实际上是通过 MASS::glmmPQL() 拟合的，并且该函数不会最大化可能性(它使用惩罚拟似然)，因此不存在可能性，也没有 AIC。

关于返回的对象包含多个模型拟合的要点仍然存在；这些函数返回混合模型形式作为一个组件，并返回 GAM 形式作为第二个组件。使用 str(my_model, max = 1) 查看这两个组件的名称，因为它们在模型拟合函数之间有所不同。