高度偏斜數據的 t 檢驗
我有一個包含數以萬計醫療費用數據觀察的數據集。該數據高度向右傾斜，並且有很多零。兩組人看起來像這樣（在這種情況下，兩個年齡組，每個年齡組 > 3000 obs）：
Min. 1st Qu. Median Mean 3rd Qu. Max. 0.0 0.0 0.0 4536.0 302.6 395300.0 Min. 1st Qu. Median Mean 3rd Qu. Max. 0.0 0.0 0.0 4964.0 423.8 721700.0如果我對此數據執行 Welch 的 t 檢驗，我會得到一個結果：
Welch Two Sample t-test data: x and y t = -0.4777, df = 3366.488, p-value = 0.6329 alternative hypothesis: true difference in means is not equal to 0 95 percent confidence interval: -2185.896 1329.358 sample estimates: mean of x mean of y 4536.186 4964.455我知道對這些數據使用 t 檢驗是不正確的，因為它非常不正常。但是，如果我對均值的差異使用置換檢驗，我會一直得到幾乎相同的 p 值（並且隨著更多的迭代它會變得更接近）。
在 R 和 permTS 中使用 perm 包和精確的 Monte Carlo
Exact Permutation Test Estimated by Monte Carlo data: x and y p-value = 0.6188 alternative hypothesis: true mean x - mean y is not equal to 0 sample estimates: mean x - mean y -428.2691 p-value estimated from 500 Monte Carlo replications 99 percent confidence interval on p-value: 0.5117552 0.7277040為什麼排列檢驗統計量如此接近 t.test 值？如果我記錄數據，那麼我會得到 0.28 的 t.test p 值，並且與置換測試相同。我認為 t 檢驗值比我在這裡得到的更垃圾。我喜歡的許多其他數據集都是如此，我想知道為什麼 t 檢驗似乎在不應該工作的情況下工作。
我在這裡擔心的是個人成本不是 iid 有許多具有非常不同的成本分佈（女性與男性、慢性病等）的人的子群體似乎違反了中心極限定理的 iid 要求，或者我不應該擔心關於那個？
t 檢驗和置換檢驗都沒有太大的能力來識別兩個這種異常偏斜分佈之間的均值差異。 因此，它們都給出了無痛 p 值，表明根本沒有意義。問題不在於他們似乎同意；那是因為他們很難發現任何差異，他們根本無法不同意！
出於某種直覺，考慮一下如果一個數據集中單個值發生變化會發生什麼。例如，假設第二個數據集中沒有出現最大值 721,700。平均值將下降大約 721700/3000，即 240。然而，平均值的差異僅為 4964-4536 = 438，甚至沒有兩倍大。這表明（儘管它沒有證明）任何對均值的比較都不會發現差異顯著。
但是，我們可以驗證 t 檢驗不適用。 讓我們生成一些與這些具有相同統計特徵的數據集。為此，我創建了混合物，其中
- 在任何情況下，數據都是零。
- 其餘數據具有對數正態分佈。
- 該分佈的參數被安排來重現觀察到的平均值和第三四分位數。
事實證明，在這些模擬中，最大值也離報告的最大值不遠。
讓我們將第一個數據集複製 10,000 次並跟踪其平均值。（當我們對第二個數據集執行此操作時，結果將幾乎相同。）這些均值的直方圖估計均值的採樣分佈。當該分佈近似正態時，t 檢驗有效；它偏離正態性的程度表明學生 t 分佈會出錯的程度。因此，作為參考，我還繪製了（紅色）正態分佈的 PDF 以適合這些結果。
我們看不到太多細節，因為有一些非常大的異常值。（這就是我提到的這種手段的敏感性的體現。）其中有 123 個——1.23%——超過 10,000 個。讓我們關注其餘部分，以便我們可以看到細節，因為這些異常值可能是由於假設的分佈的對數正態性導致的，而原始數據集不一定是這種情況。
這仍然嚴重偏斜並且明顯偏離正態近似，為問題中所述的現象提供了充分的解釋。它還讓我們了解測試可以檢測到的均值差異有多大：它必須在 3000 左右或更多才能顯得顯著。相反，如果您有大約 428 的實際差異，則可能會檢測到數倍的數據（在每組中）。 給定 50 倍的數據，我估計在 5% 的顯著性水平上檢測到這種差異的能力約為 0.4（這不好，但至少你會有機會）。
這是
R產生這些數字的代碼。# # Generate positive random values with a median of 0, given Q3, # and given mean. Make a proportion 1-e of them true zeros. # rskew <- function(n, x.mean, x.q3, e=3/8) { beta <- qnorm(1 - (1/4)/e) gamma <- 2*(log(x.q3) - log(x.mean/e)) sigma <- sqrt(beta^2 - gamma) + beta mu <- log(x.mean/e) - sigma^2/2 m <- floor(n * e) c(exp(rnorm(m, mu, sigma)), rep(0, n-m)) } # # See how closely the summary statistics are reproduced. # (The quartiles will be close; the maxima not too far off; # the means may differ a lot, though.) # set.seed(23) x <- rskew(3300, 4536, 302.6) y <- rskew(3400, 4964, 423.8) summary(x) summary(y) # # Estimate the sampling distribution of the mean. # set.seed(17) sim.x <- replicate(10^4, mean(rskew(3367, 4536, 302.6))) hist(sim.x, freq=FALSE, ylim=c(0, dnorm(0, sd=sd(sim.x)))) curve(dnorm(x, mean(sim.x), sd(sim.x)), add=TRUE, col="Red") hist(sim.x[sim.x < 10000], xlab="x", freq=FALSE) curve(dnorm(x, mean(sim.x), sd(sim.x)), add=TRUE, col="Red") # # Can a t-test detect a difference with more data? # set.seed(23) n.factor <- 50 z <- replicate(10^3, { x <- rskew(3300*n.factor, 4536, 302.6) y <- rskew(3400*n.factor, 4964, 423.8) t.test(x,y)$p.value }) hist(z) mean(z < .05) # The estimated power at a 5% significance level
