How to estimate the source of variation in R? (For two-way ANOVA)

Viewed 199

In my laboratory we used a statistical software package, but now we are using R for the statistics.

Using the statistical software, after a two-way ANOVA, a table with the results includes "Source of variation" and "% of total variation". How can I get those values from R?

Example: This data set

Temperature <- factor(c(rep("cold", times = 4),
                        rep("hot", times = 4)),
                      levels = c("cold", "hot"))

Light <- factor(rep(c(rep("blue", times = 2),
                      rep("yellow", times = 2)),
                    times = 2),
                levels = c("blue", "yellow"))

Result <- c(90.40, 85.20, 21.70, 25.30,
            75.12, 77.36, 6.11, 10.8)

Data <- data.frame(Temperature, Light, Result)

enter image description here

2-way ANOVA in R:

two_wayANOVA <- aov(data = Data,
                    formula = Result ~ Temperature * Light)
summary(two_wayANOVA)


 Df Sum Sq Mean Sq  F value  Pr(>F)    
Temperature        1    354     354   42.250 0.00289 ** 
Light              1   8723    8723 1041.366 5.5e-06 ***
Temperature:Light  1      6       6    0.725 0.44250    
Residuals          4     34       8                     
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

In the statistical software:

enter image description here

1 Answers

As @yh6 commented, one can determine the source of variation by dividing the SumSq by the total sum of "Sum Sq".

The broom package has the convenient tidy function to take the output from the aov function and make it into a data frame for easy handling.

library(broom)
output <- tidy(two_wayANOVA)

output$variationSource <- output$sumsq/sum(output$sumsq)*100
output


# A tibble: 4 x 7
  term                 df   sumsq  meansq statistic     p.value variationSource
  <chr>             <dbl>   <dbl>   <dbl>     <dbl>       <dbl>           <dbl>
1 Temperature           1  354.    354.      42.2    0.00289              3.88  
2 Light                 1 8723.   8723.    1041.     0.00000550          95.7   
3 Temperature:Light     1    6.07    6.07     0.725  0.442                0.0666
4 Residuals             4   33.5     8.38    NA      NA                   0.368 

Update
To recreate the significance column, this add this line to the script:

output$sig <- cut(output$p.value, breaks=c(0,  0.001,  0.01,  0.05,  0.1, 1), 
                                  labels=c( '***',  '**' , '*' , '.', ' ' ))
Related