如何为felm回归模型计算AIC准则？自定义函数返回空数值求助

你使用felm包构建回归模型后，自定义的AIC_adj函数返回空数值（numeric(empty)），原函数如下：

AIC_adj <- function(mod){
  # Number of observations
  n.N   <- nrow(mod$model)
  # Residuals vector
  u.hat <- residuals(mod)
  # Variance estimation
  s.sq  <- log( (sum(u.hat^2)/(n.N)))
  # Number of parameters (incl. constant) + one additional for variance estimation
  p     <-  length(coef(mod)) + 1
  
  # Note: minus sign cancels in log likelihood
  aic <- 2*p  +  n.N * (  log(2*pi) + s.sq  + 1 ) 
  
  return(aic)
}

模型系数结果：

G         D1         D2         D4         D5         D6   posratio     sd_ret       G:D1 
 0.0832205 -0.0030118 -0.0006392  0.0007287  0.0004233 -0.0013440  0.0059718  0.2344416  0.0087144 
      G:D2       G:D4       G:D5       G:D6 G:posratio   G:sd_ret 
 0.0046054 -0.0028214 -0.0033844 -0.0076368 -0.0881973  0.6204865 

问题原因 #

原函数返回空数值的核心原因是：felm对象的$model属性并非总是包含完整的原始观测数据，当模型包含固定效应、聚类分组等设定时，mod$model可能被处理为仅包含用于估计的变量子集，甚至为空，导致nrow(mod$model)返回空值，后续所有计算因此失效。

修正后的AIC计算函数 #

针对felm对象的结构特点，修改后的函数使用专门的观测数获取函数nobs()，并确保残差提取的正确性：

AIC_felm <- function(mod){
  # 获取模型实际使用的观测数量
  n.N <- nobs(mod)
  # 提取响应变量的原始残差
  u.hat <- residuals(mod, type = "response")
  # 计算残差方差的对数
  s.sq <- log(sum(u.hat^2) / n.N)
  # 参数个数：估计的回归系数 + 1个方差参数
  p <- length(coef(mod)) + 1
  
  # 基于高斯似然的AIC计算
  aic <- 2*p + n.N * (log(2*pi) + s.sq + 1)
  return(aic)
}

额外说明 #

• 对于包含固定效应的felm模型，AIC计算中仅计入显式估计的回归系数（即coef(mod)返回的结果），固定效应属于吸收的组内差异，不作为自由参数计入AIC。
• 如果模型使用了聚类标准误，AIC的计算逻辑不受影响，因为聚类标准误仅改变系数的标准误估计，不改变似然函数的结构。

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