Fisher信息量实战:用Python验证Cramér-Rao下界与MLE效率
Fisher信息量实战用Python验证Cramér-Rao下界与MLE效率1. 理解Fisher信息量的统计意义Fisher信息量是衡量概率分布中参数所携带信息多少的重要指标。想象你是一位侦探Fisher信息量就是你手中的信息放大镜——它决定了你从数据中推断参数的精确程度。对于一个参数θFisher信息量I(θ)定义为import numpy as np def fisher_information(p): 计算Bernoulli分布的Fisher信息量 return 1 / (p * (1 - p))这个简单的公式背后蕴含着深刻的统计意义当p接近0或1时Fisher信息量急剧增大意味着我们更容易准确估计极端概率在p0.5时信息量最小此时区分正负类最具挑战性Fisher信息量的核心性质与得分函数对数似然的梯度的方差直接相关决定了Cramér-Rao下界——任何无偏估计量方差的理论下限与MLE的渐进方差成反比关系注意Fisher信息量越大意味着参数估计可以达到的精度越高这就像拥有更高像素的相机能捕捉更清晰的图像。2. Cramér-Rao下界的数学本质Cramér-Rao不等式为无偏估计量的方差设定了理论下限Var(θ̂) ≥ 1/(nI(θ))这个不等式告诉我们即使是最优的无偏估计量其方差也不可能突破这个界限。让我们通过Bernoulli分布的例子来验证这一点。Bernoulli分布的Cramér-Rao下界计算参数pFisher信息量I(p)n100时的下界0.111.110.00090.34.760.00210.54.000.00253. 构建蒙特卡洛实验框架我们将设计一个完整的实验来验证MLE的效率是否达到Cramér-Rao下界。以下是实验的关键步骤def monte_carlo_experiment(p, n, num_simulations10000): 执行蒙特卡洛实验验证MLE性质 fisher_info fisher_information(p) cr_lower_bound 1 / (n * fisher_info) mle_variances [] for _ in range(num_simulations): samples np.random.binomial(1, p, sizen) p_hat np.mean(samples) mle_variances.append((p_hat - p)**2) empirical_variance np.mean(mle_variances) return { CR_lower_bound: cr_lower_bound, empirical_variance: empirical_variance, ratio: empirical_variance / cr_lower_bound }实验设计要点固定真实参数p和样本量n重复生成随机样本并计算MLE收集MLE的方差并与理论下界比较4. 结果可视化与分析让我们用不同样本量进行实验并可视化结果import matplotlib.pyplot as plt sample_sizes [10, 50, 100, 500, 1000] p 0.3 results [monte_carlo_experiment(p, n) for n in sample_sizes] plt.figure(figsize(10, 6)) plt.plot(sample_sizes, [r[CR_lower_bound] for r in results], labelCramér-Rao下界) plt.plot(sample_sizes, [r[empirical_variance] for r in results], labelMLE经验方差) plt.xlabel(样本量) plt.ylabel(方差) plt.title(MLE方差与Cramér-Rao下界比较(p0.3)) plt.legend() plt.grid(True) plt.show()典型实验结果会显示小样本时MLE方差明显高于下界随着样本量增加MLE方差逐渐接近并最终达到下界验证了MLE的渐进有效性5. 渐进正态性的实证检验MLE的另一个重要性质是渐进正态性√n(θ̂ - θ) → N(0, 1/I(θ))我们可以通过QQ图来验证这一性质from scipy import stats def check_asymptotic_normality(p, n, num_simulations5000): 检验MLE的渐进正态性 mles [] for _ in range(num_simulations): samples np.random.binomial(1, p, sizen) mles.append(np.mean(samples)) standardized (np.array(mles) - p) * np.sqrt(n * fisher_information(p)) stats.probplot(standardized, distnorm, plotplt) plt.title(fQQ图验证渐进正态性(n{n})) plt.show()当样本量足够大时QQ图上的点将近似落在一条直线上直观验证了MLE的渐进正态分布性质。6. 实际应用中的注意事项虽然理论性质优美但在实际应用中需要注意小样本问题当np或n(1-p)较小时正态近似可能不准确考虑使用Wilson区间等改进方法边界情况处理def safe_mle(samples): 处理全0或全1样本的MLE计算 p_hat np.mean(samples) n len(samples) if p_hat 0: return 1/(2*n) elif p_hat 1: return 1 - 1/(2*n) return p_hat模型误设影响当真实分布不是Bernoulli时MLE可能不再有效可通过稳健标准误等方法减轻影响7. 扩展到其他分布虽然我们以Bernoulli分布为例但这种方法可以推广到其他分布。例如对于泊松分布def poisson_fisher_information(lambd): 泊松分布的Fisher信息量 return 1 / lambd def poisson_mle(samples): 泊松分布的MLE return np.mean(samples)不同分布Fisher信息量对比分布参数Fisher信息量公式Bernoullip1/[p(1-p)]Poissonλ1/λNormal(μ已知)σ²1/(2σ⁴)Exponentialλ1/λ²理解这些差异有助于我们在不同建模场景中选择合适的分布和评估估计精度。