目的

相関係数の検定・区間推定を指定する条件で行うために必要な標本サイズを計算するなどのパワーアナリシスを行う。

使用法

power.cor.test(n=NULL, cor0=0, cor1=NULL, sig.level=0.05, power=NULL, alternative=c("two.sided", "one.sided"))

n, delta, sd, power, sig.level のどれか一つだけを NULL として指定して関数を呼び出すと，そのパラメータの値を求めることができる（sig.level を求めるときには，明示的に NULL を指定する必要がある）。

引数

n           標本サイズ
cor0        母相関係数
cor1        標本相関係数
sig.level   有意水準
power       検出力
alternative 仮説・信頼区間の種類（指定しない場合には両側検定・両側信頼区間）

ソース

インストールは，以下の 1 行をコピーし，R コンソールにペーストする
source("http://aoki2.si.gunma-u.ac.jp/R/src/power-cor-test.R", encoding="euc-jp")

# 相関係数の検定・区間推定のパワーアナリシス
power.cor.test <- function(  n=NULL,                                         # 標本サイズ
                                cor0=0,                                         # 母相関係数
                                cor1=NULL,                                      # 標本相関係数
                                sig.level=0.05,                                 # 有意水準
                                power=NULL,                                     # 検出力
                                alternative=c("two.sided", "one.sided"))        # 仮説・信頼区間の種類
{
        if (sum(sapply(list(n, cor0, cor1, power, sig.level), is.null)) != 1) {
                stop("n, delta, sd, power, sig.level のどれか一つだけを NULL にする")
        }
        alternative <- match.arg(alternative)                                        # 引数の補完
        tside <- switch(alternative, one.sided=1, two.sided=2)
        power.function <- quote(pnorm(sqrt(n-3)*abs(atanh(cor0)-atanh(cor1))-qnorm(sig.level/tside, lower.tail=FALSE)))
        if (is.null(power)) {
                power <- eval(power.function)
        }
        else if (is.null(n)) {
                n <- uniroot(function(n) eval(power.function)-power, c(4, 1e7))$root
        }
        else if (is.null(cor0)) {
                cor0 <- uniroot(function(cor0) eval(power.function)-power, c(cor1, 1))$root
        }
        else if (is.null(cor1)) {
                cor1 <- uniroot(function(cor1) eval(power.function)-power, c(cor0, 1))$root
        }
        else if (is.null(sig.level)) {
                sig.level <- uniroot(function(sig.level) eval(power.function)-power, c(1e-5, 0.99999))$root
        }
        else {
                stop("internal error")
        }
        METHOD <- "Power calculation of the one-sample correlation test."
        structure(list(n=n, cor0=cor0, cor1=cor1, sig.level=sig.level, power=power, alternative=alternative, method=METHOD), class="power.htest")
}


使用例・出力結果例

# 母相関係数=0，標本相関係数=0.4 のとき，サンプルサイズを求める
> power.cor.test(sig.level=0.05, power=0.8, cor0=0.0, cor1=0.4)

     Power calculation of the one-sample correlation test. 

              n = 46.73161
           cor0 = 0
           cor1 = 0.4
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

# 母相関係数=0.3，標本相関係数=0.4 のとき，サンプルサイズを求める
> power.cor.test(sig.level=0.05, power=0.8, cor0=0.3, cor1=0.4)

     Power calculation of the one-sample correlation test. 

              n = 605.5779
           cor0 = 0.3
           cor1 = 0.4
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

# 検出力を求めるとき
> power.cor.test(n=605, sig.level=0.05, cor0=0.3, cor1=0.4)

     Power calculation of the one-sample correlation test. 

              n = 605
           cor0 = 0.3
           cor1 = 0.4
      sig.level = 0.05
          power = 0.7996236
    alternative = two.sided

# 有意水準を求めるとき
> power.cor.test(n=605, power=0.8, cor0=0.3, sig.level=NULL, cor1=0.4)

     Power calculation of the one-sample correlation test. 

              n = 605
           cor0 = 0.3
           cor1 = 0.4
      sig.level = 0.05016064
          power = 0.8
    alternative = two.sided

 解説ページ