Thursday, February 12, 2015

py15. Beta Ditribution in Python

The Beta distribution depends on two parameters, alpha and beta. Here, a and b names are used instead. The probability is nonzero between x = 0 and x =1, and goes to zero near the two boundaries.


Here, x is a numpy array from 0 to 1, with a total number of points N. In Python, we can test for the last element as x[-1]. Both x[0] and x[-1] should be 0.


We use scipy.stats.beta for the beta distribution. If you are using IPython, type beta, dot and press tab to see beta functions. One of the beta functions is pdf which stands for the probability density function.


We also test for two strings. In Python, string tests are ASCII based, with the lower case letters higher than capital case letter. For example, 'A' is 65 (or 41 hex) and 'a' is 97 (or 61 hex).

# ex15.py
from __future__ import print_function, division
import matplotlib.pyplot as plt
import numpy as np
from scipy.stats import beta
N = 10000
x = np.linspace(0, 1, num = N)
a = 2
b = 5
y = beta.pdf(x,a,b)
plt.plot(x, y, 'b')
plt.xlabel("x")
plt.ylabel("Beta")
plt.title(r"$\alpha$ = %s,   $\beta$ = %s" % (a,b))
plt.show()
print('Total Probability:  %s' % (1/N*sum(y)))
print('Error, Total Probability:  %s' % (1-1/N*sum(y)))
m1 = 1/N*sum(x*y)
m2 = a/(a+b)
print('Approx. Mean: ', m1)
print('Mean: ', m2)
error = (m1-m2)/m2
print('Error, Mean: ', error)
print('python > R: ', 'python'>'R')
#Total Probability:  0.999899974997
#Error, Total Probability:  0.000100025002502
#Approx. Mean:  0.285685714286
#Mean:  0.285714285714
#Error, Mean:  -0.000100000000004
#python > R:  True

Output:

No comments:

Post a Comment