quasi-Poisson glm模型的Akaike criterion(AIC)计算问题咨询

代码正确性说明 #

你编写的函数核心逻辑是针对准泊松模型计算QAIC（准AIC），但存在一处方法层面的疏漏：

• 对数似然部分：直接用dpois计算泊松对数似然是合理的，准泊松的均值结构和泊松完全一致，准似然和泊松似然仅相差一个和离散参数φ相关的常数项，不会影响不同模型间的相对比较结果，这部分逻辑无误。
• 惩罚项与离散参数部分：你当前用每个模型自身的summary(model)$dispersion作为离散参数φ，是后续出现异常结果的核心原因。

AIC随自由度上升反而下降的原因 #

这种情况的成因分为两种，其中第一种属于计算逻辑错误：

• 错误成因（最常见）：QAIC的应用前提是所有参与比较的候选模型使用统一的离散参数φ，一般选择所有候选模型中参数最多、最饱和的模型估计得到的φ作为统一值。如果你每个模型单独计算自身的φ，当你增加模型自由度（新增变量）时，模型拟合效果提升会导致残差减小，估计得到的φ会随之降低，可能出现φ下降的幅度超过自由度p上升的幅度，最终让惩罚项2*p*φ不升反降，带动整体QAIC下降。
• 正常成因：如果新增变量对响应变量的解释力极强，对数似然上升的幅度远超过惩罚项的上升幅度，此时QAIC下降属于正常情况，说明新增变量确实提升了模型的拟合效率。

修正方案 #

修正后的计算逻辑和代码 #

1. 先拟合所有待比较的候选模型，选择参数最多的饱和模型计算统一的离散参数φ
2. 所有模型计算QAIC时都使用这个统一的φ，代码参考如下：

# 第一步：用最复杂的饱和模型计算统一离散参数
phi_unified <- summary(你的饱和模型名)$dispersion

# 第二步：修正后的QAIC计算函数
fqaic <- function(model, unified_phi) {
  # 计算泊松对数似然
  loglik_pois <- sum(dpois(model$y, model$fitted.values, log = TRUE))
  # 获取模型拟合的参数个数（含截距）
  param_num <- length(coef(model))
  # 计算QAIC
  qaic <- -2 * loglik_pois + 2 * param_num * unified_phi
  return(qaic)
}

# 调用示例：所有模型传入同一个统一phi
# qaic_model1 <- fqaic(model1, phi_unified)
# qaic_model2 <- fqaic(model2, phi_unified)

额外注意事项 #

• 如果你的样本量较小（样本量/参数个数 <40），建议使用带小样本修正的QAICc，在原QAIC基础上增加修正项即可：2 * param_num * (param_num + 1) / (样本量 - param_num - 1) * unified_phi
• QAIC的绝对值没有实际意义，仅需要比较不同模型之间的QAIC差值即可，一般差值小于2的模型可以认为拟合效果相当

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