統計ER

統計ソフトRの使い方を中心に、統計解析方法の解説をするブログ。ありそうでなかなか見つからないサンプルサイズ計算などニッチな方法について紹介しています。

t検定のサンプルサイズ計算はどうやる?

Toukei Consul Banner

KH Coder Consul Banner

t検定(平均値の差の検定)のサンプルサイズ計算はどうやるのだろうか。

1:1の割り付けじゃなくて、1:2とか1:3とかにするにはどうしたらいいか?

統計ソフトRでt検定のサンプルサイズを計算してみた。

t検定のサンプルサイズ計算(1:1の場合)

power.t.test()を用いる。

必要な情報は、

検出力をどうするかと、
二群の平均値の差。

標準偏差で割った標準化された差が必要だ。

検出力は慣例で80%=0.8が多い。90%=0.9もよくある。

標準化された差dが難しい。

先行研究は二群の差の標準偏差が計算されていないことがほとんどだ。

しかし、一つ朗報だ。

なんとなくの目安がある。

  • 小さな差を検出したければd=0.1~0.2
  • 中くらいな差を検出したければd=0.4~0.5
  • 大きな差を検出したければd=0.8~0.9

検出力80%、中くらいな差を検出したいとしてd=0.5とする。

power.t.test()では、

検出力をpower=で、
標準化された差をdelta=で、

指定する。

power.t.test(power=0.8,delta=0.5)

結果は一群64例必要と計算される。

二群合わせて128例が必要なサンプルサイズだ。

> power.t.test(power=0.8,delta=0.5)

     Two-sample t test power calculation 

              n = 63.76576
          delta = 0.5
             sd = 1
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

NOTE: n is number in *each* group

t検定のサンプルサイズ計算(1:1の場合)power.t.test()を使わない方法(2020年4月19日追記)

power.t.test()を使わない簡易的な方法のスクリプトを下記書籍を参照して自作した。

t.test.samplesize <- function(sig.level=0.05, power=0.8, 
delta, sd=1, alternative=c("two.sided","one.sided"))
{
 alternative <- match.arg(alternative)
 tside <- switch(alternative, one.sided=1, two.sided=2)

 Za <- qnorm(sig.level/tside, lower.tail=FALSE)
 Zb <- qnorm(power)
 n <- 2*((Za+Zb)/(delta/sd))^2

 NOTE <- "n is number in *each* group"
 METHOD <- "t-test sample size"
 structure(list(n=n,delta=delta, sd=sd, 
 sig.level=sig.level, power=power, 
 alternative=alternative, note=NOTE, method=METHOD),
 class="power.htest")
}

結果は一群63例と計算された。power.t.test()を使った時は64例と計算された。簡易的な方法でも大きな違いがない。

> t.test.samplesize(delta=0.5)

     t-test sample size 

              n = 62.79104
          delta = 0.5
             sd = 1
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

NOTE: n is number in *each* group

参考書籍

t検定のサンプルサイズ計算(割り付け比が1:1でない場合)

統計ソフトRの標準では、1:2とか1:3に割り付けるサンプルサイズ計算のプログラムは用意されていない。

自前で作ったプログラムが以下の通り。

power.t.test()で1:1のサンプルサイズを計算した後、割り付けの比(ctrl.ratio)で調整している。

power.t.test.unequal <- function(n = NULL, delta = NULL, sd = 1, sig.level = 0.05, power = NULL, 
 type = "two.sample", alternative = "two.sided", strict = FALSE, ctrl.ratio = NULL) {
 res <- power.t.test(n = n, delta = delta, sd = sd, sig.level = sig.level, power = power,
 type = type, alternative = alternative, strict = strict)
 alternative <- match.arg(alternative)
 tside <- switch(alternative, one.sided=1, two.sided=2)
 nA <- 2 * res$n / (ctrl.ratio + 1)
 nB <- 2 * res$n * ctrl.ratio / (ctrl.ratio + 1)
 Za <- qnorm(sig.level/tside, lower.tail=FALSE)
 power.u <- 1 - pnorm(Za - delta/sd * (1 / sqrt(1 / nA + 1 / nB)))
 METHOD <- paste("Two-sample", "t test power calculation")
 structure(list(n.equal = res$n, n.unequal.1 = nA, n.unequal.2 = nB, delta = delta, sd = sd, 
 sig.level = sig.level, power.equal = res$power, power.unequal = power.u, alternative = alternative,
 method = METHOD), class = "power.htest")
}

powerとdeltaは先ほどと同じにして、

1:1, 1:2, 1:3の場合を計算してみる。

まず1:1の場合。1:1は前節と同じように一群64例と計算される。

> power.t.test.unequal(power=0.8, delta=0.5, ctrl.ratio=1)

     Two-sample t test power calculation 

        n.equal = 63.76576
    n.unequal.1 = 63.76576
    n.unequal.2 = 63.76576
          delta = 0.5
             sd = 1
      sig.level = 0.05
    power.equal = 0.8
  power.unequal = 0.8060089
    alternative = two.sided

1:2は43例と86例と計算される。検出力は76%に落ちる。

> power.t.test.unequal(power=0.8, delta=0.5, ctrl.ratio=2)

     Two-sample t test power calculation 

        n.equal = 63.76576
    n.unequal.1 = 42.51051
    n.unequal.2 = 85.02102
          delta = 0.5
             sd = 1
      sig.level = 0.05
    power.equal = 0.8
  power.unequal = 0.7586038
    alternative = two.sided

1:3は32例と96例と計算される。検出力はさらに落ちて69%になる。

> power.t.test.unequal(power=0.8, delta=0.5, ctrl.ratio=3)

     Two-sample t test power calculation 

        n.equal = 63.76576
    n.unequal.1 = 31.88288
    n.unequal.2 = 95.64865
          delta = 0.5
             sd = 1
      sig.level = 0.05
    power.equal = 0.8
  power.unequal = 0.6861757
    alternative = two.sided

ここからは1:1が最も検出力が高いことがわかる。1:1の割り付けがおすすめなのは、検出力が高いからだ。

1:1割り付け出ない場合、最終的な検出力を80%に保つためには最初の検出力を上げておかねばならない。

最初の検出力を85%にしておいて、1:2割り付けのサンプルサイズを計算すると、どうなるか?

結果は、49例と98例で、最終の検出力が81%になる。

これならば、必要症例数は64例から49例に15例減り、検出力は80%を保っているため、有用かつ適切なデザイン変更だ。

ただし、合計の必要例数例数は49+98=147となり、1:1のときの合計例数128例より19例増加する。

> power.t.test.unequal(power=0.85, delta=0.5, ctrl.ratio=2)

     Two-sample t test power calculation 

        n.equal = 72.80053
    n.unequal.1 = 48.53369
    n.unequal.2 = 97.06738
          delta = 0.5
             sd = 1
      sig.level = 0.05
    power.equal = 0.85
  power.unequal = 0.8116907
    alternative = two.sided

1:3の場合は、1:1の検出力を90%にすると、1:3の検出力が80%を超える。

必要例数は43例と128例で、合計171例必要になる。1:1の時より43例増加する。

> power.t.test.unequal(power=0.90, delta=0.5, ctrl.ratio=3)

     Two-sample t test power calculation 

        n.equal = 85.03129
    n.unequal.1 = 42.51564
    n.unequal.2 = 127.5469
          delta = 0.5
             sd = 1
      sig.level = 0.05
    power.equal = 0.9
  power.unequal = 0.8060558
    alternative = two.sided

症例を集めるのが大変な場合は致し方ないが、できれば1:1にしたほうが、合計の必要なサンプルサイズは小さくなる。

t検定のサンプルサイズ計算(割り付け比に対応)power.t.test()を使わない方法(2020年4月20日追記)

1:1の場合と同様に、power.t.test()を使わない簡易的な方法のスクリプトを作成して、power.t.test()を使った場合と比較した。

power.t.test.unequal.approx <- function(
 delta = NULL, sd = 1, sig.level = 0.05, power = 0.8, 
 alternative = c("two.sided","one.sided"), 
 ctrl.ratio = NULL) {

 alternative <- match.arg(alternative)
 tside <- switch(alternative, one.sided=1, two.sided=2)

 Za0 <- qnorm(sig.level/tside, lower.tail=FALSE)
 Zb0 <- qnorm(power)
 n <- 2*((Za0+Zb0)/(delta/sd))^2

 nA <- 2 * n / (ctrl.ratio + 1)
 nB <- 2 * n * ctrl.ratio / (ctrl.ratio + 1)
 Za <- qnorm(sig.level/tside, lower.tail=FALSE)
 power.u <- 1 - pnorm(Za - delta/sd * (1 / sqrt(1 / nA + 1 / nB)))

 METHOD <- paste("Two-sample", "t test power calculation")
 structure(list(n.equal = n, n.unequal.1 = nA, 
 n.unequal.2 = nB, delta = delta, sd = sd, 
 sig.level = sig.level, power.equal = power, 
 power.unequal = power.u,  alternative = alternative,
 method = METHOD), class = "power.htest")
}

まず1:1の場合。power.t.test()の結果と比べ、各群1例少ないのは上述した通り。

> power.t.test.unequal.approx(power=0.8, delta=0.5, ctrl.ratio=1)

     Two-sample t test power calculation 

        n.equal = 62.79104
    n.unequal.1 = 62.79104
    n.unequal.2 = 62.79104
          delta = 0.5
             sd = 1
      sig.level = 0.05
    power.equal = 0.8
  power.unequal = 0.8
    alternative = two.sided

1:2の場合は、1:1の場合の検出力を85%にすると、検出力が80%を超えるのは同じ。

小さい例数の群は49から48例に1例減少、大きい例数の群は98例から96例に2例減少した。大きな違いはないと言える。

> power.t.test.unequal.approx(power=0.85, delta=0.5, ctrl.ratio=2)

     Two-sample t test power calculation 

        n.equal = 71.82718
    n.unequal.1 = 47.88479
    n.unequal.2 = 95.76957
          delta = 0.5
             sd = 1
      sig.level = 0.05
    power.equal = 0.85
  power.unequal = 0.8064989
    alternative = two.sided

1:3の場合は、1:1の場合の検出力を90%にすると、検出力が80%を超えるのはこちらも同様。

小さい例数の群は、43例から変化なし。大きい例数の群は、128例から127例に1例減少。こちらも大きな違いはないと言える。

> power.t.test.unequal.approx(power=0.90, delta=0.5, ctrl.ratio=3)

     Two-sample t test power calculation 

        n.equal = 84.05938
    n.unequal.1 = 42.02969
    n.unequal.2 = 126.0891
          delta = 0.5
             sd = 1
      sig.level = 0.05
    power.equal = 0.9
  power.unequal = 0.8015779
    alternative = two.sided

よって、簡易的な(近似的な)方法は、power.t.test()を使った方法と、実用上は差がないと言えるだろう。

t検定のサンプルサイズ計算(割り付け比に対応)power.t.test()を使わないもう少し簡単な方法(2020年9月20日追記)

EZRの「2群の平均値の比較のためのサンプルサイズの計算」を参考にして、もう少し簡単に割り付け比ありのサンプルサイズ計算ができるように改造した。

Caseの例数 n1 を計算したときに切り上げて、そこに割り付け比をかけるように変更した。

t.test.samplesize.unequal.approx2 <- function(sig.level=0.05, power=0.8, 
delta, sd=1, alternative=c("two.sided","one.sided"),ctrl.ratio=NULL)
{
 alternative <- match.arg(alternative)
 tside <- switch(alternative, one.sided=1, two.sided=2)

 Za <- qnorm(sig.level/tside, lower.tail=FALSE)
 Zb <- qnorm(power)
 n1 <- ceiling((1+1/ctrl.ratio)*((Za+Zb)/(delta/sd))^2)
 n2 <- n1*ctrl.ratio

 NOTE <- "n1 or n2 is number in *each* group"
 METHOD <- "t-test sample size"
 structure(list(n1=n1, n2=n2, "control ratio"=ctrl.ratio, delta=delta, sd=sd, 
 sig.level=sig.level, power=power, 
 alternative=alternative, note=NOTE, method=METHOD),
 class="power.htest")
}

sig.level=0.05, power=0.8, sd=1のデフォルトのままとして、delta=0.5, ctrl.ratio=1,2,3のとき、結果は以下の通り。

> t.test.samplesize.unequal.approx2(delta=0.5,ctrl.ratio=1)

     t-test sample size 

             n1 = 63
             n2 = 63
  control ratio = 1
          delta = 0.5
             sd = 1
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

NOTE: n1 or n2 is number in *each* group

> t.test.samplesize.unequal.approx2(delta=0.5,ctrl.ratio=2)

     t-test sample size 

             n1 = 48
             n2 = 96
  control ratio = 2
          delta = 0.5
             sd = 1
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

NOTE: n1 or n2 is number in *each* group

> t.test.samplesize.unequal.approx2(delta=0.5,ctrl.ratio=3)

     t-test sample size 

             n1 = 42
             n2 = 126
  control ratio = 3
          delta = 0.5
             sd = 1
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

NOTE: n1 or n2 is number in *each* group

エクセルでサンプルサイズ計算【改良版】t検定 平均値の差の検定(割り付け比にも対応)(2020年9月20日追記)

エクセルファイル販売中。よければどうぞ。

t検定 平均値の差の検定(割り付け比にも対応)サンプルサイズ計算【改良版】【エクセルでサンプルサイズ】 | HHA SHOP

使い方解説動画。こちらもよければ、ぜひ。

youtu.be

参考書籍

丹後敏郎著 無作為化比較試験 朝倉書店
3.4 割り付け数の不均衡の影響