逻辑回归 - 在 R 与 Stata 中获得 Pearson 标准化残差

机器算法验证 r 状态 残差 物流
2022-04-14 12:39:02

我正在处理涉及逻辑回归模型的任务,我需要根据其中一个预测变量绘制 pearson 标准化残差。这是基本设置:

model <- glm(outcome ~ predictor1 + predictor2, family=binomial(logit))
res <- residuals(model, "pearson")

在查看残差分布时,我发现与使用 Stata 的同事(使用 predict 和 rstandard)完全不同。它们的残差或多或少是正常的,而在我的值中存在差距(不是单个残差在 -0.05 和 1.15 之间)。这在逻辑回归的背景下确实有意义,尤其是最大预测概率不是那么高(38%)。

我想了解这里发生了什么……Stata 做了什么而 R 没有做这些残差?

2个回答

对于逻辑回归,Stata 将残差和相关数量定义为如果您将所有观测值与所有预测变量的相同值分组,计算这些观测值的成功和失败,并将逻辑回归模型拟合到结果二项式数据而不是原始伯努利数据。这是一件很有用的事情,因为(如果有多个观测值具有相同的协变量模式)结果残差的行为更像你习惯于最小二乘法的那些。

为了从 R 中获得相同的残差,我怀疑您需要对数据进行分组并将模型拟合到分组数据。但我不清楚 R 是否使用与 Stata 相同的“标准化残差”定义,因为我目前无法访问 R 文档引用的众多教科书。

这是Stata手册条目中“逻辑后估计”的“方法和公式”部分的摘录(我喜欢Stata的一件事是手册提供了所有内容的完整公式):

定义的协变量模式的观察总数。定义的协变量模式的观察中的积极响应总数。MjjYjj

个观察的 Pearson 残差定义为 ...帽子矩阵未调整对角线元素是由 ,其中是估计的参数协方差矩阵。调整后的对角线元素创建的然后是 标准化的 Pearson 残差j

rj=YjMjpjMjpj(1pj)

hUjhUj=(XVX)jjVhjhathj=Mjpj(1pj)hUj
rSjrj/1hj.

Pearson 残差是通过将每个观察的原始残差除以相应方差的平方根来获得的。这个想法是得到大约方差为 1 的东西。在你的例子中,试试这个;

set.seed(3141)
x1 <- rnorm(100)
x2 <- rnorm(100)
y <- rbinom(100, 1, 0.25)
glm1 <- glm(y~x1+x2, family=binomial)
f1 <- fitted(glm1) # the fitted probability of y=1, for each observation
plot( residuals(glm1, "pearson"), (y-f1)/sqrt(f1*(1-f1)))
abline(0,1)        # they match

出现“差距”是因为的残差位于一侧,而的残差位于另一侧。标准化残差是另一种动物;它们除以估计的残差标准差;您可以使用 R 在 R 中获得它们,但对于非线性 GLM,它在计算中使用线性近似。Y=1Y=0rstandard()

任何形式的 NB 残差在逻辑回归中往往没有太大帮助。对于独立的二元数据,唯一真正关心的是我们是否正确指定了平均值——并且在样本量适中的情况下,残差图通常无法提供评估的能力。