50% 的置信區間是否比 95% 的置信區間更可靠?
我的問題源於對 Andrew Gelman 的博客文章的評論,其中他主張使用 50% 的置信區間而不是 95% 的置信區間,儘管不是因為它們被更穩健地估計:
我更喜歡 50% 到 95% 的間隔,原因有 3 個:
- 計算穩定性,
- 更直觀的評估(50% 區間的一半應該包含真實值),
- 一種感覺,在應用中最好了解參數和預測值的位置,而不是嘗試不切實際的近乎確定性。
評論者的想法似乎是,與 50% CI 相比,如果置信區間是 95% CI,則構建置信區間的假設問題會產生更大的影響。然而,他並沒有真正解釋原因。
[…] 隨著您進入更大的間隔,您通常對模型的細節或假設變得更加敏感。例如,您永遠不會相信您已經正確識別了 99.9995% 區間。或者至少這是我的直覺。如果它是正確的,它認為應該比 95% 更好地估計 50%。或者也許是“更穩健”的估計,因為它對關於噪聲的假設不太敏感,也許?
這是真的嗎?為什麼/為什麼不?
這個答案分析了引文的含義,並提供了模擬研究的結果來說明它並幫助理解它可能想說什麼。任何人(具有基本技能)都可以輕鬆地擴展這項研究,
R
以探索其他置信區間程序和其他模型。在這項工作中出現了兩個有趣的問題。 一個是關於如何評估置信區間程序的準確性。人們對健壯性的印象取決於此。我展示了兩種不同的準確度度量,以便您可以比較它們。
另一個問題是,儘管低置信度的置信區間過程可能是穩健的,但相應的置信限可能根本不穩健。間隔往往效果很好,因為他們在一端犯的錯誤通常會抵消他們在另一端犯的錯誤。作為一個實際問題,你可以很確定你的大約一半置信區間覆蓋了它們的參數,*但實際參數可能始終位於每個區間的一個特定端附近,*具體取決於現實與模型假設的偏離程度。
健壯在統計中具有標準含義:
穩健性通常意味著對偏離圍繞潛在概率模型的假設不敏感。
(Hoaglin、Moseller 和 Tukey,了解穩健和探索性數據分析。J. Wiley (1983),第 2 頁。)
這與問題中的引用一致。要理解報價,我們仍然需要知道置信區間的預期目的。為此,讓我們回顧一下格爾曼所寫的內容。
我更喜歡 50% 到 95% 的間隔,原因有 3 個:
- 計算穩定性,
- 更直觀的評估(50% 區間的一半應該包含真實值),
- 感覺在應用程序中最好了解參數和預測值的位置,而不是嘗試不切實際的接近確定性。
由於了解預測值不是置信區間 (CI) 的目的,因此我將重點了解參數值,這正是 CI 所做的。我們稱這些為“目標”值。因此,根據定義, CI 旨在以指定的概率(其置信水平)覆蓋其目標。達到預期的覆蓋率是評估任何 CI 程序質量的最低標準。(此外,我們可能對典型的 CI 寬度感興趣。為了使帖子保持合理的長度,我將忽略這個問題。)
這些考慮促使我們研究置信區間計算可能在多大程度上誤導我們關於目標參數值的問題。 該引文可以被解讀為暗示置信度較低的 CI 可能會保留其覆蓋範圍,即使數據是由與模型不同的過程生成的。這是我們可以測試的。程序是:
- 採用包含至少一個參數的概率模型。經典的方法是從未知均值和方差的正態分佈中採樣。
- 為模型的一個或多個參數選擇 CI 過程。一個優秀的方法是根據樣本均值和样本標準差構造 CI,將後者乘以學生 t 分佈給出的因子。
- 將該程序應用於各種不同的模型——與採用的模型相差不大——以評估其在一系列置信水平上的覆蓋範圍。
例如,我就是這樣做的。我允許基礎分佈在很寬的範圍內變化,從幾乎伯努利到均勻,到正態,到指數,一直到對數正態。這些包括對稱分佈(前三個)和強烈偏斜分佈(後兩個)。對於每個分佈,我生成了 50,000 個大小為 12 的樣本。對於每個樣本,我構建了置信水平在和,它涵蓋了大多數應用程序。
現在出現了一個有趣的問題:我們應該如何衡量 CI 程序執行的好壞(或差)? 一種常見的方法只是評估實際覆蓋率和置信度之間的差異。不過,這對於高置信度水平來說可能看起來非常好。例如,如果您試圖達到 99.9% 的置信度,但您只獲得 99% 的覆蓋率,那麼原始差異僅為 0.9%。但是,這意味著您的程序未能覆蓋目標的次數是應有的十倍!出於這個原因,比較覆蓋範圍的信息更豐富的方法應該使用優勢比之類的東西。我使用 logits 的差異,這是優勢比的對數。具體來說,當期望的置信水平是實際覆蓋範圍是, 然後
很好地捕捉到了差異。當它為零時,覆蓋範圍正是預期的值。當它為負時,覆蓋率太低——這意味著 CI 過於樂觀,低估了不確定性。
那麼,問題是,當基礎模型受到干擾時,這些錯誤率如何隨置信水平而變化? 我們可以通過繪製模擬結果來回答它。 *這些圖量化了 CI 的“接近確定性”*在這個原型應用程序中的“不切實際”程度。
圖形顯示了相同的結果,但左側顯示的是 logit 標度值,而右側使用的是原始標度。Beta 發行版是 Beta(實際上是伯努利分佈)。對數正態分佈是標準正態分佈的指數。包含正態分佈是為了驗證這個 CI 過程確實達到了其預期的覆蓋範圍,並揭示了從有限模擬大小中可以預期有多少變化。(事實上,正態分佈的圖表非常接近於零,沒有顯示出明顯的偏差。)
很明顯,在 logit 尺度上,隨著置信水平的增加,覆蓋範圍變得更加分散。 不過,也有一些有趣的例外。如果我們不關心引入偏斜或長尾的模型擾動,那麼我們可以忽略指數和對數正態並專注於其餘部分。他們的行為反复無常,直到超過左右(一個logit),此時分歧已經開始。
**這項小型研究使格爾曼的主張更加具體,**並說明了他可能想到的一些現象。特別是,當我們使用置信度較低的 CI 程序時,例如,那麼即使底層模型受到強烈擾動,看起來覆蓋率仍然接近:我們認為這樣的 CI 大約有一半是正確的,另一半是不正確的。那是健壯的。相反,如果我們希望是對的,比如說,的時間,這意味著我們真的只想犯錯在那個時候,我們應該為我們的錯誤率做好準備,以防萬一世界不像我們的模型假設的那樣運作。
順便說一句,這個屬性CI 之所以成立,很大程度上是因為我們正在研究對稱置信區間。對於偏態分佈,個體置信限可能很糟糕(而且根本不穩健),但它們的錯誤通常會抵消。 通常,一條尾巴短而另一條尾巴長,導致一端覆蓋過度,另一端覆蓋不足。我相信置信限不會像相應的區間那樣穩健。
這是
R
生成圖的代碼。它很容易修改以研究其他分佈、其他置信範圍和其他 CI 程序。# # Zero-mean distributions. # distributions <- list(Beta=function(n) rbeta(n, 1/30, 1/30) - 1/2, Uniform=function(n) runif(n, -1, 1), Normal=rnorm, #Mixture=function(n) rnorm(n, -2) + rnorm(n, 2), Exponential=function(n) rexp(n) - 1, Lognormal=function(n) exp(rnorm(n, -1/2)) - 1 ) n.sample <- 12 n.sim <- 5e4 alpha.logit <- seq(0, 6, length.out=21); alpha <- signif(1 / (1 + exp(-alpha.logit)), 3) # # Normal CI. # CI <- function(x, Z=outer(c(1,-1), qt((1-alpha)/2, n.sample-1))) mean(x) + Z * sd(x) / sqrt(length(x)) # # The simulation. # #set.seed(17) alpha.s <- paste0("alpha=", alpha) sim <- lapply(distributions, function(dist) { x <- matrix(dist(n.sim*n.sample), n.sample) x.ci <- array(apply(x, 2, CI), c(2, length(alpha), n.sim), dimnames=list(Endpoint=c("Lower", "Upper"), Alpha=alpha.s, NULL)) covers <- x.ci["Lower",,] * x.ci["Upper",,] <= 0 rowMeans(covers) }) (sim) # # The plots. # logit <- function(p) log(p/(1-p)) colors <- hsv((1:length(sim)-1)/length(sim), 0.8, 0.7) par(mfrow=c(1,2)) plot(range(alpha.logit), c(-2,1), type="n", main="Confidence Interval Accuracies (Logit Scales)", cex.main=0.8, xlab="Logit(alpha)", ylab="Logit(coverage) - Logit(alpha)") abline(h=0, col="Gray", lwd=2) legend("bottomleft", names(sim), col=colors, lwd=2, bty="n", cex=0.8) for(i in 1:length(sim)) { coverage <- sim[[i]] lines(alpha.logit, logit(coverage) - alpha.logit, col=colors[i], lwd=2) } plot(range(alpha), c(-0.2, 0.05), type="n", main="Raw Confidence Interval Accuracies", cex.main=0.8, xlab="alpha", ylab="coverage-alpha") abline(h=0, col="Gray", lwd=2) legend("bottomleft", names(sim), col=colors, lwd=2, bty="n", cex=0.8) for(i in 1:length(sim)) { coverage <- sim[[i]] lines(alpha, coverage - alpha, col=colors[i], lwd=2) }