母平均の検定・区間推定のパワーアナリシス     Last modified: Jun 28, 2004

目的

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

使用法

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

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

引数

n           標本サイズ
delta       母平均との差
sd          母標準偏差
power       検出力
sig.level   有意水準
alternative 仮説・信頼区間の種類(指定しない場合には両側検定・両側信頼区間)

ソース

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

# 母平均の検定・区間推定のパワーアナリシス
power.one.sample.t.test <- function( n=NULL,                                 # 標本サイズ
                                        delta=NULL,                             # 母平均との差
                                        sd=NULL,                                # 母標準偏差
                                        sig.level=0.05,                         # 有意水準
                                        power=NULL,                             # 検出力
                                        alternative=c("two.sided", "one.sided"))# 仮説・信頼区間の種類
{
        if (sum(sapply(list(n, delta, sd, power, sig.level), is.null)) != 1) {
                stop("n, delta, sd, power, and sig.level のどれか一つだけが NULL でなければならない")
        }
        alternative <- match.arg(alternative)                                        # 省略された引数の補完
        tside <- switch(alternative, one.sided=1, two.sided=2)
        power.function <- quote(pnorm((sqrt(n)*delta/sd-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(1, 1e7))$root
        }
        else if (is.null(delta)) {
                delta <- uniroot(function(delta) eval(power.function)-power, c(-1e7, 1e7))$root
        }
        else if (is.null(sd)) {
                sd <- uniroot(function(sd) eval(power.function)-power, c(1e-7, 1e7))$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 <- "1 標本 t 検定の検出力"
        structure(list(n=n, delta=delta, sd=sd, sig.level=sig.level, power=power, alternative=alternative, method=METHOD), class="power.htest")
}


使用例・出力結果例

# 必要なサンプルサイズを求める
> power.one.sample.t.test(delta=0.5, sd=sqrt(3), power=0.8, sig=0.05, alt="two")

     Power calculation of the one-sample t-test. 

              n = 94.18656
          delta = 0.5
             sd = 1.732051
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

# 検出力を求める
> power.one.sample.t.test(delta=0.5, sd=sqrt(3), n=46, sig=0.05, alt="two")

     Power calculation of the one-sample t-test. 

              n = 46
          delta = 0.5
             sd = 1.732051
      sig.level = 0.05
          power = 0.4991726
    alternative = two.sided

# 検出できる平均値の差(delta)を求める
> power.one.sample.t.test(sd=sqrt(3), n=46, power=0.8, sig=0.05, alt="two")

     Power calculation of the one-sample t-test. 

              n = 46
          delta = 0.715461
             sd = 1.732051
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

# 検出できる標準偏差を求める
> power.one.sample.t.test(delta=0.5, n=46, power=0.8, sig=0.05, alt="two")

     Power calculation of the one-sample t-test. 

              n = 46
          delta = 0.5
             sd = 1.210445
      sig.level = 0.05
          power = 0.8
    alternative = two.sided

# 有意水準を求める
> power.one.sample.t.test(delta=0.5, sd=sqrt(3), n=46, power=0.8, alt="two", sig.level=NULL)

     Power calculation of the one-sample t-test. 

              n = 46
          delta = 0.5
             sd = 1.732051
      sig.level = 0.2642912
          power = 0.8
    alternative = two.sided

・ 解説ページ


・ 直前のページへ戻る  ・ E-mail to Shigenobu AOKI ( @si.gunma-u.ac.jp )

Made with Macintosh