为何lmtest::bptest()与check_heteroskedasticity()的BP检验结果存在差异？

在对线性模型执行异方差检验时，我使用了两个R函数实现Breusch-Pagan检验：lmtest::bptest()和performance::check_heteroskedasticity()。根据check_heteroskedasticity的文档说明，它所做的正是Breusch-Pagan（1979）的（非）恒定误差方差检验，但实际运行后两个函数给出的p值存在差异——仅check_heteroskedasticity的p值小于0.05。我想确认是否因对函数的理解或应用有误导致了这一现象，以下是复现该差异的示例：

m <- lm(mpg ~ hp + cyl, data = mtcars)

performance::check_heteroscedasticity(m)
lmtest::bptest(m, studentize = TRUE)
lmtest::bptest(m, studentize = FALSE)

运行结果：

> performance::check_heteroscedasticity(m)
Warning: Heteroscedasticity (non-constant error variance) detected (p = 0.041).

> lmtest::bptest(m, studentize = TRUE)

    studentized Breusch-Pagan test

data:  m
BP = 5.5086, df = 2, p-value = 0.06365

> lmtest::bptest(m, studentize = FALSE)

    Breusch-Pagan test

data:  m
BP = 4.4698, df = 2, p-value = 0.107

结果差异的核心原因 #

两个函数的Breusch-Pagan检验实现细节存在关键差别：

• performance::check_heteroscedasticity()默认使用学生化残差，并采用基于残差平方对原模型自变量回归的F检验计算p值。F检验在小样本（如mtcars的32个观测）下的统计表现通常比卡方检验更稳健。
• lmtest::bptest()默认使用原始残差+卡方检验；即使设置studentize=TRUE启用学生化残差，它仍然基于卡方统计量计算p值，而非F检验。

可以手动复现performance包的检验逻辑验证：

# 手动计算performance的检验结果
resid_sq <- residuals(m)^2
model_resid <- lm(resid_sq ~ hp + cyl, data = mtcars)
anova(model_resid)$Pr[1] # 输出约0.041，与check_heteroscedasticity结果一致

而lmtest::bptest(m, studentize=TRUE)的卡方统计量为5.5086，对应卡方(2)分布的p值为1-pchisq(5.5086,2)=0.06365，与输出结果匹配。

总结 #

两个函数均属于Breusch-Pagan检验的变体，差异并非使用错误，而是检验统计量的选择（F vs 卡方）和残差处理逻辑的细微差别导致：

• 若样本量较小，performance::check_heteroscedasticity()的F检验结果更可靠；
• 若需要严格匹配经典Breusch-Pagan检验的卡方统计量标准，选择lmtest::bptest()即可。

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