# Efficient sampling of a mixture model

I know there is an easy way to sample standard mixture models if the mixture components can be easily sampled. But, I have a slightly different scenario and I'm wondering if there is an analogous method (or any method at all that is nice). I'll first explain the standard scenario so there is no confusion, and then explain how my problem is different. Sorry if some of this is too remedial, I just want my question to be completely clear.

A standard mixture model has pdf

$$f_X(x) = \sum_i w_i f^i_X(x)$$

where the $w_i$'s are non-negative and sum to one and the $f^i$'s are pdfs. A pretty important mixture model is mixture of Gaussians, where the $f^i$'s are normal distributions with different means and variances, i.e.,

$$f_X(x) = \sum_i w_i \mathcal{N}(\mu_i, \sigma^2_i; x)$$

I can easily sample this distribution with two samples. I'll assume there is a latent variable that governs which of the components generates the value of the random variable. Thus, I can choose $i$ with probability $w_i$ and then sample $\mathcal{N}(\mu_i, \sigma^2_i; \cdot)$ to generate my sample $\bar{x}$.

In my case, I have a slightly different mixture model

$$f_X(x) = \sum_i w_i f^i_X(x) - \sum_j w_j f^j_X(x)$$

Basically, some of the weights are negative so the components $w_j f^j$ subtract from the density at any $x$. Assume that $f_X$ is a pdf, so it is non-negative and integrates to one. Of course, there are methods for generating samples from arbitrary distributions, but it seems that computing the cdf for this one may be hard.

What is the best way to draw a sequence of samples from this distribution?

• In the standard model, $\sum_i w_i = 1$, and $\int f_i(x) dx =1 \forall i$. What is the corresponding normalization in your model? – highBandWidth Feb 20 '11 at 5:28
• $\int f_i(x) dx = 1 ~\forall i$, but weights don't necessarily sum to one. Basically the weights are what they need to be so that f(x) is a pdf with respect to the $f_i$'s and the signs of the weighting coefficients. – RandomGuy Feb 20 '11 at 7:32
• Integrating both sides, I think you need $\sum_i w_i - \sum_j w_j = 1$. But this doesn't guarantee that total $f_X(x)>0 \forall x$. We'll just assumethat this is true. – highBandWidth Feb 20 '11 at 11:58
• Yeah, you're right. I actually have this constraint written down in my notes. No useful constraints on the weights guarantee $f_X(x) > 0$ for all $x$, but its not too troublesome due to the way the distributions arise (which is not really relevant to the sampling). – RandomGuy Feb 20 '11 at 18:33
• One can always assumes $f_X(x)$ by using the $f_X=f_X^+-f_X^-$ decomposition otherwise and adding terms to both sides of the mixture. – Xi'an Feb 3 '18 at 13:52

Since $$\int f_X^i(x) dx = 1 \forall i$$ and $$\int f_X^j(x) dx = 1 \forall j,$$ it is also required that $$\sum_i w_i - \sum_j w_j = 1.$$ Let us call $\sum_i w_i = a$, so we know that $\sum w_j = a-1$.

Here is the algorithm:

1. Pick $i$ based on weights $w_i$.
2. Sample y from $f_X^i(x)$.
3. Calculate $$p_a = \frac{\sum_i w_i f^i_X(x)- \sum_j w_j f^j_X(x)}{\sum_i w_i f_X^i(x)}$$
4. Sample $u$ from $Uniform(0,1)$.
5. If $u < p_a$, then emit y, else reject this sample and go back to 1 to try another.

This is called acceptance-rejection sampling. For more information see the wikipedia article.

• This doesn't work. Try a concrete example (say, sample $2f-g$ where $f$ is $\pm 1$ with equal probability and $g$ is always $1$). – mjqxxxx Feb 20 '11 at 13:34
• I'm not sure I understand the intuition. Perhaps you could show how the pdf of the samples generated is (or approximates) $f_X$. – RandomGuy Feb 20 '11 at 18:36
• @mjqxxxx, We need to have $2f-g \ge 0$. Your example doesn't satisfy that. – highBandWidth Feb 21 '11 at 2:12
• @RandomGuy, I have corrected my answer, hopefully the rejection sampling article should be clear. – highBandWidth Feb 21 '11 at 18:56
• @RandomGuy, take a look at this handout: web.as.uky.edu/statistics/users/viele/sta695s02/rejimp/… – highBandWidth Mar 1 '11 at 19:22