how to generate random numbers that satisfy certain distribution

[Steven D'Aprano]

On Sat, 23 Jan 2010 14:10:10 -0800, Paul Rubin wrote:

> Steven D'Aprano <steve at REMOVE-THIS-cybersource.com.au> writes:
>> The Box-Muller transform is reasonably simple, but you can't be serious
>> that it is simpler than adding twelve random numbers and subtracting
>> six!
> 
> If you want a normal distribution, using the Box-Muller transform is
> simpler because it spares you the complication of figuring out whether
> the 12-trial binomial approximation is close enough to produce reliable
> results for your specific application, which you obviously have to do if
> you are using the approximation for anything serious.

But using Box-Miller gives you the complication of figuring out whether 
you care about being thread safety, which you have to do if you're doing 
anything serious. (See the comments in the random module for the gauss 
method).


> It also involves
> writing less code than that list comprehension, since it is already
> implemented in the random module so you can just call it directly.

By that logic, the Linux kernel is simpler than a function to add one to 
the argument, because the Linux kernel has already been implemented but 
you have to write your own add_one function.

(Except in Forth, which usually comes with a word to add one to the 
number at the top of the stack.)

We can agree that, given that the random module already has two normal 
distributions, there's no real point in using the binomial approximation. 
Everything else is quibbling.

By the way, does anyone know why there is no Poisson random numbers in 
the module? The implementation is quite simple (but not as simple as the 
Linux kernel *wink*):


def poisson(lambda_=1):
    L = math.exp(-lambda_)
    k = 0
    p = 1
    while 1:
        k += 1
        p *= random.random()
        if p <= L:
            break
    return k-1


although this can be improved upon for large values of lambda_.





-- 
Steven