Skip to content

2.4 常见离散分布 ​

2.2 节把随机变量定义为样本点到实数的映射,并引入了累积分布函数 F(x),离散情形用概率质量函数 p(x)=P(X=x),连续情形用概率密度函数 f(x)。2.3 节在此基础上给出了期望 \E[X]、方差 Var(X) 与矩,用少数几个数值概括分布的位置与离散程度。这两节的工具适用于任意分布,读者可以为每个取值单独指定概率,只要非负且总和为 1。这种完全自由的赋值方式在理论上方便,在实际建模中却难以使用,因为真实数据往往来自某种可描述的重复机制。

把生成机制写下来,概率就由少数几个参数决定,分布也随之确定。这样的分布集合称为分布族(distribution family)。分布族的价值在于压缩模型:把无穷多个自由参数缩减为两三个有明确含义的数值,同时保留对数据来源的理解。日常生活中许多量都有清晰的生成机制,例如抛 n 次硬币正面出现的次数、单位时间内到达服务窗口的人数、一小时内网站的访问次数、某人一个月内的网购次数。这些量的生成过程分别是独立重复试验、稀有事件流、泊松计数与长尾计数,对应的分布族是二项分布、泊松分布与本章稍后讨论的其它计数模型。

本节按生成机制把常见的离散分布串成一条线索,而不孤立罗列。先从最基础的赋值方式出发,得到离散均匀分布与分类分布;加入独立重复试验的假设,得到伯努利分布、二项分布与多项分布;改变停止规则,得到几何分布与负二项分布;去掉有放回这一条件,得到超几何分布与泊松-二项分布;控制稀有事件的发生率,得到泊松分布;最后处理长尾与零过多这两类实际数据中常见的偏离,并给出经验分布这一非参数对象。每个分布都给出生成机制、概率质量函数、期望与方差、典型场景与参数含义,并用代码验证。下一节把同一套思路用于连续分布,本节末尾的分布关系图为 2.7 节的极限定理与 2.11 节的参数估计提供索引。

2.4.1 离散均匀分布与分类分布 ​

离散均匀分布 ​

最简单的离散分布是离散均匀分布(discrete uniform distribution)。生成机制是从有限个结果中任取一个,每个结果被选中的可能性相同。当取值集合为 {1,2,…,n} 时,概率质量函数为

p(k)=1n,k=1,2,…,n

取值集合为 {a,a+1,…,b} 时,把 n=b−a+1 代入即可。期望与方差为

\E[X]=a+b2,Var(X)=n2−112=(b−a+1)2−112

掷一颗均匀骰子是 a=1、b=6 的情形,期望为 3.5,方差约为 2.9167。离散均匀分布的参数只有取值范围的端点,含义直接:端点决定了每个结果的高度是 1/n。它在抽奖、随机抽样、随机初始化的索引选择中出现,也是逆变换采样在离散情形下的基础。

分类分布 ​

分类分布(categorical distribution)把取值从等概率推广到任意概率。设随机变量 X 取值于 k 个类别 {1,…,k},概率向量 p=(p1,…,pk) 满足 pi≥0 且 ∑ipi=1,概率质量函数为

p(X=i)=pi

离散均匀分布是 pi=1/k 的特例。分类分布是离散情形里最一般的赋值方式,任何有限取值的离散分布都可以写成它。多分类任务中神经网络最后一层经过 softmax 得到的就是一组分类分布的概率,交叉熵损失度量预测分布与真实类别之间的差异,这一点在 4.3 节讨论交叉熵时会展开。

分类分布的期望与方差取决于类别标签的数值编码。若把类别编码为一个独热向量 y∈{0,1}k,则 \E[yi]=pi,Var(yi)=pi(1−pi),且 Cov(yi,yj)=−pipj(i≠j)。若把类别直接编码为标签序号 1,…,k,则 \E[X]=∑iipi,Var(X)=∑ii2pi−(∑iipi)2。两种编码给出同一分布的不同概括量,报告结果时需要说明所用的编码。

python
import numpy as np
from scipy import stats

# 离散均匀分布:均匀骰子,取值 {1, ..., 6}
du = stats.randint(1, 7)                 # scipy 约定区间为 [low, high),即 {1, ..., 6}
print(f"PMF(k=3)   = {du.pmf(3):.4f}")   # 1/6 ≈ 0.1667
print(f"理论均值   = {du.mean():.4f}")    # (1 + 6)/2 = 3.5
print(f"理论方差   = {du.var():.4f}")     # (6^2 - 1)/12 ≈ 2.9167
print(f"模拟均值   = {du.rvs(size=100000, random_state=0).mean():.4f}")  # ≈ 3.5

# 分类分布:三类,概率 (0.5, 0.3, 0.2)
rng = np.random.default_rng(0)
probs = np.array([0.5, 0.3, 0.2])
samples = rng.choice(3, size=100000, p=probs)
print(f"各类别频率 = {np.bincount(samples, minlength=3) / 100000}")  # ≈ [0.5, 0.3, 0.2]

# 用 rv_discrete 构造任意分类分布并计算期望方差
cat = stats.rv_discrete(values=(np.arange(1, 4), probs))
print(f"分类期望 = {cat.mean():.4f}, 分类方差 = {cat.var():.4f}")     # 1.7  0.61

分类分布的方差在标签序号编码下依赖类别的排列,改变类别顺序会改变数值。分类任务中若需要与顺序无关的概括量,通常报告熵或最大类别概率,而不报告标签序号的方差。

为什么需要分布族

只给定 F 或 p 的一串数值,模型会在每个取值上多出自由度,预测新数据时无从约束。分布族用少数参数限制了可能的形状,使得从数据估计参数成为可能。参数估计的任务是把数据压缩为几个参数(2.11 节),假设检验的任务是判断这些参数是否满足某个取值(2.13 节),两条路线都以分布族为对象。

2.4.2 伯努利分布与二项分布 ​

伯努利分布 ​

伯努利分布(Bernoulli distribution)描述单次二元试验。随机变量 X 取 1 表示成功,取 0 表示失败,成功概率为 p,概率质量函数可以合并写为

p(X=x)=px(1−p)1−x,x∈{0,1}

期望与方差为

\E[X]=p,Var(X)=p(1−p)

方差在 p=0.5 时取最大值 0.25,在 p 接近 0 或 1 时趋于 0。这一形状符合直觉:试验越确定,结果的不确定性越小。伯努利分布是全部二值数据的模型,例如一次抽检是否合格、一次点击广告与否、一次投篮是否命中。它的参数 "p" 解释为长期重复下成功的频率。

二项分布 ​

把 n 次独立的伯努利试验叠加,成功次数服从二项分布(binomial distribution)。若各次试验的成功概率同为 p,则

p(X=k)=(nk)pk(1−p)n−k,k=0,1,…,n

组合数 (nk) 统计的是成功出现的位置数,pk(1−p)n−k 是某一特定排列的概率,两者相乘即得恰好 k 次成功的总概率。二项分布与伯努利分布的关系可以用随机变量之和表述:若 X1,…,Xn 独立且都服从 Bernoulli(p),则 X=∑i=1nXi 服从 Binomial(n,p)。期望与方差可以直接用这一分解推出:

\E[X]=np,Var(X)=np(1−p)

方差的推导用到独立性:独立随机变量之和的方差等于方差之和(2.3 节)。二项分布的两个参数含义清楚,n 是试验次数,p 是单次成功概率。典型场景包括投掷 n 次硬币的正面数、n 次曝光中的点击次数、n 个包裹中的破损个数。

泊松近似 ​

当 n 很大而 p 很小时,二项分布可以用一个只依赖 λ=np 的分布近似。取 λ=np 固定,把 p=λ/n 代入二项分布的概率质量函数:

p(X=k)=(nk)(λn)k(1−λn)n−k=λkk!⋅n(n−1)⋯(n−k+1)nk⋅(1−λn)n(1−λn)−k

当 n→∞ 时,第一项因式趋于 1,第二项因式趋于 e−λ,第三项因式趋于 1,于是

p(X=k)→e−λλkk!

这就是泊松极限定理。它说明试验次数很多、单次概率很小时,成功次数只由期望 np 决定。例如 n=100、p=0.05 时,二项分布的均值是 5,与 Poisson(5) 的均值相同;两者的概率质量函数在整个取值范围内接近。

正态近似与中心极限定理的预告 ​

当 n 很大且 p 不接近 0 或 1 时,二项分布的形状接近正态分布。把 X 标准化为

Z=X−npnp(1−p)

则 Z 的分布在 n 增大时接近标准正态分布 N(0,1)。这一结论是中心极限定理(central limit theorem)在伯努利序列上的特例:独立同分布随机变量之和经标准化后趋向正态分布。2.7 节将完整陈述中心极限定理,此处只需记住二项分布有两种可用近似,稀有事件用泊松近似,一般情形用正态近似。

python
import numpy as np
from scipy import stats

n, p = 100, 0.05
lam = n * p
b = stats.binom(n, p)
po = stats.poisson(lam)

print(f"二项均值 {b.mean():.3f},方差 {b.var():.3f}")     # 5.000  4.750
print(f"泊松均值 {po.mean():.3f},方差 {po.var():.3f}")   # 5.000  5.000

# 泊松近似:比较若干取值的 PMF
for k in range(4):
    print(f"k={k}: 二项 {b.pmf(k):.5f}  泊松 {po.pmf(k):.5f}")
# k=0: 二项 0.00592  泊松 0.00674
# k=1: 二项 0.03116  泊松 0.03369
# k=2: 二项 0.08118  泊松 0.08422
# k=3: 二项 0.13958  泊松 0.14037

# 正态近似:计算 P(X <= 8)
z = (8 - lam) / np.sqrt(n * p * (1 - p))
print(f"二项精确 P(X<=8) = {b.cdf(8):.4f}")           # 0.9369
print(f"正态近似        = {stats.norm.cdf(z):.4f}")   # 0.9487

正态近似在 p 接近 0 时误差较大,需要连续性修正或改用泊松近似。判断的依据是 np 与 n(1−p) 是否都不小于 5,这是常用的经验阈值。

二项分布的两种近似

泊松近似只保留 λ=np 这一个参数,适合稀有事件;正态近似保留 μ=np 与 σ2=np(1−p) 两个参数,适合大规模重复试验。选择哪一种取决于 p 的量级与 n 的大小,其中 p 是否足够小是判断泊松近似能否使用的关键。

2.4.3 多项分布、狄利克雷-多项分布与贝塔-二项分布 ​

多项分布 ​

多项分布(multinomial distribution)是二项分布向多类别的推广。进行 n 次独立试验,每次落入 k 个类别之一,类别概率为 p=(p1,…,pk)。令 Xi 表示落入第 i 类的次数,计数向量 X=(X1,…,Xk) 满足 ∑iXi=n,概率质量函数为

p(x1,…,xk)=n!x1!x2!⋯xk!∏i=1kpixi

多项系数 n!∏ixi! 统计的是把 n 次试验分配到各类别的方案数。边缘分布为 Xi∼Binomial(n,pi),因此 \E[Xi]=npi,Var(Xi)=npi(1−pi)。不同类别的计数存在负相关:

Cov(Xi,Xj)=−npipj,i≠j

负相关的来源是计数总和固定为 n:某一类计数偏高时,其它类的计数总和必然偏低。典型场景包括掷多面骰子各面出现的次数、一批用户按偏好选项分类的人数、文档中各类主题词的出现频次。

贝塔-二项分布与过散布 ​

二项分布假设单次成功概率 p 恒定。若 p 本身在个体之间波动,则把 p 视为随机变量并对其取分布,得到的计数分布会比二项分布更分散。取 p∼Beta(α,β),再取 X∣p∼Binomial(n,p),X 的边缘分布称为贝塔-二项分布(beta-binomial distribution),概率质量函数为

p(X=k)=(nk)B(k+α,n−k+β)B(α,β)

其中 B(⋅,⋅) 是贝塔函数。期望与二项分布相同:

\E[X]=nαα+β,Var(X)=nαβ(α+β+n)(α+β)2(α+β+1)

令 p=α/(α+β),方差可以改写为

Var(X)=np(1−p)⋅α+β+nα+β+1

方括号中的因子大于 1,说明贝塔-二项分布的方差系统性地大于同样均值下的二项分布。这种方差超过名义值的现象称为过散布(overdispersion),它是实际计数数据中普遍存在的偏离。当 α,β→∞ 且比值固定时,该因子趋于 1,贝塔-二项分布退化为二项分布。过散布的存在提醒模型使用者,用二项分布拟合个体间存在异质的计数会低估不确定性。

狄利克雷-多项分布 ​

把贝塔分布换成它的多元推广狄利克雷分布(Dirichlet distribution),取 p∼Dirichlet(α),再取 X∣p∼Multinomial(n,p),得到的边缘分布称为狄利克雷-多项分布(Dirichlet-multinomial distribution)。它的概率质量函数为

p(x)=n!Γ(α0)Γ(n+α0)∏i=1kΓ(xi+αi)xi!Γ(αi),α0=∑i=1kαi

边缘期望为 \E[Xi]=nαi/α0,方差为

Var(Xi)=nαiα0(1−αiα0)n+α01+α0

括号外的因子 >1,说明狄利克雷-多项分布同样允许过散布。这一分布常用于文本的主题模型(如潜在狄利克雷分配)与微生物组组成的建模,因为这些场景中各类别的比例在不同样本之间波动。

python
import numpy as np
from scipy import stats

rng = np.random.default_rng(1)

# 多项分布:n = 50,三类概率 (0.5, 0.3, 0.2)
n, p = 50, np.array([0.5, 0.3, 0.2])
counts = rng.multinomial(n, p, size=100000)
print(f"各类均值 = {counts.mean(axis=0)}")    # ≈ [25. 15. 10.]
print(f"各类方差 = {counts.var(axis=0)}")      # ≈ [12.5 10.5 8.0],与 n p (1-p) 一致

# 贝塔-二项与二项对比:n = 20,二项 p = 2/7,贝塔-二项 a = 2, b = 5
bb = stats.betabinom(20, 2, 5)
bi = stats.binom(20, 2 / 7)
print(f"二项方差      = {bi.var():.4f}")      # ≈ 4.0816
print(f"贝塔-二项方差 = {bb.var():.4f}")       # ≈ 13.7778,因子 (7+20)/(7+1) ≈ 3.375
print(f"两者均值相同  = {bi.mean():.4f}, {bb.mean():.4f}")   # 5.7143  5.7143

方差之比 ≈3.375 与理论因子吻合,验证了过散布参数的解读:α+β 越小,个体间异质越强,方差放大的倍数越大。

2.4.4 几何分布与负二项分布 ​

几何分布与无记忆性 ​

几何分布(geometric distribution)描述在独立重复的伯努利试验中,首次成功所需的试验次数。若单次成功概率为 p,则前 k−1 次失败、第 k 次成功,概率质量函数为

p(X=k)=(1−p)k−1p,k=1,2,…

期望与方差为

\E[X]=1p,Var(X)=1−pp2

几何分布具有无记忆性(memorylessness):对任意 s,t≥0,

P(X>s+t∣X>s)=P(X>t)

推导用到 P(X>m)=(1−p)m 这一简洁的生存函数:P(X>s+t∣X>s)=(1−p)s+t/(1−p)s=(1−p)t=P(X>t)。无记忆性的含义是已经失败 s 次这一信息不会改变未来还需等待多久的分布,未来与过去无关。在离散分布中,只有几何分布具备无记忆性,连续分布中的指数分布是它的对应物(2.5 节)。

负二项分布 ​

负二项分布(negative binomial distribution)描述首次成功推广为第 r 次成功所需的试验次数。若记录试验总次数,概率质量函数为

p(X=k)=(k−1r−1)pr(1−p)k−r,k=r,r+1,…

期望与方差为

\E[X]=rp,Var(X)=r(1−p)p2

负二项分布与几何分布的关系是求和:r 个独立几何随机变量之和服从负二项分布,因此 r=1 时负二项分布退化为几何分布。另一个常见参数化记录的是第 r 次成功前的失败次数,此时取值范围变为 0,1,2,…,期望为 r(1−p)/p。scipy 的 stats.nbinom(n, p) 采用后一种约定,第一个参数 n 是成功次数,第二个参数 p 是单次成功概率。

负二项分布与过散布 ​

按失败次数参数化时(记成功次数为 r,均值为 μ=r(1−p)/p),方差可以写成

Var(X)=μ+μ2r

方差与均值之比为

Var(X)\E[X]=1p=1+μr

第二项使方差大于均值,超出量由参数 r 控制:r 越小,过散布越强,r→∞ 时负二项分布趋于泊松分布。这一性质使负二项分布成为计数数据中处理过散布的常用选择,例如保险理赔次数与网站点击量的建模。

python
import numpy as np
from scipy import stats

p = 0.2
g = stats.geom(p)                       # 试验次数,取值 >= 1
print(f"几何均值 {g.mean():.3f},方差 {g.var():.3f}")   # 5.000  20.000

# 无记忆性验证
print(f"P(X>10)/P(X>5) = {g.sf(10) / g.sf(5):.6f}")      # ≈ 0.327680 = 0.8^5
print(f"P(X>5)         = {g.sf(5):.6f}")                 # 0.327680

# 负二项:scipy 约定为第 n 次成功前的失败次数
nb = stats.nbinom(5, 0.5)
po = stats.poisson(5)
print(f"负二项均值 {nb.mean():.3f},方差 {nb.var():.3f}")   # 5.000  10.000
print(f"同均值泊松方差 = {po.var():.3f}")                    # 5.000,说明负二项过散布
print(f"方差/均值 = {nb.var() / nb.mean():.3f}")             # 2.000 = 1 / p

计数约定容易混淆

几何分布与负二项分布都有两套常见计数方式:一套记录试验总次数,取值范围从 1 或 r 开始;另一套记录失败次数,取值范围从 0 开始。scipy 的 geom 使用试验次数,nbinom 使用失败次数。混用两套约定会让均值相差一个偏移量,报告结果与写代码时应先确认约定。

2.4.5 超几何分布与泊松-二项分布 ​

超几何分布 ​

前面的二项分布假设各次试验相互独立,这在有放回抽样时成立。无放回抽样时,前一次抽到的结果会改变后一次的概率,计数服从超几何分布(hypergeometric distribution)。设总体中有 N 个对象,其中 K 个具有某属性,从中无放回抽取 n 个,抽到的属性对象数 X 的概率质量函数为

p(X=k)=(Kk)(N−Kn−k)(Nn),max(0,n−(N−K))≤k≤min(n,K)

期望与方差为

\E[X]=nKN,Var(X)=nKN(1−KN)N−nN−1

方差的最后一项 N−nN−1 称为有限总体校正因子(finite population correction)。当 n=1 时校正因子为 1,当 n=N 时校正因子为 0(抽取全部对象,结果确定)。校正因子小于 1,说明无放回抽样的计数比有放回抽样更集中。当 N→∞ 且 K/N→p 固定时,校正因子趋于 1,超几何分布趋于二项分布。典型场景包括从一批产品中抽检次品数、从有限库存中抽取某类商品的数量、牌桌上的牌型计算。

泊松-二项分布 ​

泊松-二项分布(Poisson-binomial distribution)描述 n 次相互独立但成功概率不同的伯努利试验之和。设第 i 次试验的成功概率为 pi,则和 X=∑i=1nXi 的期望与方差为

\E[X]=∑i=1npi,Var(X)=∑i=1npi(1−pi)

概率质量函数的显式表达需要展开生成函数 ∏i=1n[(1−pi)+piz] 的系数,通常用递推或数值卷积计算。泊松-二项分布与二项分布的区别在于成功概率是否相同,与泊松分布的区别在于是否要求试验次数很大、单次概率很小。当 ∑ipi 较小且各 pi 都小时,泊松-二项分布可以用均值为 ∑ipi 的泊松分布近似。

python
import numpy as np
from scipy import stats

# 超几何:总体 500,其中 50 个阳性,无放回抽 10 个
hg = stats.hypergeom(M=500, n=50, N=10)     # scipy: M 总体, n 属性数, N 抽取数
bi = stats.binom(10, 0.1)
print(f"超几何均值 {hg.mean():.4f},方差 {hg.var():.4f}")    # 1.0000  0.8838
print(f"二项  均值 {bi.mean():.4f},方差 {bi.var():.4f}")   # 1.0000  0.9000
print(f"有限总体校正因子 = {hg.var() / bi.var():.4f}")        # ≈ 0.9820

# 泊松-二项:10 次独立伯努利试验之和,成功率各不相同
rng = np.random.default_rng(2)
ps = np.array([0.05, 0.08, 0.1, 0.12, 0.15, 0.1, 0.09, 0.11, 0.07, 0.13])
sims = (rng.random((200000, ps.size)) < ps).sum(axis=1)
print(f"模拟均值 = {sims.mean():.4f},理论均值 = {ps.sum():.4f}")            # 1.00  1.00
print(f"模拟方差 = {sims.var():.4f},理论方差 = {np.sum(ps * (1 - ps)):.4f}") # 0.86  0.8645

把超几何分布误当成二项分布会低估方差的风险较小(两者很接近),但在抽样比例 n/N 较大时差异明显,例如从 100 个对象中抽 50 个时校正因子降为 0.5,方差被高估一倍。

2.4.6 泊松分布与泊松过程 ​

泊松分布 ​

泊松分布(Poisson distribution)描述给定区间内稀有事件发生的次数。参数 λ 既是区间内事件数的期望,也是它的强度,概率质量函数为

p(X=k)=e−λλkk!,k=0,1,2,…

期望与方差相等:

\E[X]=λ,Var(X)=λ

这一等式是泊松分布最显著的特征,也是检验计数数据是否适合泊松模型的常用判据:若样本方差远大于样本均值,说明存在过散布,需要考虑负二项或零膨胀等替代模型。泊松分布的典型场景包括单位时间内到达服务台的顾客数、一页书上印刷错误的个数、某路段一天内发生的轻微事故数、一段音频中出现的爆音次数。

泊松极限定理 ​

上一节已在二项分布中给出泊松极限定理的推导,此处补充它的适用范围。定理成立的条件是 n→∞、p→0 且 np→λ 有限。条件中的关键点是 p 足够小,使得单次事件是稀有事件;n 足够大,使得期望 np 保持有限。定理的价值在于把二项分布的两参数问题化简为泊松分布的一参数问题,同时给出一个无需知道 n 与 p 分别取值的建模工具。当 n=1000、p=0.001 时,Binomial(1000,0.001) 与 Poisson(1) 的差异在大多数取值上小于 0.001。

泊松过程 ​

泊松分布描述固定区间上的计数,泊松过程(Poisson process)描述计数在时间上的展开。速率参数为 λ 的泊松过程满足两条性质:任意长度为 t 的区间内事件数服从 Poisson(λt);不重叠区间内的事件数相互独立。由这两条性质可以推出,相邻两次事件之间的时间间隔服从参数为 λ 的指数分布 Exponential(λ)(2.5 节),且各间隔相互独立。泊松过程因此同时给出了计数分布与间隔分布两个视角,前者用离散分布描述,后者用连续分布描述。

python
import numpy as np
from scipy import stats

# 泊松极限:n 增大、p 减小、np 固定为 5
for n in [10, 100, 1000]:
    p = 5 / n
    print(f"n={n:4d}, p={p:.4f}, 二项 PMF(3) = {stats.binom(n, p).pmf(3):.5f}")
# n=  10, p=0.5000, 二项 PMF(3) = 0.11719
# n= 100, p=0.0500, 二项 PMF(3) = 0.13958
# n=1000, p=0.0050, 二项 PMF(3) = 0.14030
print(f"泊松 PMF(3) = {stats.poisson(5).pmf(3):.5f}")   # 0.14037

# 泊松过程:速率 3 的事件流,模拟到达时刻并统计单位区间计数
rng = np.random.default_rng(3)
lamb = 3.0
gaps = rng.exponential(1 / lamb, size=20000)
arrivals = np.cumsum(gaps)
T = arrivals[-1]
counts = np.histogram(arrivals, bins=int(T), range=(0, T))[0]
print(f"单位区间计数均值 = {counts.mean():.4f}")     # ≈ 3.0
print(f"单位区间计数方差 = {counts.var():.4f}")      # ≈ 3.0,与均值接近,符合泊松特征

泊松过程的两个视角

计数视角把时间切成区间,用泊松分布描述每段的计数;间隔视角关注相邻事件的等待时间,用指数分布描述。两者由同一速率参数 λ 联系,从计数分布可以推出间隔分布,反向也成立。这一对偶关系在排队论、可靠性工程与神经放电建模中反复出现。

2.4.7 幂律分布与齐夫分布 ​

幂律分布与长尾 ​

前几节的计数分布(二项、泊松、负二项)的尾部随 k 增大而迅速衰减,属于轻尾分布。幂律分布(power-law distribution)的概率质量函数按幂次衰减:

p(X=k)∝k−α,k≥1,α>0

归一化后写成 p(X=k)=k−α/ζ(α),其中 ζ(α)=∑k=1∞k−α 是黎曼 ζ 函数。α 控制衰减速度:α 越小,尾部越厚,极端大值出现的概率越高。幂律分布的尾部衰减慢,因此称为重尾分布。当 α≤1 时 ζ(α) 发散,分布无法归一化,这样的幂律不作为合法分布使用。当 1<α≤2 时期望有限而方差发散;当 2<α≤3 时方差有限而三阶矩发散。方差不存在的含义是样本方差不随样本量增大而稳定,二值以上的起伏持续存在。

齐夫分布 ​

齐夫分布(Zipf distribution)是有限支撑上幂律的一种形式。设取值 1,2,…,N,概率质量函数为

p(X=k)=k−sHN,s,HN,s=∑i=1Ni−s

当 N→∞ 时齐夫分布趋于幂律分布。齐夫定律最初来自语言学观察:把单词按出现频率从高到低排序,第 k 位单词的频率近似与 k 成反比(s≈1)。类似的排序规律出现在城市人口、网站访问量、书籍销量与学术引用次数中。scipy 的 stats.zipf(a) 采用无限支撑形式,要求 a>1。

重尾与对数坐标 ​

识别幂律的常用方法是在双对数坐标下观察概率质量函数是否近似为直线。由 log⁡p(k)=−αlog⁡k−log⁡ζ(α) 可知,log⁡p(k) 对 log⁡k 的斜率等于 −α。这一线性关系是幂律与指数衰减的判别依据:指数衰减在对数坐标下仍为曲线,幂律在对数坐标下接近直线。需要注意,重尾分布的样本均值与样本方差在小样本下波动很大,基于均值的推断需要谨慎处理。

python
import numpy as np
from scipy import stats

z = stats.zipf(2.0)
print(f"均值 = {z.mean()}")                      # inf,因 a = 2 时期望发散
print(f"P(X=1)={z.pmf(1):.4f}")                  # 1/zeta(2) ≈ 0.6079
print(f"P(X=2)={z.pmf(2):.4f}")                  # 0.1520
print(f"P(X=10)={z.pmf(10):.4f}")                # 0.0061

# 在双对数坐标下检查线性关系
ks = np.arange(1, 500)
pmf = z.pmf(ks)
slope = np.polyfit(np.log(ks), np.log(pmf), 1)[0]
print(f"log-log 斜率 = {slope:.4f}")             # ≈ -2.0,与幂律指数对应

2.4.8 零膨胀计数模型 ​

零过多的来源 ​

实际计数数据中经常出现零的比例高于标准模型预测的现象。以泊松分布为例,Poisson(4) 的零概率约为 0.0183,但若数据中三成样本都是零,泊松模型无法解释这一比例。零过多的成因可以分成两类。第一类是结构性零(structural zero):某些个体根本不可能出现该事件,例如未暴露于某风险因素的人群永远不会发病,其计数必然为零。第二类是抽样零(sampling zero):个体有出现事件的可能,但在观测期内恰好没有出现,这类零由计数的随机性产生。

两部分的混合结构 ​

零膨胀泊松分布(zero-inflated Poisson,ZIP)用一个两部分结构描述上述机制。以概率 π 产生结构性零,以概率 1−π 进入计数部分并按泊松分布取值,概率质量函数为

p(X=0)=π+(1−π)e−λ,p(X=k)=(1−π)e−λλkk!,k≥1

期望与方差为

\E[X]=(1−π)λ,Var(X)=(1−π)λ(1+πλ)

方差大于均值,超出量为 (1−π)λ⋅πλ,说明零膨胀同时带来了过散布。把计数部分换成负二项分布,得到零膨胀负二项分布(zero-inflated negative binomial,ZINB),它在零过多之外还能容纳计数部分的过散布。ZIP 与 ZINB 常用于一个月内的网购次数、某类故障的发生次数、违章停车次数等计数数据,这些数据中不参与这类活动的人与参与了但事件未发生的人都贡献零值,两类零混在一起。

python
import numpy as np
from scipy import stats

# 零膨胀泊松:结构零比例 pi = 0.3,计数部分 lambda = 4
pi, lam = 0.3, 4.0
k = np.arange(0, 11)
pmf_zip = pi * (k == 0) + (1 - pi) * stats.poisson(lam).pmf(k)
print(f"泊松(4) 零概率 = {stats.poisson(lam).pmf(0):.4f}")   # 0.0183
print(f"ZIP 零概率     = {pmf_zip[0]:.4f}")                  # 0.3128
print(f"ZIP 均值 = {(1 - pi) * lam:.4f}")                    # 2.8000
print(f"ZIP 方差 = {(1 - pi) * lam * (1 + pi * lam):.4f}")    # 6.1600,方差约为均值的 2.2 倍

# 从 ZIP 采样并与泊松对比零的比例
rng = np.random.default_rng(5)
n = 100000
is_struct_zero = rng.random(n) < pi
counts = np.where(is_struct_zero, 0, rng.poisson(lam, size=n))
print(f"ZIP 模拟零比例 = {np.mean(counts == 0):.4f}")          # ≈ 0.3128

零膨胀模型的参数 π 与 λ 分别对应两种机制,解释时需分开讨论:π 反映不可能发生事件的群体比例,λ 反映可能发生事件群体的平均强度。把两者合并为一个均值会掩盖机制上的差异。

2.4.9 经验分布 ​

经验分布函数 ​

经验分布(empirical distribution)把观测样本本身当作分布的估计,不假定任何参数形式。给定样本 x1,…,xn,经验累积分布函数定义为

Fn(x)=1n∑i=1n1{xi≤x}

它统计不超过 x 的样本比例,是一个阶梯函数,每经过一个观测点上升 1/n。由大数定律,Fn(x) 对每个固定的 x 依概率收敛到真实分布函数 F(x);格里汶科-坎泰利定理(Glivenko-Cantelli theorem)进一步保证 Fn 在整个实数轴上一致收敛:supx|Fn(x)−F(x)|→0。经验分布因此成为非参数推断的起点。

与多项分布的关系 ​

把实轴划分成 k 个互不重叠的区间,每个样本落入其中一个区间。落入各区间的频数向量服从多项分布,参数为各区间的真实概率。经验分布正是这一多项分布的实现:它给出每个区间的相对频数。这一联系说明经验分布是多项分布的连续类比,二者的差别在于经验分布的区间可以任意细分,直至每个观测点成为独立的一类。当样本量增大时,经验分布的支撑会越来越密,逐步逼近真实分布的形态。

重抽样的基础 ​

经验分布把概率质量 1/n 均匀分配给每个观测点,因此从经验分布中抽取一个大小为 n 的新样本等价于对原样本做有放回重抽样。这一操作是自助法(bootstrap,2.10 节)的核心,用于估计统计量的抽样分布与构造置信区间。自助法不依赖分布族假设,适用于中位数、分位数、相关系数等难以写出抽样分布的统计量,代价是计算量增加与对原样本代表性的依赖。

python
import numpy as np
from scipy import stats

rng = np.random.default_rng(4)
x = rng.poisson(3, size=200)                     # 来自 Poisson(3) 的样本

def ecdf(t, data):
    return np.mean(data <= t)

print(f"经验 CDF(2) = {ecdf(2, x):.4f}")          # 与理论值接近
print(f"理论 CDF(2) = {stats.poisson(3).cdf(2):.4f}")   # 0.4232

# 经验分布的均值与方差是对样本均值与方差的估计
print(f"经验均值 = {x.mean():.4f},经验方差 = {x.var(ddof=1):.4f}")   # 接近 3 与 3

# 自助法:用经验分布重抽样,估计样本均值的标准误
boot_means = np.array([rng.choice(x, size=x.size, replace=True).mean() for _ in range(2000)])
print(f"自助法标准误 = {boot_means.std(ddof=1):.4f}")     # 接近 sqrt(3/200) ≈ 0.1225

经验分布在样本量小时支撑稀疏,尾部估计粗糙;在样本量大时计算成本随样本量线性增长。它是众多非参数方法的基础,但在需要外推到样本范围之外的场景中并不适用,因为经验分布永远给不出观测范围以外的新取值。

2.4.10 分布关系的整体图景 ​

把本节讨论的分布按生成机制排列,可以得到几条清晰的联系路径。极限关系描述一种分布在大参数极限下趋于另一种分布;混合关系描述把参数随机化后得到的更分散的分布;共轭关系描述参数估计中先验与后验保持同一分布族的情形。

关系类型分布对条件结果
极限二项 → 泊松n→∞,p→0,np→λ稀有事件计数
极限二项 → 正态n→∞,p 固定中心极限定理特例(2.7 节)
极限超几何 → 二项N→∞,K/N→p 固定有放回与无放回趋于一致
极限负二项 → 泊松过散布参数 r→∞过散布消失
极限齐夫 → 幂律N→∞无限支撑幂律
混合贝塔-二项p∼Beta(α,β) 与二项混合过散布的二项
混合狄利克雷-多项p∼Dirichlet(α) 与多项混合过散布的多项
混合零膨胀泊松退化分布与泊松混合零过多且过散布
共轭Beta 与二项先验 Beta、似然二项后验仍为 Beta(2.12 节)
共轭Dirichlet 与多项先验 Dirichlet、似然多项后验仍为 Dirichlet
分解二项 = 伯努利之和n 次独立同分布伯努利期望与方差可加
分解负二项 = 几何之和r 次独立几何r=1 时退化为几何

生成机制的共性可以概括为一句话:先在单次试验层面确定结果的可能取值与概率,再规定试验如何进行(重复次数、是否放回、停止规则),最后决定是否需要引入个体异质或结构零。单次试验的取值从两类别推广到多类别,就得到伯努利到分类、二项到多项的对应;试验条件从独立同分布放宽到独立但不同分布、再到不独立,就得到泊松-二项与超几何分布;引入参数的不确定性,就得到贝塔-二项与狄利克雷-多项分布。

2.4.11 本节小结 ​

从生成机制出发理解离散分布

本节按生成机制把离散分布组织成一条线索。离散均匀分布与分类分布给出最基础的赋值方式;伯努利分布描述单次二元试验,二项分布描述独立重复试验的成功次数,二者由求和关系相连;多项分布把类别推广到多个,贝塔-二项与狄利克雷-多项通过参数随机化引入过散布;几何分布与负二项分布来自改变停止规则,负二项分布同时是过散布计数的常用工具;超几何分布去掉有放回条件,泊松-二项分布去掉成功概率相同的条件;泊松分布刻画稀有事件计数,其期望与方差相等这一特征是检验计数模型是否合适的判据;幂律分布与齐夫分布描述长尾,方差可能不存在;零膨胀模型用两部分混合处理零过多;经验分布把样本当作分布的非参数估计,是自助法的基础。每个分布都对应一个可描述的生成过程,参数都有明确的含义,把参数的含义与数据的来源对应起来,模型的选择才有依据。

下一节把这些分布的连续类比逐一道来。离散情形中由求和定义的累积概率与期望,在连续情形中改为积分;二项分布的连续类比是贝塔分布(作为成功概率的分布),泊松分布的间隔分布是指数分布,幂律分布的连续版本是帕累托分布。连续分布还引入了离散情形没有的工具,例如逆变换采样、正态分布的经验法则与指数族的自然参数。掌握这些对象之后,2.6 节将把一元分布推广到随机向量,讨论联合分布与多元正态分布。

练习题 ​

第 1 题 概念推导 ​

设 X∼Binomial(n,p),其中 X=∑i=1nXi,各 Xi 独立同分布且服从 Bernoulli(p)。由期望与方差的可加性推导 \E[X]=np 与 Var(X)=np(1−p),并说明为什么协方差项在独立情形下消失。

参考答案

先计算单个伯努利随机变量的期望与方差。Xi 取 1 的概率为 p,取 0 的概率为 1−p,故

\E[Xi]=1⋅p+0⋅(1−p)=p,\E[Xi2]=12⋅p+02⋅(1−p)=p

由 Var(Xi)=\E[Xi2]−(\E[Xi])2=p−p2=p(1−p)。

由期望的线性性,\E[X]=∑i=1n\E[Xi]=np,这一步不要求独立性。

对和方差展开:

Var(X)=∑i=1nVar(Xi)+2∑1≤i<j≤nCov(Xi,Xj)

当各 Xi 相互独立时,Cov(Xi,Xj)=\E[XiXj]−\E[Xi]\E[Xj]=p2−p2=0,双重求和项全部消失,于是 Var(X)=np(1−p)。

独立性只用于方差的计算,期望的推导不需要独立性。若各次试验存在正相关,协方差项为正,方差会大于 np(1−p),这正是过散布的来源,也是贝塔-二项分布方差更大的原因。

第 2 题 计算推理 ​

某学习平台记录每名用户在 20 道测验题中的答对数量。用户的答对率存在个体差异,模型取单次答对概率 p∼Beta(2,5),给定 p 后答对题数服从 Binomial(20,p)。计算任意一名用户答对题数的期望与方差,并与固定 p=2/7 的二项模型对比,说明哪一个的方差更大。

参考答案

由贝塔分布的均值 \E[p]=α/(α+β)=2/7,得 \E[X]=20×2/7=5.7143。

方差用贝塔-二项分布的公式:

Var(X)=nαβ(α+β+n)(α+β)2(α+β+1)=20×2×5×(2+5+20)(2+5)2×(2+5+1)=20×27049×8=20×270392≈13.78

固定 p=2/7 的二项模型方差为

Var(X)=20×27×57=20×1049≈4.08

两个模型的均值相同,均为 5.7143,但贝塔-二项模型的方差约为二项模型的 3.375 倍。放大的倍数等于 (n+α+β)/(1+α+β)=27/8=3.375。方差更大是因为用户之间的答对率存在异质,一部分用户答对率高、一部分低,这种个体差异体现在计数上就是更大的波动。使用二项模型会低估不确定性,进而低估置信区间宽度与检验的保守性。

第 3 题 代码验证 ​

用蒙特卡洛方法验证二项分布的泊松近似:取 n=1000、p=0.003,模拟 200000 组试验,比较模拟得到的成功次数分布与 Poisson(3) 在 k=0,1,2,3,4 处的概率质量函数。再取 n=1000、p=0.3,说明此时泊松近似是否仍然合适。

参考答案
python
import numpy as np
from scipy import stats

rng = np.random.default_rng(7)
n_sim = 200000

for n, p in [(1000, 0.003), (1000, 0.3)]:
    lam = n * p
    # 用二项分布直接抽样,避免逐次伯努利试验的开销
    sims = rng.binomial(n, p, size=n_sim)
    po = stats.poisson(lam)
    print(f"n={n}, p={p}, lambda={lam:.1f}")
    for k in range(5):
        print(f"  k={k}: 模拟 {np.mean(sims == k):.5f}  泊松 {po.pmf(k):.5f}")

第一组的理论结果为:

Binomial(1000,0.003)与Poisson(3)在k=0,1,2,3,4处几乎重合

模拟频率与泊松概率的差异在千分之一以内,泊松近似成立。第二组中 n=1000、p=0.3,期望同为 300,但 Binomial(1000,0.3) 的方差为 210,而 Poisson(300) 的方差为 300,两者相差明显,泊松近似不再合适,此时应使用正态近似。泊松近似成立的条件是 p 足够小,仅检查 np 是否为常数并不充分。

第 4 题 综合应用 ​

某平台统计一批用户一个月内的网购次数,样本均值约为 1.8,样本方差约为 6.5,且约 45% 的用户下单次数为零。判断以下三种模型哪一种更合适,并说明理由:(1)Poisson(1.8);(2)负二项分布;(3)零膨胀负二项分布。

参考答案

先检查泊松分布。泊松分布要求方差等于均值,而样本方差 6.5 约为样本均值 1.8 的 3.6 倍,方差远大于均值,说明存在明显的过散布,Poisson(1.8) 不合适。

再检查零的比例。Poisson(1.8) 的零概率为 e−1.8≈0.1653,而样本中零的比例约为 0.45,远高于泊松预测值,说明数据结构中存在零过多。

负二项分布可以容纳过散布,其方差为 μ+μ2/r,通过调节 r 匹配样本方差。令 μ=1.8,解 1.8+1.82/r=6.5 得 r≈0.67。但负二项分布在匹配方差的同时,其零概率未必能达到 0.45,因此它解释了过散布,却不能直接解释零过多的具体来源。

零膨胀负二项分布同时具备两部分结构:一部分以概率 π 产生结构性零,描述不可能就诊的群体;另一部分按负二项分布描述会就诊但次数波动的群体。两个参数分别对应零过多的来源与计数部分的过散布,与样本表现出的一方差大、零比例高两个特征都吻合。三种模型中零膨胀负二项分布最合适,因为它同时匹配了方差与零比例这两个超出简单泊松模型的特征,并把零过多的成因拆分为可解释的两部分。

实际建模时还应检查是否存在过度拟合,并比较 ZIP 与 ZINB 的拟合优度:若计数部分的过散布不明显,ZIP 已经足够;若计数部分方差仍显著大于均值,则应选择 ZINB。

常见错误 ​

错误 1 · 混淆几何分布与负二项分布的两套计数约定

原因:几何分布与负二项分布既可以记录试验总次数,也可以记录失败次数,两套约定的取值范围与均值相差一个偏移量。scipy 的 geom 使用试验次数(取值从 1 开始),nbinom 使用失败次数(取值从 0 开始)。混用约定会让期望与方差的计算出现系统性偏差,进而影响后续的模型比较。

解决:确定约定后再写公式。记录试验总次数时,几何分布的均值为 1/p,负二项分布的均值为 r/p;记录失败次数时,均值分别为 (1−p)/p 与 r(1−p)/p。调用 scipy 函数前先查看文档说明参数含义,必要时用小样本手动验证均值。

错误 2 · 对强相关的计数数据套用二项分布或泊松分布

原因:二项分布要求各次试验独立且成功概率相同,泊松分布要求计数区间不重叠且事件独立。数据中存在个体异质或事件聚集时,方差会超过名义值,标准差被低估,置信区间过窄,检验的显著性被高估。

解决:先比较样本方差与样本均值。方差接近均值时可以尝试泊松模型;方差明显大于均值时改用负二项分布、贝塔-二项分布或零膨胀模型。对无放回抽样场景改用超几何分布,对成功概率不同的独立试验改用泊松-二项分布。

错误 3 · 用多项分布计算两类别计数时忽略负相关

原因:多项分布中各类计数之和固定为 n,不同类别的计数之间存在负协方差 −npipj。若把各类计数当作独立的二项变量分别处理,会遗漏这一约束,导致组合概率与联合区间估计出错。

解决:需要联合分布时直接用多项分布的概率质量函数,或由类别概率与计数总和推导。只关心单类计数时可以退化为二项分布,但涉及两类以上同时出现的问题必须保留负相关结构,不可简单相乘。

错误 4 · 把零膨胀模型的两部分参数合并解释

原因:零膨胀泊松分布由结构零与计数两部分混合而成,参数 π 与 λ 分别描述不同的机制。常见的做法是只看总体均值 \E[X]=(1−π)λ,把两个参数的效果混在一起,无法区分零过多的来源与计数的强度。

解决:解释结果时分开报告 π 与 λ。π 反映不可能发生事件的群体比例,λ 反映可能发生事件群体的平均强度。在推断中同时给出两个参数的估计与置信区间,并结合研究背景解释零的两类成因。