如何使用 R 中的 Gamma GLM 通过 AIC 进行模型选择？答案

【问题标题】：How can I do model selection by AIC with a Gamma GLM in R?如何使用 R 中的 Gamma GLM 通过 AIC 进行模型选择？
【发布时间】：2012-11-04 12:01:04
【问题描述】：

正如 glm() 的文档所解释的，glm() 返回的值的 aic 组件不是有效的 AIC：

对于高斯、Gamma 和逆高斯离散度是根据残差估计的偏差，参数的数量是系数加一。对于高斯族，MLE 使用了色散，因此这是 AIC 的有效值，但对于伽玛和逆高斯族不是。

因此需要通过其他方式获得有效的 AIC。

【问题讨论】：

标签： r glm model-comparison

【解决方案1】：

如果您想使用 step() 或 MASS::stepAIC() 模型选择函数，您可以首先通过执行以下操作确保正确计算 AIC：

GammaAIC <- function(fit){
  disp <- MASS::gamma.dispersion(fit)
  mu <- fit$fitted.values
  p <- fit$rank
  y <- fit$y
  -2 * sum(dgamma(y, 1/disp, scale = mu * disp, log = TRUE)) + 2 * p
}
GammaAICc <- function(fit){
  val <- logLik(fit)
  p <- attributes(val)$df
  n <- attributes(val)$nobs
  GammaAIC(fit) + 2 * p * (p + 1) / (n - p - 1)      
}

my_extractAIC <- function(fit, scale=0, k=2, ...){
  n <- length(fit$residuals)
  edf <- n - fit$df.residual  
  if (fit$family$family == "Gamma"){
    aic <- GammaAIC(fit)
  } else {
    aic <- fit$aic
  }
  c(edf, aic + (k - 2) * edf)
}
assignInNamespace("extractAIC.glm", my_extractAIC, ns="stats")

如果你使用glmulti包，你可以简单地通过glmulti()的crit参数指定使用上述GammaAIC()或GammaAICc()函数。

【讨论】：