Add confidence ellipse to LDA ordination plot

Use ggord to plot LDA ordination plot

Installation

devtools::install_github('fawda123/ggord')

Basic LDA ordination biplot

library(MASS)
ord <- lda(Species ~ ., iris, prior = rep(1, 3)/3)

library(ggord)
p <- ggord(ord, iris$Species)
p

A function to compute confidence ellipse

get_lda_ell <- function(ord_in, grp_in, ellipse_pro = 0.97){
    ## adapted from https://github.com/fawda123/ggord/blob/master/R/ggord.R
    require(plyr)
    axes = c('LD1', 'LD2')
    obs <- data.frame(predict(ord_in)$x[, axes])
    obs$Groups <- grp_in
    names(obs)[1:2] <- c('one', 'two')
    theta <- c(seq(-pi, pi, length = 50), seq(pi, -pi, length = 50))
    circle <- cbind(cos(theta), sin(theta))
    ell <- ddply(obs, 'Groups', function(x) {
        if(nrow(x) <= 2) {
            return(NULL)
        }
        sigma <- var(cbind(x$one, x$two))
        mu <- c(mean(x$one), mean(x$two))
        ed <- sqrt(qchisq(ellipse_pro, df = 2))
        data.frame(sweep(circle %*% chol(sigma) * ed, 2, mu, FUN = '+'))
    })
    names(ell)[2:3] <- c('one', 'two')
    ell <- ddply(ell, .(Groups), function(x) x[chull(x$one, x$two), ])
    ell
}

Using the function

library(ggplot2)
anotherEll <- get_lda_ell(ord, iris$Species, 0.97)
## Loading required package: plyr
p + geom_polygon(data = anotherEll, 
                 aes_string(color = 'Groups', group = 'Groups'),
                 lty=2, fill = NA)

comments powered by Disqus