R에서 ancova를 수행하는 방법

이 튜토리얼에서는 R에서 ANCOVA를 수행하는 방법에 대한 예를 제공합니다.

예: R의 ANCOVA

다음 변수를 사용하여 학습 기술이 시험 결과에 영향을 미치는지 여부를 테스트하기 위해 ANCOVA를 수행합니다.

• 기술적 연구 : 분석하고자 하는 독립변수
• 현재 학생 성적 : 고려하려는 공변량
• 검토점수 : 분석하고자 하는 반응변수

다음 데이터세트에는 무작위로 30명씩 세 그룹으로 나누어진 90명의 학생에 대한 정보가 포함되어 있습니다.

데이터 세트에는 각 학생이 사용한 학습 기법 (A, B 또는 C) , 학습 기법을 사용하기 시작한 당시 수업의 현재 성적, 준비를 위해 한 달 동안 학습 기법을 사용한 후 시험에서 받은 성적이 표시됩니다. 시험을 위해. 시험:

 #make this example reproducible
set.seed(10)

#create dataset
data <- data.frame(technique = rep(c("A", "B", "C"), each = 30),
                   current_grade = runif(90, 65, 95),
                   exam = c(runif(30, 80, 95), runif(30, 70, 95), runif(30, 70, 90)))

#view first six lines of dataset
head(data)

# technical current_grade exam
#1 A 80.22435 87.32759
#2 A 74.20306 90.67114
#3 A 77.80723 88.87902
#4 A 85.79306 87.75735
#5 A 67.55408 85.72442
#6 A 71.76310 92.52167

1단계: 데이터 탐색

ANCOVA 모델을 피팅하기 전에 먼저 데이터를 탐색하여 데이터를 더 잘 이해하고 결과를 왜곡할 수 있는 극단적인 이상값이 없는지 확인해야 합니다.

먼저 데이터 세트의 각 변수에 대한 요약을 표시할 수 있습니다.

 summary(data)

# technical current_grade exam      
#A:30 Min. :65.43 Min. :71.17  
# B:30 1st Qu.:71.79 1st Qu.:77.27  
# C:30 Median:77.84 Median:84.69  
# Mean:78.15 Mean:83.38  
# 3rd Qu.:83.65 3rd Qu.:89.22  
# Max. :93.84 Max. :94.76  

각 학습 기법 값( A, B , C)이 데이터에 30번 나타나는 것을 볼 수 있습니다.

또한 연구 시작 시 학생들의 현재 점수가 어떻게 분포되었는지 확인할 수 있습니다. 해당 과목의 최저점수는 65.43점, 최고점수는 93.84점, 평균점수는 78.15점이었다.

마찬가지로, 시험에서 얻은 최소 점수는 71.17점, 최대 점수는 94.76점, 평균은 83.38점임을 알 수 있습니다.

그런 다음 dplyr 패키지를 사용하여 각 학습 기술에 대한 현재 성적과 시험 결과의 평균 및 표준 편차를 쉽게 찾을 수 있습니다.

 #load dplyr
library(dplyr)

data %>%
  group_by (technical) %>%
  summarize (mean_grade = mean(current_grade),
            sd_grade = sd(current_grade),
            mean_exam = mean(exam),
            sd_exam = sd(exam))

# A tibble: 3 x 5
# technique mean_grade sd_grade mean_exam sd_exam                      
#1 A 79.0 7.00 88.5 3.88
#2 B 78.5 8.33 81.8 7.62
#3 C 76.9 8.24 79.9 5.71

각 학습방법을 활용한 학생들의 현재 성적의 평균과 표준편차는 거의 유사함을 알 수 있다.

또한 기술 B 와 C 에 비해 학습 기술 A를 사용한 학생들의 평균 시험 점수가 훨씬 더 높다는 것을 알 수 있습니다.

또한 상자 그림을 사용하여 연구 기술을 기반으로 시험 결과 분포를 시각화할 수 있습니다.

 boxplot(exam ~ technique,
data = data,
main = "Exam Score by Studying Technique",
xlab = "Studying Technique",
ylab = "Exam Score",
col = "steelblue",
border = "black"
)

마찬가지로 상자 그림을 사용하여 학습 기술을 기반으로 현재 성적 분포를 시각화할 수도 있습니다.

 boxplot(current_grade ~ technical,
data = data,
main = "Current Grade by Studying Technique",
xlab = "Studying Technique",
ylab = "Current Grade",
col = "steelblue",
border = "black"
)

2단계: 모델 가정 확인

몇 가지 기본 데이터 탐색을 수행하고 데이터에 익숙해지면 ANCOVA에 대한 다음 가정이 충족되는지 확인해야 합니다.

• 공변량과 처리는 독립적입니다 . 공변량( 현재 등급) 과 처리 (연구 기법) 가 서로 독립적인지 검증해야 합니다. 모델에 공변량 항을 추가하는 것은 공변량과 치료는 반응 변수( 시험 )에 독립적으로 작용합니다.
• 분산의 동질성 – 그룹 간의 분산이 동일한지 확인해야 합니다.

공변량과 치료가 독립적인지 확인하기 위해 현재 등급을 반응 변수로 사용하고 연구 기술을 예측 변수로 사용하여 분산 분석을 수행할 수 있습니다.

 #fit anova model
anova_model <- aov(current_grade ~ technique, data = data)
#view summary of anova model
summary(anova_model)

# Df Sum Sq Mean Sq F value Pr(>F)
#technical 2 74 37.21 0.599 0.552
#Residuals 87 5406 62.14    

p-값이 0.05보다 크므로 공변량( 현재 등급) 과 치료( 연구 기법 )가 독립적인 것으로 보입니다.

그런 다음 그룹 간 분산의 동질성이 있는지 확인하기 위해 Levene 테스트를 수행할 수 있습니다.

 #load car library to conduct Levene's Test
libary(car)

#conduct Levene's Test
leveneTest(technical exam, data = data)

#Levene's Test for Homogeneity of Variance (center = median)
# Df F value Pr(>F)    
#group 2 9.4324 0.0001961 ***
#87   

검정의 p-값은 0.0001961로, 이는 그룹 간의 분산이 동일하지 않음을 나타냅니다. 이 문제를 해결하기 위해 데이터 변환을 시도할 수 있지만 현재로서는 분산 차이에 대해 너무 걱정하지 않을 것입니다.

3단계: ANCOVA 모델 조정

다음으로, 시험 점수를 반응 변수로, 학습 기술을 예측 변수(또는 “치료”) 변수로, 현재 등급을 공변량으로 사용하여 ANCOVA 모델을 적합화할 것입니다.

이를 수행하기 위해 유형 III 제곱합을 모델에 사용하도록 지정할 수 있기 때문에 패키지에서 Anova() 함수를 사용할 것입니다. 왜냐하면 유형 I 제곱합은 순서에 따라 달라지기 때문입니다. 예측변수가 모델에 입력됩니다.

 #load car library
library(car)

#fit ANCOVA model
ancova_model <- aov(exam ~ technique + current_grade, data = data)

#view summary of model
Anova(ancova_model, type="III") 

#Answer: exam
# Sum Sq Df F value Pr(>F)    
#(Intercept) 7161.2 1 201.4621 < 2.2e-16 ***
#technical 1242.9 2 17.4830 4.255e-07 ***
#current_grade 12.3 1 0.3467 0.5576    
#Residuals 3057.0 86         

기술 에 대한 p-값이 매우 낮은 것을 볼 수 있는데, 이는 현재 성적을 통제한 후에도 학습 기술이 시험 점수에 통계적으로 유의미한 영향을 미친다는 것을 나타냅니다.

4단계: 사후 테스트

ANCOVA 결과를 통해 학습 기법이 시험 점수에 통계적으로 유의미한 영향을 미치는 것으로 나타났지만, 어떤 학습 기법이 서로 다른지 확인하려면 사후 테스트를 수행해야 합니다.

이를 위해 R의 multcomp 패키지에 있는 glht() 함수를 사용하여 다중 비교를 위한 Tukey 테스트를 수행할 수 있습니다.

 #load the multcomp library
library(multicomp)

#fit the ANCOVA model
ancova_model <- aov(exam ~ technique + current_grade, data = data)

#define the post hoc comparisons to make
postHocs <- glht(ancova_model, linfct = mcp(technique = "Tukey"))

#view a summary of the post hoc comparisons
summary(postHocs)

#Multiple Comparisons of Means: Tukey Contrasts
#
#Fit: aov(formula = exam ~ technique + current_grade, data = data)
#
#Linear Assumptions:
#Estimate Std. Error t value Pr(>|t|)    
#B - A == 0 -6.711 1.540 -4.358 0.000109 ***
#C - A == 0 -8.736 1.549 -5.640 < 1e-04 ***
#C - B == 0 -2.025 1.545 -1.311 0.393089    

#view the confidence intervals associated with the multiple comparisons
confint(postHocs)

# Simultaneous Confidence Intervals
#
#Multiple Comparisons of Means: Tukey Contrasts
#
#Fit: aov(formula = exam ~ technique + current_grade, data = data)
#
#Quantile = 2.3845
#95% family-wise confidence level
#
#Linear Assumptions:
# Estimate lwr upr     
#B - A == 0 -6.7112 -10.3832 -3.0392
#C - A == 0 -8.7364 -12.4302 -5.0426
#C - B == 0 -2.0252 -5.7091 1.6588

그 결과, A 기법을 공부하는 것과 B 기법을 공부하는 것(p-값: 0,000109), 그리고 A 기법과 C 기법을 공부하는 것 사이의 시험 결과에 통계적으로 유의미한 차이(α = 0.05)가 있음을 알 수 있습니다. (p-값: <1e-04).

또한 기술 B 와 C 간에 통계적으로 유의미한 차이(α = 0.05)가 없음 을 확인할 수 있습니다. 기술 간의 신뢰 구간도 이러한 결과를 뒷받침합니다.

따라서 우리는 수업에서 학생의 현재 성적을 통제한 후에도 학습 기법 A 를 사용하면 기법 B 와 C 에 비해 학생의 시험 성적이 통계적으로 유의하게 더 높다는 결론을 내릴 수 있습니다.