[SAS] 일원 배치 분산 분석 (One-way ANOVA, ANalysis Of VAriance) - PROC GLM, PROC ANOVA
이전 포스팅에서 음주 여부에 따른 ALT 평균에 차이가 있는지를 알아보기 위해서는 T test를 사용한다는 내용을 알아보았다. 2022.10.04 - [모평균 검정/SAS] - [SAS] 독립 표본 T검정 (Independent samples T-test) - PROC TTEST
음주 여부는 "1) 음주자", "2) 비음주자"로 나뉘는 이분형 변수다. 즉, 음주자의 ALT 평균과 비음주자의 ALT 평균을 비교한 것이었다.
만약, 두 개 이상의 범주로 나뉘는 변수에 대해 모평균의 검정을 하고자 한다면 어떻게 해야 할까? 예를 들어, 흡연 상태를 "1) 비흡연자", "2) 과거 흡연자", "3) 현재 흡연자"로 나누었고, 각 그룹의 폐기능 검사 중 하나인 FVC (Forced Vital Capacity)의 평균에 차이가 있는지 알아보고자 한다.
이때 사용할 수 있는 방법이 일원 배치 분산 분석인 One-way ANOVA (Analysis of Variance)이며 이번 포스팅에서 다뤄볼 주된 내용이다.
*실습용 데이터는 아래 링크를 클릭하면 다운로드할 수 있습니다.
2022.08.04 - [공지사항 및 소개] - 분석용 데이터 (update 22.10.11)
분석용 데이터 (update 22.10.11)
2022년 10월 11일 버전입니다. 변수는 계속하여 추가될 예정입니다. 다음 카테고리에 있는 글에서 이용된 데이터입니다. - 기술 통계 - 범주형 자료 분석 - 모평균 검정 - 반복 측정 자료 분석 - 통계
medistat.tistory.com
시작하기 위해 라이브러리를 만들고, 파일을 불러온다.
라이브러리 만드는 방법: 2022.08.05 - [통계 프로그램 사용 방법/SAS] - [SAS] 라이브러리 만들기 - LIBNAME
파일 불러오는 방법: 2022.08.05 - [통계 프로그램 사용 방법/SAS] - [SAS] 데이터 불러오기 및 저장하기 - PROC IMPORT, PROC EXPORT
*라이브러리 지정하기;
LIBNAME hong "C:/Users/User/Documents/Tistory_blog";
*파일 불러오기;
PROC IMPORT
DATAFILE="C:\Users\user\Documents\Tistory_blog\Data.xlsx"
DBMS=EXCEL
OUT=hong.df
REPLACE;
RUN;
목표: 흡연 상태에 따라 FVC 평균은 서로 다르다고 모집단 수준에서 말할 수 있는가?
전제조건 (정규성, 등분산성)
ANOVA에는 두 가지 전제조건이 필요하다. 본 포스팅의 예시에 맞추어 설명하면 다음과 같다.
1) 정규성: 현재 흡연자, 과거 흡연자, 비흡연자 별로 FVC의 분포를 보았을 때 각각의 분포는 정규성을 따른다.
2) 등분산성: 현재 흡연자, 과거 흡연자, 비흡연자 별로 FVC의 분포를 보았을 때 분산은 서로 같다.
따라서 각각을 검정해야 한다.
1) 정규성
정규성 검정에 관한 내용은 다음 링크에서 확인할 수 있다.
2022.08.12 - [기술 통계/SAS] - [SAS] 정규성 검정 - PROC UNIVARIATE
SMOK: 흡연 상태
FVC_pPRED: FVC의 백분위 예측치
PROC UNIVARIATE DATA=hong.df NORMAL PLOT;
CLASS SMOK;
VAR FVC_pPRED;
HISTOGRAM FVC_pPRED/ NORMAL (MU=EST SIGMA=EST);
RUN;
결과 - 1) 비흡연자 (SMOK=0)
표본의 수가 2,000이 되지 않으므로 Shapiro-Wilk 통계량의 p-value를 보아야 하고, 이는 0.05보다 크다. Q-Q plot, 히스토그램에서 정규성을 띠지 않는다고 할만한 근거가 없으므로 정규성을 따른다고 결론을 짓는다.
결과 - 2) 과거 흡연자 (SMOK=1)
표본의 수가 2,000이 되지 않으므로 Shapiro-Wilk 통계량의 p-value를 보아야 하고, 이는 0.05보다 크다. Q-Q plot, 히스토그램에서 정규성을 띠지 않는다고 할만한 근거가 없으므로 정규성을 따른다고 결론을 짓는다.
결과 - 3) 현재 흡연자 (SMOK=2)
표본의 수가 2,000이 되지 않으므로 Shapiro-Wilk 통계량의 p-value를 보아야 하고, 이는 0.05보다 크다. Q-Q plot, 히스토그램에서 정규성을 띠지 않는다고 할만한 근거가 없으므로 정규성을 따른다고 결론을 짓는다.
따라서 정규성 전제조건은 만족한다고 할 수 있다.
2) 등분산성
등분산성 검정에 관한 내용은 다음 링크에서 확인할 수 있다.
2022.10.04 - [모평균 검정/SAS] - [SAS] 독립 표본 T검정 (Independent samples T-test) - PROC TTEST
2022.10.05 - [모평균 검정/SAS] - [SAS] 등분산성 검정 (Homogeneity of variance) - PROC GLM
위 두 개 글에서 제시한 네 가지 방법 모두를 사용하여 등분산 검정을 시행해 보겠다.
*1) Levene;
PROC GLM DATA=hong.df;
CLASS SMOK;
MODEL FVC_pPRED=SMOK;
MEANS SMOK/ HOVTEST=LEVENE(TYPE=ABS);
RUN;
*2) O'Brien;
PROC GLM DATA=hong.df;
CLASS SMOK;
MODEL FVC_pPRED=SMOK;
MEANS SMOK/ HOVTEST=OBRIEN;
RUN;
*3) Brwon and Forsythe;
PROC GLM DATA=hong.df;
CLASS SMOK;
MODEL FVC_pPRED=SMOK;
MEANS SMOK/ HOVTEST=BF;
RUN;
*4) Bartlett;
PROC GLM DATA=hong.df;
CLASS SMOK;
MODEL FVC_pPRED=SMOK;
MEANS SMOK/ HOVTEST=BARTLETT;
RUN;
QUIT;
결과
1) Levene
2) O'Brien
3) Brown and Forsythe
4) Bartlett
어떤 방법을 사용하더라도 p-value>0.05로 등분산성이라는 전제조건을 만족한다는 것을 알 수 있다. 따라서 일원 배치 분산 분석 (One-way ANOVA)를 시행할 수 있다.
일원 배치 분산 분석 (One way ANOVA) 코드
*1)사후분석- 본페로니;
PROC GLM DATA=hong.df;
CLASS SMOK;
MODEL FVC_pPRED=SMOK;
MEANS SMOK/BON CLDIFF;
RUN;
QUIT;
*2)사후분석- 본페로니;
PROC GLM DATA=hong.df;
CLASS SMOK;
MODEL FVC_pPRED=SMOK;
LSMEANS SMOK/ADJUST=BON PDIFF;
RUN;
QUIT;
*3)사후분석- 던칸;
PROC GLM DATA=hong.df;
CLASS SMOK;
MODEL FVC_pPRED=SMOK;
MEANS SMOK/DUNCAN;
RUN;
QUIT;세 가지 방법으로 코드를 소개했지만, 모두 대동소이하다. 첫 세 줄이 ANOVA 코드로 모두 동일하다. 네 번째 줄은 사후 분석 및 출력 방법에 대한 코드다. 사후 분석에 대해서는 아래에 설명하도록 하겠다.
1) 사후 분석: 본페로니(BONferroni) 방법을 사용하고, 차이의 신뢰구간 (Confidence Limits for DIFFerence)을 나타내라.
PROC GLM DATA=hong.df; : 일반화 선형 모델 (GLM)을 시행하겠다. 데이터는 hong라이브러리의 df파일을 사용하라.
CLASS SMOK; : 흡연 상태에 따라 분산 분석 (ANOVA)를 시행하라.
MODEL FVC_pPRED=SMOK; : 흡연 상태(SMOK)에 따라 FVC % 예측치(FVC_pPRED)의 평균에 차이가 있는지 나타내라.
MEANS SMOK/BON CLDIFF; : 흡연 상태(SMOK)에 따른 FVC % 예측치(FVC_pPRED)의 평균을 보여달라. 사후 분석은 본페로니(BONferroni, BON) 방법을 사용하고, 차이의 신뢰구간 (Confidence Limits for DIFFerence, CLDIFF)을 보여달라
- CLDIFF를 사용하기 위해서는 "MEANS"로 시작하는 구문을 사용해야 한다.
2) 사후 분석: 본페로니(BONferroni) 방법을 사용하고, 차이의 p-value (P-value for DIFFerence)를 나타내라.PROC GLM DATA=hong.df; : 일반화 선형 모델 (GLM)을 시행하겠다. 데이터는 hong라이브러리의 df파일을 사용하라.CLASS SMOK; : 흡연 상태에 따라 분산 분석 (ANOVA)를 시행하라.MODEL FVC_pPRED=SMOK; : 흡연 상태(SMOK)에 따라 FVC % 예측치(FVC_pPRED)의 평균에 차이가 있는지 나타내라.
LSMEANS SMOK/ADJUST=BON PDIFF; : 흡연 상태(SMOK)에 따른 FVC % 예측치(FVC_pPRED)의 최소 제곱 평균(Least Squares MEANS, LSMEANS)을 보여달라. 사후 분석은 본페로니(BONferroni, BON) 방법을 사용하고, 차이의 p-value (P-value for DIFFerence)를 보여달라
- PDIFF를 사용하기 위해서는 "LSMEANS"로 시작하는 구문을 사용해야 한다.
3) 사후 분석: 던칸(DUNCAN) 방법을 사용하라.PROC GLM DATA=hong.df; : 일반화 선형 모델 (GLM)을 시행하겠다. 데이터는 hong라이브러리의 df파일을 사용하라.CLASS SMOK; : 흡연 상태에 따라 분산 분석 (ANOVA)를 시행하라.MODEL FVC_pPRED=SMOK; : 흡연 상태(SMOK)에 따라 FVC % 예측치(FVC_pPRED)의 평균에 차이가 있는지 나타내라.
MEANS SMOK/DUNCAN; : 흡연 상태(SMOK)에 따른 FVC % 예측치(FVC_pPRED)의 평균을 보여달라. 사후 분석은 던칸(DUNCAN) 방법을 사용하라.
결과
위 세 개의 코드에서 공통적인 부분 (ANOVA)의 결과는 다음과 같다.
P-value가 <0.0001이므로 귀무가설을 기각하고 대립 가설을 채택한다. 그렇다면 여기에서 귀무가설 및 대립 가설은 무엇이었는가?
귀무가설: $H0=$ 세 집단의 모평균은 "모두" 동일하다.
대립 가설: $H1=$ 세 집단의 모평균이 모두 동일한 것은 아니다.
우리는 대립 가설을 채택해야 하므로 "세 집단의 모평균이 모두 동일한 것은 아니다."라고 결론 내릴 것이다.
이 말은, 세 집단 중 두 개씩 골라 비교했을 때, 적어도 한 쌍에서는 차이가 난다는 것이다. 따라서 세 집단 중 두 개씩 골라 비교를 해보아야 하며, 이를 사후 분석 (post hoc analysis)이라고 한다. 세 집단에서 두 개씩 고르므로 가능한 경우의 수는 $_3 C_2=3$이다.
(1) 비흡연자 vs 과거 흡연자
(2) 비흡연자 vs 현재 흡연자
(3) 과거 흡연자 vs 현재 흡연자
사후 분석 결과 부분을 보겠다.
1) 사후 분석으로 본페로니(BONferroni) 방법을 사용하고, 차이의 신뢰구간 (Confidence Limits for DIFFerence)을 나타내라.
빨간 박스를 먼저 설명하면, 비흡연자 (SMOK=0)가 과거 흡연자 (SMOK=1)보다 FVC가 9.5328%p 더 크다. 그 차이의 95% 신뢰구간은 7.2168에서 11.8499이다. 즉 신뢰구간에 0이 포함되어있지 않으므로 유의미한 차이가 있다고 할 수 있다. 그런 경우 노란 박스의 "***" 표시를 해준다. 밑의 결과도 해석은 같다.
즉, 모든 사후 분석 결과에서 유의미한 차이가 있음을 알 수 있다.
2) 사후 분석: 본페로니(BONferroni) 방법을 사용하고, 차이의 p-value (P-value for DIFFerence)를 나타내라.
위의 방법 1)에서는 95% 신뢰구간만을 보여주고, p-value가 0.05보다 작은지만 보여줄 뿐, 정확한 p-value의 값은 알려주지 않는다. 정확한 p-value값을 알고 싶은 경우 LSMEANS구문에 PDIFF옵션을 사용하면 된다.
예를 들어 위의 빨간 박스는 그룹 1 (비흡연자)와 그룹 2(과거 흡연자) 사이의 평균의 차이에 대한 p-value를 보여준다. 위와 다르게 p-value를 보여준다는 장점은 있지만, 그 차이가 얼마인지, 신뢰구간은 얼마인지 보여주지 못하는 단점이 있다. 둘 중 어떤 방법으로 출력을 할지는 개인이 선택을 하면 된다.
3) 사후 분석: 던칸(DUNCAN) 방법을 사용하라.
이 시점에서 이런 의문이 들 수 있다.
"각각의 그룹별로 평균을 비교하면 되지, 굳이 왜 ANOVA라는 방법까지 사용하는 것인가?"
아주 논리적인 의문점이다. 하지만, 반드시 ANOVA를 사용해야 한다. 그 이유는 다음과 같다. 본 사례는 흡연 상태에 따른 조합 가능한 경우의 수가 3인데, 각각 유의성의 기준을 $\alpha=0.05$로 잡아보자. 이때 세 번의 비교에서 모두 귀무가설이 학문적인 진실인데(평균에 차이가 없다.), 세 번 모두 귀무가설을 선택할 확률은 $(1-0.95)^3 \approx 0.86$이다. (이해가 어려우면 p-value에 대한 설명 글을 읽고 오길 바란다. 2022.09.05 - [통계 이론] - [이론] p-value에 관한 고찰)
그런데, ANOVA의 귀무가설은 "모든 집단의 평균이 같다."이다. 따라서 모든 집단의 평균이 같은 것이 학문적 진실일 때, 적어도 한 번이라도 대립 가설을 선택하게 될 확률은 $1-0.86=0.14$가 된다. 즉, 유의성의 기준이 올라가게 되어, 덜 보수적인 결정을 내리게 되고, 다시 말하면 대립 가설을 잘 선택하는 쪽으로 편향되게 된다. 학문적으로는 '다중 검정 (multiple testing)을 시행하면 1종 오류가 발생할 확률이 증가하게 된다.'라고 표현한다.
따라서, 각각을 비교해보는 것이 아니라 한꺼번에 비교하는 ANOVA를 시행해야 함이 마땅하다.
여기에서 한 번 더 의문이 들 수 있다.
"사후 분석을 할 때에는 1종 오류가 발생하지 않는가?"
그렇다. 1종 오류가 발생할 확률이 있으므로, p-value의 기준을 더 엄격하게 (0.05보다 더 작게) 잡아야 한다. P-value를 보정하는 방법은 여러 가지가 나와있는데, 가장 간단하고 제일 보수적인 방법이 위에 예시에서 사용한 본페로니 방법이다. 그 외에도 많은 방법이 있고 다른 방법을 사용하기 위해서는 BON의 자리에 다음 중 하나를 입력하면 된다.
| 코드 | 내용 | 코드 | 내용 |
| BON | Bonferroni t test | SCHEFFE | Scheffe's multiple comparison |
| DUNCAN | Duncan's multiple range test | SIDAK | Pairwise t test (Sidak's inequality) |
| DUNNETT | Dunnett's two-tailed t test | SMM (혹은 GT2) | Pairwise comparisons (studentized maximum modulus, Sidak's uncorrelated t-inequality) |
| DUNNETTL | Dunnett's one-tailed t test (less) | SNK | Student-Newman-Keuls |
| DUNNELLU | Dunnett's one-tailed t test (greater) | T (혹은 LSD) | Pairwise t test |
| GABRIEL | Gabriel's multiple comparison | TUKEY | Tukey's studentized range test |
| REGWQ | Ryan-Einot-Habriel-Welsch | WALLER | Waller-Duncan k-ratio t test |
굉장히 많은 방법이 있지만, Bonferroni, Duncan, Scheffe, Tukey가 가장 많이 쓰이는 듯하다. 이 중에 Duncan, Tukey는 그룹별로 표본의 수가 같을 때만 사용해야 하므로 유의하도록 한다.
그중에 DUNCAN옵션을 사용하여 산출한 결과는 다음과 같다.
막대기가 연결되어 있으면 평균이 서로 다르다고 할 수 없는 것이고, 끊어져 있으면 서로 다르다고 할 수 있는 것이다. 위의 경우 비흡연자 (SMOK=0)은 파란 막대기가 있고, 과거 흡연자 (SMOK=1)는 빨간 막대기로 파란 막대기와 분리되어 있다. 따라서 과거 흡연자의 FVC는 비흡연자의 FVC와 유의하게 다르다고 할 수 있다. 현재 흡연자 (SMOK=2)는 초록 막대기로 파랑 및 빨간 막대기와 분리되어 있다. 따라서 현재 흡연자의 FVC는 비흡연자 및 과거 흡연자의 FVC와 유의하게 다르다고 할 수 있다. 즉 모든 그룹에서 FVC의 평균은 다르다고 할 수 있다.
어떤 사후 분석을 쓸 것인가
이 논의에 대해 정답이 따로 있는 것은 아니다. 적절한 방법을 사용하여 논문에 제시하면 되고, 어떤 것이 정답이라고 콕 집어 이야기할 수는 없다. 다만, 사후 분석 방법이 여러 가지가 있다는 것은 '사후 분석 방법에 따라 산출되는 결과가 달라질 수 있다.'는 것을 의미하고, 심지어는 '어떤 사후 분석 방법을 채택하냐에 따라 유의성 여부가 달라질 수도 있다.'는 것을 의미한다. 심지어, ANOVA에서는 유의한 결과가 나왔는데, 사후 분석을 해보니 유의한 차이를 보이는 경우가 없을 수도 있다. 따라서 어떤 사후 분석 방법 결과에 따른 결과인지 유의하여 해석할 필요가 있다.
PROC ANOVA
PROC ANOVA에 대한 설명을 듣고자 본 포스팅에 들어온 독자도 있을 것이다. 하지만 이는 다루지 않을 것이다. 왜냐하면 PROC ANOVA는 그룹별로 표본의 수가 같아야만 사용할 수 있는, 아주 제한적인 상황에서만 쓸 수 있는 코드이기 때문이다. 따라서 웬만하면 PROC GLM을 사용하길 권한다. 만약 사용하고자 한다면, PROC GLM을 PROC ANOVA로만 바꾸어 다음과 같이 사용하면 된다.
PROC ANOVA DATA=hong.df;
CLASS SMOK;
MODEL FVC_pPRED=SMOK;
MEANS SMOK;
RUN;SAS 일원 배치 분산 분석 (One-way ANOVA, ANalysis Of VAriance) 정복 완료!
작성일: 2022.10.12.
최종 수정일: 2022.11.30.
이용 프로그램: SAS v9.4
운영체제: Windows 10