반응형

[R] 일표본 윌콕슨 부호 순위 검정 (비모수 일표본 중앙값 검정: One-Sample Wilcoxon Signed Rank Test) - wilcox.test()

 어떤 분포의 평균이 특정 값인지 확인하는 방법을 이전에는 일표본 T검정 (One-Sample T test)로 시행했었다. (2022.11.03 - [모평균 검정/R] - [R] 일표본 T검정 (One-sample T-test) - t.test()) 하지만 여기에는 중요한 가정이 필요한데, 분포가 정규성을 따르는 것이다. 하지만 분포가 정규성을 따르지 않는다면 어떻게 해야 할까? 그럴 때 사용하는 것이 One-Sample Wilcoxon Signed Rank Test (일표본 윌콕슨 부호 순위 검정)이다. 

 

 이번 포스팅에서는 One-Sample Wilcoxon Signed Rank Test (일표본 윌콕슨 부호 순위 검정)에 대해서 알아볼 것이다.

 

*실습용 데이터는 아래 링크를 클릭하면 다운로드할 수 있습니다.

2022.08.04 - [공지사항 및 소개] - 분석용 데이터 (update 22.11.28)

 

분석용 데이터 (update 22.11.28)

2022년 11월 28일 버전입니다. 변수는 계속하여 추가될 예정입니다. 다음 카테고리에 있는 글에서 이용된 데이터입니다. - 기술 통계 - 범주형 자료 분석 - 모평균 검정 - 반복 측정 자료 분석 - 통계

medistat.tistory.com

 

코드를 보여드리기에 앞서 워킹 디렉토리부터 지정하겠다.

워킹 디렉토리에 관한 설명은 다음 링크된 포스트에서 볼 수 있다.

2022.08.05 - [통계 프로그램 사용 방법/R] - [R] 작업 디렉토리 (Working Directory) 지정 - getwd(), setwd()

 

setwd("C:/Users/user/Documents/Tistory_blog")

 

데이터를 불러와 df에 객체로 저장하겠다.

데이터 불러오는 방법은 다음 링크에서 볼 수 있다.

2022.08.05 - [통계 프로그램 사용 방법/R] - [R] 데이터 불러오기 : EXCEL - read_excel(), read.xlsx()

2022.08.10 - [통계 프로그램 사용 방법/R] - [R] 데이터 저장하기 : CSV 파일 - write.csv(), write_csv()

2022.08.10 - [통계 프로그램 사용 방법/R] - [R] 데이터 불러오기 : SAS file (.sas7bdat) - read.sas7bdat(), read_sas()

 

install.packages("readr")
library("readr")
df<-read_csv("Data.csv")

 

목표: 모집단에서 신경심리검사 1 원점수의 중앙값이 55이라고 할 수 있는가?

 이번 포스팅의 목적은 1000명의 데이터를 가지고, 이 1000명이 기원한 모집단에서 신경심리검사 1 원점수의 중앙값이 55이라고 할 수 있는지 판단하는 것이다.

 

전제조건 (정규성) - 코드

 만약 정규성을 따른다면 t-test를 하면 되므로 정규성 여부를 파악하도록 한다. 정규성 검정에 관한 분석 내용은 다음 글에서 살펴볼 수 있다. 

2022.08.11 - [기술 통계/R] - [R] 정규성 검정 (1) : Q-Q plot - qqnorm(), qqline()

2022.08.11 - [기술 통계/R] - [R] 정규성 검정 (2) : 히스토그램 - hist(), dnorm()

2022.08.11 - [기술 통계/R] - [R] 정규성 검정 (3) : 정량적 검정 (Shapiro-Wilk, Kolmogorov-Smirnov) - shapiro.test(), ks.test()

2022.08.12 - [기술 통계/R] - [R] 정규성 검정 (4) : 정량적 검정 (Lilliefors test) - lillie.test()

2022.08.16 - [기술 통계/R] - [R] 고급 Q-Q Plot - Van der Waerden, Rankit, Tukey, Blom

 

# 1) Q-Q plot 그리기
qqnorm(df$NP1_raw)
qqline(df$NP1_raw)
# 2) 히스토그램 그리기
hist(df$NP1_raw, prob=TRUE)
NP1_rawrange<-seq(min(df$NP1_raw),max(df$NP1_raw),length=max(max(df$NP1_raw)-min(df$NP1_raw),100))
ND<-dnorm(NP1_rawrange,mean=mean(df$NP1_raw),sd=sd(df$NP1_raw))
lines(NP1_rawrange, ND, lwd=2)
# 3) Shapiro-Wilk test 시행하기
shapiro.test(df$NP1_raw)

결과

	Shapiro-Wilk normality test

data:  df$NP1_raw
W = 0.97015, p-value = 1.81e-13

 N수가 2,000개 미만이므로 Shapiro-Wilk 통계량의 p-value를 보면 0.05 이하이며, Q-Q Plot은 대부분의 데이터가 선상에 있지 않고, 히스토그램에서도 정규성을 따르지 않는 것처럼 보인다. 따라서 일표본 T검정 (One-sample T-test)를 시행할 수 없고, 일표본 윌콕슨 부호 순위 검정 (One-Sample Wilcoxon Signed Rank Test)을 시행해야 한다.

 

일표본 윌콕슨 부호 순위 검정 (One-Sample Wilcoxon Signed Rank Test)의 귀무가설과 대립 가설

 

이번 일표본 윌콕슨 부호 순위 검정 (One-Sample Wilcoxon Signed Rank Test)의 귀무가설은 하나다.

귀무가설: $H_0=$ 모집단의 NP1_raw 중앙값은 55이다.

 

하지만 대립 가설은 3개가 된다.

1) $H_1=$ 모집단의 중앙값은 55보다 크다. (단측 검정)

2) $H_1=$ 모집단의 중앙값은 55보다 작다. (단측 검정)

3) $H_1=$ 모집단의 중앙값은 55보다 크거나 작다. (양측 검정)

 

대립 가설 별로 코드를 보고자 한다.

 

1)  $H_1=$ 모집단의 중앙값은 55보다 크다. (단측 검정)

wilcox.test(df$NP1_raw, mu=55, alternative="greater", correct=FALSE)

wilcox.test(df$NP1_raw, mu=55, alternative="greater", correct=FALSE) : Wilcoxon test를 시행하라. 검정할 변수는 df 데이터의 "NP1_raw"이고, 중앙값이 55와 같은지를 검정해라. 대립 가설은 모집단의 중앙값이 55보다 크다(greater)이고, 연속성 수정은 시행하지 마라. (연속성 수정은 다음 글을 참고하길 바란다.2022.08.30 - [통계 이론] - [이론] 연속성을 수정한 카이 제곱 검정 (Chi-squared test with Yates's correction for continuity))

 

결과

	Wilcoxon signed rank test

data:  df$NP1_raw
V = 156889, p-value = 1
alternative hypothesis: true location is greater than 55

Wilcoxon signed rank test

data:  df$NP1_raw

V = 156889, p-value = 1

alternative hypothesis: true location is greater than 55

 

P-value는 1.0000이다. 즉, 대립 가설 (모집단의 중앙값은 55보다 크다.)을 선택하는 것이 말이 안 된다는 것이다. 중앙값이 55보다 크다고 보는 것보다는 55이라고 보는 것이 더 합리적이라는 뜻이다. 

 

2)  $H_1=$ 모집단의 평균은 55보다 작다. (단측 검정)

wilcox.test(df$NP1_raw, mu=55, alternative="less", correct=FALSE)

wilcox.test(df$NP1_raw, mu=55, alternative="less", correct=FALSE) : Wilcoxon test를 시행하라. 검정할 변수는 df 데이터의 "NP1_raw"이고, 중앙값이 55와 같은지를 검정해라. 대립 가설은 '모집단의 중앙값이 55보다 작다(less)'이고, 연속성 수정은 시행하지 마라.

 

결과

	Wilcoxon signed rank test

data:  df$NP1_raw
V = 156889, p-value < 2.2e-16
alternative hypothesis: true location is less than 55

Wilcoxon signed rank test

data:  df$NP1_raw

V = 156889, p-value < 2.2e-16

alternative hypothesis: true location is less than 55

 

P-value2.2e-16보다 작다. 즉, 대립 가설 (모집단의 중앙값은 55보다 작다.)을 선택하는 것이 합당하다는 뜻이다. 중앙값이 55이라고 보는 것 보다는 55보다 작다고 보는 것이 더 합리적이라는 뜻이다. 

 

3)  $H_1=$ 모집단의 평균은 55보다 작거나 크다. (양측 검정)

wilcox.test(df$NP1_raw, mu=55, alternative="two.sided", correct=FALSE)

wilcox.test(df$NP1_raw, mu=55, correct=FALSE)

위 두 개의 코드는 같은 결과를 내준다. 왜냐하면 양측 검정 (alternative="two.sided")가 기본 값이어서 alternative 옵션을 지정해주지 않으면 alternative="two.sided"라고 인식하고 코드가 돌아가기 때문이다.

 

wilcox.test(df$NP1_raw, mu=55, alternative="two.sided", correct=FALSE) : Wilcoxon test를 시행하라. 검정할 변수는 df 데이터의 "NP1_raw"이고, 중앙값이 55와 같은지를 검정해라. 대립 가설은 '모집단의 중앙값은 55가 아니다(two.sided)'이고, 연속성 수정은 시행하지 마라.

 

결과

	Wilcoxon signed rank test

data:  df$NP1_raw
V = 156889, p-value < 2.2e-16
alternative hypothesis: true location is not equal to 55

Wilcoxon signed rank test

data:  df$NP1_raw

V = 156889, p-value < 2.2e-16

alternative hypothesis: true location is not equal to 55

 

P-value는 2.2e-16보다 작다. 즉, 대립 가설 (모집단의 중앙값은 55이 아니다.)을 선택하는 것이 합당하다는 뜻이다. 중앙값이 55이라고 보는 것보다는 55보다 작거나 크다고 보는 것이 더 합리적이라는 뜻이다. 

 

보통은 양측검정 결과를 사용하게 된다.

 

그러면 실제 중앙값은 얼마이길래 이런 결론이 나왔던 것일까?

median(df$NP1_raw)
[1] 49.96328

49.96328로 거의 50에 근사한다. 그렇다면 만약 50에 대해 검정하면 어떤 결과를 내보일까?

 

코드

wilcox.test(df$NP1_raw, mu=50, alternative="two.sided", correct=FALSE)

결과

	Wilcoxon signed rank test

data:  df$NP1_raw
V = 246178, p-value = 0.6558
alternative hypothesis: true location is not equal to 50

p-value가 0.6558>0.05로 귀무가설을 채택하는 것이 합리적이며, 중앙값이 50이 아니라고는 할 수 없는 상황이다.

 

 

중요한 가정: 대칭성

위 논의에서 빠진 정말 중요한 가정이 하나 있다. 신경심리검사 1 원점수 (NP1_raw)의 분포가 어떤 점수를 기점으로 좌우대칭이라는 것이다. 윌콕슨 부호 순위 검정은 좌우대칭일 때 그 검정력이 극대화되며, 좌우대칭이 아니면 정확도가 떨어진다. 신경심리검사 1 원점수 (NP1_raw)의 분포를 히스토그램을 통해 알아보자 (히스토그램 그리는 방법은 다음 링크를 참고하길 바란다.2022.08.11 - [기술 통계/R] - [R] 정규성 검정 (2) : 히스토그램 - hist(), dnorm())

#좀 더 자세히 보기 위해 breaks=40이라는 옵션을 넣었다.
hist(df$NP1_raw, breaks=40)

봉우리가 세 개인 Trimodal shape을 하고 있고, 중앙값인 50을 기준으로 대칭인 것처럼 보인다. 따라서 위 분석에는 문제가 없는 듯하다.

 

 

가정이 성립하지 않는다면?

우리가 궁금하고, 문제가 되는 상황은 가정이 성립하지 않는 상황이다. 실제로 이럴 때가 많기 때문이다. 신경심리검사 2 원점수인 NP2_raw데이터의 분포를 확인해보자.

hist(df$NP2_raw, breaks=80)

이는 심하게 오른쪽으로 치우쳐진 right skewed data다. 이때 중앙값에 대해 Wilcoxon signed rank test를 시행하면 어떻게 될까?

wilcox.test(df$NP2_raw, mu=median(df$NP2_raw), alternative="two.sided", correct=FALSE)
	Wilcoxon signed rank test

data:  df$NP2_raw
V = 295522, p-value = 7.213e-07
alternative hypothesis: true location is not equal to 1.004532

참 중앙값으로 "NP2_raw의 중앙값이 참 중앙값과 같니?"라고 물어보았는데 "아니요"라고 대답하는 참사가 벌어진다. 이는 분포에 대칭성이 없기 때문이다. 이럴 때 선택할 수 있는 옵션이 많지 않긴 한데, 그중 하나는 변환 (transformation)이다. 이런 right skewed data는 로그 변환 (log-transformation)을 하면 대칭성이 갖추어지는 경우가 많다. 로그 변환을 하고 히스토그램을 그려 분포를 확인해보자. (로그 변환에 관한 글은 다음 링크에서 확인할 수 있다.2022.11.25 - [통계 프로그램 사용 방법/R] - [R] 변수 계산 (산술 연산))

df$log_NP2_raw<-log(df$NP2_raw)
hist(df$log_NP2_raw, breaks=80)

충분히 대칭성을 띄고 있다고 할 수 있다.  "혹시 정규성을 띄지 않을까?"라는 생각으로 정규성을 확인해보면 다음과 같음을 알 수 있다. 정규성 검정 방법은 아래 링크들을 확인하길 바란다.

2022.08.11 - [기술 통계/R] - [R] 정규성 검정 (1) : Q-Q plot - qqnorm(), qqline()

2022.08.11 - [기술 통계/R] - [R] 정규성 검정 (2) : 히스토그램 - hist(), dnorm()

2022.08.11 - [기술 통계/R] - [R] 정규성 검정 (3) : 정량적 검정 (Shapiro-Wilk, Kolmogorov-Smirnov) - shapiro.test(), ks.test()

2022.08.12 - [기술 통계/R] - [R] 정규성 검정 (4) : 정량적 검정 (Lilliefors test) - lillie.test()

2022.08.16 - [기술 통계/R] - [R] 고급 Q-Q Plot - Van der Waerden, Rankit, Tukey, Blom

 

코드

# 1) Q-Q plot 그리기
qqnorm(df$log_NP2_raw)
qqline(df$log_NP2_raw)
# 2) 히스토그램 그리기
hist(df$log_NP2_raw, prob=TRUE)
log_NP2_rawrange<-seq(min(df$log_NP2_raw),max(df$log_NP2_raw),length=max(max(df$log_NP2_raw)-min(df$log_NP2_raw),100))
ND<-dnorm(log_NP2_rawrange,mean=mean(df$log_NP2_raw),sd=sd(df$log_NP2_raw))
lines(log_NP2_rawrange, ND, lwd=2)
# 3) Shapiro-Wilk test 시행하기
shapiro.test(df$log_NP2_raw)

결과

	Shapiro-Wilk normality test

data:  df$log_NP2_raw
W = 0.99284, p-value = 9.487e-05

히스토그램에서는 정규분포를 따르는 것 같아 보이지만, QQ plot에서 직선 위에 있지 않은 점들이 많고, Shapiro-wilk test에서 p-value가 매우 작게 나오므로 정규성을 띠지 않는다고 판단하는 것이 맞다. 따라서 일표본 t검정(One-sample t-test)이 아닌 일표본 윌콕슨 부호 순위 검정 (One sample Wilcoxon signed rank test)을 시행해야 한다. 위와 같이 log_NP2_raw의 중앙값에 대해 검정해보자.

 

코드

wilcox.test(df$log_NP2_raw, mu=median(df$log_NP2_raw), alternative="two.sided", correct=FALSE)

 

결과

	Wilcoxon signed rank test

data:  df$log_NP2_raw
V = 247173, p-value = 0.7363
alternative hypothesis: true location is not equal to 0.004521896

p-value가 0.7363으로 귀무가설을 기각할 수 없게 되어 중앙값은 0.004521896이 아니라고 할 수 없다고 결론을 내리게 된다.

 

 

그 외 다룰 내용 1: 정규성을 안 따르면 꼭 Wilcoxon test를 시행해야 하나?

 T-test의 전제조건은 표본 평균이 정규분포를 따르는 것이다. 그런데 중심 극한 정리에 따라 표본의 수 (우리 데이터에서는 n=1,000)가 커질수록 표본 평균의 분포는 정규성을 띠게 된다. 따라서 n수가 적당히 크기만 하면 표본이 정규분포를 따르는지와 관계없이 표본평균은 정규성을 따른다고 할 수 있다. 이런 이유로 정규성과 관계 없이 n수가 크기만 하면 t-test의 이용이 정당화되기도 한다. 합리적인 말이다.

 그런데 왜 정규분포를 따르는지 왜 다들 확인할까? 만약 모분포가 정규분포였다면 표본도 정규분포를 따를 확률이 높아진다. 물론 꼭 그런 건 아니지만 말이다. 따라서 표본이 정규분포를 따른다면, 모분포를 정규분포를 따른다고 할 수 있고, 그렇다면 n수가 작아도 표본 평균의 분포도 정규성을 따르니 t-test를 사용해도 되는 안전한 환경이 구축된다. 이런 이유로 정규성을 확인하게 된다.

 

그 외 다룰 내용 2: Wilcoxon signed rank test는 짝지어진 자료의 검정방법 아닌가?

 맞다. 하지만 one sample test에도 사용될 수 있다. 모수적 방법에서 대응 표본 T 검정 (paired T test)가 사실 일표본 T 검정 (One sample T test)와 같다는 것과 일맥상통하는 이야기다. (대응 표본 T검정은 다음 링크에서 확인할 수 있다.2022.11.25 - [반복 측정 자료 분석/R] - [R] 대응 표본 T검정 (Paired samples T-test) - t.test())

 

코드 정리

##워킹 디렉토리 지정
setwd("C:/Users/user/Documents/Tistory_blog")

##데이터 불러오기
install.packages("readr")
library("readr")
df<-read_csv("Data.csv")

###########NP1_raw###########

##정규성 검정
# 1) Q-Q plot 그리기
qqnorm(df$NP1_raw)
qqline(df$NP1_raw)
# 2) 히스토그램 그리기
hist(df$NP1_raw, prob=TRUE)
NP1_rawrange<-seq(min(df$NP1_raw),max(df$NP1_raw),length=max(max(df$NP1_raw)-min(df$NP1_raw),100))
ND<-dnorm(NP1_rawrange,mean=mean(df$NP1_raw),sd=sd(df$NP1_raw))
lines(NP1_rawrange, ND, lwd=2)
# 3) Shapiro-Wilk test 시행하기
shapiro.test(df$NP1_raw)

##Wilcoxon signed rank test
wilcox.test(df$NP1_raw, mu=55, alternative="greater", correct=FALSE)
wilcox.test(df$NP1_raw, mu=55, alternative="less", correct=FALSE)
wilcox.test(df$NP1_raw, mu=55, alternative="two.sided", correct=FALSE)
wilcox.test(df$NP1_raw, mu=55, correct=FALSE)

##중앙값 확인하기
median(df$NP1_raw)

##중앙값 근사치(50)으로 검정하기
wilcox.test(df$NP1_raw, mu=50, alternative="two.sided", correct=FALSE)

##NP1_raw 대칭성 확인
hist(df$NP1_raw)
hist(df$NP1_raw, breaks=40)

###########NP2_raw###########

##NP2_raw 대칭성 확인
hist(df$NP2_raw, breaks=80)

##NP2_raw 비대칭이지만 wilcoxon test 시행하기
wilcox.test(df$NP2_raw, mu=median(df$NP2_raw), alternative="two.sided", correct=FALSE)

##로그변환 후 히스토그램으로 대칭성 확인하기
df$log_NP2_raw<-log(df$NP2_raw)
hist(df$log_NP2_raw, breaks=80)

##로그변환된 NP2_raw의 정규성 확인하기
# 1) Q-Q plot 그리기
qqnorm(df$log_NP2_raw)
qqline(df$log_NP2_raw)
# 2) 히스토그램 그리기
hist(df$log_NP2_raw, prob=TRUE)
log_NP2_rawrange<-seq(min(df$log_NP2_raw),max(df$log_NP2_raw),length=max(max(df$log_NP2_raw)-min(df$log_NP2_raw),100))
ND<-dnorm(log_NP2_rawrange,mean=mean(df$log_NP2_raw),sd=sd(df$log_NP2_raw))
lines(log_NP2_rawrange, ND, lwd=2)
# 3) Shapiro-Wilk test 시행하기
shapiro.test(df$log_NP2_raw)

##로그변환된 변수로 Wilcoxon test 시행하기.
wilcox.test(df$log_NP2_raw, mu=median(df$log_NP2_raw), alternative="two.sided", correct=FALSE)

 

[R] 일표본 윌콕슨 부호 순위 검정 (비모수 일표본 중앙값 검정: One-Sample Wilcoxon Signed Rank Test) 정복 완료!

 

작성일: 2022.11.29.

최종 수정일: 2022.11.29.

이용 프로그램: R 4.2.2

RStudio v2022.07.2

RStudio 2022.07.2+576 "Spotted Wakerobin" Release

운영체제: Windows 10, Mac OS 12.6.1

반응형

+ Recent posts