이전 장에서 피셔의 유의성 검정과 네이만-피어슨의 통계적 가설검정을 살펴보았습니다. 네이만-피어슨의 검정에서 $p$ 값이 직접 사용되지는 않지만, 적어도 검정통계량이 기각역에 속하는지 여부를 판단하기 위해 사용될 수 있다는 것을 보았습니다.
$p$-value가 갖는 특징 몇가지를 알아봅시다.
아래의 내용은 Daniël Lakens의 「Improving Your Statistical Inferences」를 바탕으로 여러 문헌을 공부하며 정리한 내용입니다.
저는 통계학 전공자가 아닙니다. 잘못 서술한 내용이 있을 수 있으니 비판적으로 읽어주길 바라요.
오류를 발견하면 언제든 어디서든 꼭! 알려주세요.
글을 인용하실 때에는 이 글의 출처와 함께 원문의 서지 정보를 필히 밝혀주세요. 연구, 아이디어, 주장 등은 한 사람이 뚝딱 생각하여 만드는 것이 아니라 선행 문헌들 위에서 발전되어 온 것들입니다 (a.k.a. 거인의 어깨). 원 문헌을 저술한 저자를 명시하여 그들이 받아 마땅한 공을 밝히도록 합시다.
Null Statistical Hypothesis Testing 정리하기 (2/4)
1장. 통계적 가설검정의 개념
2장. p-value의 특징 ←
3장. p-value를 제대로 보고하기
4장. p-value와 error rate에 대한 오해와 물음
2-1. P-Value Distribution
어떤 실험 1회에서 얻은 표본 데이터로부터 $p$ 값을 얻을 수 있습니다. 같은 실험을 셀 수 없을만큼 매우 많이 시행하여 얻은 $p$ 값은 어떨까요? 실험 하나에서 얻은 $p$ 값이 같은 실험을 무한히 많이 시행했을 때에도 비슷하게 유지될까요? 영가설이 실제로 거짓인 실험을 무한히 많이 시행했을 때 낮은 $p$-value가 많이 나올까요? 또는 영가설이 실제로 참인 실험을 무한히 많이 시행했을 때 높은 $p$-value가 많이 나올까요?
우리는 이를 시뮬레이션을 통해 알아볼 수 있습니다. 영가설이 참인 데이터를 만들어 놓고, 그 데이터를 통해 계산된 $p$ 값이 어떤 경향성을 갖는지 확인할 수 있는 것이지요.
아래에서는 시뮬레이션을 통해 각각의 상황에서 $p$ 값이 갖는 분포를 확인하겠습니다.
- 영가설 ($H_0$)이 실제로 거짓인 경우: 낮은 p가 나오는 빈도가 높습니다.
시뮬레이션 조건은 다음과 같습니다:
- 100,000번 반복.
- $H_0:\mu_1 = \mu_2 = 100$, $\sigma_1 = \sigma_2 = 15$, $n_1 = n_2 = 140$
- Student’s $t$-test, power $= 0.80$
# Adapted from Lakens (2024), section 1.4
p <- numeric(100000) # store all simulated *p*-values
# When H_0 FALSE
for (i in 1:100000) { # for each simulated experiment
x <- rnorm(n = 140, mean = 100, sd = 15) #
y <- rnorm(n = 140, mean = 105, sd = 15) # H_0 FALSE
p[i] <- t.test(x, y, var.equal = TRUE)$p.value # store the *p*-value
}
(sum(p < 0.05) / 100000) # compute power
hist(p, breaks = 20) # plot a histogram
위 히스토그램에서 확인할 수 있듯, 영가설이 실제로 거짓인 경우 낮은 $p$ 값이 많이 나온다는 점을 알 수 있습니다.
- 영가설 ($H_0$)이 실제로 참인 경우: p값이 나오는 빈도가 동일합니다.
시뮬레이션 조건은 다음과 같습니다:
- 100,000번 반복.
- $H_0:\mu_1 = \mu_2 = 100$, $\sigma_1 = \sigma_2 = 15$, $n_1 = n_2 = 140$
-
Student’s $t$-test, $\alpha = 0.05$
즉, 영가설이 실제로 참일 때는 어떤 p값이든 동일한 확률로 나올 수 있습니다.
- 영가설이 실제로 참일 때 이를 잘못 기각하는 오류가 나올 확률, 즉 1종오류율은 $\alpha = 0.05$입니다.
-
위 시뮬레이션의 빨간 점선은 $y=5000$을 나타내는 선입니다. 실제로 $100,000$번 중 $5,000$번 즉 전체의 $0.05\%$에서 $p$가 통계적으로 유의하지 않은 결과를 확인할 수 있습니다.즉, $p<0.05$인 경우의 수는 전체의 $0.05\%$ 이하인 것이죠.
- 예컨대 $0.04 < p \le 0.05$일 확률은 $0.89 < p \le 0.90$일 확률과 같다는 말이에요!
- 영가설이 참일 때 $p$ 값은 균등분포uniform distribution를 따릅니다.
- 단, 이러한 분포는 Pearson’s $\chi^2$ test나 Fisher’s exact test와 같이 검정통계량이 discrete한 경우에는 적용되지 않는 것으로 알려져 있습니다 (Wang et al., 2019; Murdoch et al., 2008).
2-2. Lindley’s Paradox
$p$ 값에 관한 또 다른 성질은 린들리 역설Lindley’s Paradox과 관련이 있습니다. 린들리 역설은 제프리 역설(Jeffreys’s paradox)이나 제프리-린들리 역설 (Jeffreys-Lindley paradox)로 불리기도 합니다. 린들리 역설에 대한 개념은 Spanos (2013)의 글을, 빈도주의frequentist 검정에서 린들리 역설이 이해되는 방식에 관해서는 Maier & Lakens (2022)의 글을 주로 참고하여 정리하였습니다.
먼저 린들리가 제시한 ‘Large $n$ problem’에서 시작합시다 (Lindley, 1957). 평균 $\theta \in (-\infty,\infty)$와 분산 $\sigma^2 > 0$를 갖는 independent & identically distributed simple normal model (이하 NIIDNormal Independent Identically Distributed model) $X_k$, 즉 \(X_k \sim NIID(\theta, \sigma^2), k= 1, 2,\ldots, n, \ldots\)에 대해, large $n$ problem은 다음과 같습니다:
Large $n$ problem
$\alpha$를 유의수준으로 하는 빈도주의 검정은 어떤 점 영가설point null hypothesis (e.g. $H_0:\theta = \theta_0$)도 기각되도록 하는 충분히 큰 표본 크기 $n$을 항상 가진다는 오류에 취약하다.
린들리는 large $n$ problem이 베이지안 가설검정Bayesian hypothesis test의 결과와 충돌하는 역설적인 상황을 다음과 같이 설명합니다.
Jeffreys-Lindley paradox
어떤 prior에 대해서, $n\rightarrow\infty$이면 ($H_0$에 대한 posterior probability) $\rightarrow 1$이다.
즉, $n$이 충분히 클 때 빈도주의 접근에서 기각되는 영가설이 베이지안 접근에서 기각되지 않을 수 있는 역설적 상황이 벌어지는 것이다.
동일한 영가설이 빈도주의 검정에서 기각되는데 베이지안 검정에서는 기각되지 않는다니!
한편 린들리 역설은 빈도주의 검정에서 다음과 같이 이해될 수도 있습니다1 (Maier & Lakens, 2022):
Jeffreys-Lindley paradox in Frequentist Test
검정력이 증가할수록 효과가 있을 때보다 없을 때 $p < \alpha = 0.05$인 $p$-value가 나올 가능성이 높아지는 역설적인 현상을 말한다.
- 검정력이 높은 상황 (e.g. 표본 크기 $n$이 클 때)에 대부분의 $p$ 값이 $\alpha$ 값 이하에 집중되며, 작은 $p$ 값 (e.g. $p<0.05$)은 효과가 있을 때 보다도 없을 때 관찰될 가능성이 더 높습니다.
- $p$ 값이 $\alpha$뿐만 아니라 검정력과도 함께 해석되어야 하는 이유일 것 같습니다.
- 다음은 제가 작성한 간단한 시뮬레이션 결과를 바탕으로 그린 그래프입니다:
각각 검정력이 $50\%$, $80\%$, $99\%$이도록 한 $H_1$이 참인 데이터와 $H_0$가 참인 데이터, 총 4개 데이터에서 Student $t$-test의 결과로 얻은 $p$ 값의 분포가 어떻게 다른지 보고자 하였습니다. 표본 크기는 $d = 0.5$와 각각의 검정력을 만족하도록 정하였고요. 각 시행을 $100,000$번 반복하여 검정을 $100,000$번 진행한 후 $p$ 값을 얻었습니다. 얻은 $p$ 값의 히스토그램을 ggplot2 package의 geom_density() 함수를 통해 추정한 density로 나타냈습니다. geom_density()는 stat::density()를 사용하는데, stat::density()는 기본값으로 Gaussian kernel density estimation을 통해 추정합니다. 사용한 코드를 간단히 알아보겠습니다.
먼저 시뮬레이션 결과를 p_null, p_99, p_80, p_50에 저장하고, 이를 p_power에 tibble 타입으로 저장합니다.
코드 보기
# alpha = 0.05, H0 TRUE
p_null <- numeric(100000)
for (i in 1:100000) { # for each simulated experiment
x <- rnorm(n = 100, mean = 0, sd = 1) # Simulate data
y <- rnorm(n = 100, mean = 0, sd = 1) # Simulate data
p_null[i] <- t.test(x, y, var.equal = TRUE)$p.value # store the *p*-value
}
(sum(p_null < 0.05) / 100000) # compute alpha
# d = 0.5, 99% power, H1 TRUE
p_99 <- numeric(100000)
for (i in 1:100000) { # for each simulated experiment
x <- rnorm(n = 150, mean = 0, sd = 1) # Simulate data
y <- rnorm(n = 150, mean = 0.5, sd = 1) # Simulate data
p_99[i] <- t.test(x, y, var.equal = TRUE)$p.value # store the *p*-value
}
(sum(p_99 < 0.05) / 100000) # compute power
# d = 0.5, 80% power, H1 TRUE
p_80 <- numeric(100000)
for (i in 1:100000) { # for each simulated experiment
x <- rnorm(n = 64, mean = 0, sd = 1) # Simulate data
y <- rnorm(n = 64, mean = 0.5, sd = 1) # Simulate data
p_80[i] <- t.test(x, y, var.equal = TRUE)$p.value # store the *p*-value
}
(sum(p_80 < 0.05) / 100000) # compute power
# d = 0.5, 50% power, H1 TRUE
p_50 <- numeric(100000)
for (i in 1:100000) { # for each simulated experiment
x <- rnorm(n = 32, mean = 0, sd = 1) # Simulate data
y <- rnorm(n = 32, mean = 0.5, sd = 1) # Simulate data
p_50[i] <- t.test(x, y, var.equal = TRUE)$p.value # store the *p*-value
}
(sum(p_50 < 0.05) / 100000) # compute power
p_power <- tibble(p_null = p_null, p_50 = p_50, p_80 = p_80,
p_99 = p_99, blank = numeric(100000))이후 결과를 ggplot2 패키지를 통해 그립니다.
코드 보기
cols <- c(p_99 = "black", p_80 = "purple", p_50 = "skyblue",p_null = "red")
p_density <- ggplot(p_power) +
geom_density(aes(x = p_99, color = "p_99"), linewidth = 0.7, bw = bw_99) +
geom_density(aes(x = p_80, color = "p_80"), linewidth = 0.7, bw = bw_80) +
geom_density(aes(x = p_50, color = "p_50"), linewidth = 0.7, bw = bw_50) +
geom_density(aes(x = p_null, color = "p_null"), linewidth = 0.7, bw = bw_null,
linetype = "dashed") +
...이때 geom_density()의 bw 파라미터는 bandwidth를 의미하는데, bandwidth는 Scott의 방법
으로 정하였습니다2.
코드 보기
# Bandwidth parameter in Gaussian KDE: Scott's bandwidth
bw_50 <- 1.06 * sd(p_50) * length(p_50)^(-0.2)
bw_80 <- 1.06 * sd(p_80) * length(p_80)^(-0.2)
bw_99 <- 1.06 * sd(p_99) * length(p_99)^(-0.2)
bw_null <- 1.06 * sd(p_null) * length(p_null)^(-0.2)결론
Long-term 시뮬레이션을 통해 $p$ 값이 지니는 두 가지 특성을 알아보았습니다.
- 영가설이 실제로 거짓일 때는 낮은 $p$ 값이 더 많이 나오지만, 영가설이 실제로 참일 때에는 $p$가 어떤 값이든 비슷한 확률로 나올 수 있습니다.
- 검정력이 매우 큰 상황에서 영가설이 거짓일 때보다 참일 때 작은 $p$가 나올 가능성이 높아질 수 있습니다.
Reference
- Lindley, D. V. (1957). A Statistical Paradox. Biometrika, 44(1/2), 187-192. https://doi.org/10.2307/2333251
- Maier, M., & Lakens, D. (2022). Justify Your Alpha: A Primer on Two Practical Approaches. Advances in Methods and Practices in Psychological Science, 5(2), 25152459221080396. https://doi.org/10.1177/25152459221080396
- Murdoch, D. J., Tsai, Y.-L., & Adcock, J. (2008). P-Values are Random Variables. The American Statistician, 62(3), 242-245. https://doi.org/10.1198/000313008X332421
- Spanos, A. (2013). Who Should Be Afraid of the Jeffreys-Lindley Paradox? Philosophy of Science, 80(1), 73-93. https://doi.org/10.1086/668875
- Wagenmakers, E.-J., & Ly, A. (2023). History and nature of the Jeffreys–Lindley paradox. Archive for History of Exact Sciences, 77(1), 25-72. https://doi.org/10.1007/s00407-022-00298-3
- Wang, B., Zhou, Z., Wang, H., M. Tu, X., & Feng, C. (2019). The p-value and model specification in statistics. General Psychiatry, 32(3), e100081. https://doi.org/10.1136/gpsych-2019-100081
Leave a comment