Simulace normálního rozdělení

import numpy as np
import matplotlib.pyplot as plt

# pdf normálního rozdělení N(mi,sigma) pro kontrolu
def Gaussian(x,mu,sigma):
    return 1/(np.sqrt(2*np.pi)*sigma)*np.exp(-(x-mu)**2/(2*sigma))

výběr a inicializace generátoru U(0,1)

rng=np.random.default_rng() #defaultní generátor PCG64
print(rng)

Generator(PCG64)

pomocí centrální limitní věty

$r_i \in U(0,1)$
$x=\left (\sum_{i=1}^{12} r_i\right )-6$
$x \in N(0,1)$

N=100000 #počet vygenerovaných náhodných proměnných z N(0,1)
Nsum=12 #počet náhodných proměnných z U(0,1), které budeme sčítat
x=np.zeros(N)
for j in range(N):
    for i in range(Nsum):
        x[j]+=rng.random()    
    x[j]=x[j]-6    

fig,ax=plt.subplots()
ax.hist(x,bins=100,density=True)
xp=np.linspace(-5,5,1000)
ax.set_xlim(-5,5)
ax.set_xlabel("x")
ax.set_ylabel("f(x)")
ax.set_title("simulace pomocí centrální limitní věty")
ax.plot(xp,Gaussian(xp,0,1),c='red')
plt.show()

Pomocí transformace rovnoměrného rozdělení

$r_1,r_2 \in U(0,1)$
$x_1=\sqrt{-2\text{ln}r_1}\cos(2\pi r_2)$
$x_2=\sqrt{-2\text{ln}r_1}\sin(2\pi r_2)$
$x_1,x_2 \in N(0,1)$

x=np.zeros(N)
for i in range(0,N,2):
    r1=rng.random()
    r2=rng.random()
    x[i]=np.sqrt(-2*np.log(r1))*np.cos(2*np.pi*r2)
    x[i+1]=np.sqrt(-2*np.log(r1))*np.sin(2*np.pi*r2)

fig,ax=plt.subplots()
ax.hist(x,bins=100,density=True)
xp=np.linspace(-5,5,1000)
ax.set_xlim(-5,5)
ax.set_xlabel("x")
ax.set_ylabel("f(x)")
ax.set_title("Simulace pomocí transformace rovnoměrného rozdělení")
ax.plot(xp,Gaussian(xp,0,1),c='red')
plt.show()

pomocí funkce numpy.random.normal()

x=rng.normal(0,1,N)

fig,ax=plt.subplots()
ax.hist(x,bins=100,density=True)
xp=np.linspace(-5,5,1000)
ax.set_xlim(-5,5)
ax.set_xlabel("x")
ax.set_ylabel("f(x)")
ax.set_title("Simulace pomocí funkce numpy.random.normal()")
ax.plot(xp,Gaussian(xp,0,1),c='red')
plt.show()