10
$\begingroup$

I want to simulate a sample from a mixture normal distribution such that

$$p\times\mathcal{N}(\mu_1,\sigma_1^2) + (1-p)\times\mathcal{N}(\mu_2,\sigma_2^2) $$

is restricted to the interval $[0,1]$ instead of $\mathbb{R}$. This means I want to simulate a truncated mixture of normal distributions.

I know that there are some algorithm to simulate a truncated normal (i.e from this question) and corresponding package in R to do this. But how can I simulate a truncated mixture normal? Is it the same if I simulate two truncated normal of $\mathcal{N}(\mu_1,\sigma_1^2)$ and $\mathcal{N}(\mu_2,\sigma_2^2$) to make a truncated mixture normal?

Improve this question
$\endgroup$
3
  • 3
    $\begingroup$ If it's on the unit interval, why not use betas instead of normals? For $\alpha=\beta>1$, the distribution is symmetric and unimodal and bounded on the unit interval. $\endgroup$
    – Sycorax ♦
    Commented Mar 23, 2015 at 14:45
  • 2
    $\begingroup$ If you don’t need your simulations to be very fast, you can do it using rejection sampling: (1) sample $x$ from the mixture of two normals, (2) if $x$ is not in $[0,1]$, go back to step 1, (3) output $x$. (but user777 is right, do you have a good reason to chose this distribution instead of a mixture of betas?) $\endgroup$
    – Elvis
    Commented Mar 23, 2015 at 14:53
  • 1
    $\begingroup$ @user777 a Truncated Gaussian mixture has a different distribution from a Beta distribution and cannot be swapped out just because you can enforce symmetry and the same support. $\endgroup$
    – mjnichol
    Commented May 11, 2015 at 21:15

1 Answer 1

Reset to default
16
$\begingroup$

Simulation from a truncated normal is easily done if you have access to a proper normal quantile function. For instance, in R, simulating $$ \mathcal{N}_a^b(\mu,\sigma^2)$$where $a$ and $b$ denote the lower and upper bounds can be done by inverting the cdf $$\dfrac{\Phi(\sigma^{-1}\{x-\mu\})-\Phi(\sigma^{-1}\{a-\mu\})}{\Phi(\sigma^{-1}\{b-\mu\})-\Phi(\sigma^{-1}\{a-\mu\})} $$ e.g., in R

x = mu + sigma * qnorm( pnorm(a,mu,sigma) + 
     runif(1)*(pnorm(b,mu,sigma) - pnorm(a,mu,sigma)) )

Otherwise, I developed a truncated normal accept-reject algorithm twenty years ago.

If we consider the truncated mixture problem, with density $$ f(x;\theta) \propto \left\{p\varphi(x;\mu_1,\sigma_1)+(1-p)\varphi(x;\mu_2,\sigma_2)\right\}\mathbb{I}_{[a,b]}(x) $$ it is a mixture of truncated normal distributions but with different weights: $$ f(x;\theta) \propto p\left\{\Phi(\sigma_1^{-1}\{b-\mu_1\})-\Phi(\sigma_1^{-1}\{a-\mu_1\}) \right\}\dfrac{\sigma_1^{-1}\phi(\sigma_1^{-1}\{x-\mu_1\})}{\Phi(\sigma_1^{-1}\{b-\mu_1\})-\Phi(\sigma_1^{-1}\{a-\mu_1\})} \\[15pt] +(1-p)\left\{\Phi(\sigma_2^{-1}\{b-\mu_2\})-\Phi(\sigma_2^{-1}\{a-\mu_2\}) \right\}\dfrac{\sigma_2^{-1}\phi(\sigma_2^{-1}\{x-\mu_2\})}{\Phi(\sigma_2^{-1}\{b-\mu_2\})-\Phi(\sigma_1^{-1}\{a-\mu_2\})} $$ Therefore, to simulate from a truncated normal mixture, it is sufficient to take $$x=\begin{cases} x_1\sim\mathcal{N}_a^b(\mu_1,\sigma_1^2) &\text{with probability }\\ &\qquad p\left\{\Phi(\sigma_1^{-1}\{b-\mu_1\})-\Phi(\sigma_1^{-1}\{a-\mu_1\}) \right\}\big/\mathfrak{s}\\ x_2\sim\mathcal{N}_a^b(\mu_2,\sigma_2^2) &\text{with probability }\\ &\qquad(1-p)\left\{\Phi(\sigma_2^{-1}\{b-\mu_2\})-\Phi(\sigma_2^{-1}\{a-\mu_2\}) \right\}\big/\mathfrak{s} \end{cases} $$ where \begin{align} \mathfrak{s}=&p\left\{\Phi(\sigma_1^{-1}\{b-\mu_1\})-\Phi(\sigma_1^{-1}\{a-\mu_1\}) \right\}+ \\ &(1-p)\left\{\Phi(\sigma_2^{-1}\{b-\mu_2\})-\Phi(\sigma_2^{-1}\{a-\mu_2\}) \right\} \end{align}

Improve this answer
$\endgroup$
11
  • 1
    $\begingroup$ Ah! I think I see the issue. It is because the entire distribution is being truncated, not each distribution separately. If each sub-distribution of the mixture were individually truncated before being added into the mixture then we would be able to simply sample from the distribution according to the relative weights of each sub-distribution, right? $\endgroup$
    – mjnichol
    Commented May 11, 2015 at 22:03
  • 1
    $\begingroup$ @mjnichol It is a mixture but with different weights than $p$ and $1-p$. $\endgroup$
    – Xi'an
    Commented May 12, 2015 at 5:21
  • 2
    $\begingroup$ @mjnichol: in that case, you would have $$p\mathcal{N}_a^b(\mu_1,\sigma_1^2)+(1-p)\mathcal{N}_a^b(\mu_2,\sigma_2^2)$$so yes indeed this would work. $\endgroup$
    – Xi'an
    Commented May 12, 2015 at 17:50
  • 1
    $\begingroup$ Just reflect on the fact that$$p/\sigma_1+(1-p)/\sigma_2\ne 1$$ $\endgroup$
    – Xi'an
    Commented Apr 4 at 9:53
  • 1
    $\begingroup$ The $\sigma_i^{-1}$ are part of the densities, not "multiplying" these densities. $\endgroup$
    – Xi'an
    Commented Apr 5 at 8:01

Your Answer

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