I am using the SIBER package to develop a plot with three groups, I would like to color each group according to a legend. The SIBER pdf says to use the clr command, however, when using this command, this colors each of the boxes for all three x axis groups the same color. The pdf says this is basically an hdr boxplot. I am curious if anyone knows how to color each box according to my legend? It seems like the reason this is difficult is because each x axis group has a box for 3 different CI, therefore requiring three different shades of the same color for each group. the code is below.

ylab = expression("Standard Ellipse Area " ('\u2030' ^2) ),
bty = "L",
las = 1,
main = "",
prn = TRUE,# this prints the CI
clr = # this is where I am not sure what to do to color each of the boxes according to my legend

#so this is using the iris set. You would need rjags (which you have to download
#from the internet and then install package) and the SIBER package. Thanks!
library(ellipse)
library(rjags)
library(SIBER)
data(iris)
str(iris)
irisSIBERdat_1 <- as.data.frame(cbind(iris$Sepal.Length,iris$Sepal.Width, iris$Species, rep(1,length.out = length(iris$Sepal.Width))))
colnames(irisSIBERdat_1) = c("iso1", "iso2", "group", "community")
my.siber.data_1 <- createSiberObject(irisSIBERdat_1)
group.ML_1 <- groupMetricsML(my.siber.data_1)
print(group.ML_1)
#options for running jags
parms <- list()
parms$n.iter <- 2 * 10^4 # number of iterations to run the model for
parms$n.burnin <- 1 * 10^3 # discard the first set of values
parms$n.thin <- 10 # thin the posterior by this many
parms$n.chains <- 2 # run this many chains
#define the priors
priors <- list()
priors$R <- 1 * diag(2)
priors$k <- 2
priors$tau.mu <- 1.0E-3
#fit the ellipses which uses an Inverse Wishart prior
#on the covariance matrix Sigma, and a vague normal prior on the
#means. Fitting is via the JAGS method.
ellipses.posterior <- siberMVN(my.siber.data_1, parms, priors)
SEA.B <- siberEllipses(ellipses.posterior)
#below are the colors I want for group one through three
clr_1<- c("royalblue1","dodgerblue3", "deepskyblue3")
clr_2<- c("tomato","salmon", "lightcoral")
clr_3<- c("limegreen","green1", "mediumspringgreen")
siberDensityPlot(SEA.B, xticklabels = c("group 1", "group 2", "group 3"),
ylab = expression("Standard Ellipse Area " ('\u2030' ^2) ),
bty = "L",
las = 1,
main = "",
prn = TRUE,# this prints the CI
clr = # this is where I am not sure what to do to color each of the groups different shades according to my legend on a separate graph
)

For illustration, I choose the confidence intervals to be displayed as 99, 75, 50, 25 (see probs parameter). I gave each of them a different colour. As you can see, the colors correspond to the same confidence interval in each boxplot, and this help comparison between groups. Coloring the different groups in different colors would not make sense as you loose the visual aid for comparing the same confidence interval across different groups.