手动计算逻辑回归中的logLik

Question

我使用一个名为GMMAT的R软件包（函数：glmmkin()）运行混合模型逻辑回归，并用遗传关系矩阵调整模型。模型输出包括以下内容（摘自用户手册）：

• theta：离散参数估计[1]和方差分量参数估计[2]
• coefficients：固定效应参数估计（包括截距）。
• linear.predictors：线性预测值。
• fitted.values：原始尺度上的拟合均值。
• Y：与最终工作向量长度相等的向量。
• P：维数等于样本大小的投影矩阵。
• residuals：原始尺度上的残差。不受离散参数重新缩放。
• cov：固定效应（包括截距）的协方差矩阵。
• converged：收敛的逻辑指示器。

我试图获取对数似然以计算方差解释程度。我的第一反应是拆开logLik.glm函数以手动计算AIC，但我卡在了尝试计算AIC时。我使用了这里的答案。

我用stats::glm()运行了逻辑回归进行了健全性检查，在那里model1$aic为4013.232，但是使用我找到的Stack Overflow答案，我得到了30613.03。

我的问题是：有人知道如何使用上述输出手动计算R中的逻辑回归对数似然？

- mkv8

Answer 1

这里没有统计分析，只是我从查看glm.fit中看到的解决方案。这仅适用于在拟合模型时未指定权重的情况（或者如果您这样做了，则需要在模型对象中包含这些权重）。

get_logLik <- function(s_model, family = binomial(logit)) {
    n <- length(s_model$y)
    wt <- rep(1, n) # or s_model$prior_weights if field exists
    deviance <- sum(family$dev.resids(s_model$y, s_model$fitted.values, wt))
    mod_rank <- sum(!is.na(s_model$coefficients)) # or s_model$rank if field exists

    aic <- family$aic(s_model$y, rep(1, n), s_model$fitted.values, wt, deviance) + 2 * mod_rank
    log_lik <- mod_rank - aic/2
    return(log_lik)
}

例如...

model <- glm(vs ~ mpg, mtcars, family = binomial(logit))
logLik(model)
# 'log Lik.' -12.76667 (df=2)

sparse_model <- model[c("theta", "coefficients", "linear.predictors", "fitted.values", "y", "P", "residuals", "cov", "converged")]
get_logLik(sparse_model)
#[1] -12.76667

- Chrisss