วิธีดำเนินการวิเคราะห์ความแปรปรวนแบบสองทางใน r

การวิเคราะห์ ความแปรปรวนแบบสองทาง (“การวิเคราะห์ความแปรปรวน”) ใช้เพื่อพิจารณาว่ามีความแตกต่างที่มีนัยสำคัญทางสถิติระหว่างค่าเฉลี่ยของกลุ่มอิสระสามกลุ่มขึ้นไปที่ถูกแบ่งออกเป็นสองปัจจัยหรือไม่

บทช่วยสอนนี้จะอธิบายวิธีดำเนินการวิเคราะห์ความแปรปรวนแบบสองทางใน R

ตัวอย่าง: การวิเคราะห์ความแปรปรวนแบบสองทางใน R

สมมติว่าเราต้องการพิจารณาว่าความเข้มข้นของการออกกำลังกายและเพศส่งผลต่อการลดน้ำหนักหรือไม่ ในกรณีนี้ ปัจจัยสองประการที่เรากำลังดูคือ การออกกำลังกาย และ เพศ และตัวแปรการตอบสนองคือ การลดน้ำหนัก โดย มีหน่วยวัดเป็นปอนด์

เราสามารถใช้การวิเคราะห์ความแปรปรวนแบบสองทางเพื่อตรวจสอบว่าการออกกำลังกายและเพศส่งผลต่อการลดน้ำหนักหรือไม่ และเพื่อพิจารณาว่ามีปฏิสัมพันธ์ระหว่างการออกกำลังกายและเพศต่อการลดน้ำหนักหรือไม่

เรากำลังคัดเลือกผู้ชาย 30 คนและผู้หญิง 30 คนเพื่อเข้าร่วมการทดลอง โดยสุ่มมอบหมายให้คนละ 10 คนทำตามโปรแกรมที่ไม่ออกกำลังกาย ออกกำลังกายเบา ๆ หรือออกกำลังกายหนัก ๆ เป็นเวลาหนึ่งเดือน

รหัสต่อไปนี้จะสร้าง data frame ที่เราจะใช้งานด้วย:

 #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

สำรวจข้อมูล

ก่อนที่จะปรับแบบจำลองการวิเคราะห์ความแปรปรวนแบบสองทางให้เหมาะสม เราสามารถเข้าใจข้อมูลได้ดีขึ้นโดยการค้นหาค่าเฉลี่ยและค่าเบี่ยงเบนมาตรฐานของการลดน้ำหนักสำหรับแต่ละกลุ่มการรักษาทั้งหกกลุ่มโดยใช้แพ็คเกจ 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 

นอกจากนี้เรายังสามารถสร้าง boxplot สำหรับแต่ละกลุ่มการรักษาทั้งหกกลุ่มเพื่อให้เห็นภาพการกระจายตัวของการลดน้ำหนักสำหรับแต่ละกลุ่ม:

 #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
)

เราจะเห็นได้ทันทีว่าทั้ง 2 กลุ่มที่ออกกำลังกาย หนักๆ ดูเหมือนจะมีค่าการลดน้ำหนักที่สูงกว่า เรายังเห็นได้ว่าผู้ชายมักจะมีค่าการลดน้ำหนักที่สูงกว่าผู้หญิงทั้งในกลุ่มออกกำลังกาย ที่เข้มข้น และ เบา

ต่อไป เราจะปรับแบบจำลองการวิเคราะห์ความแปรปรวนแบบสองทางให้พอดีกับข้อมูลของเรา เพื่อดูว่าความแตกต่างทางสายตาเหล่านี้มีนัยสำคัญทางสถิติจริง ๆ หรือไม่

ติดตั้งโมเดล ANOVA แบบสองทาง

ไวยากรณ์ทั่วไปสำหรับการปรับโมเดล ANOVA แบบสองทางให้เหมาะสมใน R คือ:

aov (ตัวแปรการตอบสนอง ~predictor_variable1 *predictor_variable2, data = ชุดข้อมูล)

โปรดทราบว่าเครื่องหมาย * ระหว่างตัวแปรทำนายทั้งสองบ่งชี้ว่าเราต้องการทดสอบเอฟเฟกต์การโต้ตอบระหว่างตัวแปรทำนายทั้งสองด้วย

ในตัวอย่างของเรา เราสามารถใช้โค้ดต่อไปนี้เพื่อให้พอดีกับโมเดล ANOVA แบบสองทาง โดยใช้ 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

การตรวจสอบสมมติฐานแบบจำลอง

ก่อนที่เราจะไปไกลกว่านั้น เราต้องตรวจสอบว่าเป็นไปตามสมมติฐานของแบบจำลองของเรา เพื่อให้ผลลัพธ์ของแบบจำลองของเราเชื่อถือได้ โดยเฉพาะอย่างยิ่ง การวิเคราะห์ความแปรปรวนแบบสองทางจะถือว่า:

1. ความเป็นอิสระ – การสังเกตของแต่ละกลุ่มจะต้องเป็นอิสระจากกัน เนื่องจากเราใช้ การออกแบบแบบสุ่ม จึงควรเป็นไปตามสมมติฐานนี้ ดังนั้นเราจึงไม่ต้องกังวลกับเรื่องนี้มากเกินไป

2. ความปกติ – ตัวแปรตามควรมีการแจกแจงแบบปกติโดยประมาณสำหรับการรวมกันของกลุ่มของทั้งสองปัจจัยแต่ละกลุ่ม

วิธีหนึ่งในการทดสอบสมมติฐานนี้คือการสร้างฮิสโตแกรมของส่วนที่เหลือของแบบจำลอง หากส่วนที่เหลือมีการกระจายตามปกติโดยประมาณ ก็ควรจะเป็นไปตามสมมติฐานนี้

 #define model residuals
reside <- model$residuals

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

ส่วนที่เหลือจะมีการกระจายตามปกติโดยประมาณ ดังนั้นเราจึงสรุปได้ว่าเป็นไปตามสมมติฐานภาวะปกติ

3. ความแปรปรวนที่เท่ากัน – ความแปรปรวนสำหรับแต่ละกลุ่มจะเท่ากันหรือเท่ากันโดยประมาณ

วิธีหนึ่งในการตรวจสอบสมมติฐานนี้คือทำการทดสอบ Levene เพื่อหาความเท่าเทียมกันของความแปรปรวนโดยใช้แพ็คเกจ รถ :

 #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 เราจึงสามารถสรุปได้ว่าเป็นไปตามสมมติฐานของเราในเรื่องความเท่าเทียมกันของความแปรปรวนระหว่างกลุ่ม

วิเคราะห์ความแตกต่างในการรักษา

เมื่อเราตรวจสอบแล้วว่าเป็นไปตามสมมติฐานของแบบจำลองแล้ว เราก็สามารถทำการ ทดสอบหลังการทดสอบ เพื่อระบุได้อย่างแน่ชัดว่ากลุ่มการรักษาใดที่แตกต่างกัน

สำหรับการทดสอบหลังการทดสอบ เราจะใช้ฟังก์ชัน 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-value: 0.990364)

นอกจากนี้เรายังสามารถเห็นภาพช่วงความเชื่อมั่น 95% ที่เกิดจากการทดสอบ Tukey โดยใช้ฟังก์ชัน 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) ทำการทดสอบ HSD ของ Post-hoc Tukey

สำหรับผู้ชาย โปรแกรมการออกกำลังกาย เข้มข้น ส่งผลให้น้ำหนักลดลงมากกว่าโปรแกรม เบาๆ (p < 0.0001) หรือ ไม่ออกกำลังกายเลย (p < 0.0001) อย่างมีนัยสำคัญ นอกจากนี้ ในผู้ชาย การรับประทานอาหาร เบาๆ ส่งผลให้น้ำหนักลดลงมากกว่า การไม่ออกกำลังกายอย่าง มีนัยสำคัญ (p < 0.0001)

สำหรับผู้หญิง โปรแกรมการออกกำลังกาย เข้มข้น ส่งผลให้น้ำหนักลดลงมากกว่าโปรแกรม เบาๆ (p < 0.0001) หรือ ไม่ออกกำลังกายเลย (p < 0.0001) อย่างมีนัยสำคัญ

การตรวจสอบภาวะปกติและการทดสอบ Levene ดำเนินการเพื่อตรวจสอบว่าเป็นไปตามสมมติฐานของ ANOVA