Bookshelf icon indicating copy to clipboard operation
Bookshelf copied to clipboard

Add normal distribution in the random module

Open aksiome opened this issue 1 year ago • 2 comments

aksiome avatar Dec 01 '24 11:12 aksiome

I think this issue should be splitted in two issues:

  • Fast normal distrib (described below)
  • Random number in any distribution (a bit more complexe... and costly, probably not what you want for the next release)

A fast method to get a random number that follows a normal distribution is to cumulate several random numbers that follows a uniform distribution.

def random_normal(mu, sigma):
    N = 10 # Number of uniform r.v. to cumulate
    r = 0
    for _ in range(N):
        r += np.random.random()*2-1 # Cumulate uniform r.v. in [-1,1]
    return mu + sigma * r / np.sqrt(N/3)

With N> 10 or above, the obtained distributions are very satisfying: Image

VForiel avatar May 03 '25 17:05 VForiel

To generate a random number in any distribution, I suggest to use a MCMC method such as the well known Metropolis algorithm.

The idea consist in placing a walker randomly in the desired area in wich we want to generate a random number. We then generate randomly a new position for the walker in this area. If the distribution is lower at the new position than at the old one, we accept the move only if a random number between 0 and 1 is lower than the ratio of the distribution at these two positions.

def metropolis(dist, a, b):
    x = np.random.random() * (b - a) + a
    N = 10
    for _ in range(N):
        x2 = np.random.random() * (b - a) + a
        if np.random.random() < dist(x2) / dist(x):
            x = x2
    return x

Image

VForiel avatar May 03 '25 17:05 VForiel