兩個樣本 t 檢驗的功效

October 15, 2011

我試圖了解兩個獨立樣本 t 檢驗的功率計算（不假設方差相等，所以我使用了 Satterthwaite）。

這是我發現有助於理解該過程的圖表：

所以我假設給定以下關於兩個總體的信息並給出樣本量：
mu1<-5
mu2<-6
sd1<-3
sd2<-2
n1<-20
n2<-20
我可以計算與具有 0.05 上尾概率相關的空值下的臨界值：
df<-(((sd1^2/n1)+(sd2^2/n2)^2)^2) / ( ((sd1^2/n1)^2)/(n1-1) + ((sd2^2/n2)^2)/(n2-1)  )
CV<- qt(0.95,df) #equals 1.730018
然後計算備擇假設（在這種情況下，我學到的是“非中心 t 分佈”）。我使用非中心分佈和上面找到的臨界值計算了上圖中的 beta。這是R中的完整腳本：
#under alternative
mu1<-5
mu2<-6
sd1<-3
sd2<-2
n1<-20
n2<-20


#Under null
Sp<-sqrt(((n1-1)*sd1^2+(n2-1)*sd2^2)/(n1+n2-2))
df<-(((sd1^2/n1)+(sd2^2/n2)^2)^2) / ( ((sd1^2/n1)^2)/(n1-1) + ((sd2^2/n2)^2)/(n2-1)  )
CV<- qt(0.95,df)


#under alternative
diff<-mu1-mu2
t<-(diff)/sqrt((sd1^2/n1)+ (sd2^2/n2))
ncp<-(diff/sqrt((sd1^2/n1)+(sd2^2/n2)))


#power
1-pt(t, df, ncp)
這給出了 0.4935132 的功率值。

這是正確的方法嗎？我發現如果我使用其他功率計算軟件（例如 SAS，我認為我已經設置了與下面的問題等效的軟件），我會得到另一個答案（來自 SAS，它是 0.33）。

SAS 代碼：
proc power;
     twosamplemeans test=diff_satt
        meandiff = 1
        groupstddevs = 3 | 2
        groupweights = (1 1)
        ntotal = 40
        power = .
       sides=1;
  run;
最終，我希望獲得一種理解，使我能夠查看更複雜程序的模擬。

編輯：我發現了我的錯誤。本來應該

1-pt(CV, df, ncp) 不是 1-pt(t, df, ncp)

您已經接近了，但需要進行一些小的更改：

均值的真實差異通常被視為，而不是相反。

G*電源使用作為自由度- 在這種情況下的分佈（不同的方差，相同的組大小），按照 Cohen 的建議，如此處所述

SAS 可能會使用韋爾奇公式或薩特斯韋特公式來計算 df 給定的不等方差（在您引用的此 pdf中找到） - 結果中只有 2 個有效數字無法分辨（見下文）

如n1, n2, mu1, mu2, sd1, sd2您的問題中所定義：
> alpha   <- 0.05
> dfGP    <- n1+n2 - 2                     # degrees of freedom (used by G*Power)
> cvGP    <- qt(1-alpha, dfGP)             # crit. value for one-sided test (under the null)
> muDiff  <- mu2-mu1                       # true difference in means
> sigDiff <- sqrt((sd1^2/n1) + (sd2^2/n2)) # true SD for difference in empirical means
> ncp     <- muDiff / sigDiff              # noncentrality parameter (under alternative)
> 1-pt(cvGP, dfGP, ncp)                    # power
[1] 0.3348385
這與G*Power的結果相匹配，這是一個很好的解決這些問題的程序。它還顯示 df、臨界值和 ncp，因此您可以單獨檢查所有這些計算。

編輯：使用 Satterthwaite 的公式或 Welch 的公式變化不大（仍然是 0.33*）：
# Satterthwaite's formula
> var1  <- sd1^2
> var2  <- sd2^2
> num   <- (var1/n1 + var2/n2)^2
> denST <- var1^2/((n1-1)*n1^2) + var2^2/((n2-1)*n2^2)
> (dfST <- num/denST)
[1] 33.10309

> cvST <- qt(1-alpha, dfST)
> 1-pt(cvST, dfST, ncp)
[1] 0.3336495

# Welch's formula
> denW <- var1^2/((n1+1)*n1^2) + var2^2/((n2+1)*n2^2)
> (dfW <- (num/denW) - 2)
[1] 34.58763

> cvW   <- qt(1-alpha, dfW)
> 1-pt(cvW, dfW, ncp)
[1] 0.3340453
（注意我把一些變量名稍微改成了t, df,diff也是內置函數的名字，還要注意你的代碼的分子df是錯誤的，它有一個錯位^2，一個^2太多了，應該是((sd1^2/n1) + (sd2^2/n2))^2）

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

comments powered by Disqus

兩個樣本 t 檢驗的功效

相關問答