COMMunication
RESearch.NET

User Tools

Site Tools

You are here: index » r » twoway_anova
Trace:
r:twoway_anova

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision		Previous revision
r:twoway_anova [2017/11/06 08:57] – [3] hkimscilr:twoway_anova [2018/10/19 08:39] (current) – [e.g. 5] hkimscil
Line 1: Line 1:
 +data sample: {{:r:diet.csv}}
 +
 ====== Twoway ANOVA ====== ====== Twoway ANOVA ======
 recap of oneway ANOVA recap of oneway ANOVA
Line 78: Line 80:
 45 Carrier 3    Office 5 15.04 45 Carrier 3    Office 5 15.04
  
-> library(ggplot2) +> library("ggpubr") 
-> ggplot(delivery.df, aes(Time, Destination, colour = Service)) + geom_point() +> ggline(delivery.df, x = "Destination", y = "Time", color = "Service", add = c("mean_se", "dotplot"))
-> +
 </code> </code>
  
Line 217: Line 218:
 [{{:r:interaction effect eg04.png|Interaction effect example 4}}] [{{:r:interaction effect eg04.png|Interaction effect example 4}}]
  
-===== 5 ===== +===== e.g., 5 ===== 
-download [[http://dl.dropbox.com/u/10246536/Web/RTutorialSeries/dataset_anova_twoWay_comparisons.csv|dataset_anova_twoWay_interactions.csv]] +download {{:r:dataset_anova_twoWay_comparisons.csv|dataset_anova_twoWay_interactions.csv}} 
  
-<code>data <- read.csv("dataset_ANOVA_TwoWayComparisons.csv")+<code>stressdata <- read.csv("http://commres.net/wiki/_media/r/dataset_anova_twoway_comparisons.csv")
 > #display the data > #display the data
-> data+> stressdata
    Treatment   Age StressReduction    Treatment   Age StressReduction
 1     mental young              10 1     mental young              10
Line 252: Line 253:
 27   medical   old               0</code> 27   medical   old               0</code>
  
-<code>> aov(d$StressReduction~d$Treatment*d$Age, d) +<code> 
-Call: +> stressdata <- read.csv("http://commres.net/wiki/_media/r/dataset_anova_twoway_comparisons.csv") 
-   aov(formula = d$StressReduction ~ d$Treatment * d$Age, data = d) +>  
- +> a.mod <- aov(StressReduction~Treatment*Age, data=stressdata)
-Terms: +
-                d$Treatment d$Age d$Treatment:d$Age Residuals +
-Sum of Squares           18   162                 0        18 +
-Deg. of Freedom           2     2                 4        18 +
- +
-Residual standard error: 1 +
-Estimated effects may be unbalanced +
-> a.mod <- aov(d$StressReduction~d$Treatment*d$Age, d)+
 > summary(a.mod) > summary(a.mod)
-                  Df Sum Sq Mean Sq F value  Pr(>F)     +              Df Sum Sq Mean Sq F value  Pr(>F)     
-d$Treatment        2     18       9       9 0.00195 **  +Treatment      2     18       9       9 0.00195 **  
-d$Age              2    162      81      81   1e-09 *** +Age            2    162      81      81   1e-09 *** 
-d$Treatment:d$Age  4      0       0       0 1.00000     +Treatment:Age  4      0       0       0 1.00000     
-Residuals         18     18       1                    +Residuals     18     18       1                    
 --- ---
 Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
-></code> +>  
 +</code> 
  
-<code>> pairwise.t.test(d$StressReduction, d$Treatment, p.adj = "none")+<code>> TukeyHSD(a.mod, which="Treatment") 
 +  Tukey multiple comparisons of means 
 +    95% family-wise confidence level
  
- Pairwise comparisons using t tests with pooled SD +Fit: aov(formula = StressReduction ~ Treatment * Age, data = stressdata)
  
-data:  d$StressReduction and d$Treatment +$Treatment 
 +                 diff        lwr       upr     p adj 
 +mental-medical      2  0.7968988 3.2031012 0.0013531 
 +physical-medical    1 -0.2031012 2.2031012 0.1135025 
 +physical-mental    -1 -2.2031012 0.2031012 0.1135025
  
-         medical mental +>  
-mental   0.13    -      +</code>  
-physical 0.45    0.45  +<code> 
 +> TukeyHSD(a.mod, which="Age") 
 +  Tukey multiple comparisons of means 
 +    95% family-wise confidence level
  
-P value adjustment method: none</code>  +Fit: aov(formula = StressReduction ~ Treatment * Age, data = stressdata)
-<code>> pairwise.t.test(d$StressReduction, d$Age, p.adj = "none")+
  
- Pairwise comparisons using t tests with pooled SD +$Age 
 +          diff       lwr       upr    p adj 
 +old-mid     -3 -4.203101 -1.796899 1.54e-05 
 +young-mid    3  1.796899  4.203101 1.54e-05 
 +young-old    6  4.796899  7.203101 0.00e+00 
 +> 
 +</code> 
  
-data:  d$StressReduction and d$Age  
- 
-      mid     old     
-old   2.5e-05 -       
-young 2.5e-05 2.3e-10 
- 
-P value adjustment method: none  
->  
-</code> 
  
 <WRAP info>위는 아래의 linear model 을 이용하여도 가능. 사실 모든 ANOVA 테스트는 linear model이기도 함(lm) <WRAP info>위는 아래의 linear model 을 이용하여도 가능. 사실 모든 ANOVA 테스트는 linear model이기도 함(lm)
Line 317: Line 316:
 </code> </code>
  
-<code>> pairwise.t.test(d$StressReduction, d$Treatment, p.adj = "none") 
  
- Pairwise comparisons using t tests with pooled SD +===== e.g., 6 =====
  
-data:  d$StressReduction and d$Treatment +<code>tdata <- ToothGrowth 
 +tdata 
 +str(tdata) 
 +tdata$dose <- factor(tdata$dose) 
 +# or a nicer way 
 +tdata$dose <- factor(tdata$dose, levels=c(0.5, 1, 2), labels=c("dhalf", "d1", "d2"))
  
-         medical mental +tdata
-mental   0.13    -      +
-physical 0.45    0.45  +
  
-P value adjustment method: none</code>  +table(tdata$supp, tdata$dose) 
-<code>> pairwise.t.test(d$StressReduction, d$Age, p.adj = "none")+</code>
  
- Pairwise comparisons using t tests with pooled SD +==== Visualize the data ====
  
-data:  d$StressReduction and d$Age +<code>install.packages("ggpubr")
  
-      mid     old     +library("ggpubr") 
-old   2.5e-05 -       +ggboxplot(tdata, x = "dose", y = "len", color = "supp", 
-young 2.5e-05 2.3e-10+          palette = c("red", "blue")) 
 +ggline(tdata, x = "dose", y = "len", color = "supp", 
 +       add = c("mean_se", "dotplot"), 
 +       palette = c("red", "blue")) 
 + 
 +# prettier 
 +boxplot(len ~ supp * dose, data=tdata, frame = FALSE,  
 +      col = c("#00AFBB", "#E7B800"), ylab="Tooth Length") 
 + 
 +# or  
 +interaction.plot(x.factor = tdata$dose, trace.factor =  
 +      tdata$supp, response = tdata$len, fun = mean,  
 +      type = "b", legend = TRUE,  
 +      xlab = "Dose", ylab="Tooth Length", 
 +      pch=c(1,19), col = c("#00AFBB", "#E7B800"))
  
-P value adjustment method: none  
->  
 </code> </code>
  
 +==== Compute ANOVA test ====
  
 +<code>res.aov2 <- aov(len ~ supp + dose, data = tdata)
 +summary(res.aov2)
 +            Df Sum Sq Mean Sq F value   Pr(>F)    
 +supp         1  205.4   205.4   14.02 0.000429 ***
 +dose         2 2426.4  1213.2   82.81  < 2e-16 ***
 +Residuals   56  820.4    14.7                     
 +---
 +Signif. codes:  
 +0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1</code>
 +
 +
 +<code>res.aov3 <- aov(len ~ supp * dose, data = tdata)
 +res.aov3 <- aov(len ~ supp + dose + supp:dose, data = tdata)
 +summary(res.aov3)
 +            Df Sum Sq Mean Sq F value   Pr(>F)    
 +supp         1  205.4   205.4  15.572 0.000231 ***
 +dose         2 2426.4  1213.2  92.000  < 2e-16 ***
 +supp:dose    2  108.3    54.2   4.107 0.021860 *  
 +Residuals   54  712.1    13.2                     
 +---
 +Signif. codes:  
 +0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
 +</code>
 +
 +From the ANOVA results, you can conclude the following, based on the p-values and a significance level of 0.05:
 +
 +  * the p-value of supp is 0.000429 (significant), which indicates that the levels of supp are associated with significant different tooth length.
 +  * the p-value of dose is < 2e-16 (significant), which indicates that the levels of dose are associated with significant different tooth length.
 +  * the p-value for the interaction between supp*dose is 0.02 (significant), which indicates that the relationships between dose and tooth length depends on the supp method.
 +
 +==== pair-wise test ====
 +
 +<code>> TukeyHSD(res.aov3, which = "dose")
 +  Tukey multiple comparisons of means
 +    95% family-wise confidence level
 +
 +Fit: aov(formula = len ~ supp * dose, data = tdata)
 +
 +$dose
 +        diff       lwr       upr   p adj
 +d1-dh  9.130  6.362488 11.897512 0.0e+00
 +d2-dh 15.495 12.727488 18.262512 0.0e+00
 +d2-d1  6.365  3.597488  9.132512 2.7e-06
 +</code>
 +
 +<code>pairwise.t.test(tdata$len, tdata$dose, p.adjust.method = "BH")
 +</code>
 +==== Homogeneity ====
 +
 +<code>plot(res.aov3, 1)</code>
 +
 +<code>install.packages("car")
 +library(car)
 +leveneTest(len ~ supp*dose, data = tdata)
 +
 +Levene's Test for Homogeneity of Variance (center = median)
 +      Df F value Pr(>F)
 +group  5  1.7086 0.1484
 +      54               
 +</code>
 +p value is greater than .05, which indicates that there is no evidence of significant difference between variances of groups.
 +
 +==== check normality ====
 +<code>plot(res.aov3, 2)
 +</code>
 +The residuals are plotted against the quantiles of normal distribution (the straight line).
 +
 +==== unbalance design ====
 +<code>> library(car)
 +> my_anova <- aov(len ~ supp * dose, data = tdata)
 +> 
 +> Anova(my_anova, type="III"
 ++ )
 +Anova Table (Type III tests)
 +
 +Response: len
 +             Sum Sq Df F value    Pr(>F)    
 +(Intercept) 1750.33  1 132.730 3.603e-16 ***
 +supp         137.81  1  10.450  0.002092 ** 
 +dose         885.26  2  33.565 3.363e-10 ***
 +supp:dose    108.32  2   4.107  0.021860 *  
 +Residuals    712.11 54                      
 +---
 +Signif. codes:  
 +0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
 +</code>
r/twoway_anova.1509928056.txt.gz · Last modified: by hkimscil
Donate Powered by PHP Valid HTML5 Valid CSS Driven by DokuWiki