ねこすたっと

ねこの気持ちと統計について悩む筆者の備忘録的ページ。

サンプルサイズ計算:2群の連続アウトカム(正規分布)を比較する [R]

自分用のリファレンスとして、サンプルサイズ計算に関する記事をシリーズで書いています。 なるべく体裁を統一するために、以下のように決めています。

  • 群を示す添字について:
    • c, C = 対照群
    • t, T = 介入群*1
    • e, E = 要因群
    • 添字なし = 全体
  • 使用する記号:
    • p = 反応割合
    • μ = 平均値
    • σ = 標準偏差
    • Φ = 割り付け比 Nt/Nc(デフォルトは1)
    • α = αエラー(デフォルトは0.05)
    • β = βエラー(デフォルトは0.2)
    • z = 標準正規分布の累積分布関数

想定シナリオ

疾患Dに対する治療として、従来薬Cと新規薬Tのどちらが優れているだろうか?
ランダム化比較試験で患者を1:1に割り付けて、疾患活動度スコアを比較したいと考えている。
これまでの知見から、治療開始後のスコアは従来薬Cでは60点、新規薬Tでは50点となることが見込まれる。 両側有意水準 5%, 検出力 80%として、必要なサンプルサイズはどれくらいだろうか?

このシナリオでサンプルサイズを計算するのに必要な条件を整理してみると、以下のようになる。

  • 帰無仮説  H_0: \mu_T = \mu_C
  • 対立仮説  H_A: \mu_T \neq \mu_C
  • 有意水準:α = 0.05
  • 検出力:1-β = 0.80
  • 割り付け比:φ =  N_T / N_C = 1
  • 従来薬C群における治療後活動度スコアの平均値(想定値): \hat{\mu}_C = 60
  • 新規薬T群における治療後活動度スコアの平均値(想定値): \hat{\mu}_T = 50
  • 治療後活動度スコアの標準偏差(想定値): \hat{\sigma} = 10(両群共通)

新規薬Tの効果は、2群の平均値の差  \delta = \mu_T - \mu_S で推定され、効果サイズと呼ばれる。 効果サイズの大小は、スコアのバラツキ具合、つまり標準偏差  \sigma(standard deviation, SD)と比べて相対的に判断する必要があるので、

 
\begin{aligned}
\Delta = \frac{\delta}{\sigma}
\end{aligned}

で計算される「標準効果サイズ(別名:Cohen's d) 」が用いられる。
標準効果サイズは現実的には0.1-1.0の範囲で設定され、目安は以下のとおり。

  • Δ=0.2:効果が小さい
  • Δ=0.5:効果は中等度
  • Δ=0.8:効果が大きい

設定する効果サイズが小さいほど、必要なサンプルサイズは大きくなる。小さい魚をすくおうとすると、細かい目の網が必要になるイメージ。 ただし出来るだけ先行研究をもとにして効果サイズを見積もる方がよい。

方法0:手計算

なるべく正確に計算しようとすると、以下の2点を考慮して計算式を選ばなければならない。

  • 非心t分布(non-central t [NCT] distribution)*2
    対立仮説のもとでの効果サイズの分布は、本来NCT分布に従うので、検出力はNCT分布をもとに計算する必要があるが、自由度(平たく言えばサンプルサイズ)が増えるほど正規分布に近づく。
  • 不等分散:
    標準偏差σが2群で等しいと仮定できないの場合は、そのことを考慮した補正が必要(Satterhwaiteの補正など。他にもあるのかは調べてません。)

上記の2点を考慮するととても複雑な式になってしまうので、「非心t分布は正規分布に近似」かつ「等分散」の条件の下でサンプルサイズを計算するには以下の式を使う。

 
\begin{aligned}
N = \frac{(1+\phi)^2}{\phi}\frac{(z_{1-\alpha/2} + z_{1-\beta})^2}{\Delta ^2} + \frac{z_{1-\alpha/2}^2}{2}
\end{aligned}

末尾の項は「分散既知のz検定」と「分散未知のStudentのt検定」の違いを補う項。有意水準α=0.05とすると、この項目はおおよそ2になるので、分散が未知の場合はサンプルサイズは2例ほど余計に必要になる。

上の式に、

  •  \Delta = \frac{(\mu_T - \mu_C)}{\sigma} = 10/10 = 1.0
  • φ = 1
  •  z_{1-\alpha/2} = 1.96
  •  z_{ 1 - \beta} = 0.8416

を代入すると、N = 33.3...となるので、サンプルサイズは両群で34例(1群で17例)必要なことが分かる。

方法1:pwr.t.test( )を使って計算する

まずはpwrパッケージの関数を使ってみる。この関数でサンプルサイズを計算するときは以下の引数を指定する。

  • d:Cohenのd(標準化効果サイズ)。ここでは(60-50)/10=1を指定。
  • sig.level:有意水準
  • power:検出力
  • type:"two.sample"(2群比較)、"one.sample"(1群の平均を既知の値と比較)、"paired"(対応のある比較)から選ぶ。ここでは"two.sample"を指定。
  • alternative:対立仮説の種類。ここでは"two.sided"(両側検定)を指定。
library(pwr)
pwr.t.test(d = (60-50)/10,
           type = "two.sample",alternative="two.sided",
           sig.level=0.05,
           power=0.8)
     Two-sample t test power calculation 

              n = 16.71472
              d = 1
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

NOTE: n is number in *each* group

1群あたり17例、両群で34例必要という結果になった。

方法2:TwoSampleMean.Equality( )を使って計算する

次にTrialSizeパッケージの関数を使ってみる。
この関数では対立仮説(alternative hypothesis)は関数名に含まれているので、引数として指定する必要はない。*3

指定する引数は以下のとおり。

  • alpha:αエラー
  • beta:βエラー
  • sigma:2群をまとめた標準偏差(pooled SD)
  • k:割り付け比
  • margin:両群の平均の差。標準化していない効果サイズ。ここでは60-50=10を指定。
library(TrialSize)
TwoSampleMean.Equality(alpha = 0.05,
                       beta = 0.2,
                       sigma = 10,
                       k=1,
                       margin = 10)
[1] 15.69776

1群あたり16例、両群で32例必要という結果。方法0との微妙な違いは許容範囲か。

方法3:power.prop.test( )を使って計算する

このシナリオなら基本関数を使っても計算できる。
サンプルサイズ計算のときに指定する主な引数は以下のとおり。

  • delta:両群の平均の差
  • sd:標準偏差
  • sig.level:αエラー
  • power:1-βエラー
  • type:"two.sample"(2群比較)、"one.sample"(1群の平均を既知の値と比較)、"paired"(対応のある比較)から選ぶ。ここでは"two.sample"を指定。
  • alternative:対立仮説。"one.sided"(片側)か"two.sided"(両側)を指定。
power.t.test(delta = 10,
             sd = 10,
             sig.level = 0.05,
             power = 0.8,
             type = "two.sample",
             alternative = "two.sided")
     Two-sample t test power calculation 

              n = 16.71477
          delta = 10
             sd = 10
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

NOTE: n is number in *each* group

1群あたり17例、両群で34例必要という結果になった。

おわりに

  • 今回から連続値アウトカムのサンプルサイズ計算です。あと3回分で一旦お休みします。

参考資料

  • 奥村先生の効果サイズに関する記事です。

oku.edu.mie-u.ac.jp

  • サンプルサイズ計算の本なのにあまり眠くならない。手計算はこの本を参考にしました。

*1:interventionのIを使うと1と区別しにくいのでtrial

*2:t分布を一般化したもの。一般的には中心が0ではなく、左右対称でない。

*3:"equality"は"equivalent"(同等性)と混同しやすいので注意。この研究の仮説は「優越性(superiority)」なんだけど、このパッケージで"superiority"は非劣性試験と対応させた「十分な優越性(substantially superiority)」の意味で使われているので、これも注意。