# Simulate the values for 10000 random variables using .rvs.
data= stats.norm.rvs(loc = 0, scale = 1, size =1000) # loc specifies the mean while scale the std

print("Mean: %g " %np.mean(data))
print("SD: %g" %np.std(data, ddof=1))

plt.figure()
plt.hist(data, bins = 20, color = '#F50057', normed = True, label = 'Gaussian(0,1)')
plt.xlim(-5,5)
pylab.legend(loc='upper right')
plt.title("Gaussian Random Variables")
plt.show()