From: jpenderg@nyx10.cs.du.edu (jon pendergrass) Newsgroups: sci.math Subject: Generating pseudo-random numbers Date: 13 Jan 1995 17:42:33 -0700 I am looking for a method of generating random numbers with a probability density function f(x), where f(x) is some appropriate function. I've tried the "rejection" method from 'Numerical Recipes' but it requires a lot of work to change to some new p.d.f. Also, the technique must be easily implemented in Fortran or C. Thanks in advance. jsp@crc.com ============================================================================== 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 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 ============================================================================== From: martin@balboa.eng.uci.edu (Martin Dowd) Newsgroups: sci.math Subject: Re: Generating pseudo-random numbers Date: 17 Jan 1995 21:40:58 GMT In article <3fh7ou$ptd@mp.cs.niu.edu>, Dave Rusin wrote: >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. There's a good discussion of this method in Knuth's Art of Computer Programming, Volume 2 I think. He gives a really fast version. He also gives a neat special method, called the power method, for the normal distribution. - Martin Dowd ============================================================================== From: jpenderg@nyx.cs.du.edu (jon pendergrass) Newsgroups: sci.math Subject: Generating pseudo-random numbers Date: 18 Jan 1995 20:14:01 -0700 Thanks to Mssrs. Rusin and Dowd. I'll look up the reference.