検定の知識整理(3) 検定力について1

有意水準・効果量・標本数・検定力はカルテットである

他の2つの条件が同じであれば、以下の3つの性質がある。

有意水準を大きくすると検定力も大きくなる。有意水準を小さくすると検定力も小さくなる。
効果量が大きければ検定力も大きい。効果量が小さければ検定力も小さい。
標本数を大きくすると検定力も大きくなる。標本数を小さくすると検定力も小さくなる。

第2種の誤り(偽陰性)を十分に小さくした(検定力を大きくした)状態で、出来るだけ少ない標本数で有意性を検出することが検定の理想(らしい)。

t検定の検出力を行ってみる。関数は、次の通り。

power.t.test(n=NULL,delta=NULL,sd=1,sig.level=0.05,type=c("two.sample","one.sample","paired"),alternative = c("two.sided","one.sided"),power=NULL,strict=FALSE)

ここで、各引数は次の意味を持つ。

n：群毎の標本数

delta：真の平均値差(μ実験―μ対称)。sd=1ならdelta=効果量。

sd：2群に共通の標準偏差。

sig.level：有意水準(第1種の誤りの確率)。既定値は5%。

power：検定力 1-β

type：t検定の型。対応なしの2群なら"two.sample", 1標本なら"one.sample",対応ありなら"paired"。期待値は対応なしの2群。

alternative：両側検定か、片側検定か。両側ならtwo.sided。片側ならone.sided。既定値は両側検定。

strict：理論値。Tなら両側の場合で厳密な計算をする。Fなら簡便解。

早速、何通りかで使ってみる。対応なしの2群、両側検定とする。

効果量を固定したときの有意水準・検定力の変化に対する標本数の変化を考える

次のコードを実行。

res <- NULL
result <- NULL
k <- 1
alpha <- c(0.01,0.05,0.10,0.50,0.80,0.90)

for(i in alpha) {
  res[[k]] <- matrix(0,100,2)
  for(j in 1:99) {
    try(result[[k]] <- power.t.test(n=NULL,delta=0.5,sd=1,sig.level=i,
                           type="two.sample",alternative="two.sided",
                           power=j/100,strict=FALSE))
    try(res[[k]][j,] <- c(j/100,round(result[[k]]$n)+1))
  }
  k <- k + 1
}

for(i in 1:k) {
  res[[i]] <- subset(res[[i]],res[[i]][,1]>0)
  if (i == 1) {
    plot(res[[i]][,1],res[[i]][,2],type="l",
         xlab="Power",ylab="Number of samples",
         main="Relation between power and number of samples in unpaired t-test")
  } else {
    lines(res[[i]][,1],res[[i]][,2],col=i)
  }
}

legend("topleft",legend=alpha,col=1:k,lty=1)

結果は以下の通り。

f:id:nakhirot:20131010165208p:plain