Tests of probability distributions

#!/usr/bin/env python
# (C) 2000 Huaiyu Zhu <hzhu@users.sourceforge.net>.  Licence: GPL
# $Id: test_probs.py,v 1.7 2000/08/06 07:21:17 hzhu Exp $
"""
Statistical distributions
"""

from MatPy.Stats.utils import histplot
if __name__ == "__main__":	waittime = None
else:						waittime = 1

i = 0
def test(g, dist, x, __name__):
	"For dist: generate rand, draw pdf, cdf and cdfc on points x"
	global i
	
	print dist.__doc__
	print dist.rand((2,5))
	g.title("The %s distribution" % dist.__class__.__name__)

	n = 400
	y = dist.rand((n,1))
	
	histplot(y, x=30, waittime=1, g=g, scaled=1)
	g.holdon
	curves = [dist.pdf(x), dist.cdf(x), dist.cdfc(x)]
	g.plot(curves, xs=[x]*3, names=["pdf", "cdf", "cdfc"])

	g.hardcopy("probs_%s.ps" %i)
	i = i + 1
	g.holdoff
	wait(waittime)

#------------------------------------------------------------------

from MatPy.gplot import Gplot, wait
g = Gplot()

from MatPy.Stats.distribs import binomial, negbinomial, poisson
from MatPy.Stats.distribs import beta, gamma, normal, chi2, t, F

from MatPy.Matrix import r_range
x = r_range(100)/10.

#------------------------------------------------------------------

print "-"*30, "Testing distributions in [-Inf, Inf]"
for dist in [normal(), t(2)]:
	test(g, dist, x-5, __name__)

print "-"*30, "Testing distributions in [0, Inf]"
for dist in [gamma(2,1), chi2(4), F(5,9)]:
	test(g, dist, x+1e-100, __name__)


"""
binomial,
negbinomial,
poisson,
beta,
gamma,
normal,
chi2,
t,
F,
"""

Result obtained with
>>> from MatPy.tests import test_probs

------------------------------ Testing distributions in [-Inf, Inf]
Normal distribution:
	pdf =  1/(b*sqrt(pi)) * exp(-((x-a)/b)^2 / 2)
	
[-0.946  1.57   1.06   0.253  0.351
 -0.219 -1.43   1.07   0.582  1.73  ]
<MatPy.gplot.Gplot instance at 0x8199e04>
Student t distribution (n>0):
	pdf =  gam((n+1)/2) / (sqrt(n*pi)* gam(n/2)) * (1+x^2/n)^(-(n+1)/2)
	
[-8.9   -0.61   0.397 -0.263 -1    
 -0.851 -0.223  0.0603  1.71  -0.419 ]
<MatPy.gplot.Gplot instance at 0x8199e04>
------------------------------ Testing distributions in [0, Inf]
Gamma distribution (a>0, b>0):
	pdf = b/gamma(a) * (b*x)^(a-1) * exp(-b*x)
	
[ 0.948  1.42   0.311  3.95   2.33 
  1.04   0.762  1.99   2.3    2.09  ]
<MatPy.gplot.Gplot instance at 0x8199e04>
Chi square distribution with n degrees of freedom (n>0):
	pdf_t = 1/(2* gam(n/2)) * (t/2)^(n/2-1) * exp(-t/2)
	= pdfg(x, n/2, 1/2)
	
[ 5.01   5.46   2.91   2.99   4.75 
  1.29   1.22   1.74   13     5.41  ]
<MatPy.gplot.Gplot instance at 0x8199e04>
F distribution (n1>0, n2>0):
	(also known as Snedcor's density or the	variance ratio density)
	pdf = density of (u1/n1)/(u2/n2),
	where u1 and u2 are Chi variates with n1 and n2 degrees of freedom.
	
[ 2.28   0.158  0.474  0.718  1.42 
  0.421  3.45   1.76   0.627  1.52  ]
<MatPy.gplot.Gplot instance at 0x8199e04>