# [R] Mixed Beta Disrubutions

Bob Aronoff bobaronoff at gmail.com
Wed Jan 6 19:19:29 CET 2016

```I am working to understand the same issues with my datasets.
Adapting Dr Zeileis' posting I have written a humble function.
It creates a two cluster beta mix on a vector of data.
It seems to be working well on my datasets.
You are welcome to try on yours.

regards,
Bob

_____

bi.modal.beta<-function (vals){

# vals is a vector data points in range (0,1).  There should not be any 0,1, or NA in vals

library(betareg)
library(flexmix)

d <- data.frame(y = vals)

#create bimodal beta model
m <- betamix(y ~ 1 | 1, data = d, k = 2)

#extract beta parameters (mean and precision) with anti-link functions
mu <- plogis(coef(m)[,1])
phi <- exp(coef(m)[,2])

#convert mean & precision to alpha & beta
a <- mu * phi
b <- (1 - mu) * phi

#report parameters
print("alpha:")
print( a)
print("beta:")
print(b)

#plot in order to inspect result
ys <- seq(0, 1, by = 0.01)
p <- prior(m\$flexmix)
# p is equivalent to the percentage of each distribution

hist(d\$y, breaks = 0:25/25, freq = FALSE,
main = "Bimodal Betamix", xlab = "values")
lines(ys, p * dbeta(ys, shape1 = a, shape2 = b) , lwd = 2, col="red")
lines(ys, p * dbeta(ys, shape1 = a, shape2 = b) , lwd = 2,col="blue")
}

____

```