Three ANOVAs

Intro to the Experiment

This experiment showcased the lifespan of female mice who were assigned to differing differing diets. The treatment groups, or diets, included NP, N/N85, N/R50, R/R50, N/R50 lopro, and N/R40.

library(tidyverse)
library(openintro)
library(lmPerm)
library(DescTools)
library(Sleuth3)
library(grid)

data("case0501")
names(case0501)

## [1] "Lifetime" "Diet"

Part A

ggplot(data=case0501, aes(x=Diet, y=Lifetime, fill=Diet)) +
  geom_boxplot(outlier.color="red", outlier.size=2) +
  labs(title="Lifetime of Mice") +
  theme_classic()

summary(
  ANOVA.model <- aov(Lifetime ~ Diet, data=case0501)
)

##              Df Sum Sq Mean Sq F value Pr(>F)    
## Diet          5  12734  2546.8    57.1 <2e-16 ***
## Residuals   343  15297    44.6                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

PostHocTest(ANOVA.model, method="scheffe", conf.level=NA, ordered=FALSE)

## 
##   Posthoc multiple comparisons of means: Scheffe Test 
## 
## $Diet
##       N/N85   N/R40   N/R50   NP      R/R50 
## N/R40 < 2e-16 -       -       -       -     
## N/R50 1.1e-11 0.3289  -       -       -     
## NP    0.0063  < 2e-16 < 2e-16 -       -     
## R/R50 9.7e-12 0.6644  0.9986  < 2e-16 -     
## lopro 1.6e-05 0.0022  0.4440  1.4e-15 0.2695

After the ANOVA analysis and Posthoc test, we see that R/R50 and N/R50 are equally effective treatments.

Part B

GraphNormality <- function(model) {
  residual.data <- data.frame(e=model$residuals)
  H             <- shapiro.test(residual.data$e)
  ggplot(residual.data, aes(sample=e)) + 
    stat_qq() +
    geom_abline(color="blue", intercept=mean(residual.data$e), slope=sd(residual.data$e)) +
    labs(title="Are the Residuals Normally Distributed ", 
         subtitle = paste("Shapiro-Wilks test p-value =", signif(H$p.value,5)) ) +
    theme_classic()
}

GraphHomogeneity <- function(response, predictor, dataset=NULL) {
  H <- bartlett.test(response~predictor)    
  ggplot(data=dataset, aes(x=predictor, y=response, fill=predictor)) +
    geom_boxplot(outlier.color = "red", outlier.size=3) +
    geom_jitter(width=0.2) +
    labs( title="Homogenous Variances",
         subtitle=paste("Bartlett's test p-value =", signif(H$p.value,5)) ) +    
    theme_classic()
}

GraphNormality(ANOVA.model)

GraphHomogeneity(case0501$Lifetime, case0501$Diet, dataset=case0501)

Anova assumptions were violated because the graphs do not appear to be normally distributed or homogenous.

Part C

summary(
  perm.model <- aovp(Lifetime ~ Diet, data=case0501)
)

## [1] "Settings:  unique SS "

## Component 1 :
##              Df R Sum Sq R Mean Sq Iter  Pr(Prob)    
## Diet          5    12734    2546.8 5000 < 2.2e-16 ***
## Residuals   343    15297      44.6                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

class(perm.model) 

## [1] "aovp" "aov"  "lmp"  "lm"

PostHocTest(perm.model, method="scheffe", conf.level=NA, ordered=FALSE)

## 
##   Posthoc multiple comparisons of means: Scheffe Test 
## 
## $Diet
##       N/N85   N/R40   N/R50   NP      R/R50 
## N/R40 < 2e-16 -       -       -       -     
## N/R50 1.1e-11 0.3289  -       -       -     
## NP    0.0063  < 2e-16 < 2e-16 -       -     
## R/R50 9.7e-12 0.6644  0.9986  < 2e-16 -     
## lopro 1.6e-05 0.0022  0.4440  1.4e-15 0.2695

( kw.test <- kruskal.test(Lifetime ~ Diet, data=case0501) )

## 
##  Kruskal-Wallis rank sum test
## 
## data:  Lifetime by Diet
## Kruskal-Wallis chi-squared = 159.01, df = 5, p-value < 2.2e-16

class(kw.test)

## [1] "htest"

The results of the alternate nonparametric tests were in fact consistent with the previous ANOVA test findings.

There does appear to be clear “best” and “worst” diets, but that being said there are also diets that are considered to be “just as good” or “just as bad”.

Source Reference

R. Weindruch, R. L. Walford, S Fligiel, and D. Guthrie, “The Retardation of Aging in Mice by Dietary Restriction: Longevity, Cancer, Immunity, and Lifetime Energy Intake,” Journal of Nutrition 116(4) (1986):641-54).

…