This script reproduces all of the analysis and graphs for the MANOVA of the Iris data in the paper and also includes other analyses not described there. It is set up as an R script that can be “compiled” to HTML, Word, or PDF using knitr::knit(). This is most convenient within R Studio via the File -> Compile Notebook option.

Load packages and the data

library(heplots)
library(car)    # actually, loaded by heplots
data(iris)

Initial scatterplots and data ellipses

op <- par(mfcol=c(1,2), mar=c(5,4,1,1)+.1)
scatterplot(Sepal.Width ~ Sepal.Length | Species, data=iris,
        ellipse=TRUE, levels=0.68, smoother=NULL, reg.line=FALSE, grid=FALSE,
        legend.coords=list(x=7, y=4.4), col=c("red", "blue", "black"))  

scatterplot(Sepal.Width ~ Sepal.Length | Species, data=iris,
        ellipse=TRUE, levels=0.68, smoother=NULL, grid=FALSE,
        reg.line=FALSE, cex=0, 
        legend.plot=FALSE, col=c("red", "blue", "black"))

par(op)

Using the covEllipse function

Uncentered and centered, first two variables

covEllipses(iris[,1:4], iris$Species, 
    fill=c(rep(FALSE,3), TRUE))

covEllipses(iris[,1:4], iris$Species, center=TRUE, 
    fill=c(rep(FALSE,3), TRUE), fill.alpha=.1, label.pos=c(1:3,0))

All pairs when more than two are specified

covEllipses(iris[,1:4], iris$Species, 
    fill=c(rep(FALSE,3), TRUE), variables=1:4, 
    fill.alpha=.1)

covEllipses(iris[,1:4], iris$Species, center=TRUE,
    fill=c(rep(FALSE,3), TRUE), variables=1:4, 
    label.pos=c(1:3,0), fill.alpha=.1)

view in PCA space

NB: scale.=FALSE by default

iris.pca <- prcomp(iris[,1:4])
summary(iris.pca)
## Importance of components:
##                          PC1    PC2    PC3     PC4
## Standard deviation     2.056 0.4926 0.2797 0.15439
## Proportion of Variance 0.925 0.0531 0.0171 0.00521
## Cumulative Proportion  0.925 0.9777 0.9948 1.00000
op <- par(mfcol=c(1,2), mar=c(5,4,1,1)+.1)
covEllipses(iris.pca$x, iris$Species, 
    fill=c(rep(FALSE,3), TRUE), 
    label.pos=1:4, fill.alpha=.1, asp=1)

covEllipses(iris.pca$x, iris$Species, 
    fill=c(rep(FALSE,3), TRUE), center=TRUE,
    label.pos=1:4, fill.alpha=.1, asp=1)

par(op)

# all variables
covEllipses(iris.pca$x, iris$Species, 
    fill=c(rep(FALSE,3), TRUE), variables=1:4, 
    label.pos=1:4, fill.alpha=.1)

covEllipses(iris.pca$x, iris$Species, center=TRUE,  
    fill=c(rep(FALSE,3), TRUE), variables=1:4, 
    label.pos=1:4, fill.alpha=.1)