我的意思是：标准误差的 Huber/White/sandwich 估计量。在我看来，在连续结果的情况下，标准误差的稳健估计相当简单，因为每个观察的残差方差是作为回归的平方（估计）残差计算的。但我无法弄清楚这如何适用于二元结果：在估计每个观察的概率之后，我如何估计该观察的残差方差，因为它必然有 0 或 1 的结果，对应于logit标度，到或者? 特别是，我正在考虑我有连续预测变量的情况（所以我的一组回归变量的每个值只出现一次）。
PS：我读到了一些关于使用稳健标准误差进行逻辑回归的批评，因为如果方差的估计有偏差，那么参数估计本身也是有偏差的（在二项式情况下，给定的平均值和方差是相关的）。但是，我不怀疑稳健的标准误差在这里是否有意义，而只是想知道它们是如何计算的。
如果您使用所谓的工作残差并且回归量和残差适当加权，则该定义完全类似。工作权重是用于广义线性模型（包括逻辑回归）的常用迭代加权最小二乘 (IWLS) 算法的一部分。
主要思想是使用对数似然相对于平均值或相应的线性预测变量的导数作为偏差的度量。或者，您可以使用所谓的分数或估计函数作为偏差的度量，即对数似然相对于线性预测变量的系数的导数。所有这些想法都在我们sandwich的 R 包中使用，以实现不同三明治估计器的面向对象的实现。有关详细信息，请参阅 Zeileis（2006，三明治估计器的面向对象计算，doi:10.18637/jss.v016.i09）。
R 中使用逻辑回归进行劳动力市场参与的工作说明是：
data("SwissLabor", package = "AER")
m <- glm(participation ~ ., data = SwissLabor, family = binomial)
可以使用我们的sandwich包计算经典的 Eicker/Huber/White 三明治协方差矩阵：
library("sandwich")
vcov_sandwich <- sandwich(m)
您可以手动轻松提取工作残差和回归量，并使用工作权重适当加权：
r <- residuals(m, "working") * sqrt(weights(m, "working"))
X <- model.matrix(m) * sqrt(weights(m, "working"))
然后，协方差矩阵的公式与线性回归模型中的公式相同
vcov_byhand <- solve(t(X) %*% X) %*% t(X) %*% diag(r^2) %*% X %*% solve(t(X) %*% X)
all.equal(vcov_sandwich, vcov_byhand)
## [1] TRUE
两个警告：（1）当然，不应该像这样计算三明治，但可以使用更有效和数值更稳定的矩阵计算。(2) 对于具有独立二元响应的数据，这些三明治协方差对任何事物都不是鲁棒的。它们是一致的，没关系，但没有办法在不错误指定模型方程的情况下错误指定可能性。