手动计算回归杠杆统计量及R实现验证问题咨询

杠杆值（又称帽子值）的核心计算逻辑来自多元线性回归的帽子矩阵对角线元素，公式为h_ii = diag(X(X^T X)^{-1}X^T)，其中X为回归的设计矩阵。你可以按照以下步骤实现手动计算并验证结果：

实现代码 #

# 1. 构造你提供的数据集
df <- data.frame(
  response = c(-0.1999981, -0.2796748, -0.3590789, -0.4381798, -0.5169470, -0.5953499),
  x = c(2.000000, 1.997601, 1.994407, 1.990421, 1.985645, 1.980083),
  y = c(-4.794927e-09, -3.995733e-02, -7.981885e-02, -1.195688e-01, -1.591913e-01, -1.986710e-01)
)

# 2. 按你的回归方程（无截距项）手动计算杠杆值
# 构造无截距的设计矩阵
X <- as.matrix(df[, c("x", "y")])
# 计算(X转置乘X)的逆矩阵
XTX_inv <- solve(t(X) %*% X)
# 生成帽子矩阵
hat_matrix <- X %*% XTX_inv %*% t(X)
# 提取对角线元素即为各观测的杠杆值
manual_leverage <- diag(hat_matrix)

# 3. 与lm内置函数结果对比验证
fit <- lm(response ~ 0 + x + y, data = df)
lm_leverage <- hatvalues(fit)

# 输出两组结果对比
print("手动计算杠杆值：")
print(manual_leverage)
print("lm函数输出杠杆值：")
print(lm_leverage)

结果说明 #

• 运行代码后可看到手动计算结果与R内置函数输出完全一致，你可以用该结果核对你之前的计算值
• 如果你使用的是带截距项的回归模型，只需将设计矩阵修改为cbind(1, df[,c("x","y")])，同时将lm公式改为response ~ x + y即可
• 高杠杆点的常用判断阈值为2*p/n，其中p为模型参数个数（无截距时为自变量个数，带截距时为自变量个数+1），n为观测样本量，超过阈值的观测可判定为高杠杆点

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