From: rusin@washington.math.niu.edu (Dave Rusin)
Newsgroups: sci.math
Subject: Re: Generating pseudo-random numbers
Date: 17 Jan 1995 20:02:38 GMT

In article <3f76lp$1kb@nyx10.cs.du.edu>,
jon pendergrass <jpenderg@nyx10.cs.du.edu> wrote:
>I am looking for a method of generating random numbers with a probability 
>density function f(x), where f(x) is some appropriate function. 

How well do you know your function? If you can compute its integral
C(x)=\int_0^x f(t)dt (that's the cumulative distribution function) then
you could let  g be the inverse of that function and then generate the
random numbers by taking  x = g(z)  where the inputs  z  are randomly
generated in [0,1] with uniform density. Indeed, given an interval
(a,b),  the probability that the number  x  so generated lies in the
interval is the probability that the chosen  z  lies in the interval
(g^-1(a), g^-1(b) ), i.e., (C(a), C(b)). Since  z  is generated randomly
with a uniform distribution, this probability is just the length
(C(b)-C(a)) of the last interval. For  b  close to  a, this is roughly
C'(a).(b-a), i.e., f(a).(b-a). I assume this is what you mean when you
say "x must be generated randomly with a probability density function  f".

Observe that it is not necessary to be able to write  g  in any closed
form. Once z  is randomly generated, you can compute  x=g(z)  by
solving the equation  z=C(x)  for  x.  This can be done with Newton's
method, for example, which is particularly well suited here since this
method calls for evaluations of  C', which is just the original  f  you
have. Of course, there are the usual conditions to be concerned with but
Newton's method is typically quite fast (and provably so).

dave
