ANOVA and Multiple Comparision
백승훈
2023-10-04
Data 불러오기
getwd()## [1] "C:/data"library(readxl)
data_anova <- read_excel("data_anova.xlsx")데이터 변환
ANOVA 수행을 위한 일반표 형태로 변환한다.
data_anova <- as.data.frame(data_anova)
data_anova_a <- c()
for (i in 1:6){
data_anova_a <- append(data_anova_a,data_anova[,i])
i <- i+1
}
group <- rep(1:6, each=300)
group <- as.factor(group)
# combining into data.frame
group_df <- data.frame(data_anova_a, group)
sapply(group_df, class)## data_anova_a group
## "numeric" "factor"데이터 시각화 및 주요 기술통계
데이터의 분포를 박스 플롯으로 시각화하고 주요 기술통계를 본다.
# boxplot
boxplot(data_anova_a ~ group,
main = "Boxplot of evaluation for each question",
xlab = "Factor Levels : questions",
ylab = "evauation")
# descriptive statistics by group
tapply(data_anova_a, group, summary)## $`1`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 58.79 66.17 69.67 69.67 73.16 83.35
##
## $`2`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 58.22 70.78 75.06 74.74 78.36 91.31
##
## $`3`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 54.33 65.00 68.85 68.46 71.95 80.73
##
## $`4`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 49.84 61.58 65.25 64.92 68.48 77.39
##
## $`5`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 61.51 70.05 73.35 73.24 76.47 88.90
##
## $`6`
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 60.08 67.94 71.28 71.24 74.58 82.90ANOVA 수행
제안별 수용도에서 하나 이상의 차이가 있는지 확인하기 위해 ANOVA를 수행한다.
# one-wayANOVA
result <- aov(data_anova_a ~ group, data = group_df)
summary(result)## Df Sum Sq Mean Sq F value Pr(>F)
## group 5 18602 3720 149.4 <2e-16 ***
## Residuals 1794 44664 25
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# anova(result)
model.tables(result) # 전체 평균 대비 각 그룹 수용도의 차이## Tables of effects
##
## group
## group
## 1 2 3 4 5 6
## -0.710 4.362 -1.921 -5.462 2.865 0.866일원 분산분석 결과 최소한 하나 이상의 오퍼에 대한 수용도 차이가 있다고 보인다. (p<0.05이므로)
오퍼별 평균 대비 상대적 값을 본 결과, 2,5번 오퍼에 대한 수용도가 높은 반면 3,4번에 대한 수용도는 낮은 것으로 보인다.
정규성 검정
잔차의 정규성을 확인하기 위해서 aov() 함수의 결과에서 residuals를 뽑아서 정규성에 대해 가설검정을 해본다. p-값> 0.05이므로 잔차의 정규성을 기각하지 못했다.
# plot(result)
shapiro.test(result$residuals)##
## Shapiro-Wilk normality test
##
## data: result$residuals
## W = 0.99843, p-value = 0.09167등분산성 검정
등분산성을 가설검정할 때 요인의 수준이 2 개인 경우에는 두 모집단의 분산의 차이 검정을 F-분포를 사용하여 수행할 수 있다. 이 경우에는 var.test() 함수를 이용하여 가설검정을 수행하였다. 그러나 요인의 수준이 3 개 이상인 경우에는 bartlet.test()을 사용한다. p-값>0.05이므로 등분산 가정을 기각하지 못했다.
# Bartlett test to test the null hypothesis of equal group variances: 등분산 검정
bartlett.test(data_anova_a ~ group, data = group_df)##
## Bartlett test of homogeneity of variances
##
## data: data_anova_a by group
## Bartlett's K-squared = 8.6907, df = 5, p-value = 0.1221다중 비교
ANOVA 분석을 통해 평균이 같다는 귀무가설이 기각되면, 어는 수준 간의 차이에 의해 귀무가설이 기각되었는지 살펴보아야 한다. 이를 살펴보는 방법 중 하나가 가능한 모든 두 수준 간에 평균의 차이가 있었는지를 가설검정해 보는 것이다. 이는 두 모집단 간의 모평균 비교와 같아진다. 다만 수준이 많은 경우 단순히 두 수준 간의 평균의 차이를 가설검정하는 것은 여러 개의 가설을 동시적으로 검정하는 것이 되어서 기각역 또는 p-값의 조정이 필요하다. R의 기본 기능에 포함되어 있는 TukeyHSD() 함수는 이러한 수준별 쌍 비교를 수행해 준다.
##------------------------------------------------------------------
## multiple comparison
## - Tukey's HSD(honestly significant difference) test : TukeyHSD()
##------------------------------------------------------------------
TukeyHSD(result)## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = data_anova_a ~ group, data = group_df)
##
## $group
## diff lwr upr p adj
## 2-1 5.071821 3.9095902 6.23405206 0.0000000
## 3-1 -1.211064 -2.3732952 -0.04883339 0.0353955
## 4-1 -4.751949 -5.9141802 -3.58971834 0.0000000
## 5-1 3.574484 2.4122533 4.73671518 0.0000000
## 6-1 1.575414 0.4131828 2.73764465 0.0015912
## 3-2 -6.282885 -7.4451164 -5.12065452 0.0000000
## 4-2 -9.823770 -10.9860013 -8.66153948 0.0000000
## 5-2 -1.497337 -2.6595678 -0.33510595 0.0033332
## 6-2 -3.496407 -4.6586383 -2.33417648 0.0000000
## 4-3 -3.540885 -4.7031159 -2.37865403 0.0000000
## 5-3 4.785549 3.6233176 5.94777949 0.0000000
## 6-3 2.786478 1.6242471 3.94870897 0.0000000
## 5-4 8.326434 7.1642026 9.48866445 0.0000000
## 6-4 6.327363 5.1651321 7.48959392 0.0000000
## 6-5 -1.999071 -3.1613015 -0.83683960 0.0000149위의 모든 비교에서 유의미한 차이가 있다고 보이므로 앞서 봤던 평균 순서 그대로 차이가 있다고 판단된다.
\[2>5>6>1>3>4\]
##------------------------------------------------------------------
## multiple comparison
## - Scheffe test : scheffe.test()
##------------------------------------------------------------------
# install.packages('agricolae')
library(agricolae)
scheffe.test(result,"group",alpha = 0.05,group=TRUE,console=TRUE)##
## Study: result ~ "group"
##
## Scheffe Test for data_anova_a
##
## Mean Square Error : 24.89632
##
## group, means
##
## data_anova_a std r se Min Max Q25 Q50
## 1 69.66722 4.773053 300 0.2880759 58.79182 83.35249 66.16948 69.67084
## 2 74.73904 5.358928 300 0.2880759 58.21940 91.30598 70.77921 75.06022
## 3 68.45616 4.941499 300 0.2880759 54.33159 80.72676 65.00495 68.85282
## 4 64.91527 5.265326 300 0.2880759 49.83550 77.38503 61.58378 65.24671
## 5 73.24171 4.692675 300 0.2880759 61.51371 88.89632 70.04562 73.35434
## 6 71.24264 4.869756 300 0.2880759 60.07674 82.90091 67.93614 71.27999
## Q75
## 1 73.16403
## 2 78.35956
## 3 71.95412
## 4 68.47830
## 5 76.46883
## 6 74.57515
##
## Alpha: 0.05 ; DF Error: 1794
## Critical Value of F: 2.219085
##
## Minimum Significant Difference: 1.357044
##
## Means with the same letter are not significantly different.
##
## data_anova_a groups
## 2 74.73904 a
## 5 73.24171 b
## 6 71.24264 c
## 1 69.66722 d
## 3 68.45616 d
## 4 64.91527 e또 다른 다중 비교 방법인 scheffe test 결과도 거의 비슷하게 나왔지만 1과 3은 차이가 없는 것으로 판단했다.
\[2>5>6>1=3>4\]