I have some values with unknown joint distribution, but I am assuming that the marginal distributions are two-part Normal Mixtures. I am modelling the dependency between the distributions via vine-copulas and pairwise copula constructions.

What I want to do is simulate new values from these distributions, taking into consideration the dependency between them.

What I did:

  • Estimated the parameters of the marginals using R (bayesmix)
  • Used an empirical distribution function to get uniform values from the sample ones, so that I could construct the copula (I could also have used the estimated distributions for the transform)
  • Constructed a vine copula and simulated values from it, meaning I now have uniform values for all variables (VineCopula package)

What I still need to do:

  • Use the generated uniform values that possess the dependency information to get the actual values from the distribution. The problem here comes from the fact that I do not know the inverse CDF of the Normal Mixture and how to do this in R. If it was any of the standard distributions, for which I could calculate the inverse (or the inverse is already implemented), there would be no problem.

So my question is, how can I do this? Is there a way to do this?

I would prefer answers with both theory and R, but will be perfectly satisfied with either.

  • $\begingroup$ I do not know if this will work or not, you can have a try. You already know how to use inverse in R, via qnorm. So how about $w_1$*qnorm()_1 + $w_2$*qnorm()_2 $\endgroup$
    – Alice
    Commented Jul 22, 2020 at 10:53
  • $\begingroup$ I meant just multiply the weight by the inverse of each margins. $\endgroup$
    – Alice
    Commented Jul 22, 2020 at 10:54
  • 2
    $\begingroup$ @Alice: this is not correct, the quantile function of the mixture is not the mixture of the quantile functions. $\endgroup$
    – Xi'an
    Commented Jul 22, 2020 at 11:10

1 Answer 1


The inverse cdf approach means solving (in $x^i$) an equation of the form $(i=1,2)$ $$\omega_1^i \Phi(\{x^i-\mu_1^i\}/\sigma_1^i)+(1-\omega_1^i) \Phi(\{x^i-\mu_2^i\}/\sigma_2^i) = u^i$$ Since this equation has no analytical solution, it need be solved by an numerical resolution. For instance, here is a raw R rendering of this resolution


when $u^i=0.3$, with solution

[1] -1.740754
[1] -5.121608e-06

Note that both components do not need to be simulated this way. More precisely, $X^1$ can be generated from the corresponding Normal mixture, then transformed into $U^1$ by the mixture cdf, then $U^2$ can be generated from the copula conditional distribution, and $X^2$ derived by the above numerical inversion.

  • 2
    $\begingroup$ +1. It's no more difficult to provide a general solution for all Normal mixtures, as in qmix <- Vectorize(function(q, mu, sigma, weights, ...) { f <- function(x) sum(weights * pnorm(x, mu, sigma)) - q; ifelse(0 < q & 1-q > 0, uniroot(f, range(qnorm(q, mu, sigma)), ...)$root, qnorm(q)) }, "q") As an example of its use (and a decent test) you can plot the quantile function, as in curve(qmix(q, mu=c(1,-1), sigma=c(1/0.2,1/0.7), weights=c(0.2,0.8)), xname="q", ylab="x", n=1001). Your example uses q = 0.3. This will work for inverting Gaussian KDEs, too. $\endgroup$
    – whuber
    Commented Jul 22, 2020 at 14:17
  • 1
    $\begingroup$ (+1) Thank you for the generalisation and the comment about KDE inverse cdf generation, although I would not venture that path when generating them! $\endgroup$
    – Xi'an
    Commented Jul 23, 2020 at 17:50

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.