别再死记公式了！用Python从零推导极大似然估计，5分钟搞懂核心思想

统计学课本里那些密密麻麻的数学推导总是让人望而生畏？今天我们用Python代码和可视化图形，带你绕过复杂公式，直击 极大似然估计 （Maximum Likelihood Estimation, MLE）的核心思想。无需微积分基础，只要会写几行Python代码，你就能理解这个支撑现代机器学习的统计基石。

想象你是一位考古学家，在遗址中发现10枚古币，其中7枚是正面朝上。你会如何推测这批古币的"真实"正面概率？这就是MLE要解决的经典问题—— 通过观察到的数据，反推最可能产生这些数据的模型参数 。下面我们用一个正态分布的例子，分三步实现这个思想：

1. 模拟实验环境：生成观测数据

任何统计推断都始于数据。我们先设定一个"真实"的正态分布（均值μ=5，标准差σ=1），然后从中抽样生成模拟观测数据：

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)
true_mu = 5  # 真实均值
true_sigma = 1  # 真实标准差
sample_size = 100  # 样本量
observations = np.random.normal(true_mu, true_sigma, sample_size)

用直方图可视化这些数据点，可以看到它们围绕5左右波动：

plt.hist(observations, bins=20, density=True, alpha=0.6)
plt.title('模拟观测数据分布')
plt.xlabel('数值')
plt.ylabel('频率')
plt.show()

注意：现实中我们只有观测数据，不知道真实的μ和σ。这里生成数据只是为了模拟实验场景。

2. 构建似然函数：量化参数可能性

似然函数L(θ)衡量的是： 在给定参数θ下，观察到当前数据的概率 。对于正态分布，其概率密度函数(PDF)为：

$$ f(x|\mu,\sigma) = \frac{1}{\sigma\sqrt{2\pi}} e^{-\frac{(x-\mu)^2}{2\sigma^2}} $$

因为各数据点独立，整体似然是所有数据点概率的乘积：

def likelihood(mu, sigma, data):
    # 计算单个数据点的概率密度
    single_prob = (1/(sigma * np.sqrt(2*np.pi))) * np.exp(-0.5*((data-mu)/sigma)**2)
    # 返回所有数据点的联合概率（乘积）
    return np.prod(single_prob)

但直接计算乘积会遇到 数值下溢 问题（大量小于1的数相乘结果趋近于0）。因此实践中使用 对数似然 （log-likelihood），将乘积转为求和：

def log_likelihood(mu, sigma, data):
    n = len(data)
    const = -n * np.log(sigma * np.sqrt(2*np.pi))
    sum_sq = -np.sum((data - mu)**2) / (2 * sigma**2)
    return const + sum_sq

3. 寻找最大似然：可视化参数空间

现在我们在μ的可能取值范围内（比如0到10），计算每个μ对应的对数似然值：

mu_candidates = np.linspace(0, 10, 100)
log_likelihoods = [log_likelihood(mu, true_sigma, observations) for mu in mu_candidates]

plt.plot(mu_candidates, log_likelihoods)
plt.axvline(x=true_mu, color='r', linestyle='--', label='真实均值')
plt.xlabel('候选μ值')
plt.ylabel('对数似然值')
plt.title('对数似然函数曲线')
plt.legend()
plt.show()

你会看到曲线在μ=5附近达到峰值——这正是MLE的核心理念： 选择使观测数据出现概率最大的参数值 。用代码找出精确的最大值点：

max_idx = np.argmax(log_likelihoods)
mle_mu = mu_candidates[max_idx]
print(f"极大似然估计的μ值: {mle_mu:.2f}")

运行结果通常会接近5（如4.85），与真实值非常接近。样本量越大，估计越准确。

4. 扩展到多参数估计：联合优化μ和σ

现实中σ也未知。我们可以扩展似然函数，同时优化两个参数：

from mpl_toolkits.mplot3d import Axes3D

# 生成参数网格
mu_grid = np.linspace(4, 6, 50)
sigma_grid = np.linspace(0.5, 1.5, 50)
M, S = np.meshgrid(mu_grid, sigma_grid)

# 计算每个(μ,σ)组合的对数似然
LL = np.array([[log_likelihood(m, s, observations) 
               for m in mu_grid] for s in sigma_grid])

# 3D可视化
fig = plt.figure(figsize=(10, 7))
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(M, S, LL, cmap='viridis')
ax.set_xlabel('μ')
ax.set_ylabel('σ')
ax.set_zlabel('对数似然')
plt.title('双参数对数似然函数曲面')
plt.show()

通过 scipy.optimize 可以自动找到最大值点：

from scipy.optimize import minimize

def neg_log_likelihood(params, data):
    mu, sigma = params
    return -log_likelihood(mu, sigma, data)  # 最小化负对数似然

initial_guess = [4, 1]  # 初始猜测
result = minimize(neg_log_likelihood, initial_guess, args=(observations,))
mle_mu, mle_sigma = result.x
print(f"MLE估计结果: μ={mle_mu:.2f}, σ={mle_sigma:.2f}")

5. 从数学到实践：MLE的应用陷阱

虽然MLE理论优美，但实际应用中需要注意：

• 样本量敏感性 ：小样本下估计可能偏差较大
• 模型误设风险 ：如果假设的分布模型与真实数据生成过程不符，结果将失真
• 计算复杂度 ：对于复杂模型，似然函数可能难以优化

一个改进方案是 正则化 ，通过在似然函数中加入惩罚项防止过拟合：

def regularized_log_likelihood(params, data, lambda_reg=0.1):
    mu, sigma = params
    # 原始对数似然 + L2正则化项（惩罚大的σ值）
    return log_likelihood(mu, sigma, data) - lambda_reg * sigma**2

最后分享一个实用技巧：在Python中，许多统计模型（如 statsmodels 库）已经内置了MLE实现。理解原理后，你可以直接调用：

import statsmodels.api as sm
model = sm.MLE(observations, model_function)
results = model.fit()
print(results.summary())