python模拟常见分布及做相关假设检验

1. 常见分布

按随机变量类型为离散还是连续（取值个数为有限还是无限），可以将事件分布分为离散型分布与连续型分布。

常见离散型分布有：

• 二项分布（伯努利分布）
• 泊松分布
• 几何分布
• 负二项分布
• 超几何分布

常见连续分布：

• 均匀分布
• 指数分布
• 正态分布

1.1 生成特定分布的随机变量

python生成随机分布主要使用Numpy和Scipy的stats，代码如下：

• 离散随机变量

import numpy as np
s=np.random

#二项分布，参数为试验次数与每次成功概率
s.binomial(10,0.6,10)
array([8, 4, 4, 7, 4, 5, 5, 7, 7, 6])

#二项分布，参数为试验次数与每次成功概率
s.binomial(10,0.6,10)
array([8, 4, 4, 7, 4, 5, 5, 7, 7, 6])

#几何分布，参数为成功概率
s.geometric(0.6,(10,))
array([1, 1, 1, 2, 3, 1, 2, 1, 1, 1])

#返回成功3次，每次成功概率为0.6，前会失败的次数
s.negative_binomial(3,0.6,10)
array([1, 1, 2, 0, 1, 6, 4, 0, 3, 1])

#参数为关注事物的数量，不关注事物的数量，采样的数量，实验次数
#返回为样本中包含关注事物的个数
s.hypergeometric(10,5,2,10)
array([1, 0, 1, 2, 2, 0, 1, 2, 2, 2])

随着λ越来越大，泊松分布越来越对称，越来越接近正态分布。

• 连续随机变量

#正态分布，均值，标准差，试验次数
s.normal(0,1,10)
s.standard_normal(10)
array([-1.3756928 , -0.69395138, -0.05536103, -0.49320461,  0.53539012,
       -0.65562987, -0.60314098,  0.53543226,  0.23511523,  0.1169497 ])

#指数分布,单位时间接待3个客户，那么下一个客户需要的时间为多少，这里试验10次。注意该方法参数为指数分布参数λ的倒数
s.exponential(1/3,10)
array([0.04577549, 0.04436587, 0.0040138 , 0.0359573 , 0.65450092,
       0.58551778, 0.02194068, 0.01608489, 0.01789454, 0.04151969])

1.2 计算统计分布的分布函数

使用概率来描述变量的分布情况。对于离散随机变量，使用概率质量函数（PMF）描述其分布。对于连续型随机变量，使用概率密度函数（PDF）和累计函数，描述其分布。

l=range(10)
#二项分布的概率
stats.binom.pmf(l,10,0.6)

#array([1.04857600e-04, 1.57286400e-03, 1.06168320e-02, 4.24673280e-02,
#       1.11476736e-01, 2.00658125e-01, 2.50822656e-01, 2.14990848e-01,
#       1.20932352e-01, 4.03107840e-02])

#泊松分布概率
stats.poisson.pmf(l,mu=3)

#array([0.04978707, 0.14936121, 0.22404181, 0.22404181, 0.16803136,
#      0.10081881, 0.05040941, 0.02160403, 0.00810151, 0.0027005 ])

二项分布概率质量函数：

泊松分布概率质量函数：

# 计算均匀分布U(0,1)的PDF
x = np.linspace(1,3,10)
#loc起始位置，scale区间长度
stats.uniform.pdf(x,loc=1, scale=2)
#array([0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5])

# 计算正态分布N(0,1)的PDF
x = np.linspace(-3,3,10)
stats.norm.pdf(x,loc=0, scale=1)

#array([0.00443185, 0.02622189, 0.09947714, 0.24197072, 0.37738323,
#       0.37738323, 0.24197072, 0.09947714, 0.02622189, 0.00443185])

# 计算指数分布E(1)的PDF
x = np.linspace(0,10,10)
#scale为指数分布参数λ
stats.expon.pdf(x,loc=0,scale=1)

#array([1.00000000e+00, 3.29192988e-01, 1.08368023e-01, 3.56739933e-02,
#       1.17436285e-02, 3.86592014e-03, 1.27263380e-03, 4.18942123e-04,
#       1.37912809e-04, 4.53999298e-05])

正态分布概率密度函数：

指数函数概率密度函数：

随着λ变大，等待时间变化：

2.假设检验

scipy的stats模块也可以方便地做假设检验

2.1 正态假设检验

Shapiro-Wilk test 做正态假设检验

from scipy.stats import shapiro
data_nonnormal = np.random.exponential(size=100)
data_normal = np.random.normal(size=100)

def normal_judge(data):
    stat, p = shapiro(data)
    if p > 0.05:
        return 'stat={:.3f}, p = {:.3f}, probably gaussian'.format(stat,p)
    else:
        return 'stat={:.3f}, p = {:.3f}, probably not gaussian'.format(stat,p)

# output
normal_judge(data_nonnormal)
'stat=0.827, p = 0.000, probably not gaussian'

2.2 独立性检验（卡方检验）

• 卡方检验
  目的：检验两组类别变量是相关的还是独立的

  H0: 两个样本是独立的

  H1: 两组样本不是独立的

from scipy.stats import chi2_contingency
table = [[10, 20, 30],[6,  9,  17]]
#stat统计量，p显著度，dof自由度，expected期望边界
stat, p, dof, expected = chi2_contingency(table)
print('stat=%.3f, p=%.3f' % (stat, p))
if p > 0.05:
	print('Probably independent')
else:
	print('Probably dependent')
	
stat=0.272, p=0.873
Probably independent

expected
array([[10.43478261, 18.91304348, 30.65217391],
       [ 5.56521739, 10.08695652, 16.34782609]])

2.3 均值检验

• T-test
  目的：检验两个独立样本集的均值是否具有显著差异(进而判断分布是否相同，因为这里默认方差相同)

  H0: 均值是相等的

  H1: 均值是不等的

from scipy.stats import ttest_ind
data1 = np.random.normal(size=10)
data2 = np.random.normal(size=10)
#这里默认总体方差相同，判断两者均值是否相同
stat, p = ttest_ind(data1, data2)
print('stat=%.3f, p=%.3f' % (stat, p))
if p > 0.05:
	print('Probably the same mean value')
else:
	print('Probably different mean value')

stat=-0.761, p=0.457
Probably the same mean value

• ANOVA
  目的：与t-test类似，ANOVA可以检验两组及以上独立样本集的均值是否具有显著差异

  H0: 均值是相等的

  H1: 均值是不等的

from scipy.stats import f_oneway
data1 = np.random.normal(size=10)
data2 = np.random.normal(size=10)
data3 = np.random.normal(size=10)
stat, p = f_oneway(data1, data2, data3)
print('stat=%.3f, p=%.3f' % (stat, p))
if p > 0.05:
	print('Probably the same distribution')
else:
	print('Probably different distributions')
 
# output
# stat=0.189, p=0.829
# Probably the same distribution

2.4 同分布检验

• Mann-Whitney U Test
  目的：检验两个样本集的分布是否相同（直接判断分布是否一致）

  H0: 两个样本集的分布相同

  H1: 两个样本集的分布不同

from scipy.stats import mannwhitneyu
data1 = [0.873, 2.817, 0.121, -0.945, -0.055, -1.436, 0.360, -1.478, -1.637, -1.869]
data2 = [1.142, -0.432, -0.938, -0.729, -0.846, -0.157, 0.500, 1.183, -1.075, -0.169]
stat, p = mannwhitneyu(data1, data2)
print('stat=%.3f, p=%.3f' % (stat, p))
if p > 0.05:
	print('Probably the same distribution')
else:
	print('Probably different distributions')

# output
# stat=40.000, p=0.236
# Probably the same distribution

参考

本内容均来自DataWhale概率统计组队学习。