lý thuyết ngày 4: Phân tích ANOVA

Nhập dữ liệu

A = c(8, 9, 11, 4, 7, 8, 5)
B = c(7, 17, 10, 14, 12, 24, 11, 22)
C = c(28, 21, 26, 11, 24, 19)
D = c(26, 16, 13, 12, 9, 10, 11, 17, 15)
x = c(A, B, C, D)
group = c(rep("A", 7), rep("B", 8), rep("C", 6),
rep("D", 9))
data = data.frame(x, group)
data

##     x group
## 1   8     A
## 2   9     A
## 3  11     A
## 4   4     A
## 5   7     A
## 6   8     A
## 7   5     A
## 8   7     B
## 9  17     B
## 10 10     B
## 11 14     B
## 12 12     B
## 13 24     B
## 14 11     B
## 15 22     B
## 16 28     C
## 17 21     C
## 18 26     C
## 19 11     C
## 20 24     C
## 21 19     C
## 22 26     D
## 23 16     D
## 24 13     D
## 25 12     D
## 26  9     D
## 27 10     D
## 28 11     D
## 29 17     D
## 30 15     D

Phân tích ANOVA

av = aov(x ~ group)
summary(av)

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## group        3  642.3  214.09   8.197 0.000528 ***
## Residuals   26  679.1   26.12                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Có sự khác biệt giữa các nhóm

Phân tích hậu định (posthoc analysis)

• Có sự khác biệt giữa các nhóm
• Nhưng giữa nhóm nào? # Lựa chọn phương pháp nào # Phương pháp Bonferroni
• Ngưỡng 𝛼 toàn diện (global) vd: 𝛼 = 0.05
• nếu có k test, thì adjusted alpha = alpha/k # Phương pháp nào cho ra kết quả với khoảng tin cậy ngắn nhất là phương pháp tối ưu nhất # R code - Tukey’s Method

av = aov(x ~ group)
TukeyHSD(av)

##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = x ~ group)
## 
## $group
##           diff          lwr        upr     p adj
## B-A  7.1964286  -0.05969765 14.4525548 0.0525014
## C-A 14.0714286   6.27132726 21.8715299 0.0002134
## D-A  6.9047619  -0.16073856 13.9702624 0.0571911
## C-B  6.8750000  -0.69675602 14.4467560 0.0850381
## D-B -0.2916667  -7.10424368  6.5209103 0.9994049
## D-C -7.1666667 -14.55594392  0.2226106 0.0597131

Tuy nhiên nếu trình bày như trên rất rối

Nên trình bày băng biểu đồ

plot(TukeyHSD(av), las=2 , col="brown")

# Tóm lược # • Phân tích phương sai là phương pháp tách phương sai theo từng nguồn - Có nhiều phương pháp phân tích hậu định (posthoc) - Phương pháp với khoảng tin cậy ngắn nhất là tối ưu nhấtà TukeyHSD? # R codes cho biểu đồ false positive

k = 1:100
alpha0.05 = 1-(1-0.05)^k
alpha0.01 = 1-(1-0.01)^k
alpha0.001 = 1-(1-0.001)^k
k = c(1:100, 1:100, 1:100)
alpha = c(alpha0.05, alpha0.01, alpha0.001)
P.value = c(rep(0.05, 100), rep(0.01, 100), rep(0.001, 100))
dat = data.frame(alpha, k, P.value)

library(ggplot2)

p = ggplot(data=dat, aes(x=k, y=alpha, col=factor(P.value))) + geom_line()

p

p + scale_y_continuous(breaks=seq(0, 1, 0.05)) + scale_x_continuous(breaks=seq(0, 100, 10)) + theme(legend.position="top") + labs(x="Số lần kiểm định giả thuyết", y="Xác suất phát hiện ít nhất 1 kết quả significant")

R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

summary(cars)

##      speed           dist       
##  Min.   : 4.0   Min.   :  2.00  
##  1st Qu.:12.0   1st Qu.: 26.00  
##  Median :15.0   Median : 36.00  
##  Mean   :15.4   Mean   : 42.98  
##  3rd Qu.:19.0   3rd Qu.: 56.00  
##  Max.   :25.0   Max.   :120.00

Including Plots

You can also embed plots, for example:

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.