CODEWITHSUNDEEP
pythonnumpyrandomscipyscientific-computing
aph
aph

Reputation: 1855

drawing gaussian random variables using scipy erfinv

I would like to draw npts random variables distributed as a gaussian with mean mu and dispersion sigma. I know how to do this in Numpy:

x = np.random.normal(loc=mu, scale=sigma, size=npts)
print(np.std(x), np.mean(x))
0.1998, 0.3997

This should also be possible to do using scipy.special.erfinv via inverse transforms, beginning from a uniform distribution:

u = np.random.uniform(0, 1, npts)

However, I can't figure out how to get the scaling correct. Has anyone done this before?

Upvotes: 0

Views: 719

Answers (2)

user2285236
user2285236

Reputation:

Try this:

mean = 100
sigma = 7
x = mean + 2**0.5 * sigma * erfinv(np.random.uniform(size=10**5) * 2 - 1)
x.mean(), x.std()
Out: (99.965915366042381, 7.0062395839075107)

enter image description here

The conversion from erf to normal distribution is from John D. Cook's blog.

Upvotes: 6

sascha
sascha

Reputation: 33522

The best approach would probably be an own implementation of highly-optimized gaussian-samplers like:

But here is something simple (which i saw a long time ago here). This will be less efficient than the approaches above (as it's using the quantile function / percent point function which has no closed-form representation for the normal-distribution and will be approximated):

import numpy as np
from scipy.stats.distributions import norm

uniform_samples = np.random.rand(100000)
gaussian_samples = norm(loc=0, scale=1).ppf(uniform_samples)

Upvotes: 2

Related Questions