Как выполнить односторонний дисперсионный анализ в r

Однофакторный дисперсионный анализ используется для определения наличия или отсутствия статистически значимой разницы между средними значениями трех или более независимых групп.

Этот тип теста называется однофакторным дисперсионным анализом, поскольку мы анализируем влияние предикторной переменной на переменную отклика.

Примечание . Если бы нас вместо этого интересовало влияние двух переменных-предсказателей на переменную отклика, мы могли бы выполнить двусторонний дисперсионный анализ .

Как выполнить односторонний дисперсионный анализ в R

В следующем примере показано, как выполнить однофакторный дисперсионный анализ в R.

Фон

Предположим, мы хотим определить, оказывают ли три разные программы упражнений различное влияние на потерю веса. Переменной-предиктором, которую мы изучаем, является программа упражнений, а переменной-ответом является потеря веса, измеряемая в фунтах.

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

Мы набираем 90 человек для участия в эксперименте, в котором случайным образом назначаем 30 человек следовать программе А, программе Б или программе С в течение месяца.

Следующий код создает фрейм данных, с которым мы будем работать:

 #make this example reproducible
set.seed(0)

#create data frame
data <- data.frame(program = rep(c("A", "B", "C"), each = 30),
                   weight_loss = c(runif(30, 0, 3),
                                   runif(30, 0, 5),
                                   runif(30, 1, 7)))

#view first six rows of data frame
head(data)

# program weight_loss
#1 A 2.6900916
#2 A 0.7965260
#3 A 1.1163717
#4 A 1.7185601
#5 A 2.7246234
#6 A 0.6050458

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

Изучите данные

Прежде чем даже использовать однофакторную модель ANOVA, мы можем лучше понять данные, найдя среднее и стандартное отклонение потери веса для каждой из трех программ с помощью пакета dplyr :

 #load dplyr package
library (dplyr)

#find mean and standard deviation of weight loss for each treatment group
data %>%
  group_by (program) %>%
  summarize (mean = mean(weight_loss),
            sd = sd(weight_loss))

# A tibble: 3 x 3
# program mean sd
#      
#1 A 1.58 0.905
#2 B 2.56 1.24 
#3 C 4.13 1.57  

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

 #create boxplots
boxplot(weight_loss ~ program,
data = data,
main = "Weight Loss Distribution by Program",
xlab = "Program",
ylab = "Weight Loss",
col = "steelblue",
border = "black")

Из этих коробчатых диаграмм мы видим, что средняя потеря веса является самой высокой у участников программы C, а средняя потеря веса является самой низкой у участников программы A.

Мы также можем видеть, что стандартное отклонение («длина» коробчатой диаграммы) для потери веса немного выше в программе C по сравнению с двумя другими программами.

Далее мы применим однофакторную модель ANOVA к нашим данным, чтобы увидеть, действительно ли эти визуальные различия статистически значимы.

Однофакторная аппроксимация модели ANOVA

Общий синтаксис для подбора однофакторной модели ANOVA в R:

aov (переменная ответа ~ переменная_предиктора, данные = набор данных)

В нашем примере мы можем использовать следующий код для соответствия однофакторной модели ANOVA, используя Weight_loss в качестве переменной ответа и программу в качестве переменной-предиктора. Затем мы можем использовать функцию summary() для отображения результата нашей модели:

 #fit the one-way ANOVA model
model <- aov(weight_loss ~ program, data = data)

#view the model output
summary(model)

# Df Sum Sq Mean Sq F value Pr(>F)    
#program 2 98.93 49.46 30.83 7.55e-11 ***
#Residuals 87 139.57 1.60                     
#---
#Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Из результатов модели мы видим, что программа переменных-предикторов статистически значима на уровне значимости 0,05.

Другими словами, существует статистически значимая разница между средней потерей веса в результате трех программ.

Проверка предположений модели

Прежде чем идти дальше, нам необходимо убедиться, что предположения нашей модели выполняются, чтобы результаты нашей модели были надежными. В частности, односторонний дисперсионный анализ предполагает:

1. Независимость – наблюдения каждой группы должны быть независимы друг от друга. Поскольку мы использовали рандомизированный план (то есть мы распределяли участников по программам упражнений случайным образом), это предположение должно соблюдаться, чтобы нам не приходилось слишком беспокоиться об этом.

2. Нормальность – зависимая переменная должна иметь примерно нормальное распределение для каждого уровня предикторной переменной.

3. Равная дисперсия – дисперсии для каждой группы равны или примерно равны.

Один из способов проверить предположения о нормальности и равной дисперсии — использовать функциюplot() , которая создает четыре графика проверки модели. В частности, нас особенно интересуют следующие два сюжета:

• Остатки против. подобрано – этот график показывает взаимосвязь между остатками и подобранными значениями. Мы можем использовать этот график, чтобы примерно оценить, является ли дисперсия между группами примерно одинаковой или нет.
• График QQ – на этом графике отображаются стандартизованные остатки в сравнении с теоретическими квантилями. Мы можем использовать этот график, чтобы примерно оценить, выполняется ли предположение о нормальности.

Следующий код можно использовать для создания графиков проверки модели:

 plot(model)

График QQ, приведенный выше, позволяет нам проверить предположение о нормальности. В идеале стандартизированные остатки должны располагаться вдоль прямой диагональной линии графика. Однако на графике выше мы видим, что остатки немного отклоняются от линии в сторону начала и конца. Это указывает на то, что наше предположение о нормальности может быть нарушено.

Остатки против. Скорректированный график выше позволяет нам проверить наше предположение о равных дисперсиях. В идеале мы хотели бы, чтобы остатки распределялись одинаково для каждого уровня подобранных значений.

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

Чтобы формально проверить равенство дисперсий, мы могли бы выполнить тест Левена, используя пакет car :

 #load car package
library (car)

#conduct Levene's Test for equality of variances
leveneTest(weight_loss ~ program, data = data)

#Levene's Test for Homogeneity of Variance (center = median)
# Df F value Pr(>F)  
#group 2 4.1716 0.01862 *
#87                  
#---
#Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

P-значение теста составляет 0,01862 . Если мы используем уровень значимости 0,05, мы отвергнем нулевую гипотезу о том, что дисперсии одинаковы во всех трех программах. Однако если мы используем уровень значимости 0,01, мы не отвергнем нулевую гипотезу.

Хотя мы могли бы попытаться преобразовать данные, чтобы гарантировать выполнение наших предположений о нормальности и равенстве дисперсий, сейчас мы не будем слишком об этом беспокоиться.

Анализ различий в лечении

После того как мы убедились, что предположения модели выполняются (или обоснованно выполняются), мы можем провести апостериорный тест , чтобы точно определить, какие группы лечения отличаются друг от друга.

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

 #perform Tukey's Test for multiple comparisons
TukeyHSD(model, conf.level=.95) 

#Tukey multiple comparisons of means
# 95% family-wise confidence level
#
#Fit: aov(formula = weight_loss ~ program, data = data)
#
#$program
# diff lwr upr p adj
#BA 0.9777414 0.1979466 1.757536 0.0100545
#CA 2.5454024 1.7656076 3.325197 0.0000000
#CB 1.5676610 0.7878662 2.347456 0.0000199

Значение p указывает, существует ли статистически значимая разница между каждой программой. Результаты показывают, что существует статистически значимая разница между средней потерей веса каждой программы на уровне значимости 0,05.

Мы также можем визуализировать 95% доверительные интервалы, полученные в результате теста Тьюки, с помощью функцииplot(TukeyHSD()) в R:

 #create confidence interval for each comparison
plot(TukeyHSD(model, conf.level=.95), las = 2)

Результаты доверительных интервалов согласуются с результатами проверки гипотез.

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

Это соответствует тому, что все значения p для наших тестов гипотез составляют менее 0,05.

Отчет о результатах однофакторного дисперсионного анализа

Наконец, мы можем сообщить о результатах однофакторного дисперсионного анализа таким образом, чтобы суммировать результаты:

Для изучения влияния программы упражнений был проведен однофакторный дисперсионный анализ.   на похудение (измеряется в фунтах). Была статистически значимая разница между эффектами трех программ на потерю веса (F(2, 87) = 30,83, p = 7,55e-11). Апостериорно были проведены тесты Тьюки HSD.

Средняя потеря веса участников программы С значительно превышает среднюю потерю веса участников программы Б (р <0,0001).

Средняя потеря веса участников программы С значительно превышает среднюю потерю веса участников программы А (р <0,0001).

Кроме того, средняя потеря веса участников программы Б была значительно больше, чем средняя потеря веса участников программы А (р = 0,01).

Дополнительные ресурсы

В следующих руководствах представлена дополнительная информация об однофакторном дисперсионном анализе:

Введение в однофакторный дисперсионный анализ
Руководство по использованию апостериорного тестирования с помощью ANOVA
Полное руководство: Как сообщить о результатах ANOVA