Как провести ancova в r

В этом руководстве представлен пример выполнения ANCOVA в R.

Пример: ANCOVA в R

Мы проведем ANCOVA, чтобы проверить, влияет ли методика обучения на результаты экзамена, используя следующие переменные:

• Техническое исследование : независимая переменная, которую мы хотим проанализировать.
• Текущая оценка учащегося : ковариата, которую мы хотим принять во внимание.
• Оценка по обзору : переменные ответа, которые мы хотим проанализировать.

Следующий набор данных содержит информацию о 90 студентах, случайно разделенных на три группы по 30 человек.

Набор данных показывает технику обучения, которую использовал каждый учащийся (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.

Мы также можем визуализировать распределение результатов экзамена в зависимости от методики обучения, используя коробчатые диаграммы :

 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, поэтому ковариата ( текущая степень) и лечение ( метод исследования ) кажутся независимыми.

Затем, чтобы проверить однородность дисперсии между группами, мы можем выполнить тест Левена:

 #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, используя оценку экзамена в качестве переменной ответа, технику обучения в качестве предикторной переменной (или «лечения») и текущую оценку в качестве ковариаты.

Мы будем использовать функцию Anova() в пакете, потому что для этого просто чтобы иметь возможность указать, что мы хотим использовать сумму квадратов типа III для модели, поскольку сумма квадратов типа 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 показали нам, что методика обучения оказала статистически значимое влияние на результаты экзаменов, нам необходимо провести апостериорное тестирование , чтобы определить, какие методы обучения отличаются друг от друга.

Для этого мы можем использовать функцию glht() в пакете multcomp в R, чтобы выполнить тест Тьюки для множественных сравнений:

 #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

Из результата мы видим, что существует статистически значимая разница (при α = 0,05) в результатах экзамена между изучением техники A и изучением техники B (значение p: 0,000109), а также между техникой A и техникой C. (значение p: <1e-04).

Мы также можем видеть, что нет статистически значимой разницы (при α = 0,05) между методами B и C. Доверительные интервалы между методами также подтверждают эти результаты.

Таким образом, можно сделать вывод, что использование методики обучения А приводит к статистически значимо более высокой оценке учащихся на экзамене по сравнению с методиками В и С , даже после контроля текущей оценки учащегося в классе.