The Inversion Method (Generating Random Sequences IV)

In the first post of this series, I discussed how to generate permutations of sequences using the Fisher-Yates method and I explained (although indirectly) how a linear congruential generator works. In a second post, I explained how to generate 2D points uniformly and randomly distributed a triangle, discussing the method of rejection. In a third post, I’ve discussed how to generate points on a sphere.

All these methods have something in common: they are based on the uniform (pseudo)random generator, and they map uniform numbers onto a shape (or move numbers around, in the first case). What if we need another density function than uniform?

We discussed the rejection method here, but in a somewhat cursory fashion. It’s worth discussing (I may in a future post) but its iterative nature is not always suitable, and it’s finicky to get right. No, in this post I’m going to discuss the inversion method that may be easier to get right when the density function considered is playing nice. Let met explain.

Let f(x) be a density function that gives the value for P(X=x) for a random variable X. It is a function such that 0 \leqslant f(x) \leqslant 1 for all x \in R (for some range R, which may or may not be the reals \mathbb{R}) and such that

\displaystyle\int_R f(x) dx = 1

The cumulative density function, F (that is a notation taken by many author where the lowercase letter is the density and the uppercase letter the cumulative density), which is the function yielding P(X\leqslant x), is given by:

\displaystyle F(x)=\int_{L}^x f(t) dt

where L is the lower bound (or infimum, or -\infty, whatever delimits the range R on the left). The function F(x) is monotone non-decreasing, and it may be invertible, that is, F^{-1}(F(x)) may exist and return x. If it is invertible, then we can use it to generate a random number drawn with density f(x).

The key observation is that 0 \leqslant F(x) \leqslant 1 for all x, and that, therefore, using 0\leqslant u \leqslant 1 will allow us to access the whole range R through F^{-1}(u). Drawing 0\leqslant u \leqslant 1 uniformly will make F^{-1}(u) distributed as f(x)!

OK, maybe an illustration will make that clear(er):

On the graph, we have f and F, which are our density and cumulative density functions, respectively. (let us say that indeed the curve F is the integral of the curve f.) We see that three uniformly spaced numbers on the y axis map, through F^{-1}(u) to very different places on the x axis and that the spaces between the points on the x are spaced non-uniformly. The space between a' and b' is such that F(b')-F(a')=b-a!

Let us now see how it would work with a real density function, say, the exponential distribution (the example in [2]).

We have

f(x)=\lambda e^{-\lambda x},

for a distribution parameter \lambda (the “rate parameter”, related to the expectation). The cumulative distribution function is

F(x)=1- e^{-\lambda x}.

We compute the inverse:

y=1- e^{-\lambda x}

-y+1=e^{-\lambda x}

1-y=e^{-\lambda x}

\ln(1-y)=-\lambda x


We can now feed a function F^{-1} with a uniform generator on [0,1].

* *

In Numerical Recipes [2], we find a rather lengthy, but paradoxically terse, discussion on all kind of generators. In [2], the discussion is minimalist, but outlines a few methods correctly. For some distributions, there are better methods than inverses. Especially when the inverse is analytically troublesome. For example, the ziggurat algorithm, a rejection-based algorithm applicable to monotone-decreasing densities (possibly dealing with symmetries), seems to be a runner-up against more traditional methods.

[1] W. H. Press, S. A. Teukolsky, W. T. Vetterling, B. P. Flanner — Numerical Recipes in C: The Art of Scientific Computing — 3rd ed, Cambridge University Press, 2007

[2] S. M. Ross — Introduction to Probability and Statistics for Engineers and Scientists — 3rd ed, Academic Press, 2004

2 Responses to The Inversion Method (Generating Random Sequences IV)

  1. […] I’ve discussed permutations, generating uniform random points in a triangle, on a sphere, the inversion method, even recycling expensive random bits. In this entry, I will use the rejection method to generate […]

  2. nicole says:

    I really like your plot, but for clarity, the CDF is generally described as the cumulative distribution function rather than density, which is reserved for the PDF.

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: