上一节我们讨论的都是随机事件,某一个随机事件可能包含若干个随机试验样本空间中的随机结果,如果对于每一个可能的实验结果都关联一个特定的值,这样就形成了一个随机变量。
例如抛一个骰子,将抛出的骰子的值作为随机变量的值;足球比赛,将某一只球队进球的个数作为随机变量的值;抛一根标枪,抛出的距离作为随机变量的值;今年一年的降水量作为随机变量等等。
随机变量的取值并不是连续的,而是有限个数值,或者是可以计数的无限个数值,这样的随机变量被称为离散随机变量。回顾一下上面提出的四个例子,投一个骰子,将抛出的骰子的值作为随机变量的值,这些值只可能是 1,2,3,4,5,6,所以是离散随机变量分布;足球比赛,将某一只球队进球的个数作为随机变量的值,这个进球数可能是无限个,只不过数值很大的时候概率很低而已,这也是离散随机变量分布。但是对于抛一根标枪,抛出的距离作为随机变量的值和今年一年的降水量作为随机变量这些是无法计数的,被称为连续随机变量。
对于离散随机变量,搞清楚每个值的概率也是很重要的,将随机变量的每个值映射到其概率上,这就是概率质量函数(PMF).记随机变量为 X,PMF 为 P(X=k) 。接下来,我们就来详细说明几种典型的随机变量以及其概率质量函数。
我们有如下几个例子:
以上这些例子中,都可以理解为在同样的条件下重复地、相互独立地进行的一种随机试验,其特点是该随机试验只有两种可能:发生或者不发生。我们假设该项试验独立重复地进行了 n 次,那么就称这一系列重复独立的随机试验为 n 重伯努利试验。n 重伯努利试验结果的分布就是二项分布
二项分布的 PMF为: P(X = k) = C_n^kp^k(1-p)^{n-k}
根据 PMF 推导期望与方差,假设伯努利实验的随机变量只有两个值 0(不发生),1(发生),那么 n 次试验的期望为: E[X] = 1C_n^1p(1-p)^{n-1} + 2C_n^2p^2(1-p)^{n-2} + ... + (n-1)C_n^{n-1}p^{n-1}(1-p) + nC_n^np^{n}
根据 kC_n^k = nC_{n-1}^{k-1} ,则有
E[X] = nC_{n-1}^0p(1-p)^{n-1} + nC_{n-1}^1p^2(1-p)^{n-2} + ... + nC_{n-1}^{n-2}p^{n-1}(1-p) + nC_n^np^{n}
= np(C_{n-1}^0(1-p)^{n-1} + C_{n-1}^1p^1(1-p)^{n-2} + ... + C_{n-1}^{n-2}p^{n-2}(1-p) + C_n^np^{n-1}) = np\sum_{k=0}^{n-1}C_{n-1}^{n-1-k}p^k(1-p)^{n-1-k}
根据 (p+q)^n = C_n0pnq^0 + C_n1p{n-1}q^{1} +… + C_n{n-1}p{1}q^{n-1} + C_nnp0q^n = \sum_{k=0}{n}C_{n}{n-k}p{n-k}qk ,则有: E[X]=np(1−p+p)n−1=np
n 次试验的方差为:
D[X] = E(X - E[X])^2 = E[X^2 - 2xE[X] + E[X]^2] = E[X^2] - 2E[X]E[X] + E[X]^2 = E[X^2] - E[X]^2 其中 E[X^2] 为 E[X^2] = 1^2C_n^1p(1-p)^{n-1} + 2^2C_n^2p^2(1-p)^{n-2} + ... + (n-1)^2C_n^{n-1}p^{n-1}(1-p) + n^2C_n^np^{n} = \sum_{k=0}^{n}k^2C_n^kp^k(1-p)^{n-k} 根据 k^2C_n^k = knC_{n-1}^{k-1} = (k - 1 + 1)nC_{n-1}^{k-1} = n(k - 1)C_{n-1}^{k-1} + nC_{n-1}^{k-1} = n(n-1)C_{n-2}^{k-2} + nC_{n-1}^{k-1} ,则有:
E[X^2] = n(n-1)\sum_{k=0}^{n}C_{n-2}^{k-2}p^k(1-p)^{n-k} + n\sum_{k=0}^{n}C_{n-1}^{k-1}p^k(1-p)^{n-k}
设 k+j=n,则有: E[X^2] = n(n-1)\sum_{j=0}^{n}C_{n-2}^{n-j-2}p^{n-j}(1-p)^{j} + n\sum_{j=0}^{n}C_{n-1}^{n-j-1}p^{n-j}(1-p)^{j} =n(n-1)p^2\sum_{j=0}^{n}C_{n-2}^{n-j-2}p^{n-j-2}(1-p)^{j} + np\sum_{j=0}^{n}C_{n-1}^{n-j-1}p^{n-j-1}(1-p)^{j} =n(n-1)p^2(p+1-p)^{n-2} + np(p+1-p)^{n-1} = n(n-1)p^2 + np 则: E[X^2] - E[X]^2 = n(n-1)p^2 + np - n^2p^2 = n^2p^2 - np^2 + np - n^2p^2 = np - np^2 = np(1-p)
接下来我们通过一个小程序来直观感受下二项分布:
from scipy.stats import binom
import matplotlib.pyplot as plt
import numpy as np
import seaborn
#我个人比较喜欢这个主题
seaborn.set_style("whitegrid")
#使用内置库 binom 模拟 n = 10, p = 0.3 的情况,随机 100000 次
binom_sim = data = binom.rvs(n=10, p=0.3, size=100000)
#计算本次平均值以及标准差
print ("Mean: %g" % np.mean(binom_sim))
print ("Std: %g" % np.std(binom_sim, ddof=1))
#绘图,利用归一化的直方图,这样展示的很接近概率质量函数
plt.hist(binom_sim, bins=10, density=True, stacked=True)
plt.show()
可以看到输出为:
可以看到平均值很接近期望,标准差也是很接近我们用公式算出来的方差开根。
当 n 比较大, p 比较小的时候,二项分布可以近似为 泊松分布。
由此可以推导出,泊松分布的期望方差为: E[X] = np = \lambda D[X] = np(1-p) = np * 1 = np = \lambda
泊松分布的 PMF 为: P(X = k) = \lim_{n \to \infty, p \to 0}C_n^kp^k(1-p)^{n-k} = \lim_{n \to \infty, p \to 0} \frac{n(n-1)...(n-k+1)}{k!}p^k(1-p)^{n-k} = \lim_{n \to \infty, p \to 0} \frac{n(n-1)...(n-k+1)}{k!}(\frac{\lambda}{n})^k(1-\frac{\lambda}{n})^{n-k} 分别来看每一部分: \lim_{n \to \infty, p \to 0} n(n-1)...(n-k+1)= n^k 由于 n 趋近于无穷,所以 \lim_{n \to \infty}(1+\frac{x}{n})^n=e^x 所以有: \lim_{n \to \infty, p \to 0}(1-\lambda/n)^n = e^{-\lambda} 最后: \lim_{n \to \infty, p \to 0}(1-\lambda/n)^k = 1^k = 1 合并带入: \lim_{n \to \infty, p \to 0} \frac{n(n-1)...(n-k+1)}{k!}(\frac{\lambda}{n})^k(1-\frac{\lambda}{n})^{n-k} n→∞,p→0limk!n(n−1)...(n−k+1)(nλ)k(1−nλ)n−k = \lim_{n \to \infty, p \to 0}\frac{n^k}{k!}\frac{\lambda^k}{n^k}e^{-\lambda} =n→∞,p→0limk!nknkλke−λ = \frac{\lambda^k}{k!}e^{-\lambda}
为何要近似为泊松分布呢?在生活中我们会根据历史数据来预测结果,同时有很多事件可以抽象为泊松分布,例如:
对于这种,推测概率很难,但是可以通过历史数据描述其期望的,我们一般通过抽象为泊松分布来计算它的先验概率。
接下来我们继续通过程序模拟:
from scipy.stats import poisson
import matplotlib.pyplot as plt
import numpy as np
import seaborn
#我个人比较喜欢这个主题
seaborn.set_style("whitegrid")
#使用内置库 poisson 模拟 λ=5 的情况,随机 10000 次
poisson_sim = poisson.rvs(mu=5, size=10000)
#计算本次平均值以及标准差
print ("Mean: %g" % np.mean(poisson_sim))
print ("Std: %g" % np.std(poisson_sim, ddof=1))
#绘图,利用归一化的直方图,这样展示的很接近概率质量函数
plt.hist(poisson_sim, density=True, stacked=True)
plt.gca().axes.set_xticks(range(0, 11))
plt.show()
输出为:
类比二项分布的例子,我们稍微做下修改:
这些随机变量的分布就是几何分布,其 PMF 为:
P(X=k)=p(1−p)k−1
根据 PMF 推导期望与方差,假设伯努利实验的随机变量只有两个值 0(不发生),1(发生),那么期望为:
E[X] = p + 2p(1-p) + 3p(1-p)^2 + ... + np(1-p)^{n-1}
(1-p)E[X] = p(1-p) + 2p(1-p)^2 + 3p(1-p)^3 + ... + np(1-p)^{n}
E[X]-(1-p)E[X] = p + p(1-p) + p(1-p)^2 + ... + p(1-p)^{n-1} - np(1-p)^{n} E[X
pE[X] = p\frac{1-(1-p)^n}{1-(1-p)} - np(1-p)^n
E[X] = \frac{1-(1-p)^n}{1-(1-p)} - n(1-p)^n
当 n 趋近于无穷的时候:
E[X] = \frac{1}{p}
方差为:
D[X] = E[X^2] - E[X]^2
E[X^2] = p + 2^2p(1-p) + 3^2p(1-p)^2 + ... + n^2p(1-p)^{n-1}
\frac{E[X^2]}{p} = 1 + 2^2(1-p) + 3^2(1-p)^2 + ... + n^2(1-p)^{n-1}
假设 q = 1-p
,则有:
\frac{E[X^2]}{1-q} = 1 + 2^2q + 3^2q^2 + ... + n^2q^{n-1}
根据求导公式,则有:
\frac{E[X^2]}{1-q} = (q + 2q^2 + 3q^3 + ... + nq^{n})'
根据推导期望的结果,可以知道:
q + 2q^2 + 3q^3 + ... + nq^{n} = \frac{E(X)(q)}{1-q} = \frac{q}{(1-q)^2}
根据求导公式 [f(x)/g(x)]'=[f'(x)g(x)-f(x)g'(x)]/[g(x)]^2 则有:
\frac{E[X^2]}{1-q} = (\frac{q}{(1-q)^2})' = \frac{(1-q)^2 - ((1-q)^2)'q}{(1-q)^4} = \frac{1-2q+q^2-(-2q + 2q^2)}{(1-q)^4} = \frac{1-q^2}{(1-q)^4} = \frac{1+q}{(1-q)^3}
则:
E[X^2] = \frac{1+q}{(1-q)^2} = \frac{2-p}{p^2}
则:
D[X] = E[X^2] - E[X]^2 = \frac{2-p}{p^2} - \frac{1}{p^2} = \frac{1-p}{p^2}
接下来我们继续通过程序模拟:
from scipy.stats import geom
import matplotlib.pyplot as plt
import numpy as np
import seaborn
#我个人比较喜欢这个主题
seaborn.set_style("whitegrid")
#使用内置库 geom 模拟 p = 0.5 的情况,随机 100000 次
geom_rv = geom(p=0.5)
geom_sim = geom_rv.rvs(size=100000)
#计算本次平均值以及标准差
print ("Mean: %g" % np.mean(geom_sim))
print ("Std: %g" % np.std(geom_sim, ddof=1))
#绘图,利用归一化的直方图,这样展示的很接近概率质量函数
plt.hist(geom_sim, bins=20, density=True, stacked=True)
#设置要显示的坐标值
plt.gca().axes.set_xticks(range(1,11))
plt.show()
根据公式计算出的期望为 \frac{1}{p}=2 ,与模拟的接近。根据公式计算出方差为 \frac{1-p}{p^2}=2 ,标准差约为 1.41421 也和模拟的接近。