統計ER

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

反復測定分散分析 群分けあり Split Plot Design のサンプルサイズ計算

Toukei Consul Banner

KH Coder Consul Banner

反復測定分散分析 群分けあり Split Plot Designと呼ばれるデザインの分析におけるサンプルサイズ計算は、いろいろ探してもなかなか出てこない。

ようやく方法を見つけて、やっとのことで読み解いたので、共有したい。

統計ソフトRの関数の形になっているので、コピペで使ってもらえれば幸いである。

反復測定分散分析 群分けあり Split Plot Designとは?

同じ対象者のある指標を反復して測定したときのデータを分析したいことがある。

対象者をいくつかの群に振り分けて比較したいことが多い。

ここでは、この解析方法を「反復測定分散分析 群分けあり」と呼んでいる。

Split Plot Designとも呼ばれる方法だ。

参考となる過去記事はこちら。

toukeier.hatenablog.com

反復測定分散分析 群分けあり Split Plot Design のサンプルサイズ計算はどうやる?

反復測定分散分析 群分けあり Split Plot Design のサンプルサイズ計算はどうやるか?

分散分析のサンプルサイズ計算関数 pwr.anova.test()の一部を抜き出す

まずは、分散分析のサンプルサイズ計算スクリプトを pwr パッケージの pwr.anova.test() から抜き出してみる。

pwr.anova.test.n <- function(k, f, sig.level=0.05, power=0.8)
{
    p.body <- quote({
        lambda <- k*n*f^2
        pf(qf(sig.level, k-1, (n-1)*k, lower=FALSE), k-1, (n-1)*k, lambda, lower=FALSE)
    })

    n <- uniroot(function(n) eval(p.body)-power, c(2+1e-10, 1e+09))$root

    NOTE <- "n is number in each group"
    METHOD <- "Balanced one-way analysis of variance power calculation"
    structure (list(k=k, n=n, f=f, sig.level=sig.level, power=power, note=NOTE, method=METHOD), 
        class="power.htest")
}

例えば、4グループの比較で、effect size (f)が中程度で0.25だったとすると、一群45例必要と計算される。

> pwr.anova.test.n(k=4,f=0.25)

     Balanced one-way analysis of variance power calculation 

              k = 4
              n = 44.59927
              f = 0.25
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

これをもとに、反復測定分散分析 群分けあり Split Plot Design のサンプルサイズ計算のスクリプトを作成してみた。

反復測定分散分析 群分けあり Split Plot Design のサンプルサイズ計算関数を作成した

比較する群の数をp、反復して測定する回数をq、群間比較のeffect sizeをf.A.CR、対象者の測定値における自己相関をrhoとすると、以下のような関数になる。

pwr.split.plot.test.n <- function(p, q, f.A.CR, rho, sig.level=0.05, power=0.8)
{
    p.body <- quote({
        f.A <- f.A.CR/sqrt(1+(q-1)*rho)
        lambda <- p*q*n*f.A^2
        pf(qf(sig.level, p-1, (n-1)*p, lower=FALSE), p-1, (n-1)*p, lambda, lower=FALSE)
    })

    n <- uniroot(function(n) eval(p.body)-power, c(2+1e-10, 1e+09))$root

    NOTE <- "n is number in each group"
    METHOD <- "Split-plot design analysis of variance power calculation"
    structure (list(p=p, q=q, n=n, f.A.CR=f.A.CR, rho=rho,
        sig.level=sig.level, power=power, note=NOTE, method=METHOD), 
        class="power.htest")
}

群の数を4、反復回数も4、f.A.CRを0.25、rhoを0.5とすると、一群29例と計算される。

> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.25, rho=0.5)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 28.2526
         f.A.CR = 0.25
            rho = 0.5
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

rhoが1(完璧な自己相関)の場合、一群45例と計算される。この場合は通常の反復測定しない分散分析と同じになる。完璧に自己相関するなら1回の測定でも同じということだ。

> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.25, rho=1)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 44.59927
         f.A.CR = 0.25
            rho = 1
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

rhoが0(全くの自己相関なし)の場合、一群12例と計算される。これは通常あり得ない。同じ対象者は必ず少しは相関する。例えば、血圧が高めの人は何度測っても高めのはず。低めの人はいつも低め。これは感覚的にわかるだろう。

> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.25, rho=0)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 11.92611
         f.A.CR = 0.25
            rho = 0
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

反復測定分散分析 群分けあり Split Plot Design のサンプルサイズ計算をする際に大事なこと

反復測定分散分析 群分けあり Split Plot Design のサンプルサイズ計算をする際に決めるのが大変なのは、f.A.CRとあるeffect sizeとrhoという対象者内の自己相関だろう。

f.A.CRは、いわゆる群間の標準偏差と、群内の標準偏差の比である。パイロットスタディや先行研究などがあれば参考にできるが、そうでない場合は、困ってしまう。

行動科学分野では、以下のような目安があるので、それに従うこともできる。

  • Small Effect: f=0.1
  • Medium Effect: f=0.25
  • Large Effect: f=0.4

自己相関 rho のほうはもっと困ってしまう。これもパイロットスタディや先行研究のデータがある場合は、自己相関を計算してみてもいいかもしれない。

利用可能なデータがない場合は、rho=c(0.3, 0.5, 0.7)などいくつか振って計算してみるという方法もあるかもしれない。

以下、effect sizeを振ってみた計算結果、rhoを振ってみた計算結果を列挙してみる。

# effect size=0.1
> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.1, rho=0.5)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 171.3325
         f.A.CR = 0.1
            rho = 0.5
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

# effect size=0.25
> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.25, rho=0.5)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 28.2526
         f.A.CR = 0.25
            rho = 0.5
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

# effect size=0.4
> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.4, rho=0.5)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 11.67164
         f.A.CR = 0.4
            rho = 0.5
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

# rho=0.3
> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.25, rho=0.3)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 21.71697
         f.A.CR = 0.25
            rho = 0.3
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

# rho=0.5
> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.25, rho=0.5)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 28.2526
         f.A.CR = 0.25
            rho = 0.5
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

# rho=0.7
> pwr.split.plot.test.n(p=4, q=4, f.A.CR=0.25, rho=0.7)

     Split-plot design analysis of variance power calculation 

              p = 4
              q = 4
              n = 34.79044
         f.A.CR = 0.25
            rho = 0.7
      sig.level = 0.05
          power = 0.8

NOTE: n is number in each group

まとめ

反復測定分散分析 群分けあり Split Plot Design のためのサンプルサイズ計算を紹介した。

pwr パッケージのpwr.anova.test()の一部を抜き出して、修正したものである。

計算に必要はeffect sizeと自己相関を決めねばならず、それが悩ましいポイントだ。

何らかの前提・仮定をおいて、進めていくしかないと思う。

必ずしも先行研究が見つかるわけでもなく、パイロットスタディを行っている時間もないことも多いからである。

参考文献

Bradley, D.R., Russell, R.L. Some cautions regarding statistical power in split-plot designs. Behavior Research Methods, Instruments, & Computers 30, 462–477 (1998).

Some cautions regarding statistical power in split-plot designs | SpringerLink