9.3 Welch t-test의 검정력 구하기

\(X_1 \sim N(10, 2^2)\) 과 \(X_2 \sim N(12, 2^2)\) 으로부터 각각 10개의 난수를 발생시켜 t 검정을 하는 것은 다음과 같다.

x1 = rnorm(10, 10, 2)
x2 = rnorm(10, 12, 2)
t.test(x1, x2)

    Welch Two Sample t-test

data:  x1 and x2
t = -1.6245, df = 13.539, p-value = 0.1273
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -2.9748010  0.4153159
sample estimates:
mean of x mean of y 
 10.28960  11.56934 

만약 위의 것을 Sim.t() 라는 함수를 만들어 수행한다면 다음과 같이 할 수 있다.

Sim.t = function(mu1, mu2, sig1, sig2, n1, n2)
{
  x1 = rnorm(n1, mu1, sig1)
  x2 = rnorm(n2, mu2, sig2)
  t.test(x1, x2)$p.value
}
Sim.t(10, 12, 2, 2, 10, 10)

[1] 0.0006798321

R에서 기본적으로 제공하는 함수 power.t.test()를 이용하여 각 군당 10명, 두 군의 평균의 차이는 2, 표준편차는 2로 동일한 경우에 대한 검정력은 다음과 같이 구할 수 있다.

power.t.test(n=10, delta=2, sd=2)

     Two-sample t test power calculation 

              n = 10
          delta = 2
             sd = 2
      sig.level = 0.05
          power = 0.5619846
    alternative = two.sided

NOTE: n is number in *each* group

위 함수는 각 군의 표본 크기가 다르거나, 표준편차가 다른 경우에는 사용할 수 없다. 따라서, 다음과 같은 simulation 함수를 작성하여 좀 더 일반적인 경우에도 구할 수 있다.

Power.t = function(mu1, mu2, sig1, sig2, n1, n2, N=2000)
{
  Res = rep(NA, N)
  for (i in 1:N) {
    x1 = rnorm(n1, mu1, sig1)
    x2 = rnorm(n2, mu2, sig2)
    Res[i] = t.test(x1, x2)$p.value
  }
  return(sum(Res < 0.05)/N)
}

만약 power.t.test()를 이용하여 풀었던 문제를 우리가 작성한 Power.t()를 이용하여 푼다면 다음과 같다.

Power.t(10, 12, 2, 2, 10, 10)

[1] 0.546

두 군의 크기가 다르거나 표준편차가 다른 경우에도 다음과 같이 구할 수 있다.

Power.t(10, 12, 2, 3, 30, 30)

[1] 0.8415

지금까지는 귀무가설(\(H_0\))이 \(\mu_1 = \mu_2\)이고, 따라서, 대립가설이 \(\mu_1 \neq \mu_2\)인 경우이다. 이 경우 두 평균의 차이가 얼마인지에 대한 설정이 없으므로, 이것만으로는 검정력(power)을 구할 수 없다. 이 때 가정하는 두 군의 차이값에 대한 미지의 진실을 가정하고 이에 따라 검정력을 구하게 된다. 이 때 두 값의 차이에 대한 가정을 ’최소검출자’라고도 한다.

만약 귀무가설이 \(\mu_1 - \mu_2 < 2\) (즉, 대립가설이 \(\mu_1 - \mu_2 \geq 2\)이다.) 가 된다면 위의 R script를 어떻게 바꾸어야 하는지 생각해보라.