Generating a normal distn

Stolen from here:

http://telliott99.blogspot.com/2010/11/normal-distribution-construction-in.html

Here is an example with the normal distribution that will seem trivial after the t-distribution (here).

The basic form of the normal is exp {-x2/2}. We define that as a Python function f(x), vectorize it, and construct an array X of discrete points from -10 to +10 with interval dx = 0.001. We apply the vectorized function to the array to get the relative densities. So that we obtain the correct area under the curve, we multiply the height (value of f(x)) of each piece by its width, dx. When we sum up all the pieces, the total is equal to √2π as seen in the printout. We divide by this value to normalize the distribution so that its total area is equal to 1 and it becomes a pdf.

The form of the normal that includes a term for the standard deviation is normal(x). Everything is as before, except we substitute z for x, and at the end we find that the normalizing constant is 1/σ√2π. We plot it to have something pretty to look at.

% matplotlib inline
from __future__ import division
import numpy as np
import pandas as pd
import math
import matplotlib.pyplot as plt

@np.vectorize
def f(x):
return math.e**(-0.5*(x**2))

@np.vectorize
def normal(x,mu=0,sigma=1):
z = (x-mu)/sigma
return math.e**(-0.5*(z**2))

dx = 0.001
X = np.arange(-10,10+dx,dx)
pdf = f(X)
pdf *= dx
print round(sum(pdf),3),
S = math.sqrt(2*math.pi)
print round(S,3),
pdf /= S

2.507 2.507

sigma = 2
pdfn = normal(X,3,sigma)
pdfn *= dx
print round(sum(pdfn),3),
Sn = math.sqrt(2*math.pi) * sigma
print round(Sn,3),
pdfn /= Sn

5.012 5.013

plt.plot(X,pdf,color='r',lw=4)
plt.plot(X,pdfn,color='k',lw=4)
ax = plt.axes()
ax.set_xlim(-6,6)
m = max(pdf)
ymin,ymax = -m/100,m*1.1
ax.set_ylim(ymin,ymax)
plt.text(3,0.75*ymax,s='$\sigma = 1$', color='r',fontsize=24)
plt.text(3,0.65*ymax,s='$\sigma = 2$', color='k',fontsize=24)