# 教材位置：第2章，第2.8节；图 2-2。
# 单位、变量与形状：values 为 (4000,n)，means 为 (4000,)；原等待时间用统一时间单位，标准
# 化量无量纲。
# 学习目标：用偏斜的指数分布演示中心极限定理，理解均值的标准化和标准误。
# 核心关系：独立同分布、方差有限且非零时，sqrt(n)(样本均值-μ)/σ 的分布随 n 增大
# 趋向标准正态分布；本例指数分布的尺度 θ=2，均值和标准差均为 θ。
# 数据与输出：对 n=1、5、50 各重复 4000 次抽样，打印理论标准误和均值的样本标准差，
# 保存三幅标准化均值直方图与标准正态密度的对照图 e04_clt.png。
# 运行方式：安装 NumPy、SciPy、Matplotlib 后执行 python e04_clt.py；图片写入当前工作目录
# 。

import numpy as np
# stats 提供理论分布的概率密度函数，用于绘制参照曲线。
from scipy import stats
import matplotlib.pyplot as plt
# 在创建第一幅图之前配置可跨系统选择的中文字体。
from plot_style import configure_chinese_plot

configure_chinese_plot()

# 生成器在后面的循环中持续推进；不会在每个样本量下重新设置种子。
rng = np.random.default_rng(0)
# theta 是尺度参数，不是速率参数；对应的指数分布速率为 1/theta。
# repetitions 是重复实验次数，决定直方图的模拟精细程度；n 决定每次均值使用的样本量。
theta, repetitions = 2.0, 4000
# 生成一行三列子图，并统一坐标尺度，以便直接比较分布形状的变化。
fig, axes = plt.subplots(1, 3, figsize=(6.3, 2.8), sharex=True, sharey=True)
# 标准正态参照曲线的横坐标网格；300 是画线点数，与模拟样本量无关。
z = np.linspace(-4, 4, 300)
# zip 将三个坐标轴分别与三个 n 配对；大 n 下的标准化均值通常更接近正态。
for ax, n in zip(axes, [1, 5, 50]):
    # values 形状为 (4000, n)：每行是一轮实验，每列是一轮中的一个独立观测。
    # NumPy 的 exponential 第一个参数为尺度，因此这里 E[X]=theta、Var(X)=theta²。
    values = rng.exponential(theta, size=(repetitions, n))
    # 沿 axis=1 对每行取平均，将 (4000, n) 压缩成 (4000,) 的均值数组。
    means = values.mean(axis=1)
    # 减去总体均值，再除以理论标准误 theta/sqrt(n)；不是除以原观测的标准差 theta。
    # 本例使用已知总体参数，所以对应中心极限定理中的标准化量，而不是 t 统计量。
    standardized = (means - theta) / (theta / np.sqrt(n))
    # density=True 按“频数/(总数×组宽)”归一化，使所有柱形总面积为 1。
    # 柱高是密度而非点概率，区间概率要看该区间覆盖的面积；bins 控制分组数量。
    ax.hist(standardized, bins=40, density=True, alpha=0.6)
    # norm.pdf 给出标准正态密度，不是累积概率；连续随机变量在单点的概率为 0。
    ax.plot(z, stats.norm.pdf(z), color="C1")
    # xlim 只限制显示范围，不会丢弃上述模拟数据，也不会重新计算标准化结果。
    ax.set(title=f"n={n}", xlabel="标准化样本均值", xlim=(-4, 4))
    # 对 4000 个均值求样本标准差，估计“均值的抽样标准差”，应接近 理论标准误。
    # ddof=1 使用 repetitions-1 作分母；有限次模拟下二者不必完全相等。
    print(n, "理论标准误", theta / np.sqrt(n), "均值的样本标准差", means.std(ddof=1))
    # 检查每轮恰好产出一个均值，防止将求平均的轴写错而改变实验含义。
    assert means.shape == (repetitions,)
# 中心极限定理 描述的是均值经标准化后的分布；它不声称原始指数观测变成正态。
# 增加重复次数会使直方图更平滑，却不能代替增大 n 所带来的正态近似改善。
axes[0].set_ylabel("概率密度")
fig.tight_layout()
# 在当前工作目录保存图片后关闭图形，便于脚本在无交互环境下连续运行。
fig.savefig("e04_clt.png", dpi=160)
plt.close(fig)