Assuming that all the mixture parameters are known, how can one sample from a mixture of $\texttt{Gamma}(\alpha,\beta)$ distributions: $$\theta \sim \pi \texttt{Gamma}(\alpha_1,\beta_1)+(1-\pi)\texttt{Gamma}(\alpha_2,\beta_2)$$ Is it done by sampling a whole bunch of $\theta$s from either of the distributions and putting them in a pool and at the end computing their corresponding probability according to the mixture, and then sampling proportional to the accumulated computed probabilities?

  • $\begingroup$ r solution Generating random variables from a mixture of Normal distributions $\endgroup$
    – user230053
    Oct 12, 2019 at 19:48
  • $\begingroup$ How can I apply the EM algorithm to the simulated data (a mixture of two Gamma distributions) in order to estimate the final parameters like p, beta, etc.? alphas are assumed to be fixed and known. $\endgroup$
    – Mostafa
    Mar 16, 2021 at 22:08

1 Answer 1


Actually, it is much simpler. Assuming that you have $K$-component mixture of gamma distributions the algorithm to draw sample of size $N$ is as follows:

Repeat $N$ times:
1. draw $k$ from categorical distribution parametrized by vector $\pi$,
2. draw single value from $k$-th gamma distribution parametrized by $\alpha_k$, $\beta_k$.

You can use similar algorithm to draw samples from mixture of any distributions. Notice that this follows exactly from the definition of mixture distribution:

mixture distribution is the probability distribution of a random variable that is derived from a collection of other random variables as follows: first, a random variable is selected by chance from the collection according to given probabilities of selection, and then the value of the selected random variable is realized.

If you know R, this translates to the following example:

# density
dmixgamma <- function(x, pi, alpha, beta) {
  k <- length(pi)
  n <- length(x)
  rowSums(vapply(1:k, function(i) pi[i] * dgamma(x, alpha[i], beta[i]), numeric(n)))

# random generation
rmixgamma <- function(n, pi, alpha, beta) {
  k <- sample.int(length(pi), n, replace = TRUE, prob = pi)
  rgamma(n, alpha[k], beta[k])


pi <- c(4/10, 6/10)
alpha <- c(20, 15)
beta <- c(10, 25)

hist(rmixgamma(1e5, pi, alpha, beta), 100, freq = FALSE)
xx <- seq(0, 10, by = 0.001)
lines(xx, dmixgamma(xx, pi, alpha, beta), col = "red")

enter image description here

  • $\begingroup$ Ahhh why didn't I thought about it this way. Sweet! $\endgroup$ Aug 2, 2016 at 10:07
  • $\begingroup$ A secondary question: is there stat package for c++ that you would recommend? $\endgroup$ Aug 2, 2016 at 10:08
  • $\begingroup$ @user3639557 sorry but I'm using C++ only accompanied with R via Rcpp, but this is probably not what you are asking ;) $\endgroup$
    – Tim
    Aug 2, 2016 at 10:11
  • 2
    $\begingroup$ @user3639557 but if you wanted Rcpp example of sampling from mixture, then here you are: github.com/twolodzko/extraDistr/blob/master/src/… (if you accidentally spot any bugs, posting comments and rising issues are welcomed) $\endgroup$
    – Tim
    Aug 2, 2016 at 10:15

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.