矩母函数 M(s) 实战:5大常见分布推导与Python代码验证

📅 2026/7/6 22:07:55 👁️ 阅读次数 📝 编程学习
矩母函数 M(s) 实战:5大常见分布推导与Python代码验证

矩母函数 M(s) 实战:5大常见分布推导与Python代码验证

在概率论与统计学的工具箱中,矩母函数(Moment Generating Function, MGF)是一种强大却常被低估的数学工具。不同于直接处理概率密度函数或累积分布函数,MGF通过生成各阶矩的方式,为随机变量的分析提供了统一框架。本文将聚焦五种常见概率分布——泊松、二项、指数、正态和均匀分布,不仅展示其矩母函数的数学推导过程,更通过Python代码实现符号计算与数值验证,让抽象理论落地为可执行的实践。

1. 矩母函数核心概念与Python基础配置

矩母函数M(s)定义为随机变量X的期望值E[e^(sX)],这个看似简单的指数变换蕴含着随机变量的全部矩信息。MGF之所以强大,源于三个关键特性:

  1. 唯一性定理:在适当条件下,MGF与概率分布一一对应
  2. 矩生成特性:n阶导数在s=0处的值给出第n阶矩
  3. 卷积简化:独立随机变量和的MGF等于各MGF的乘积

注意:MGF并非总是存在,当积分或级数不收敛时(如柯西分布),需要改用特征函数

配置Python符号计算环境:

import numpy as np import matplotlib.pyplot as plt from sympy import symbols, exp, summation, oo, integrate, simplify, latex from sympy.stats import density, E, P, Poisson, Binomial, Exponential, Normal, Uniform s, x, λ, n, p, μ, σ, a, b = symbols('s x λ n p μ σ a b', real=True)

2. 泊松分布:离散事件的动态建模

泊松分布是描述单位时间内稀有事件发生次数的经典模型,其概率质量函数为: [ P(X=k) = \frac{e^{-λ}λ^k}{k!}, \quad k=0,1,2,... ]

推导过程: [ M(s) = E[e^{sX}] = \sum_{k=0}^\infty e^{sk} \frac{e^{-λ}λ^k}{k!} = e^{-λ} \sum_{k=0}^\infty \frac{(λe^s)^k}{k!} = e^{λ(e^s - 1)} ]

Python验证代码:

def poisson_mgf(s_val, λ_val): X = Poisson('X', λ_val) mgf = E(exp(s*X)) return mgf.doit().subs({s: s_val, λ: λ_val}) # 数值验证 λ_example = 2.5 s_vals = np.linspace(-1, 0.5, 100) mgf_vals = [float(poisson_mgf(s, λ_example)) for s in s_vals] theoretical = [np.exp(λ_example*(np.exp(s)-1)) for s in s_vals] plt.plot(s_vals, mgf_vals, 'b-', label='SymPy计算') plt.plot(s_vals, theoretical, 'r--', label='理论值') plt.title(f'泊松分布MGF验证 (λ={λ_example})') plt.legend(); plt.show()

矩计算应用

# 计算期望和方差 mgf = exp(λ*(exp(s)-1)) first_moment = mgf.diff(s).subs(s, 0) # 输出: λ second_moment = mgf.diff(s, 2).subs(s, 0) # 输出: λ**2 + λ variance = second_moment - first_moment**2 # 输出: λ

3. 二项分布:独立试验的聚合分析

二项分布描述n次独立伯努利试验的成功次数,其PMF为: [ P(X=k) = \binom{n}{k} p^k (1-p)^{n-k} ]

MGF推导: [ M(s) = \sum_{k=0}^n e^{sk} \binom{n}{k} p^k (1-p)^{n-k} = (pe^s + 1 - p)^n ]

Python实现:

def binomial_mgf(s_val, n_val, p_val): X = Binomial('X', n_val, p_val) return E(exp(s*X)).doit().subs({s: s_val, n: n_val, p: p_val}) # 参数设置 n_example, p_example = 10, 0.4 s_range = np.linspace(-2, 1, 100) binomial_theory = [(p_example*np.exp(s) + 1 - p_example)**n_example for s in s_range] # 可视化对比 plt.figure(figsize=(10,5)) plt.plot(s_range, [binomial_mgf(s, n_example, p_example) for s in s_range], 'bo') plt.plot(s_range, binomial_theory, 'r-') plt.title('二项分布MGF验证'); plt.show()

复合分布案例

# 独立二项变量和的MGF n1, n2, p = 5, 8, 0.3 mgf_product = binomial_mgf(s, n1, p) * binomial_mgf(s, n2, p) simplify(mgf_product) # 输出: (0.3*exp(s) + 0.7)**13 → 符合Bin(n1+n2,p)

4. 指数分布:无记忆性的连续时间建模

指数分布的概率密度函数: [ f(x) = λe^{-λx}, \quad x ≥ 0 ]

MGF推导: [ M(s) = \int_0^\infty e^{sx} λe^{-λx} dx = λ \int_0^\infty e^{-(λ-s)x} dx = \frac{λ}{λ-s} \quad (s < λ) ]

Python验证:

def exponential_mgf(s_val, λ_val): X = Exponential('X', λ_val) mgf = E(exp(s*X)) return mgf.doit().subs({s: s_val, λ: λ_val}) # 收敛域验证 λ_set = 1.5 s_critical = λ_set # MGF在s=λ处发散 print(f"s=λ时的计算结果: {exponential_mgf(λ_set-0.0001, λ_set)}") # 接近+∞ print(f"s>λ时的积分: {integrate(λ*exp((s-λ)*x), (x, 0, oo)).subs({λ: λ_set, s: λ_set+0.1})}") # 发散

可靠性工程应用

# 并联系统MTTF计算 λ1, λ2 = 0.5, 0.8 mgf1 = λ1/(λ1 - s) mgf2 = λ2/(λ2 - s) system_mgf = 1 - (1 - mgf1)*(1 - mgf2) # 并联系统可靠度 mttf = -system_mgf.diff(s).subs(s, 0) # 平均失效时间: 1/λ1 + 1/λ2 - 1/(λ1+λ2)

5. 正态分布:高斯模型的矩分析

标准正态分布N(0,1)的MGF推导: [ M(s) = \int_{-\infty}^\infty e^{sx} \frac{1}{\sqrt{2π}} e^{-x^2/2} dx = e^{s^2/2} ]

一般正态分布N(μ,σ²)的MGF: [ M(s) = e^{μs + \frac{1}{2}σ^2s^2} ]

Python符号推导:

def normal_mgf(s_val, μ_val, σ_val): X = Normal('X', μ_val, σ_val) return E(exp(s*X)).doit().subs({s: s_val, μ: μ_val, σ: σ_val}) # 中心极限定理演示 def clt_demo(sample_size=30): samples = np.random.exponential(scale=1, size=(1000, sample_size)) sample_means = samples.mean(axis=1) s_vals = np.linspace(-0.5, 0.5, 20) empirical_mgf = [np.mean(np.exp(s*sample_means)) for s in s_vals] theoretical = [normal_mgf(s, 1, 1/np.sqrt(sample_size)) for s in s_vals] plt.plot(s_vals, empirical_mgf, 'bo', label='样本均值MGF') plt.plot(s_vals, theoretical, 'r-', label='N(1,1/√n)理论值') plt.legend(); plt.show()

6. 均匀分布:有限支撑集的典型代表

连续均匀分布U(a,b)的MGF: [ M(s) = \frac{e^{bs} - e^{as}}{s(b-a)} ]

离散均匀分布的MGF推导: [ M(s) = \frac{1}{n} \sum_{k=1}^n e^{s k} = \frac{e^s (e^{ns} - 1)}{n (e^s - 1)} ]

Python实现对比:

def continuous_uniform_mgf(s_val, a_val, b_val): X = Uniform('X', a_val, b_val) return E(exp(s*X)).doit().subs({s: s_val, a: a_val, b: b_val}) # 特殊情形验证 a_test, b_test = 0, 1 s_special = 0 # 直接计算会得到0/0不定式 lhopital = exp(b_test*s)*(b_test) - exp(a_test*s)*(a_test) / (b_test - a_test) print(f"s→0极限值: {lhopital.subs(s, 0)}") # 应为1,符合M(0)=1的性质 # 离散均匀分布案例 def discrete_uniform_mgf(s_val, n_val): k = symbols('k', integer=True) return summation(exp(s*k)/n_val, (k, 1, n_val)).subs({s: s_val, n: n_val})

7. 综合应用:矩母函数的工程实践

金融风险建模案例

# 投资组合损失分布建模 λ1, λ2 = 0.3, 0.7 # 两种风险事件发生率 severity1 = 10000 # 风险事件1的损失金额 severity2 = Normal('S2', 5000, 1500) # 风险事件2的损失分布 # 复合过程MGF计算 mgf_poisson1 = exp(λ1*(exp(s*severity1) - 1)) mgf_poisson2 = exp(λ2*(E(exp(s*severity2)) - 1)) total_mgf = mgf_poisson1 * mgf_poisson2 # 计算VaR(95%)的近似解 from scipy.optimize import fsolve def tail_prob(q): return lambda x: 1 - np.real(np.exp(-x*q)*total_mgf.subs(s, x)) - 0.95 var_95 = fsolve(tail_prob(1), 0.01)[0] # 数值反变换

假设检验中的MGF应用

# 似然比检验统计量的MGF def likelihood_ratio_mgf(s_val, n_val, θ0, θ1): mgf_under_h0 = (θ0/θ1)**(n_val*s) * ( (1-θ0)/(1-θ1) )**(n_val*s) return mgf_under_h0.subs({s: s_val, n: n_val}) # 计算检验势函数 def power_function(critical_val, sample_size, true_θ): s_star = symbols('s*') equation = likelihood_ratio_mgf(s_star, sample_size, 0.5, true_θ) - exp(s_star*critical_val) beta = solve(equation, s_star)[0] return 1 - exp(-beta*critical_val)*likelihood_ratio_mgf(beta, sample_size, 0.5, true_θ)