The inversion method

Inversion is a general method for simulating random variables. It makes use of the fact that the transformation $X=F^{-1}(U)$ gives a random variable $X$ with distribution function $F$ provided the inverse function $F^{-1}$ exists. This is a simple consequence of the change of variables formula, this time with $g(U)=F^{-1}(U)$. Since $g^{-1}(x) = F(x)$, the density of $X$ becomes $\frac{d}{dx}F(x)= f(x),$ which is the probability density corresponding to the distribution function $F$. The method is illustrated in Figure 1.2.

Figure 1.2: Generation of a simulated value $x$ from the distribution function $F$ by the method of inversion

Example: Note that applying this procedure in the exponential case, for which $F(x)=1-e^{-x}, x>0$, gives $F^{-1}(U)=-\log(1-U)$ rather than $-\log(U)$ as given previously. However there is no real inconsistency here since $U$ and $1-U$ have the same distribution.

Example: Another simple example is the Weibull distribution which has distribution function

$$F(x; \alpha, \beta)= 1-\exp(-(x/\beta)^{\alpha}), \, x > 0,\alpha > 0, \beta > 0.$$ (1.5)

This is simulated by $x= \beta [-\log(1-u)]^{1/ \alpha}$ or equivalently $\beta [-\log(u)]^{1/ \alpha}$. The Weibull distribution is used extensively in reliability theory as a model of the distribution of breaking strengths in material testing.

Example: Similarly, the logistic distribution with distribution function

$$F(x;\alpha,\beta)=\frac{1}{1+e^{\alpha-\beta x}},$$ (1.6)

is simulated by $x=[\alpha - \log((1-u)/u)]/\beta$.

Example: Finally the Pareto distribution with probability density function

$$f(x;\lambda) = \lambda x^{-(\lambda +1)}\;\;$$   for$$\;\; x>1 \;\;$$   and$$\;\; \lambda > 2,$$

has distribution function

$$F(x,\lambda) = \int_1^x f(y;\lambda) dy = 1- x^{-\lambda},$$

so that

$$x= F^{-1}(u) = (1-u)^{-1/\lambda}.$$

To use the inversion method, the inverse function $F^{-1}$ either has to be available explicitly, as in the exponential, Weibull, logistic and Pareto cases, or has to be computable in a reasonable amount of time. The equality $x=F^{-1}(u)$ is equivalent to $u=F(x)$, so that finding $x$ for given $u$ is equivalent to finding a root of the equation $F(x)-u=0$. When $F$ is strictly monotone, there is only one root and standard numerical root-finding algorithms can be used, provided of course that $F(x)$ itself is easy to evaluate. If it is required to sample repeatedly from the same distribution, it may be worthwhile devoting some time to the development of an accurate approximation to $F^{-1}$ beforehand. For the standard normal distribution $N(0,1)$, a numerical approximation to $F^{-1}$ is given in Abramovitz and Stegun [1]. However inversion is usually not the method of choice for the normal distribution since there are other simple methods such as Algorithm 1.1.

The main advantage of the inversion method is that it is generally easy to verify that a computer algorithm which uses it is written correctly. In this sense the method is efficient, i.e. in saving the time of the programmer. However, there are usually competing methods which will run faster at the expense of mathematical and algorithmic complexity. We proceed to discuss some of these.