引言
频率学派把参数看作固定的未知常数,用重复采样的频率解释概率;贝叶斯学派则把参数本身看作随机变量,用先验分布表达在看到数据前的信念,用似然吸收数据证据,得到后验分布。后验不只是给出一个点估计,而是给出参数的完整不确定性刻画——这在样本少、需要决策风险量化的场景里价值巨大。
代价是后验往往没有解析解,必须借助采样或近似推断。本文从贝叶斯定理出发,讲清共轭先验的便利、MCMC 与变分推断两大推断路线,以及概率图模型如何用图结构编码条件独立性,最后用 PyMC 完成端到端实战。
前置:概率分布与假设检验见 /ml-model-evaluation/;因果图与 do 算子的区别见 /ml-causal-inference/。
目录
- 1. 贝叶斯定理与三要素
- 2. 共轭先验
- 3. MLE、MAP 与后验预测
- 4. 概率图模型
- 5. 蒙特卡洛与 MCMC
- 6. 变分推断
- 7. PyMC 实战
- 8. 贝叶斯 A/B 测试与决策
- 9. 常见坑与选型
- 10. 总结
1. 贝叶斯定理与三要素
贝叶斯定理把联合概率的两种分解联系起来:
P(θ | D) = P(D | θ) · P(θ) / P(D)
四个部分各有名字:
| 名称 | 符号 | 含义 |
|---|---|---|
| 先验 | P(θ) | 看到数据前对参数的信念 |
| 似然 | P(D | θ) | 给定参数时观测到数据的概率 |
| 后验 | P(θ | D) | 看到数据后更新的信念 |
| 证据 | P(D) | 归一化常数,∫P(D|θ)P(θ)dθ |
分母 P(D) 与 θ 无关,常写作正比形式:
后验 ∝ 似然 × 先验
import numpy as np
from scipy import stats
# 抛硬币:先验 Beta(1,1)(均匀),观测 7 次正面 3 次反面
prior_a, prior_b = 1, 1
heads, tails = 7, 3
post_a, post_b = prior_a + heads, prior_b + tails # Beta 共轭更新
x = np.linspace(0, 1, 100)
print("后验均值:", round(post_a / (post_a + post_b), 4))
print("后验 95% 区间:", stats.beta.ppf([0.025, 0.975], post_a, post_b).round(3))
2. 共轭先验
当先验与后验属于同一分布族时,称先验为共轭先验(Conjugate Prior),此时后验有闭式解,无需采样。这是贝叶斯推断里最省事的一类情形。
| 似然 | 共轭先验 | 后验更新 |
|---|---|---|
| Bernoulli / Binomial | Beta(α, β) | α+成功次数, β+失败次数 |
| Poisson | Gamma(α, β) | α+Σx, β+n |
| Normal(方差已知) | Normal(μ₀, σ₀²) | 精度加权平均 |
| Normal(均值已知) | Inverse-Gamma | 形状、尺度累加 |
| Categorical / Multinomial | Dirichlet(α) | α_k + 类别计数 |
2.1 Beta-Binomial 的直觉
Beta(α, β) 可以理解成「先验里已经看过 α-1 次成功、β-1 次失败」。观测数据直接加到参数上,先验越弱(α, β 越小),数据影响力越大。
def beta_update(alpha, beta, successes, failures):
return alpha + successes, beta + failures
# 强先验 Beta(50,50),观测 7/3
print(beta_update(50, 50, 7, 3)) # (57, 53) —— 数据被先验拉回
# 弱先验 Beta(1,1)
print(beta_update(1, 1, 7, 3)) # (8, 4) —— 数据主导
2.2 Dirichlet-Multinomial
多类别场景下 Dirichlet 是 Beta 的推广:
import numpy as np
alpha = np.array([1.0, 1.0, 1.0]) # 三类先验
counts = np.array([10, 5, 2]) # 观测计数
posterior = alpha + counts
print("后验参数:", posterior)
print("后验均值:", (posterior / posterior.sum()).round(3))
3. MLE、MAP 与后验预测
3.1 三种点估计的关系
| 方法 | 目标 | 与先验关系 |
|---|---|---|
| MLE | max P(D | θ) | 不使用先验 |
| MAP | max P(θ | D) ∝ P(D|θ)P(θ) | 用先验作为正则 |
| 后验均值 | E[θ | D] | 用整个后验 |
MAP 估计等价于「似然 + 先验正则项」。例如正态先验下 MAP 就是 L2 正则(岭回归),拉普拉斯先验下就是 L1 正则(Lasso)——正则化本质上是引入了先验信念。
from scipy.optimize import minimize
def neg_log_posterior(theta, data, prior_mu=0.0, prior_sigma=1.0):
n = len(data)
ll = -0.5 * np.sum((data - theta) ** 2) # 高斯似然
lp = -0.5 * ((theta - prior_mu) / prior_sigma) ** 2 # 高斯先验
return -(ll + lp)
data = np.array([1.2, 0.8, 1.5, 0.9, 1.1])
res = minimize(neg_log_posterior, x0=[0.0], args=(data,))
print("MAP 估计:", round(res.x[0], 4))
print("MLE(样本均值):", round(data.mean(), 4))
先验均值 0 把 MAP 估计往 0 拉,因此 MAP 小于 MLE——这就是收缩(Shrinkage)。
3.2 后验预测分布
贝叶斯推断的最终目标是预测,而不是估计参数。后验预测分布把参数不确定性积分掉:
P(x_new | D) = ∫ P(x_new | θ) P(θ | D) dθ
比「用点估计代入」更稳健,因为它对参数的不确定性做了平均,避免过拟合。
4. 概率图模型
概率图模型(Probabilistic Graphical Model, PGM)用图结构紧凑地表示多变量联合分布,节点是随机变量,边是依赖关系。
4.1 有向图:贝叶斯网络
贝叶斯网络用有向无环图(DAG)表示因果/依赖方向,联合分布按父节点分解:
P(X1, ..., Xn) = Π P(Xi | Parents(Xi))
例如「下雨 → 洒水器 → 草地湿」的网络:
from pgmpy.models import DiscreteBayesianNetwork
from pgmpy.factors.discrete import TabularCPD
model = DiscreteBayesianNetwork([("Rain", "Wet"), ("Sprinkler", "Wet")])
cpd_rain = TabularCPD("Rain", 2, [[0.8], [0.2]])
cpd_spr = TabularCPD("Sprinkler", 2, [[0.5], [0.5]])
cpd_wet = TabularCPD("Wet", 2,
[[1.0, 0.1, 0.1, 0.01], [0.0, 0.9, 0.9, 0.99]],
evidence=["Rain", "Sprinkler"], evidence_card=[2, 2])
model.add_cpds(cpd_rain, cpd_spr, cpd_wet)
print("模型合法:", model.check_model())
4.2 条件独立与 d-分离
图结构编码了条件独立关系。d-分离(d-separation) 给出判定规则:若给定观测集 Z 后 X 与 Y 之间所有路径都被阻断,则 X ⊥ Y | Z。三种基本结构:
链式 X → M → Y 给定 M 则阻断
分叉 X ← M → Y 给定 M 则阻断
对撞 X → M ← Y 给定 M 反而打通(对撞偏差)
对撞结构是贝叶斯网络最容易出错的地方:两个本独立的变量,在给定共同后代后变得相关。
4.3 无向图与朴素贝叶斯
马尔可夫随机场(MRF)用无向图表示对称依赖,用势函数定义联合分布。朴素贝叶斯则是最简单的贝叶斯网络:假设特征在类别给定下条件独立:
P(y | x1..xn) ∝ P(y) Π P(xi | y)
from sklearn.naive_bayes import GaussianNB
from sklearn.datasets import load_iris
from sklearn.model_selection import cross_val_score
X, y = load_iris(return_X_y=True)
print("朴素贝叶斯 CV:", cross_val_score(GaussianNB(), X, y, cv=5).mean().round(4))
5. 蒙特卡洛与 MCMC
绝大多数真实模型的后验无解析解,需要采样。蒙特卡洛的核心思想是:只要能从未知分布采样,就能用样本均值近似期望。
5.1 蒙特卡洛积分
# 用采样估计 E[X] 与 P(X > 2),X ~ 标准正态
samples = np.random.standard_normal(100000)
print("均值估计:", round(samples.mean(), 4))
print("P(X>2) 估计:", round((samples > 2).mean(), 4))
print("P(X>2) 理论:", round(1 - stats.norm.cdf(2), 4))
5.2 Metropolis-Hastings
当无法直接采样时,构造马尔可夫链,使其平稳分布恰好是目标后验。MH 算法:
1. 从提议分布 q(θ' | θ) 采样候选 θ'
2. 计算接受率 α = min(1, [p(θ')q(θ|θ')] / [p(θ)q(θ'|θ)])
3. 以概率 α 接受,否则留在原地
def metropolis_hastings(log_pdf, init, n_samples=20000, step=0.5):
samples = np.zeros(n_samples)
current = init
for i in range(n_samples):
proposal = current + np.random.normal(0, step)
log_ratio = log_pdf(proposal) - log_pdf(current)
if np.log(np.random.uniform()) < log_ratio:
current = proposal
samples[i] = current
return samples
# 目标:标准正态(未归一化)
log_target = lambda x: -0.5 * x ** 2
chain = metropolis_hastings(log_target, 0.0)
print("后验均值:", round(chain[5000:].mean(), 4)) # 丢弃预热
5.3 Gibbs 与 HMC/NUTS
- Gibbs 采样:逐个变量、用其全条件分布采样,适合共轭结构多的高维模型。
- HMC(Hamiltonian Monte Carlo):借梯度信息在参数空间「弹跳」,在高维下比随机游走的 MH 高效得多。
- NUTS(No-U-Turn Sampler):HMC 的自适应变体,自动调步长与轨迹长度,是 PyMC/Stan 的默认采样器。
5.4 收敛诊断
MCMC 必须检查收敛,不能直接信样本:
| 指标 | 含义 | 合格阈值 |
|---|---|---|
| R-hat | 链间与链内方差比 | < 1.01 |
| ESS | 有效样本量 | 越大越好,> 400 |
| Trace 图 | 链的轨迹 | 应像毛毛虫,不漂移 |
| 自相关 | 样本相关性 | 迅速衰减 |
6. 变分推断
采样在高维大数据下仍慢。变分推断(Variational Inference, VI) 用一个简单分布 q(θ) 近似后验,把推断变成优化:
目标:最小化 KL(q(θ) ‖ p(θ|D))
等价于:最大化 ELBO = E_q[log p(D, θ)] - E_q[log q(θ)]
ELBO 分解为「期望似然」减去「KL 散度」:前者鼓励拟合数据,后者鼓励 q 接近先验。
import numpy as np
def elbo(mu_q, log_sigma_q, data, prior_sigma=1.0):
sigma_q = np.exp(log_sigma_q)
# 期望对数似然(高斯模型)
exp_ll = -0.5 * np.sum((data - mu_q) ** 2) - len(data) * log_sigma_q
# KL(q || prior)
kl = 0.5 * (sigma_q**2 / prior_sigma**2 + mu_q**2 / prior_sigma**2
- 1 - 2 * log_sigma_q + 2 * np.log(prior_sigma))
return exp_ll - kl
data = np.array([1.2, 0.8, 1.5, 0.9, 1.1])
print("ELBO@mu=1.1:", round(elbo(1.1, np.log(0.4), data), 3))
VI 比 MCMC 快几个数量级,适合大规模,代价是近似可能低估后验方差。常用方法有均场近似、随机变分推断、黑盒变分推断(BBVI)。
7. PyMC 实战
用 PyMC 做贝叶斯线性回归,把「广告投入预测销量」写成概率模型。
import numpy as np
import pymc as pm
rng = np.random.default_rng(0)
n = 100
x = rng.normal(0, 1, n)
y = 2.5 * x + 1.0 + rng.normal(0, 0.5, n)
with pm.Model() as model:
intercept = pm.Normal("intercept", mu=0, sigma=5)
slope = pm.Normal("slope", mu=0, sigma=5)
sigma = pm.HalfNormal("sigma", sigma=1)
mu = intercept + slope * x
pm.Normal("y", mu=mu, sigma=sigma, observed=y)
trace = pm.sample(draws=1000, tune=1000, chains=4, cores=4,
random_seed=1, progressbar=False)
print(pm.summary(trace, var_names=["intercept", "slope", "sigma"]))
输出会给出每个参数的后验均值、标准差、HDI(最高密度区间)与 R-hat。相比最小二乘的点估计,贝叶斯回归额外告诉你斜率有多确定——HDI 是否跨过 0 直接回答了「这个效应是否可信」。
7.1 后验预测检查
with model:
post_pred = pm.sample_posterior_predictive(trace, random_seed=1,
progressbar=False)
pred_mean = post_pred.posterior_predictive["y"].mean(dim=["chain", "draw"])
print("预测 RMSE:", round(np.sqrt(((pred_mean - y) ** 2).mean()), 4))
8. 贝叶斯 A/B 测试与决策
贝叶斯方法在 A/B 测试中的优势:直接给出「B 优于 A 的概率」,无需 p 值与多重比较校正,且支持序贯查看。
from scipy import stats
import numpy as np
# A: 100/1000, B: 130/1000
a_post = stats.beta(1 + 100, 1 + 900)
b_post = stats.beta(1 + 130, 1 + 870)
samples = np.random.default_rng(0).beta
sa = a_post.rvs(200000, random_state=0)
sb = b_post.rvs(200000, random_state=0)
print("P(B > A):", round((sb > sa).mean(), 4))
print("预期提升:", round((sb - sa).mean() / sa.mean() * 100, 2), "%")
print("提升 95% 区间:", np.percentile(sb - sa, [2.5, 97.5]).round(4))
8.1 损失函数决策
光看「B 更好」不够,还要看选错的期望损失:
loss_choose_b = np.maximum(sa - sb, 0).mean()
loss_choose_a = np.maximum(sb - sa, 0).mean()
print("选 B 的期望损失:", round(loss_choose_b, 5))
print("选 A 的期望损失:", round(loss_choose_a, 5))
当「选 A 的期望损失」低于业务可接受阈值时即可停止实验上线 B。
8.2 与贝叶斯优化的联系
贝叶斯优化用高斯过程代理目标函数,再用采集函数(EI/UCB)选下一个评估点——它的「后验更新」正是本文的贝叶斯推断框架,具体实现见 /ml-automl-hpo/。
9. 常见坑与选型
| 现象 | 根因 | 处理 |
|---|---|---|
| R-hat 远大于 1.01 | 链未收敛 | 增加 tune、重参数化、检查模型 |
| 采样极慢 | 高维、共轭少 | 换 NUTS、中心化参数、用 VI |
| 后验被先验主导 | 先验太强或数据太少 | 用弱信息先验,或收集更多数据 |
| 后验区间过窄 | VI 低估方差 | 用全秩 VI 或换 MCMC 验证 |
| 结果不可复现 | 随机种子未固定 | 固定 random_seed,报告多链结果 |
9.1 采样还是变分
数据量小、维度低、要精确后验?
├── 是 → MCMC(NUTS)
└── 否
├── 大数据、要快速迭代 → 变分推断 / 随机 VI
├── 有共轭结构 → Gibbs / 闭式解
└── 只是要个点估计 + 正则 → MAP(等价于带正则的 MLE)
10. 总结
10.1 核心要点
- 贝叶斯推断 = 先验 + 似然 → 后验,后验给出完整的不确定性刻画。
- 共轭先验让后验有闭式解;MAP 等价于带先验正则的 MLE。
- 概率图模型用图结构编码条件独立,对撞结构是最易踩的坑。
- 无解析解时用 MCMC(精确但慢)或 变分推断(快但近似)。
10.2 何时选贝叶斯
当样本稀少、需要不确定性量化、需要序贯更新(如在线 A/B)、或需要把领域知识写进先验时,贝叶斯方法明显优于频率方法。反之,纯预测任务、样本充足时,频率学派方法往往更简单够用。想进一步理解采样与优化的数值方法,可延伸阅读 Python 专题 ,若关注算法与工程实现深度,可参考 AI/ML 专题 。
继续阅读
探索更多技术文章
浏览归档,发现更多关于系统设计、工具链和工程实践的内容。