Inverse transform sampling - CDF is not invertible

Question

Suppose the cumulative distribution function $F$ is given but not invertible to use the inverse transform sampling technique (to compute $X=F^{-1}(Y)$). Do we have other alternative methods? I would appreciate to know the name of all possible methods...

Answer 1

The inverse cdf method operates even when the cdf is not invertible, using the generalised inverse$$F^-(u)=\sup\{x;\ F(x)\le u\}\tag{1}$$ which is always defined for $u\in(0,1)$.

When there is no solution in $x$ to the equation $$F(x)=u$$ it means that $F$ has a jump between a value less than $u$, $u-\epsilon$ and a value more than $u$, $u+\eta$. Hence the distribution has a point mass with mass $\eta+\epsilon$ at a point $x_0$, with $$F(x_0^{-})=u-\epsilon\quad\text{and}\quad F(x_0)=u+\eta$$ In that case (1) leads to$$F^-(u)=x_0$$

Similarly, if the equation $$F(x)=u$$ has an infinite number of solutions, say $x\in [a,b)$, it means that the cdf is constant over this interval and hence that all values in (a,b) have zero probability to occur. In that case, $$F^-(u)=b$$ [which is of course a convention since the event $U=F(a)$ has probability zero to occur.

Answer 2

In low dimensions a good alternative is to use rejection sampling from the pdf $f_X$ (in high dimensions this becomes very inefficient).

Say $f_X$ is your pdf for some random variable $X$, which you want to sample from in the interval $I=[x_\mathrm{min}, x_\mathrm{max}]$. Then you can draw samples $x_i$ uniformly from $I$ and accept/reject them with the probability $f_X(x_i)$, i.e. you draw another uniformly distributed random number $u_i\in[0, \max(f_X(x)|x \in I)]$ and if $u_i \lt f_X(x_i)$ you accept that sample point otherwise you reject it.

Answer 3

Methods can be quite different depending on the distribution you want to simulate. There many good books on simulation methods and lots of info about simulation on the Internet. Methods often exploit relationships among distributions (such as Poisson, exponential, gamma, chi-squared, F, beta---see discussion in comments). Ultimately, almost all current computer simulation ultimately uses standard uniform output from a pseudo-random generator.

Sometimes, methods depend on the current state of technology. For example, the normal CDF is one of the many CDFs that cannot be expressed in closed form.

• The first method of simulating normal variates in the 1950's, when computation beyond basic arithmetic was expensive, was to use $Z = \sum_{i=1}^{12} U_i - 6 \stackrel{aprx}{\sim}\mathsf{Norm}(0,1),$ by the CLT, where $U_i \stackrel{iid}{\sim}\mathsf{Unif}(0,1).$
• Subsequently, the Box-Muller method was used to obtain two independent standard normal random variables from to independent standard uniform ones.
• Currently, it is common to use a very accurate rational approximation to the standard normal quantile function to get one standard normal deviate from one standard uniform deviate. I believe the function runif in R uses Michael Wichura's rational approximation, which is accurate up double-precision computer representation.