在R中使用sandwich包后，如何基于含稳健标准误的模型做预测？

一、直接用predict()做点预测是否可行？ #

完全可行。原因是：

• 稳健标准误（如sandwich估计器）仅修正系数的标准误，用于改进假设检验和系数置信区间的准确性，但不会改变模型的系数点估计（系数仍然是glm输出的MLE估计值）。
• predict(model, newdata = ..., type = "response")的点预测值完全基于系数点估计计算，因此和使用普通标准误的模型结果一致，不需要做任何调整。

二、如何在预测中纳入稳健标准误？ #

稳健标准误的作用是量化预测的不确定性，即构建稳健的预测区间。具体步骤如下：

核心思路 #

预测的不确定性来自两部分：系数估计的方差（用稳健方差矩阵修正），以及模型的随机误差（泊松模型的离散性）。我们需要先计算线性预测尺度（log尺度）的稳健方差，再转换到响应尺度（相对风险尺度）。

具体代码实现 #

1. 准备新数据（示例）

set.seed(123)
# 模拟与原数据分布匹配的新数据starwars2
starwars2 <- data.frame(
  height = rnorm(10, mean = mean(data$height, na.rm = T), sd = sd(data$height, na.rm = T)),
  mass = rnorm(10, mean = mean(data$mass, na.rm = T), sd = sd(data$mass, na.rm = T)),
  birth_year = rnorm(10, mean = mean(data$birth_year, na.rm = T), sd = sd(data$birth_year, na.rm = T))
)

2. 计算稳健预测区间

library(sandwich)
library(lmtest)

# 获取新数据的模型矩阵（包含截距项，与原模型公式一致）
X_new <- model.matrix(formula(model), data = starwars2)

# 提取稳健系数方差-协方差矩阵（HC1是流行病学中常用的稳健类型，可替换为HC3等）
vcov_robust <- vcovHC(model, type = "HC1")

# 计算log尺度的线性预测值
lin_pred <- predict(model, newdata = starwars2, type = "link")

# 计算每个新观测的线性预测方差（基于稳健协方差矩阵）
lin_pred_var <- apply(X_new, 1, function(x) t(x) %*% vcov_robust %*% x)

# 转换为响应尺度的预测区间（95%置信水平）
# 注意：先在log尺度计算区间再指数转换，比直接在响应尺度计算更准确
lower_ci <- exp(lin_pred - 1.96 * sqrt(lin_pred_var))
upper_ci <- exp(lin_pred + 1.96 * sqrt(lin_pred_var))

# 整合点预测与稳健区间
pred_results <- data.frame(
  点预测值 = predict(model, newdata = starwars2, type = "response"),
  95%稳健区间下限 = lower_ci,
  95%稳健区间上限 = upper_ci
)

print(pred_results)

3. 考虑泊松离散性的完整预测区间

如果需要同时考虑系数估计不确定性和泊松模型的随机离散误差，可以在方差中加入泊松的均值方差项（泊松方差等于均值）：

# 计算总方差：系数估计方差转换到响应尺度 + 泊松离散方差
response_var <- exp(2 * lin_pred) * lin_pred_var
total_var <- response_var + pred_results$点预测值

# 构建包含离散误差的95%预测区间（下界设为0避免负数）
lower_ci_full <- pmax(pred_results$点预测值 - 1.96 * sqrt(total_var), 0)
upper_ci_full <- pred_results$点预测值 + 1.96 * sqrt(total_var)

# 更新结果
pred_results$包含离散误差的95%下限 <- lower_ci_full
pred_results$包含离散误差的95%上限 <- upper_ci_full

print(pred_results)

补充说明 #

• 对于你用泊松模型拟合二分类结局的场景，稳健标准误恰好能处理这类模型常见的过度离散问题（实际方差大于泊松假设的均值方差），因此使用稳健预测区间比普通区间更可靠。
• 所有计算的核心是用sandwich包的稳健协方差矩阵替换普通的vcov(model)，从而得到更准确的不确定性估计。

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