这篇文章的契机是专业英语课上老师布置的一个关于验证中心极限定理的作业。通过这次作业主要进行了以下几点的学习或者复习，在文中都有提到：

Python类的继承

几种常用的概率的PDF和简写

服从某种概率的随机数生成方法

1. 常见的概率分布

常见的概率分布分为三类，分别是：continuous distribution(连续分布)，discrete distribution(离散分布)，mixed discrete/continuous distribution(混合分布)。这里主要介绍7种常用的概率分布。

1.1 常见离散分布Bernoulli distribution(伯努利分布)

0-1分布的随机变量X的取值只能是0或1，它的分布规律如下表格所示。常记为X~0-1(p)或X~B(1,p)。其概率密度函数如下：

binomial distribution(二项分布)

独立重复n次伯努利实验，若事件的发生的概率为p，那么n次实验中事件发生k次的概率如公式所示。那么我们就说随机变量X服从二项分布，记为X~B(n,p)。

geometrical distribution(几何分布)

独立重复伯努利实验，直到第k次事件才发生的概率如下公式所示。那么我们称随机变量X服从几何分布，记为X~G(p)。

poisson distribution(泊松分布)

在一个单位内(时间，面积，空间)等等，某个时间发生K次的概率，例如1次，2次等等，这些共同组成了泊松分布。记为X~π(p)。

1.2 常见连续分布normal distribution(正态分布)

正态分布在统计学上十分重要，经常用在自然和社会科学来代表一个不明的随机变量。记做x~N(u,sigma^2)

exponential distribution(指数分布)

指数分布可以用来表示独立随机事件发生的时间间隔，比如旅客进入机场的时间间隔、打进客服中心电话的时间间隔、中文维基百科新条目出现的时间间隔等等。记为x~Exponential(λ)

uniform distribution(均匀分布)

在a-b之间具有相同的发生概率，在其他的地方不存在发生的可能性，称为均匀分布。记作x~U(a,b)

2. Python程序实现

本次全部使用python编写，包括使用python生成符合分布的随机数也由自己编写。

2.1 Python程序结构设计

此次使用python进行编程实现，首先将这次的项目分成两个部分，一个是定义了各个分布类的头文件，一个是用来实际操作的源文件。而由于分布类都有类似的操作，因此首先定义一个描述这些分布的基类，其他的分布类作为子类从这个基类中继承方法。现在设计的基类方法(通用)如下：

产生一个符合该概率分布的值

产生n个符合该概率分布的值

生成多组包含n个符合该概率分布的值的数组

绘制这个概率分布的基本图形

而每一个子类(特定的概率函数)需要分别定义下面两个函数：

初始化函数，为概率的参数赋值

CDF的反函数，因为每种概型的PDF都不同。具体作用会在下面说

本次用到了numpy，marplotlib，random三个库，大家可以自行导入，命名规则就是常用的规则，下面是基类的声明：

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48class :

def __init__(self, name):

self.name = name

# 定义CDF反函数，空着

def cumulative_function(self):

return 1

# 产生一个随机数，子类可直接调用

def generate_one_num(self):

seed = random.random()

output = self.cumulative_function(seed)

return output

# 产生多个随机数，子类可直接调用

def generate_numbers(self, num):

nums = []

while(num):

nums.append(self.generate_one_num())

num -= 1

return nums

# 产生多组随机数，子类可直接调用

def generate_clusters(self, sample_num, group_num):

'''产生多个服从分布的数'''

arr = np.array(self.generate_numbers(sample_num))

for i in range(group_num - 1):

arr = np.row_stack((arr, np.array(self.generate_numbers(sample_num))))

return arr

# 绘图，子类可直接调用

def figure(self, num_list, num_agg, f_style = 'b*'):

bar_xstep = np.array([0] * num_agg)

num_list = np.array(num_list)

step = (max(num_list) - min(num_list)) / num_agg

min_num = min(num_list)

# 统计各个步长内的数据数

x_lable = []

for i in range(num_agg):

temp = len(num_list[(min_num + step * i < num_list) & (num_list <= min_num + (step * (i + 1)))])

bar_xstep[i] = temp

x_lable.append('%.2f'%(min_num + step * i))

plt.xlabel('Number')

plt.ylabel('Probability')

plt.title('The Figure of %s'%self.name)

plt.bar(np.array(x_lable), bar_xstep, align = 'edge')

plt.show()

由于连续的概率分布函数CDF的反函数是一个连续的函数，方便定义，因此选择指数，正态，均匀分布这三种进行编写。

2.2 产生符合指定分布的随机数

产生符合分布的随机数，首先需要我们能够产生较为均匀的随机数作为种子(seed)，种子的值域应该在(0, 1)之间(因为CDF的值域是0-1)。这些种子和不同的概率结合就可以计算出符合某种分布的随机数。那么我们首先如何获得种子呢？在Python的random.random()方法中已经定义了这种随机数产生器，基于梅森螺旋法，随机度较高，且在0-1之间，因此可以直接使用。

其次使我们的核心问题，如何差生符合某种分布的随机数呢？下面将写。

step1 定义概率密度函数的参数

step2 对概率密度函数f(x)积分获得累积分布函数F(x)

step3 求累积分布函数的反函数

step4 将种子输入累积分布函数的反函数，求解值

在这里我们以指数分布为例(简单只写大于0的)，第一步的参数我们最后再赋，以方便理解：

概率密度函数：

$$

f(x) = lambda e^{-lambda x}

$$

积分以后，它的累积分布函数如下：

$$

F(x) = 1-e^{-lambda x}

$$

再求解它的反函数：

$$

F^{-1}(y) = -frac{In(1-y)}{lambda}

$$

再将我们产生的一个随机数作为输入就好了，例如：

$$

F^{-1}(0.5) = -frac{In(1-0.5)}{lambda}

$$

最后lambda由我们自己定义就可以求解出来了。最后在代码中就是这样呈现(请无视函数名，因为这是CDF的反函数)：

1

2

3def cumulative_function(self, x):

P = - np.log(1 - x) / self.lamda

return P

均匀分布也可以使用类似的方法求解，而正态分布需要先理解该如何使用程序进行积分…笔者没有做出来，直接调用了numpy里的函数numpy.random.normal()。

最后三种概型的子类定义如下(都是继承了Base父类，方法通用)：

指数分布：

1

2

3

4

5

6

7

8class C_Exponient(Base):

def __init__(self, lamda):

self.lamda = lamda

self.name = 'Exponiential Distribution'

def cumulative_function(self, x):

P = - np.log(1 - x) / self.lamda

return P

均匀分布：

1

2

3

4

5

6

7

8

9

10

11

12

13

14class C_Uniform(Base):

def __init__(self, a, b):

self.a = a

self.b = b

self.name = 'Uniform Distribution'

def cumulative_function(self, x):

if x < 0:

P = -np.inf

elif 0 <= x <= 1:

P = (self.b - self.a) * x + self.a

else:

P = np.inf

return P

正态分布：

1

2

3

4

5

6

7

8

9class C_Normal(Base):

def __init__(self, mu, sigma):

self.sigma = sigma

self.mu = mu

self.name = 'Normal Distribution'

def cumulative_function(self, x):

P = np.random.normal(self.mu, self.sigma, 1).tolist()[0]

return P

在控制台文件中已经调试成功了！

3. 验证中心极限定理

首先需要解释一下什么是中心极限定理。

the central limit theorem (CLT) establishes that, in some situations, when independent random variables are added, their properly normalized sum tends toward a normal distribution (informally a “bell curve”) even if the original variables themselves are not normally distributed.

中心极限定理说明，在适当的条件下，大量相互独立随机变量的均值经适当标准化后依分布收敛于正态分布。

摘自维基百科。详细可见传送门，还有个比较有用的 概率分布表。

在这里，我们分别做一些简单的对比。分别进行下面三组实验：

产生10组数据，每组1000个样本数据。求每组内的1000个数据的均值和方差。再看这10组均值的分布和方差的分布。

产生100组数据，每组1000个样本数据。求每组内的1000个数据的均值和方差。再看这100组均值的分布和方差的分布。

产生10000组数据，每组1000个样本数据。求每组内的1000个数据的均值和方差。再看这10000组均值的分布和方差的分布。

下面是我的控制台程序，每次只需要改概型的类即可。

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19from statistic import *

SAMPLE_NUM = 1000

GROUP_NUM = 1000

BAR_STEP = 10

def main():

# 一组样本的分布

Dtype = C_Exponient(1.5)

sample_list = Dtype.generate_numbers(SAMPLE_NUM)

Dtype.figure(sample_list, BAR_STEP)

# 多组样本的均值方法分布

sample_cluster = Dtype.generate_clusters(SAMPLE_NUM, GROUP_NUM)

mean, var = np.mean(sample_cluster, axis = 1), np.var(sample_cluster, axis = 1)

Dtype.figure(mean, BAR_STEP)

Dtype.figure(var, BAR_STEP)

if __name__ == '__main__':

main()

最后是运行结果的分析，可以看到三种不同的分布随着组数的增加，最后的结果如下：

这是三种分布产生的随机数：

首先是不同组数的指数分布的均值方差分布对比:

其次是不同组数的正态分布的均值方差分布对比:

其次是不同组数的均匀分布的均值方差分布对比: