ANOVA

Getting Started

To demonstrate the various types of sums of squares, we’ll create a data frame called df_disease taken from the SAS documentation. The corresponding data can be found here.

The Model

For this example, we’re testing for a significant difference in stem_length using ANOVA. In R, we’re using lm() to run the ANOVA, and then using broom::glance() and broom::tidy() to view the results in a table format.

lm_model <- lm(y ~ drug + disease + drug*disease, df_disease)

The glance function gives us a summary of the model diagnostic values.

lm_model %>% 
  glance() %>% 
  pivot_longer(everything())

# A tibble: 12 × 2
   name               value
   <chr>              <dbl>
 1 r.squared        0.456  
 2 adj.r.squared    0.326  
 3 sigma           10.5    
 4 statistic        3.51   
 5 p.value          0.00130
 6 df              11      
 7 logLik        -212.     
 8 AIC            450.     
 9 BIC            477.     
10 deviance      5081.     
11 df.residual     46      
12 nobs            58      

The tidy function gives a summary of the model results.

lm_model %>% tidy()

# A tibble: 12 × 5
   term           estimate std.error statistic      p.value
   <chr>             <dbl>     <dbl>     <dbl>        <dbl>
 1 (Intercept)      29.3        4.29    6.84   0.0000000160
 2 drug2            -1.33       6.36   -0.210  0.835       
 3 drug3           -13          7.43   -1.75   0.0869      
 4 drug4           -15.7        6.36   -2.47   0.0172      
 5 disease2         -1.08       6.78   -0.160  0.874       
 6 disease3         -8.93       6.36   -1.40   0.167       
 7 drug2:disease2    6.58       9.78    0.673  0.504       
 8 drug3:disease2  -10.9       10.2    -1.06   0.295       
 9 drug4:disease2    0.317      9.30    0.0340 0.973       
10 drug2:disease3   -0.900      9.00   -0.100  0.921       
11 drug3:disease3    1.10      10.2     0.107  0.915       
12 drug4:disease3    9.53       9.20    1.04   0.306       

The Results

You’ll see that R print the individual results for each level of the drug and disease interaction. We can get the combined F table in R using the anova() function on the model object.

lm_model %>% 
  anova() %>% 
  tidy() %>% 
  kable()

term	df	sumsq	meansq	statistic	p.value
drug	3	3133.2385	1044.4128	9.455761	0.0000558
disease	2	418.8337	209.4169	1.895990	0.1617201
drug:disease	6	707.2663	117.8777	1.067225	0.3958458
Residuals	46	5080.8167	110.4525	NA	NA

We can add a Total row, by using add_row and calculating the sum of the degrees of freedom and sum of squares.

lm_model %>%
  anova() %>%
  tidy() %>%
  add_row(term = "Total", df = sum(.$df), sumsq = sum(.$sumsq)) %>% 
  kable()

term	df	sumsq	meansq	statistic	p.value
drug	3	3133.2385	1044.4128	9.455761	0.0000558
disease	2	418.8337	209.4169	1.895990	0.1617201
drug:disease	6	707.2663	117.8777	1.067225	0.3958458
Residuals	46	5080.8167	110.4525	NA	NA
Total	57	9340.1552	NA	NA	NA

Sums of Squares Tables

Unfortunately, it is not easy to get the various types of sums of squares calculations in using functions from base R. However, the rstatix package offers a solution to produce these various sums of squares tables. For each type, you supply the original dataset and model to the. anova_test function, then specify the ttype and se detailed = TRUE.

Type I

df_disease %>% 
  rstatix::anova_test(
    y ~ drug + disease + drug*disease, 
    type = 1, 
    detailed = TRUE) %>% 
  rstatix::get_anova_table() %>% 
  kable()

Effect	DFn	DFd	SSn	SSd	F	p	p<.05	ges
drug	3	46	3133.239	5080.817	9.456	5.58e-05	*	0.381
disease	2	46	418.834	5080.817	1.896	1.62e-01		0.076
drug:disease	6	46	707.266	5080.817	1.067	3.96e-01		0.122

Type II

df_disease %>% 
  rstatix::anova_test(
    y ~ drug + disease + drug*disease, 
    type = 2, 
    detailed = TRUE) %>% 
  rstatix::get_anova_table() %>% 
  kable()

Effect	SSn	SSd	DFn	DFd	F	p	p<.05	ges
drug	3063.433	5080.817	3	46	9.245	6.75e-05	*	0.376
disease	418.834	5080.817	2	46	1.896	1.62e-01		0.076
drug:disease	707.266	5080.817	6	46	1.067	3.96e-01		0.122

Type III

df_disease %>% 
  rstatix::anova_test(
    y ~ drug + disease + drug*disease, 
    type = 3, 
    detailed = TRUE) %>% 
  rstatix::get_anova_table() %>% 
  kable()

Effect	SSn	SSd	DFn	DFd	F	p	p<.05	ges
(Intercept)	20037.613	5080.817	1	46	181.414	0.00e+00	*	0.798
drug	2997.472	5080.817	3	46	9.046	8.09e-05	*	0.371
disease	415.873	5080.817	2	46	1.883	1.64e-01		0.076
drug:disease	707.266	5080.817	6	46	1.067	3.96e-01		0.122

Type IV

In R there is no equivalent operation to the Type IV sums of squares calculation in SAS.