Python实现12种概率分布

目录

Beta 分布 2

Beta 分布（连续型） 3

高斯分布/正态分布（连续型） 3

标准正态分布 4

范畴分布 6

多项式分布（离散型） 7

二项分布（离散型） 8

伯努利分布 9

指数分布 10

狄利克雷分布（连续型） 11

Gamma 分布（连续型） 13

学生t-分布 14

Beta 分布

import numpy as np

from matplotlib import pyplot as plt

def gamma_function(n):

cal = 1

for i in range(2, n):

cal *= i

return cal

def beta(x, a, b):

gamma = gamma_function(a + b) / \

(gamma_function(a) * gamma_function(b))

y = gamma * (x ** (a - 1)) * ((1 - x) ** (b - 1))

return x, y, np.mean(y), np.std(y)

for ls in [(1, 3), (5, 1), (2, 2), (2, 5)]:

a, b = ls[0], ls[1]

# x in [0, 1], trial is 1/0.001 = 1000

x = np.arange(0, 1, 0.001, dtype=np.float)

x, y, u, s = beta(x, a=a, b=b)

plt.plot(x, y, label=r'$\mu=%.2f,\ \sigma=%.2f,'

r'\ \alpha=%d,\ \beta=%d$' % (u, s, a, b))

plt.legend()

plt.savefig('graph/beta.png')

β分布图

Beta 分布（连续型）

import numpy as np

from matplotlib import pyplot as plt

def gamma_function(n):

cal = 1

for i in range(2, n):

cal *= i

return cal

def chi_squared(x, k):

c = 1 / (2 ** (k/2)) * gamma_function(k//2)

y = c * (x ** (k/2 - 1)) * np.exp(-x /2)

return x, y, np.mean(y), np.std(y)

for k in [2, 3, 4, 6]:

x = np.arange(0, 10, 0.01, dtype=np.float)

x, y, _, _ = chi_squared(x, k)

plt.plot(x, y, label=r'$k=%d$' % (k))

plt.legend()

plt.savefig('graph/chi-squared.png')

plt.show()

Beta 分布（连续型）图

高斯分布/正态分布（连续型）

import numpy as np

from matplotlib import pyplot as plt

def gaussian(x, n):

u = x.mean()

s = x.std()

# divide [x.min(), x.max()] by n

x = np.linspace(x.min(), x.max(), n)

a = ((x - u) ** 2) / (2 * (s ** 2))

y = 1 / (s * np.sqrt(2 * np.pi)) * np.exp(-a)

return x, y, x.mean(), x.std()

x = np.arange(-100, 100) # define range of x

x, y, u, s = gaussian(x, 10000)

plt.plot(x, y, label=r'$\mu=%.2f,\ \sigma=%.2f$' % (u, s))

plt.legend()

plt.savefig('graph/gaussian.png')

plt.show()

高斯分布（连续型）图

标准正态分布

import numpy as np

from matplotlib import pyplot as plt

def normal(x, n):

u = x.mean()

s = x.std()

# normalization

x = (x - u) / s

# divide [x.min(), x.max()] by n

x = np.linspace(x.min(), x.max(), n)

a = ((x - 0) ** 2) / (2 * (1 ** 2))

y = 1 / (s * np.sqrt(2 * np.pi)) * np.exp(-a)

return x, y, x.mean(), x.std()

x = np.arange(-100, 100) # define range of x

x, y, u, s = normal(x, 10000)

plt.plot(x, y, label=r'$\mu=%.2f,\ \sigma=%.2f$' % (u, s))

plt.legend()

plt.savefig('graph/normal.png')

plt.show()

标准正态分布图

范畴分布

#Multi-Bernoulli 分布又称为范畴分布（Categorical distribution）

import random

import numpy as np

from matplotlib import pyplot as plt

def categorical(p, k):

return p[k]

n_experiment = 100

p = [0.2, 0.1, 0.7]

x = np.arange(n_experiment)

y = []

for _ in range(n_experiment):

pick = categorical(p, k=random.randint(0, len(p) - 1))

y.append(pick)

u, s = np.mean(y), np.std(y)

plt.scatter(x, y, label=r'$\mu=%.2f,\ \sigma=%.2f$' % (u, s))

plt.legend()

plt.savefig('graph/categorical.png')

plt.show()

范畴分布图

多项式分布（离散型）

import numpy as np

from matplotlib import pyplot as plt

import operator as op

from functools import reduce

def factorial(n):

return reduce(op.mul, range(1, n + 1), 1)

def const(n, a, b, c):

"""

return n! / a! b! c!, where a+b+c == n

"""

assert a + b + c == n

numer = factorial(n)

denom = factorial(a) * factorial(b) * factorial(c)

return numer / denom

def multinomial(n):

"""

:param x : list, sum(x) should be `n`

:param n : number of trial

:param p: list, sum(p) should be `1`

"""

# get all a,b,c where a+b+c == n, a<b<c

ls = []

for i in range(1, n + 1):

for j in range(i, n + 1):

for k in range(j, n + 1):

if i + j + k == n:

ls.append([i, j, k])

y = [const(n, l[0], l[1], l[2]) for l in ls]

x = np.arange(len(y))

return x, y, np.mean(y), np.std(y)

for n_experiment in [20, 21, 22]:

x, y, u, s = multinomial(n_experiment)

plt.scatter(x, y, label=r'$trial=%d$' % (n_experiment))

plt.legend()

plt.savefig('graph/multinomial.png')

plt.show()

多项式分布（离散型）图

二项分布（离散型）

import numpy as np

from matplotlib import pyplot as plt

import operator as op

from functools import reduce

def const(n, r):

r = min(r, n-r)

numer = reduce(op.mul, range(n, n-r, -1), 1)

denom = reduce(op.mul, range(1, r+1), 1)

return numer / denom

def binomial(n, p):

q = 1 - p

y = [const(n, k) * (p ** k) * (q ** (n-k)) for k in range(n)]

return y, np.mean(y), np.std(y)

for ls in [(0.5, 20), (0.7, 40), (0.5, 40)]:

p, n_experiment = ls[0], ls[1]

x = np.arange(n_experiment)

y, u, s = binomial(n_experiment, p)

plt.scatter(x, y, label=r'$\mu=%.2f,\ \sigma=%.2f$' % (u, s))

plt.legend()

plt.savefig('graph/binomial.png')

plt.show()

二项分布（离散型）图

伯努利分布

import random

import numpy as np

from matplotlib import pyplot as plt

def bernoulli(p, k):

return p if k else 1 - p

n_experiment = 100

p = 0.6

x = np.arange(n_experiment)

y = []

for _ in range(n_experiment):

pick = bernoulli(p, k=bool(random.getrandbits(1)))

y.append(pick)

u, s = np.mean(y), np.std(y)

plt.scatter(x, y, label=r'$\mu=%.2f,\ \sigma=%.2f$' % (u, s))

plt.legend()

plt.savefig('graph/bernoulli.png')

plt.show()

伯努利分布图

指数分布

import numpy as np

from matplotlib import pyplot as plt

def exponential(x, lamb):

y = lamb * np.exp(-lamb * x)

return x, y, np.mean(y), np.std(y)

for lamb in [0.5, 1, 1.5]:

x = np.arange(0, 20, 0.01, dtype=np.float)

x, y, u, s = exponential(x, lamb=lamb)

plt.plot(x, y, label=r'$mu=%.2f, sigma=%.2f,'

r' lambda=%d$' % (u, s, lamb))

plt.legend()

plt.savefig('graph/exponential.png')

plt.show()

指数分布图

狄利克雷分布（连续型）

from random import randint

import numpy as np

from matplotlib import pyplot as plt

def normalization(x, s):

"""

:return: normalizated list, where sum(x) == s

"""

return [(i * s) / sum(x) for i in x]

def sampling():

return normalization([randint(1, 100),

randint(1, 100), randint(1, 100)], s=1)

def gamma_function(n):

cal = 1

for i in range(2, n):

cal *= i

return cal

def beta_function(alpha):

"""

:param alpha: list, len(alpha) is k

:return:

"""

numerator = 1

for a in alpha:

numerator *= gamma_function(a)

denominator = gamma_function(sum(alpha))

return numerator / denominator

def dirichlet(x, a, n):

"""

:param x: list of [x[1,...,K], x[1,...,K], ...], shape is (n_trial, K)

:param a: list of coefficient, a_i > 0

:param n: number of trial

:return:

"""

c = (1 / beta_function(a))

y = [c * (xn[0] ** (a[0] - 1)) * (xn[1] ** (a[1] - 1))

* (xn[2] ** (a[2] - 1)) for xn in x]

x = np.arange(n)

return x, y, np.mean(y), np.std(y)

n_experiment = 1200

for ls in [(6, 2, 2), (3, 7, 5), (6, 2, 6), (2, 3, 4)]:

alpha = list(ls)

# random samping [x[1,...,K], x[1,...,K], ...], shape is (n_trial, K)

# each sum of row should be one.

x = [sampling() for _ in range(1, n_experiment + 1)]

x, y, u, s = dirichlet(x, alpha, n=n_experiment)

plt.plot(x, y, label=r'$\alpha=(%d,%d,%d)$' % (ls[0], ls[1], ls[2]))

plt.legend()

plt.savefig('graph/dirichlet.png')

plt.show()

狄利克雷分布（连续型）图

Gamma 分布（连续型）

import numpy as np

from matplotlib import pyplot as plt

def gamma_function(n):

cal = 1

for i in range(2, n):

cal *= i

return cal

def gamma(x, a, b):

c = (b ** a) / gamma_function(a)

y = c * (x ** (a - 1)) * np.exp(-b * x)

return x, y, np.mean(y), np.std(y)

for ls in [(1, 1), (2, 1), (3, 1), (2, 2)]:

a, b = ls[0], ls[1]

x = np.arange(0, 20, 0.01, dtype=np.float)

x, y, u, s = gamma(x, a=a, b=b)

plt.plot(x, y, label=r'$\mu=%.2f,\ \sigma=%.2f,'

r'\ \alpha=%d,\ \beta=%d$' % (u, s, a, b))

plt.legend()

plt.savefig('graph/gamma.png')

plt.show()

Gamma 分布（连续型）图

学生t-分布

import numpy as np

from matplotlib import pyplot as plt

def gamma_function(n):

cal = 1

for i in range(2, n):

cal *= i

return cal

def student_t(x, freedom, n):

# divide [x.min(), x.max()] by n

x = np.linspace(x.min(), x.max(), n)

c = gamma_function((freedom + 1) // 2) \

/ np.sqrt(freedom * np.pi) * gamma_function(freedom // 2)

y = c * (1 + x**2 / freedom) ** (-((freedom + 1) / 2))

return x, y, np.mean(y), np.std(y)

for freedom in [1, 2, 5]:

x = np.arange(-10, 10) # define range of x

x, y, _, _ = student_t(x, freedom=freedom, n=10000)

plt.plot(x, y, label=r'$v=%d$' % (freedom))

plt.legend()

plt.savefig('graph/student_t.png')

plt.show()

学生t-分布图