BIOINFORMATICS WITH PARK-KLEIS

지난 글에서는 여러번 비교를 할 때 발생할 수 있는 1종 오류 가능성의 증가에 대해 살펴보았다.

2023.02.01 - [STATISTICS] - 다중비교 (Multiple comparisons test) [1]

여러 번 비교를 하게 되면, 1종 오류의 가능성이 증가하게 되는데, 그럼 이를 조정할 필요가 있다.

이걸 예를 들어 다시 간략하게 설명해 보겠다.

A그룹과 B그룹, C그룹이 있고, 세 그룹의 평균을 비교한다고 하자.

그럼,

귀무가설 H0 : A = B = C

대립가설 H1 : H0가 아니다.

대립가설이 성립하려면 ① A ≠ B, ② A ≠ C, ③ B ≠ C 이 셋 중 한 경우라도 만족하면 된다. 

이게 지난 글에서 예시로 들었던 세 개의 문제가 있는 두 과목을 선택하는 것이었다.

두 교과목에서 다루는 세 개의 문제의 난이도가 같다고 할 때,

 - 한 교과목 (아랍어 과목) 은 세 개의 문제를 모두 풀어야 학점을 받을 수 있고,

 - 다른 교과목 (중국어 과목) 은 세 개의 문제 중 한 개라도 풀면 학점을 받을 수 있다고 하면,

어떤 교과목을 선택할 것인가?

웬만한 사람들은 세 개의 문제 중 한 개라도 풀면 학점을 받을 수 있는 교과목인 중국어 과목을 선택할 것이다.

다중비교도 이와 같다.

세 개의 문제를 A, B, C 라고 생각해보면,

① A ≠ B, ② A ≠ C, ③ B ≠ C

이 셋 중 한 경우라도 만족하면 귀무가설이 기각된다는 문제가 생긴다.

만약 두 교과목을 공평하게 만들고자 한다면, 아랍어 과목을 중국어 과목보다 더 쉽게 만들거나, 아니면 중국어 과목을 아랍어 과목보다 훨씬 어렵게 만들어야 한다.

제 1종의 오류 ɑ (실제로 귀무가설이 맞는데, 이를 기각하는 문제)는 중요한 문제이고, 이를 최소화 할 필요가 있다.

이를 해결하기 위한 방법은 여러가지가 있는데, 많이 사용되는 방법 두 가지를 소개하고자 한다.

1. FWER control

먼저 첫 번째 방법은

FWER (Family Wise Error Rate) 을 컨트롤 하는 Bonferroni correction이다.

Bonferroni correction은 매우 직관적이고 간단한 방법이다.

ɑ를 검정하려는 가설의 수로 나누어 주면 끝!

FWER를 ɑ(=0.05)로 유지하기 위해, 만약 가설의 수가 5개라면, ɑ(=0.05)를 5로 나누어주면 된다.

그럼 ① A ≠ B, ② A ≠ C, ③ B ≠ C 이 세 가지 가설 검정에서 p-value가 0.05/5=0.01보다 작아야지만 귀무가설을 기각할 수 있게 된다.

만약 가설의 수가 10개라면, 0.05/10=0.005보다 작아야지만 귀무가설을 기각할 수 있다.

Bonferroni correction의 단점은 매우 conservative 하다는 것인데 이는 power가 떨어진다는 말과 일맥상통한다. 만약 유전자데이터처럼 매우 많은 가설을 검정한다면, 기각을 할 수 있는 p-value가 매우 작아지게 된다.

2. FDR control

두 번째 방법은

FDR (False Discovery Rate) 을 컨트롤 하는 Benjamini-Hochberg procedure 이다.

이 방법은 Bonferroni correction 보다는 덜 conservative한 방법이다.

생각보다 간단하니 다음 과정을 천천히 따라해보면 이해가 될 것이다.

검정을 하려는 m개의 가설이 있다면, m개의 p-value를 가장 작은 것부터, 큰 순서대로 나열한다.

가장 작은 p-value를 p(1)이라고 하고, 그 다음 작은 p-value는 p(2), 이런 순서대로 가장 큰 p(m)까지

p(1), p(2), p(3), ..., p(m) 로 나열한다.

FDR을 ɑ 미만으로 컨트롤 하기 위해서 p(k) < k*ɑ/m을 만족하는 가장 큰 k를 찾는 것이다.

예를 들어, ɑ를 가장 많이 사용되는 유의수준인 0.05라고 하자. (ɑ=0.05)

검정하고자 하는 가설이 다음과 같이 10개 (m=10) 이고,

$$  H_{0A}, H_{0B}, H_{0C}, H_{0D},H_{0E},H_{0F},H_{0G},H_{0H},H_{0I}, H_{0J} $$

각 가설의 p-value는 다음과 같다.

$$  P_{A}=0.034, P_{B}=0.047, P_{C}=0.0001, P_{D}=0.0004, P_{E}=0.028, $$

$$ P_{F}=0.0095, P_{G}=0.003, P_{H}=0.00005, P_{I}=0.0019, P_{J}=0.051 $$

이를 p-value가 작은 것부터 나열하면,

k=1, p = 0.00005

k=2, p = 0.0001

k=3, p = 0.0004

k=4, p = 0.0019

k=5, p = 0.003

k=6, p = 0.0095

k=7, p = 0.028

k=8, p = 0.034

k=9, p = 0.047

k=10, p=0.051 

이와 같이 정리할 수 있다.

p(k) < k*ɑ/m을 만족하는 가장 큰 k 찾는 것이라 했으므로,

p-value가 가장 크게 나온 Hr 가설에 대입해보면, cutoff = k*ɑ/m = k*0.05/10 = k*0.005이다.

k를 가설개수(=m)인 10부터 시작한다.

1. k*ɑ/m = 10*0.005=0.05 이기 때문에 cutoff가 0.05,

    k=10일 때의 귀무가설 하에서 p-value는 p(10) = 0.051이므로, 해당 가설은 기각할 수 없다.

2. k*ɑ/m = 9*0.005=0.045 이기 때문에 cutoff가 0.045,

    k=9일 때의 귀무가설 하에서 p-value는 p(9) = 0.047이므로, 해당 가설은  기각할 수 없다.

3. k*ɑ/m = 8*0.005=0.04 이기 때문에 cutoff가 0.04,

    k=8일 때의 귀무가설 하에서 p-value는 p(8) = 0.034이고, 이는 cutoff인 0.04보다 작으므로, 기각 가능하다.

따라서 가장 큰 k는 8 이고, 그 밑으로는 계산할 필요가 없이 전부 기각 가능하다고 결론 내릴 수 있다.

⇨ p-value가 0.051, 0.047인 가설 J와 가설 B을 제외한 나머지 가설들은 귀무가설을 기각할 수 있다.

위 과정을 R 로 구현하면 다음과 같다.

위 R 결과에서도 보듯이 앞에서부터 8개 가설만이 adjusted p-value가 0.05보다 작음을 알 수 있다.

그럼 다중 검정은 여기에서 끝!