In AndrewLJackson/SIBER: Stable Isotope Bayesian Ellipses in R

# remove previously loaded items from the current environment and remove previous graphics.
rm(list=ls())
graphics.off()

# Here, I set the seed each time so that the results are comparable. 
# This is useful as it means that anyone that runs your code, *should*
# get the same results as you, although random number generators change 
# from time to time.
set.seed(1)

library(SIBER)

# load in the included demonstration dataset
data("demo.siber.data")
#
# create the siber object
siber.example <- createSiberObject(demo.siber.data)

par(mfrow=c(1,1))

community.hulls.args <- list(col = 1, lty = 1, lwd = 1)
group.ellipses.args  <- list(n = 100, p.interval = 0.95, lty = 1, lwd = 2)
group.hull.args      <- list(lty = 2, col = "grey20")

# this time we will make the points a bit smaller by 
# cex = 0.5
plotSiberObject(siber.example,
                  ax.pad = 2, 
                  hulls = F, community.hulls.args, 
                  ellipses = F, group.ellipses.args,
                  group.hulls = F, group.hull.args,
                  bty = "L",
                  iso.order = c(1,2),
                  xlab=expression({delta}^13*C~'permille'),
                  ylab=expression({delta}^15*N~'permille'),
                  cex = 0.5
                  )



addEllipse(mu = siber.example$ML.mu[[1]][,,1],
           sigma = siber.example$ML.cov[[1]][,,1],p.interval = 0.95)
addEllipse(mu = siber.example$ML.mu[[1]][,,2],
           sigma = siber.example$ML.cov[[1]][,,2],p.interval = 0.95)
addEllipse(mu = siber.example$ML.mu[[1]][,,3],
           sigma = siber.example$ML.cov[[1]][,,3],p.interval = 0.95)
addEllipse(mu = siber.example$ML.mu[[2]][,,1],
           sigma = siber.example$ML.cov[[2]][,,1],p.interval = 0.95)
addEllipse(mu = siber.example$ML.mu[[2]][,,2],
           sigma = siber.example$ML.cov[[2]][,,2],p.interval = 0.95)
addEllipse(mu = siber.example$ML.mu[[2]][,,3],
           sigma = siber.example$ML.cov[[2]][,,3],p.interval = 0.95)




# You can add more ellipses by directly calling plot.group.ellipses()
# Add an additional p.interval % prediction ellilpse
plotGroupEllipses(siber.example, n = 100, p.interval = 0.95,
                    lty = 2, lwd = 2)