為什麼 R 中的 lm 和 biglm 為相同的數據給出不同的 p 值？

這是一個小例子：

MyDf<-data.frame(x=c(1,2,3,4), y=c(1.2, .7, -.5, -3))

現在有了base::lm：

> lm(y~x, data=MyDf) %>% summary

Call:
lm(formula = y ~ x, data = MyDf)

Residuals:
   1     2     3     4 
-0.47  0.41  0.59 -0.53 

Coefficients:
           Estimate Std. Error t value Pr(>|t|)  
(Intercept)   3.0500     0.8738   3.491   0.0732 .
x            -1.3800     0.3191  -4.325   0.0495 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.7134 on 2 degrees of freedom
Multiple R-squared:  0.9034,    Adjusted R-squared:  0.8551 
F-statistic: 18.71 on 1 and 2 DF,  p-value: 0.04952

現在，biglm從biglm包中嘗試同樣的事情：

XX<-biglm(y~x, data=MyDf) 
print(summary(XX), digits=5)

Large data regression model: biglm(y ~ x, data = MyDf)
Sample size =  4 
            Coef     (95%      CI)      SE       p
(Intercept)  3.05  1.30243  4.79757 0.87378 0.00048
x           -1.38 -2.01812 -0.74188 0.31906 0.00002

請注意，我們需要printanddigits來查看 p 值。係數和標準誤相同，但 p 值非常不同。為什麼會這樣？

要查看哪些 p 值是正確的（如果有的話），讓我們對原假設為真的模擬數據重複計算。在當前設置中，計算是對 (x,y) 數據的最小二乘擬合，零假設是斜率為零。在這個問題中有四個 x 值 1、2、3、4，估計誤差在 0.7 左右，所以讓我們將其合併到模擬中。

這是設置，編寫為每個人都可以理解，即使是那些不熟悉R.

beta <- c(intercept=0, slope=0)
sigma <- 0.7
x <- 1:4
y.expected <-  beta["intercept"] + beta["slope"] * x

模擬生成獨立的錯誤，將它們添加到 中y.expected，調用lm以進行擬合，併summary計算 p 值。儘管這效率低下，但它正在測試使用的實際代碼。 我們仍然可以在一秒鐘內進行數千次迭代：

n.sim <- 1e3
set.seed(17)
data.simulated <- matrix(rnorm(n.sim*length(y.expected), y.expected, sigma), ncol=n.sim)
slope.p.value <- function(e) coef(summary(lm(y.expected + e ~ x)))["x", "Pr(>|t|)"]
p.values <- apply(data.simulated, 2, slope.p.value)

正確計算的 p 值將在和當原假設為真時。這些 p 值的直方圖將使我們能夠直觀地檢查這一點——它看起來是否大致水平——並且均勻性的卡方檢驗將允許進行更正式的評估。這是直方圖：

h <- hist(p.values, breaks=seq(0, 1, length.out=20))

而且，對於那些可能認為這不夠統一的人，這裡是卡方檢驗：

chisq.test(h$counts)

X 平方 = 13.042，df = 18，p 值 = 0.7891

此測試中的大 p 值表明這些結果與預期的均勻性一致。換句話說，lm是正確的。

那麼，p 值的差異從何而來？讓我們檢查一下可能被調用來計算 p 值的公式。在任何情況下，測試統計量將是

等於估計係數之間的差異和假設的（和正確的值），表示為係數估計的標準誤差的倍數。在問題中，這些值是

對於截距估計和

用於斜率估計。通常這些將與學生進行比較自由度參數為的分佈（數據量）減（估計係數的數量）。讓我們計算它的截距：

pt(-abs(3.05/0.87378), 4-2) * 2

[1] 0.0732

（這個計算乘以左尾學生概率由因為這是一個測試反對雙向選擇.) 它與lm輸出一致。

另一種計算將使用標準正態分佈來近似學生分配。讓我們看看它產生了什麼：

pnorm(-abs(3.05/0.87378)) * 2

[1] 0.000482

果然：biglm假設零分佈統計是標準正態。這是多大的錯誤？biglm使用代替重新運行前面的模擬lm會給出這個 p 值的直方圖：

這些 p 值中幾乎有 18% 小於，“顯著性”的標準閾值。這是一個巨大的錯誤。

---

我們可以從這個小小的調查中學到的一些教訓是：

1. 不要對小數據集使用從漸近分析（如標準正態分佈）得出的近似值。
2. 了解您的軟件。

引用自：https://stats.stackexchange.com/questions/255520