用Python代码可视化:从正态分布到卡方分布,理解统计学的核心桥梁

统计学中有两个经典分布——正态分布和卡方分布,它们不仅是理论基石,更是数据分析的实用工具。但对于初学者来说,公式推导往往让人望而生畏。本文将用Python代码和可视化手段,带你直观理解这两个分布的关系。

1. 正态分布:自然界中的"标准答案"

正态分布被称为"上帝创造的曲线",在自然界中无处不在。从人类身高到测量误差,从考试成绩到股票收益率,都能看到它的身影。

用Python生成正态分布数据只需几行代码:

import numpy as np
import matplotlib.pyplot as plt

# 生成10000个标准正态分布随机数
data = np.random.normal(0, 1, 10000)

# 绘制直方图
plt.hist(data, bins=50, density=True, alpha=0.6, color='g')

# 绘制理论PDF曲线
x = np.linspace(-4, 4, 100)
plt.plot(x, 1/np.sqrt(2*np.pi)*np.exp(-x**2/2), 'r', linewidth=2)

plt.title('标准正态分布')
plt.show()

这段代码展示了:

  • np.random.normal 生成正态分布随机数
  • 直方图展示实际数据分布
  • 红色曲线是理论概率密度函数(PDF)

关键参数解析

参数 含义 标准正态分布取值
μ 均值(分布中心) 0
σ 标准差(分布宽度) 1

2. 卡方分布:正态分布的平方和

卡方分布是统计学中最重要的分布之一,尤其在假设检验中广泛应用。它的定义简单却深刻:k个独立标准正态随机变量的平方和服从自由度为k的卡方分布。

让我们用代码验证这个定义:

from scipy.stats import chi2

k = 3  # 自由度
samples = 10000

# 方法1:直接使用卡方分布生成
chi_data = chi2(k).rvs(samples)

# 方法2:通过正态分布构造
normal_data = np.random.normal(0, 1, (samples, k))
constructed_chi = np.sum(normal_data**2, axis=1)

# 比较两种方法
plt.hist(chi_data, bins=50, density=True, alpha=0.5, label='直接生成')
plt.hist(constructed_chi, bins=50, density=True, alpha=0.5, label='构造生成')
plt.legend()
plt.title(f'自由度为{k}的卡方分布')
plt.show()

运行这段代码,你会发现两种方法生成的分布几乎完全重合,验证了卡方分布的定义。

3. 自由度k如何影响卡方分布

卡方分布的形状高度依赖自由度k。随着k增大,分布会发生有趣的变化:

fig, axes = plt.subplots(2, 2, figsize=(10, 8))
k_values = [1, 5, 10, 30]

for k, ax in zip(k_values, axes.ravel()):
    data = chi2(k).rvs(10000)
    ax.hist(data, bins=50, density=True, alpha=0.6)
    ax.set_title(f'k={k}')
    ax.set_xlim(0, 50)
    
plt.tight_layout()
plt.show()

观察这个可视化,你会发现:

  • 当k=1时,分布高度右偏
  • 随着k增大,分布逐渐对称
  • 当k=30时,已接近对称的钟形曲线

4. 卡方分布趋近正态分布的数学原理

这种现象不是巧合,而是有严格的数学基础——中心极限定理。虽然卡方分布本身就是多个随机变量的和,但只有当k足够大时,和的分布才会趋近正态。

让我们用代码量化这个逼近过程:

from scipy.stats import norm

k_values = [1, 2, 5, 10, 20, 50]
plt.figure(figsize=(10, 6))

for k in k_values:
    # 标准化卡方变量
    mean = k
    std = np.sqrt(2*k)
    data = (chi2(k).rvs(10000) - mean) / std
    
    # 绘制KDE曲线
    sns.kdeplot(data, label=f'k={k}')

# 标准正态分布作为参照
x = np.linspace(-4, 4, 100)
plt.plot(x, norm.pdf(x), 'k--', label='标准正态')
plt.legend()
plt.title('标准化卡方分布随k增大的变化')
plt.show()

这个可视化清晰地展示了:

  1. 当k较小时,分布与正态分布差异明显
  2. 随着k增大,曲线逐渐向标准正态靠拢
  3. 当k=50时,两者已几乎重合

收敛速度分析

自由度k 与正态分布的KS检验统计量
10 0.042
30 0.015
50 0.008
100 0.004

5. 实际应用:卡方检验示例

理解这些分布关系后,我们来看一个实际应用场景——卡方拟合优度检验。这是验证观测数据是否符合某种理论分布的常用方法。

observed = np.array([25, 30, 45])  # 观测频数
expected = np.array([30, 30, 40])  # 期望频数

# 计算卡方统计量
chi2_stat = np.sum((observed - expected)**2 / expected)
p_value = 1 - chi2(2).cdf(chi2_stat)  # 自由度为类别数-1

print(f"卡方统计量: {chi2_stat:.2f}")
print(f"P值: {p_value:.4f}")

在这个例子中:

  • 自由度为2(3个类别-1)
  • 计算出的统计量服从卡方分布
  • P值帮助我们判断差异是否显著

6. 进阶可视化:动态观察分布变化

为了更直观地理解,我们可以创建动态可视化:

from matplotlib.animation import FuncAnimation

fig, ax = plt.subplots(figsize=(8, 5))
k_values = range(1, 51)
line, = ax.plot([], [], 'r-')
ax.set_xlim(0, 80)
ax.set_ylim(0, 0.2)

def update(k):
    ax.clear()
    data = chi2(k).rvs(10000)
    ax.hist(data, bins=50, density=True, alpha=0.6)
    
    # 叠加正态分布曲线
    x = np.linspace(0, 80, 100)
    ax.plot(x, norm.pdf(x, k, np.sqrt(2*k)), 'r-')
    ax.set_title(f'自由度为{k}的卡方分布 vs 正态分布')
    return ax

ani = FuncAnimation(fig, update, frames=k_values, interval=200)
plt.close()

这段代码生成动画,可以清晰看到随着k增大,卡方分布如何逐渐接近对应的正态分布。在实际项目中,这种动态可视化往往比静态图表更能揭示统计规律。

更多推荐