Як виконати двосторонній дисперсійний аналіз у r

Двосторонній дисперсійний аналіз («дисперсійний аналіз») використовується, щоб визначити, чи існує статистично значуща різниця між середніми значеннями трьох або більше незалежних груп, які були розділені на два фактори.

У цьому посібнику пояснюється, як виконати двосторонній дисперсійний аналіз у R.

Приклад: двосторонній дисперсійний аналіз R

Скажімо, ми хочемо визначити, чи впливають інтенсивність вправ і стать на втрату ваги. У цьому випадку два фактори, які ми розглядаємо, — це фізичні вправи та стать , а змінна відповіді — це втрата ваги, виміряна у фунтах.

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

Ми набираємо 30 чоловіків і 30 жінок для участі в експерименті, у якому ми випадковим чином призначаємо 10 із кожного для виконання програми без фізичних вправ, з легкими вправами або з інтенсивними вправами протягом місяця.

Наступний код створює фрейм даних, з яким ми будемо працювати:

 #make this example reproducible
set.seed(10)

#create data frame
data <- data.frame(gender = rep(c("Male", "Female"), each = 30),
                   exercise = rep(c("None", "Light", "Intense"), each = 10, times = 2),
                   weight_loss = c(runif(10, -3, 3), runif(10, 0, 5), runif(10, 5, 9),
                                   runif(10, -4, 2), runif(10, 0, 3), runif(10, 3, 8)))

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

# gender exercise weight_loss
#1 Male None 0.04486922
#2 Male None -1.15938896
#3 Male None -0.43855400
#4 Male None 1.15861249
#5 Male None -2.48918419
#6 Male None -1.64738030

#see how many participants are in each group
table(data$gender, data$exercise)

# Intense Light None
# Female 10 10 10
# Male 10 10 10

Вивчіть дані

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

 #load dplyr package
library(dplyr)

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

# A tibble: 6 x 4
# Groups: gender [2]
# gender exercise means sd
#          
#1 Female Intense 5.31 1.02 
#2 Female Light 0.920 0.835
#3 Female None -0.501 1.77 
#4 Male Intense 7.37 0.928
#5 Male Light 2.13 1.22 
#6 Male None -0.698 1.12 

Ми також можемо створити коробковий графік для кожної з шести груп лікування, щоб візуалізувати розподіл втрати ваги для кожної групи:

 #set margins so that axis labels on boxplot don't get cut off
by(mar=c(8, 4.1, 4.1, 2.1))

#create boxplots
boxplot(weight_loss ~ gender:exercise,
data = data,
main = "Weight Loss Distribution by Group",
xlab = "Group",
ylab = "Weight Loss",
col = "steelblue",
border = "black", 
las = 2 #make x-axis labels perpendicular
)

Ми відразу бачимо, що дві групи, які брали участь у інтенсивних вправах, здається, мають вищі значення втрати ваги. Ми також бачимо, що чоловіки, як правило, мають вищі показники втрати ваги, ніж жінки, як у групах інтенсивних , так і в групах легких вправ.

Далі ми адаптуємо модель двостороннього дисперсійного аналізу до наших даних, щоб побачити, чи справді ці візуальні відмінності є статистично значущими.

Підгонка двосторонньої моделі ANOVA

Загальний синтаксис для підгонки двосторонньої моделі ANOVA в R такий:

aov(змінна відповіді ~predictor_variable1 *predictor_variable2, дані = набір даних)

Зауважте, що * між двома змінними предиктора вказує на те, що ми також хочемо перевірити ефект взаємодії між двома змінними предиктора.

У нашому прикладі ми можемо використати наведений нижче код, щоб підібрати модель двостороннього дисперсійного аналізу, використовуючи weight_loss як змінну відповіді та стать і фізичні вправи як дві змінні предиктора.

Потім ми можемо використати функцію summary() для відображення результату нашої моделі:

 #fit the two-way ANOVA model
model <- aov(weight_loss ~ gender * exercise, data = data)

#view the model output
summary(model)

# Df Sum Sq Mean Sq F value Pr(>F)    
#gender 1 15.8 15.80 11.197 0.0015 ** 
#exercise 2 505.6 252.78 179.087 <2e-16 ***
#gender:exercise 2 13.0 6.51 4.615 0.0141 *  
#Residuals 54 76.2 1.41                   
#---
#Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

З результатів моделі ми бачимо, що стать , фізичні вправи та взаємодія між двома змінними є статистично значущими на рівні значущості 0,05.

Перевірка припущень моделі

Перш ніж йти далі, нам потрібно перевірити, чи виконуються припущення нашої моделі, щоб результати нашої моделі були надійними. Зокрема, двосторонній ANOVA передбачає:

1. Незалежність – спостереження кожної групи мають бути незалежними одне від одного. Оскільки ми використовували рандомізований дизайн , це припущення має виконуватися, тож нам не потрібно про це надто турбуватися.

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

Одним із способів перевірити це припущення є створення гістограми залишків моделі. Якщо залишки розподілені приблизно нормально, це припущення має виконуватися.

 #define model residuals
reside <- model$residuals

#create histogram of residuals
hist(resid, main = "Histogram of Residuals", xlab = "Residuals", col = "steelblue")

Залишки розподілені приблизно нормально, тому ми можемо припустити, що припущення нормальності виконується.

3. Рівна дисперсія – дисперсії для кожної групи рівні або приблизно рівні.

Одним із способів перевірити це припущення є виконання тесту Левена на рівність дисперсій за допомогою комплекту автомобіля :

 #load car package
library(car)

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

#Levene's Test for Homogeneity of Variance (center = median)
# Df F value Pr(>F)
#group 5 1.8547 0.1177
#54  

Оскільки p-значення тесту перевищує наш рівень значущості 0,05, ми можемо припустити, що наше припущення про рівність дисперсій між групами виконується.

Проаналізуйте відмінності лікування

Після того, як ми переконалися, що припущення моделі виконуються, ми можемо провести пост-тест, щоб точно визначити, які групи лікування відрізняються одна від одної.

Для нашого post hoc тесту ми використаємо функцію TukeyHSD() для виконання тесту Tukey для кількох порівнянь:

 #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 ~ gender * exercise, data = data)
#
#$gender
# diff lwr upr p adj
#Male-Female 1.026456 0.4114451 1.641467 0.0014967
#
#$exercise
# diff lwr upr p adj
#Light-Intense -4.813064 -5.718493 -3.907635 0.0e+00
#None-Intense -6.938966 -7.844395 -6.033537 0.0e+00
#None-Light -2.125902 -3.031331 -1.220473 1.8e-06
#
#$`gender:exercise`
# diff lwr upr p adj
#Male:Intense-Female:Intense 2.0628297 0.4930588 3.63260067 0.0036746
#Female:Light-Female:Intense -4.3883563 -5.9581272 -2.81858535 0.0000000
#Male:Light-Female:Intense -3.1749419 -4.7447128 -1.60517092 0.0000027
#Female:None-Female:Intense -5.8091131 -7.3788841 -4.23934219 0.0000000
#Male:None-Female:Intense -6.0059891 -7.5757600 -4.43621813 0.0000000
#Female:Light-Male:Intense -6.4511860 -8.0209570 -4.88141508 0.0000000
#Male:Light-Male:Intense -5.2377716 -6.8075425 -3.66800066 0.0000000
#Female:None-Male:Intense -7.8719429 -9.4417138 -6.30217192 0.0000000
#Male:None-Male:Intense -8.0688188 -9.6385897 -6.49904786 0.0000000
#Male:Light-Female:Light 1.2134144 -0.3563565 2.78318536 0.2185439
#Female:None-Female:Light -1.4207568 -2.9905278 0.14901410 0.0974193
#Male:None-Female:Light -1.6176328 -3.1874037 -0.04786184 0.0398106
#Female:None-Male:Light -2.6341713 -4.2039422 -1.06440032 0.0001050
#Male:None-Male:Light -2.8310472 -4.4008181 -1.26127627 0.0000284
#Male:None-Female:None -0.1968759 -1.7666469 1.37289500 0.9990364

P-значення вказує на те, чи є статистично значуща різниця між кожною групою.

Наприклад, в останньому рядку вище ми бачимо, що група чоловіків, які не займалися фізичними вправами, не відчула статистично значущої різниці у втраті ваги порівняно з групою жінок, які не займалися фізичними вправами (p-значення: 0,990364).

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

 #set axis margins so labels don't get cut off
by(mar=c(4.1, 13, 4.1, 2.1))

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

Звіт про результати двостороннього дисперсійного аналізу

Нарешті, ми можемо повідомити про результати двостороннього дисперсійного аналізу таким чином, щоб узагальнити результати:

Двосторонній дисперсійний аналіз проводився для вивчення впливу статі ( чоловіча, жіноча) і програми вправ (немає, легка, інтенсивна) на втрату ваги (вимірюється у фунтах). Була статистично значуща взаємодія між впливом статі та фізичних вправ на втрату ваги (F(2, 54) = 4,615, p = 0,0141). Були проведені додаткові тести Tukey HSD.

Для чоловіків інтенсивна програма тренувань призвела до значно більшої втрати ваги, ніж легка програма (p < 0,0001) або програма без тренувань (p < 0,0001). Крім того, у чоловіків легка дієта призвела до значно більшої втрати ваги, ніж відсутність режиму фізичних вправ (p <0,0001).

Для жінок інтенсивна програма тренувань призвела до значно більшої втрати ваги, ніж легка програма (p < 0,0001) або програма без тренувань (p < 0,0001).

Були проведені перевірки нормальності та тест Левена, щоб переконатися, що припущення ANOVA були виконані.