misturaGama <- function(pis, as, bs, y, TO = 12, YLIM = 0.30, ...){
    stopifnot(length(pis) == length(as))
    stopifnot(length(pis) == length(bs))
    logCs <- as*log(bs) - lgamma(as)
    n <- length(y)
    sY <- sum(y)
    as.post <- as + sY
    bs.post <- bs + n
    logCs.post <- as.post * log(bs.post) - lgamma(as.post)
    ws <- pis * exp(logCs - logCs.post)
    pis.post <- ws/sum(ws)
    f <- function(theta){
        sapply(theta, function(x) sum(pis * dgamma(x, as, bs)))
    }
    f.post <- function(theta){
        sapply(theta, function(x) sum(pis.post * dgamma(x, as.post, bs.post)))
    }
    curve(f.post, 0, TO, ylim=c(0, YLIM),
          xlab = expression(theta), ylab=expression(paste("f(",theta,"|y)")))
    curve(dgamma(x, sY, n), add = T,  0, TO, col = "blue")
    curve(f, 0, TO, add = T, col = "red")
    return(list(priori = data.frame(pesos = pis, alpha=as, beta = bs),
                post = data.frame(pesos = pis.post, alpha=as.post, beta = bs.post),
                dados = c(n = n, somaY = sY)
           ))
}

## Dois componentes e amostra de tamanho 1 
misturaGama(c(0.5, 0.5), as = c(4, 12), bs = c(2,2), y = 3, YLIM = 0.35)
misturaGama(c(0.5, 0.5), as = c(4, 12), bs = c(2,2), y = 4, YLIM = 0.35)
misturaGama(c(0.5, 0.5), as = c(4, 12), bs = c(2,2), y = 5, YLIM = 0.35)

## amostras de tamanho 2
misturaGama(c(0.5, 0.5), as = c(4, 12), bs = c(2,2), y = c(3, 3), YLIM = 0.35)
misturaGama(c(0.5, 0.5), as = c(4, 12), bs = c(2,2), y = c(4, 4), YLIM = 0.35)
misturaGama(c(0.5, 0.5), as = c(4, 12), bs = c(2,2), y = c(5, 5), YLIM = 0.35)

## alterando pesos e parâmetros
misturaGama(c(0.35, 0.65), as = c(4, 12), bs = c(2,2), y = 4)
misturaGama(c(0.35, 0.65), as = c(4, 18), bs = c(2,3), y = 4)
misturaGama(c(0.5, 0.5), as = c(4, 18), bs = c(2,3), y = 3)
misturaGama(c(0.5, 0.5), as = c(4, 18), bs = c(2,3), y = 4)

## ammostra de tamanho 3
misturaGama(c(0.5, 0.5), as = c(4, 18), bs = c(2,3), y = c(3, 4, 5), YLIM = 0.40)
misturaGama(c(0.35, 0.65), as = c(4, 12), bs = c(2,2), y = c(3, 4, 5), YLIM = 0.40)
misturaGama(c(0.35, 0.65), as = c(4, 12), bs = c(2,2), y = c(3, 3, 4), YLIM = 0.40)

## três componentes e amostra de tamanho 1
misturaGama(c(0.3, 0.3, 0.4), as = c(4, 12, 20), bs = c(2,2,2), y = 3)
misturaGama(c(0.5, 0.3, 0.2), as = c(4, 12, 20), bs = c(2,2,2), y = 3, YLIM = 0.32)