概率密度函数的参数估计：最大似然估计、贝叶斯估计与贝叶斯学习

问题描述

\qquad 在模式分类问题中，贝叶斯判决的核心是通过贝叶斯公式计算后验概率：

$p({\omega }_i|\mathbf{x})=\frac{p(\mathbf{x}|{\omega }_i)p({\omega }_i)}{\sum _{j=1}^{c}p(\mathbf{x}|{\omega }_j)p({\omega }_j)}$

\qquad 其中，$p({\omega }_i)$ 为先验概率，$p(\mathbf{x}|{\omega }_i)$ 为类条件概率密度。

例: 假设单位某个人数较多部门的“男性职工(${\omega }_1$类)”人数为 $n_1$，“女性职工(${\omega }_2$类)”人数为 $n_2$，用身高 $x_1$ 和体重 $x_2$ 作为特征向量 $\mathbf{x}=[x_1,x_2{]}^T$。
　
先验概率：　　　　$p({\omega }_1)=\frac{n_1}{n_1+n_2}$ 和 $p({\omega }_2)=\frac{n_2}{n_1+n_2}$
　
类条件概率密度：　$\{\begin{array}{ll}p(\mathbf{x}|{\omega }_1)& \text{男性职工的(身高,体重)分布}\\ p(\mathbf{x}|{\omega }_2)& \text{女性职工的(身高,体重)分布}\end{array}$
　
贝叶斯判决：将新样本 ${\mathbf{x}}^{\ast }$ 判为最大后验概率所对应的类别 $\underset{j}{\mathrm{argmax}}p({\omega }_j|{\mathbf{x}}^{\ast })$（“最小错误概率准则”）
考虑该单位其他部门的某个人，通过测量其身高、体重组成特征向量 ${\mathbf{x}}^{\ast }=[x_1^{\ast },x_2^{\ast }{]}^T$ 来判断是男性(${\omega }_1$类)还是女性(${\omega }_2$类)。若 $p({\omega }_1|{\mathbf{x}}^{\ast })>p({\omega }_2|{\mathbf{x}}^{\ast })$，则认为是男性；否则，认为是女性。

\qquad 显然，如何有效地估计“类条件概率密度函数” $p(\mathbf{x}|{\omega }_i)$ 最为关键，常用的估计方法有：参数方法（最大似然估计，贝叶斯估计，递归贝叶斯学习）、非参数方法（$\text{Parzen}$窗，$\text{kNN}$）以及朴素贝叶斯等方法。

在上例中，若男性职工(${\omega }_1$类)样本集为 D 1 = { x 1 , ⋯   , x n 1 } \mathcal D_1=\{\boldsymbol x_1,\cdots,\boldsymbol x_{n_1}\} D1​={x1​,⋯,xn1​​}
　　　　　　女性职工(${\omega }_2$类)样本集为 D 2 = { x n 1 + 1 , ⋯   , x n 1 + n 2 } \mathcal D_2=\{\boldsymbol x_{{n_1+1}},\cdots,\boldsymbol x_{n_1+n_2}\} D2​={xn1​+1​,⋯,xn1​+n2​​}　
通过已观测到的样本来计算“类条件概率密度”就是要计算 $p(\mathbf{x}|{\omega }_1,{\mathcal{D}}_1)$ 和 $p(\mathbf{x}|{\omega }_2,{\mathcal{D}}_2)$

\qquad 不失一般性，假设针对某个特定的类 $\omega$ 有 $n$ 个观测数据 D = { x 1 , ⋯   , x n } \mathcal D=\{\boldsymbol x_1, \cdots, \boldsymbol x_n\} D={x1​,⋯,xn​}。“类条件概率密度函数估计”问题就变成了，根据 $\omega$ 类的观测数据集 $\mathcal{D}$ 来估计“条件概率密度函数” $p(\mathbf{x}|\mathcal{D})$。

$(1)$ 若观测数据 $\mathcal{D}$ 也固定不变，就可以把“条件概率密度函数” $p(\mathbf{x}|\mathcal{D})$ 简记为“概率密度函数” $p(\mathbf{x})$。

$(2)$ 若概率密度函数的参数形式 $p(\mathbf{x}|\mathbf{\theta })$ 已知、而参数 $\mathbf{\theta }$ 未知，估计“概率密度函数”的问题就转变为对 $p(\mathbf{x}|\mathbf{\theta })$ 的参数 $\mathbf{\theta }$ 的估计问题。

在上例中，若假设男性和女性的类条件概率密度函数 $p(\mathbf{x}|{\omega }_i,{\mathcal{D}}_i,{\mathbf{\theta }}_i)$ 都服从某个正态分布：
此时，$p(\mathbf{x}|{\omega }_1,{\mathcal{D}}_1,{\mathbf{\theta }}_1)\sim \mathcal{N}({\mathbf{\mu }}_1,{\Sigma }_1)$ ，其中参数 ${\mathbf{\theta }}_1=({\mathbf{\mu }}_1,{\Sigma }_1)$
　　　$p(\mathbf{x}|{\omega }_2,{\mathcal{D}}_2,{\mathbf{\theta }}_2)\sim \mathcal{N}({\mathbf{\mu }}_2,{\Sigma }_2)$ ，其中参数 ${\mathbf{\theta }}_2=({\mathbf{\mu }}_2,{\Sigma }_2)$
　
已知某个特定的类 $\omega$ 服从正态分布 $\mathcal{N}(\mathbf{\mu },\Sigma )$，只是参数 $\mathbf{\theta }=(\mathbf{\mu },\Sigma )$ 未知，针对该类的观测数据集固定为 $\mathcal{D}$，那么估计该类的“类条件概率密度函数” $p(\mathbf{x}|\omega ,\mathcal{D},\mathbf{\theta })$ 就可以简记为估计 $p(\mathbf{x}|\mathbf{\theta })\sim \mathcal{N}(\mathbf{\mu },\Sigma )$ 的参数 $\mathbf{\theta }=(\mathbf{\mu },\Sigma )$ 的值。

\qquad

1. 最大似然估计

\qquad 已知观测数据集 D = { x 1 , ⋯   , x n } \mathcal D=\{\boldsymbol x_1, \cdots, \boldsymbol x_n\} D={x1​,⋯,xn​}，同时假设概率密度函数的参数形式 $p(\mathbf{x}|\mathbf{\theta })$ 是已知的、而参数 $\mathbf{\theta }$ 未知。基于观测数据集 $\mathcal{D}$ 对概率密度函数 $p(\mathbf{x}|\mathbf{\theta })$ 的参数 $\mathbf{\theta }$ 进行最大似然估计时，参数 $\mathbf{\theta }$ 是一个确定而未知的量。

$(1)$ 由于观测数据 D = { x 1 , ⋯   , x n } \mathcal D=\{\boldsymbol x_1, \cdots, \boldsymbol x_n\} D={x1​,⋯,xn​} 是独立抽取的，可定义观测数据集 $\mathcal{D}$ 的似然函数：

$p(\mathcal{D}|\mathbf{\theta })=p({\mathbf{x}}_1,\cdots ,{\mathbf{x}}_n|\mathbf{\theta })=\prod _{i=1}^{n}p({\mathbf{x}}_i|\mathbf{\theta })$

$(2)$ 当参数 $\mathbf{\theta }$ 取不同的值，显然 $p(\mathcal{D}|\mathbf{\theta })$ 的值也会有不同的变化。

\qquad 　　若参数 $\mathbf{\theta }=\hat{\mathbf{\theta }}$ 时，$p(\mathcal{D}|\hat{\mathbf{\theta }})=\underset{\mathbf{\theta }}{\max }p(\mathcal{D}|\mathbf{\theta })$，那么 $\hat{\mathbf{\theta }}$ 就是概率密度函数参数估计的最大似然解，由观测数据集 $\mathcal{D}$ 估计出的概率密度函数就为 $p(\mathbf{x}|\hat{\mathbf{\theta }})$。

\qquad
最大似然估计的求解

$(1)$ 定义对数似然函数：

$\ell (\mathbf{\theta })=\ln p(\mathcal{D}|\mathbf{\theta })$

$(2)$ 最大似然解：

$\hat{\mathbf{\theta }}=\underset{\mathbf{\theta }}{\mathrm{argmax}}p(\mathcal{D}|\mathbf{\theta })$

\quad\quad\quad 或者说 $p(\mathcal{D}|\hat{\mathbf{\theta }})=\underset{\mathbf{\theta }}{\max }p(\mathcal{D}|\mathbf{\theta })$

$(3)$ 通过计算 ${\nabla }_{\mathbf{\theta }}\ell (\mathbf{\theta })=0$（对每个 ${\theta }_i$ 求偏导）求得 $\hat{\mathbf{\theta }}$ 的值：

${\nabla }_{\mathbf{\theta }}\ell (\mathbf{\theta })=\sum _{i=1}^n{\nabla }_{\mathbf{\theta }}\ln p({\mathbf{x}}_i|\mathbf{\theta })=0$

\quad\quad\quad 或者说 $\frac{\partial \ell (\mathbf{\theta })}{\partial {\theta }_i}=0,\forall i$

\quad 参数 $\mathbf{\theta }$ 的最大似然解 $\hat{\mathbf{\theta }}$ 实际上是一个“完全由观测样本决定的统计量”，也就是说：

$\hat{\mathbf{\theta }}=g({\mathbf{x}}_1,\cdots ,{\mathbf{x}}_n)$

“线性回归参数模型”的例子，可参考《线性模型参数求解的最大似然估计》

例1：以一维情况为例，假设概率密度函数的参数形式 $p(x|\theta )\sim \mathcal{N}(\mu ,{\sigma }^2)$，已知方差 ${\sigma }^2=18^2$，未知参数 $\theta =\mu$（真实值为 $20$）。

from scipy.stats import multivariate_normal
from scipy.stats import norm
import numpy as np
import matplotlib.pyplot as plt
mean = 20
std = 18
var = std**2
num = 100
x = np.linspace(-50, 100, num*2, endpoint=False)
y = norm.pdf(x, mean, std)
y1 = multivariate_normal.pdf(x, mean, var)
plt.figure(1, figsize=(8,4))
plt.plot(x, y), plt.xlabel('x'), plt.ylabel('y')
plt.title('$\mu=20, \sigma=18^2$')
# 生成数据集 D，共有num个样本点
rx = mean + std * np.random.randn(num)
ry = np.zeros_like(rx)
plt.plot(rx, ry, 'r.', markersize='3')
# MLE (假设均值未知，方差已知)
mu = np.linspace(mean-10, mean+10, 100)
likelihood = mu.copy()
# 计算数据集 D 对不同均值的似然值
for i in range(len(mu)):
    likelihood[i] = np.prod(norm.pdf(rx, mu[i], std))
mean_ml = mu[np.argmax(likelihood)]
print(mean_ml)
plt.figure(2, figsize=(8,4))
plt.plot(mu, likelihood)# 显示不同均值时的似然函数值，横坐标为未知参数（均值）
plt.title('$\hat\mu=$'+str(mean_ml))
plt.show()

运行结果：

Matlab版本可参考《正态分布的最大似然估计》

\qquad

2. 贝叶斯参数估计

\qquad 参数方法假设概率密度函数的参数形式 $p(\mathbf{x}|\mathbf{\theta })$ 已知，参数向量 $\mathbf{\theta }$ 未知。

$(1)$ “贝叶斯参数估计”方法将参数向量 $\mathbf{\theta }$ 看成是一个随机向量

\qquad 　　在观测任何训练样本之前，我们对参数向量 $\mathbf{\theta }$ 的全部了解只有其先验概率 $p(\mathbf{\theta })$

$(2)$ 根据先验概率 $p(\mathbf{\theta })$ 抽取出一个参数向量 $\mathbf{\theta }=\hat{\mathbf{\theta }}$，从而确定了概率密度函数 $p(\mathbf{x}|\hat{\mathbf{\theta }})$

\qquad 　　若观测数据 D = { x 1 , ⋯   , x n } \mathcal D=\{\boldsymbol x_1, \cdots, \boldsymbol x_n\} D={x1​,⋯,xn​} 是由这个概率密度函数 $p(\mathbf{x}|\hat{\mathbf{\theta }})$ 经过 $n$ 次独立抽取的
\qquad 　　就可以求出观测数据集 $\mathcal{D}$ 关于参数向量 $\hat{\mathbf{\theta }}$ 的似然值 $p(\mathcal{D}|\hat{\mathbf{\theta }})=\prod _{i=1}^{n}p({\mathbf{x}}_i|\hat{\mathbf{\theta }})$

$(3)$ 应用贝叶斯公式，可求出当前观测数据集 $\mathcal{D}$ 关于参数向量 $\mathbf{\theta }=\hat{\mathbf{\theta }}$ 的后验概率

$p(\mathbf{\theta }|\mathcal{D})=\frac{p(\mathcal{D}|\mathbf{\theta })p(\mathbf{\theta })}{\int p(\mathcal{D}|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }}$

\qquad 　　若 $p(\hat{\mathbf{\theta }}|\mathcal{D})=\max p(\mathbf{\theta }|\mathcal{D})$，且 $p(\mathbf{\theta }|\mathcal{D})$ 在 $\hat{\mathbf{\theta }}$ 附近形成一个显著的尖峰
\qquad 　　就可以认为概率密度函数 $p(\mathbf{x})=p(\mathbf{x}|\mathcal{D})\approx p(\mathbf{x}|\hat{\mathbf{\theta }})$

$\begin{array}{rl}p(\mathbf{x}|\mathcal{D})& =\int p(\mathbf{x},\mathbf{\theta }|\mathcal{D})d\mathbf{\theta }\\ & =\int p(\mathbf{x}|\mathbf{\theta },\mathcal{D})p(\mathbf{\theta }|\mathcal{D})d\mathbf{\theta }\\ & =\int p(\mathbf{x}|\mathbf{\theta })p(\mathbf{\theta }|\mathcal{D})d\mathbf{\theta }\text{训练样本集 }\mathcal{D}\text{与测试样本}\mathbf{x}\text{ 的选取是独立抽取的}\end{array}$
　
\qquad 若 $p(\mathbf{\theta }|\mathcal{D})$ 正好在 $\hat{\mathbf{\theta }}$ 附近形成一个显著的尖峰，$p(\hat{\mathbf{\theta }}|\mathcal{D})$ 的值比其他位置的 $p(\mathbf{\theta }|\mathcal{D})$ 都大很多，比如 $p(\mathbf{\theta }|\mathcal{D})\to 0,\forall \mathbf{\theta }\ne \hat{\mathbf{\theta }}$，此时，显然有 $p(\mathbf{x}|\mathcal{D})\approx p(\mathbf{x}|\hat{\mathbf{\theta }})$
　
\qquad 最理想的情况是 $p(\mathbf{\theta }|\mathcal{D})=\delta (\mathbf{\theta }-\hat{\mathbf{\theta }})$，此时 $p(\mathbf{x}|\mathcal{D})=p(\mathbf{x}|\hat{\mathbf{\theta }})$

\qquad
例2：以一维情况为例，假设概率密度函数的参数形式 $p(x|\theta )\sim \mathcal{N}(\mu ,{\sigma }^2)$，且方差 ${\sigma }^2$ 已知，未知参数 $\theta =\mu$。

$(1)$ 贝叶斯参数估计方法认为 $\theta$ 是一个随机变量，并假设 $p(\theta )\sim \mathcal{N}({\mu }_0,{\sigma }_0^2)$，其中 ${\mu }_0,{\sigma }_0^2$ 是已知的

$(2)$ 根据 $p(\theta )\sim \mathcal{N}({\mu }_0,{\sigma }_0^2)$ 抽取一个具体的 $\mu$ 值（比如取 $\mu =\hat{\mu }$），也就确定了概率密度函数 $p(x|\theta )\sim \mathcal{N}(\hat{\mu },{\sigma }^2)$，并认为观测数据 D = { x 1 , ⋯   , x n } \mathcal D=\{x_1, \cdots, x_n\} D={x1​,⋯,xn​} 是根据该概率分布独立地抽取出的 $n$ 个样本，从而计算出观测数据集 $\mathcal{D}$ 关于参数向量 $\hat{\mu }$ 的似然值 $p(\mathcal{D}|\hat{\mu })=\prod _{i=1}^{n}p(x_i|\hat{\mu })$

$(3)$ 由贝叶斯公式，可知 $p(\mu |\mathcal{D})$ 同样服从正态分布

$p(\mu |\mathcal{D})=\frac{p(\mathcal{D}|\mu )p(\mu )}{\int p(\mathcal{D}|\mu )p(\mu )d\mu }=\cdots =\frac{1}{\sqrt{2\pi }{\sigma }_n}\exp \left\{-\frac{1}{2}{\left(\frac{\mu -{\mu }_n}{{\sigma }_n}\right)}^2\right\}$

\qquad 其中，$\{\begin{array}{l}\frac{1}{{\sigma }_n^2}=\frac{n}{{\sigma }^2}+\frac{1}{{\sigma }_0^2}\\ \\ \frac{{\mu }_n}{{\sigma }_n^2}=\frac{n}{{\sigma }_n^2}{\hat{\mu }}_n+\frac{{\mu }_0}{{\sigma }_0^2}\end{array}$，${\hat{\mu }}_n=\frac{1}{n}\sum _{k=1}^n{x}_k$ 为样本均值

\qquad 从而得到：$\{\begin{array}{l}{\mu }_n=\left(\frac{n{\sigma }_0^2}{n{\sigma }_0^2+{\sigma }^2}\right){\hat{\mu }}_n+\frac{{\mu }_0{\sigma }^2}{n{\sigma }_0^2+{\sigma }^2}\\ \\ {\sigma }_n^2=\frac{{\sigma }_0^2{\sigma }^2}{n{\sigma }_0^2+{\sigma }^2}\end{array}$

详细推导过程见《Pattern Classification》3.4节公式（26）~（35）
贝叶斯参数估计通常需要满足共轭性$\text{(conjugacy)}$的条件，也就是“后验概率（正⽐于先验和似然函数的乘积）与先验概率有着相同的函数形式”，可参考《Pattern Recognition and Machine Learning》第2章

取自于《Pattern Classification》Fig 3.2（左图为二维示意图，右图为三维示意图）
　
可以看出，观测样本数 $n$ 越大，越使得 $p(\mu |\mathcal{D})$ 形成一个显著的尖峰，也就表示未知参数 $\mu$ 的值更多可能出现 ${\mu }^{\ast }={\mathrm{argmax}}_{\mu }p(\mu |\mathcal{D})$ 的附近（当$n\to \infty$时，${\sigma }_n^2\to 0$，尖峰越显著）。
显然，${\mu }_n\to {\hat{\mu }}_n(n\to \infty )$ 趋于样本均值，${\sigma }_n^2=\frac{{\sigma }_0^2{\sigma }^2}{n{\sigma }_0^2+{\sigma }^2}\to 0(n\to \infty )$

\qquad

3. 顺序(sequential)贝叶斯学习

\qquad 贝叶斯学习，是一个“增量学习$\text{(incremental learning)}$”过程。其主要特点是，学习过程随着观测数据的不断获取而不断进行下去。

$(1)$ 若 $n-1$ 时刻观测样本集为 D n − 1 = { x 1 , ⋯   , x n − 1 } \mathcal D^{n-1}=\{\boldsymbol x_1, \cdots, \boldsymbol x_{n-1}\} Dn−1={x1​,⋯,xn−1​}，随着第 $n$ 时刻观测到数据 ${\mathbf{x}}_n$，观测样本集更新为 D n = D n − 1 ∪ { x n } = { x 1 , ⋯   , x n − 1 , x n } \mathcal D^{n}=\mathcal D^{n-1}\cup\{\boldsymbol x_n\}=\{\boldsymbol x_1, \cdots, \boldsymbol x_{n-1},\boldsymbol x_n\} Dn=Dn−1∪{xn​}={x1​,⋯,xn−1​,xn​}

$(2)$ 由于观测样本集 $\mathcal{D}$ 一直在更新，$p(\mathbf{\theta }|\mathcal{D})$ 的值也一直在更新，直到 $p(\mathbf{\theta }|{\mathcal{D}}^n)$ 随着 $n\to \infty$ 而形成一个显著的尖峰

\qquad
贝叶斯学习的实现步骤如下：

$(1)$ 第 $0$ 时刻

\qquad 观测数据集 ${\mathcal{D}}^0=\varnothing$，定义 $p(\mathbf{\theta }|{\mathcal{D}}^0)=p(\mathbf{\theta })$，$p(\mathbf{\theta })$ 是未知参数向量的先验概率。

$(2)$ 第 $1$ 时刻

\qquad 从 $p(\mathbf{\theta }|{\mathcal{D}}^0)$ 中抽取一个 $\mathbf{\theta }$ 值用于确定概率密度函数 $p(\mathbf{x}|\mathbf{\theta })$，同时观测第一个数据 ${\mathbf{x}}_1$，观测数据集更新为 D 1 = { x 1 } \mathcal D^{1}=\{\boldsymbol x_1\} D1={x1​}，然后计算似然函数 $p({\mathbf{x}}_1|\mathbf{\theta })$，最后使用贝叶斯公式计算后验概率：

$p(\mathbf{\theta }|{\mathcal{D}}^1)=\frac{p({\mathbf{x}}_1|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^0)}{\int p({\mathbf{x}}_1|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^0)d\mathbf{\theta }}$，此处 $p(\mathbf{\theta }|{\mathcal{D}}^0)=p(\mathbf{\theta })$

与上节中描述的贝叶斯参数估计不同的是，该时刻的后验概率 $p(\mathbf{\theta }|{\mathcal{D}}^1)$ 将作为下一时刻贝叶斯公式中的先验概率，也就是先验概率由 $p(\mathbf{\theta }|{\mathcal{D}}^0)=p(\mathbf{\theta })$ 更新为 $p(\mathbf{\theta }|{\mathcal{D}}^1)$

$(3)$ 第 $2$ 时刻

\qquad 从 $p(\mathbf{\theta }|{\mathcal{D}}^1)$ 中抽取一个 $\mathbf{\theta }$ 值用于确定概率密度函数 $p(\mathbf{x}|\mathbf{\theta })$，同时观测第二个数据 ${\mathbf{x}}_2$，观测数据集更新为 D 2 = { x 1 , x 2 } \mathcal D^{2}=\{\boldsymbol x_1,\boldsymbol x_2\} D2={x1​,x2​}，然后计算数据 ${\mathbf{x}}_2$ 的似然函数 $p({\mathbf{x}}_2|\mathbf{\theta })$，最后使用贝叶斯公式计算后验概率：

$p(\mathbf{\theta }|{\mathcal{D}}^2)=\frac{p({\mathbf{x}}_2|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^1)}{\int p({\mathbf{x}}_2|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^1)d\mathbf{\theta }}$

该时刻的后验概率 $p(\mathbf{\theta }|{\mathcal{D}}^2)$ 将作为下一时刻贝叶斯公式中的先验概率，也就是先验概率由 $p(\mathbf{\theta }|{\mathcal{D}}^1)$ 更新为 $p(\mathbf{\theta }|{\mathcal{D}}^2)$

$\cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots$

$(4)$ 第 $n$ 时刻

\qquad 从 $p(\mathbf{\theta }|{\mathcal{D}}^{n-1})$ 中抽取一个 $\mathbf{\theta }$ 值用于确定概率密度函数 $p(\mathbf{x}|\mathbf{\theta })$，同时观测第 $n$ 个数据 ${\mathbf{x}}_n$，观测数据集更新为 D n = { x 1 , ⋯   , x n } \mathcal D^{n}=\{\boldsymbol x_1,\cdots,\boldsymbol x_n\} Dn={x1​,⋯,xn​}，然后计算数据 ${\mathbf{x}}_n$ 的似然函数 $p({\mathbf{x}}_n|\mathbf{\theta })$，最后使用贝叶斯公式计算后验概率：

$p(\mathbf{\theta }|{\mathcal{D}}^n)=\frac{p({\mathbf{x}}_n|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^{n-1})}{\int p({\mathbf{x}}_n|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^{n-1})d\mathbf{\theta }}$

该时刻的后验概率 $p(\mathbf{\theta }|{\mathcal{D}}^n)$ 将作为下一时刻贝叶斯公式中的先验概率，也就是先验概率由 $p(\mathbf{\theta }|{\mathcal{D}}^{n-1})$ 更新为 $p(\mathbf{\theta }|{\mathcal{D}}^n)$

由于样本抽取是独立的，因此 $p({\mathcal{D}}^n|\mathbf{\theta })=p({\mathbf{x}}_n|\mathbf{\theta })p({\mathcal{D}}^{n-1}|\mathbf{\theta })$
　
可得到： $p(\mathbf{\theta }|{\mathcal{D}}^n)=\frac{p({\mathcal{D}}^n|\mathbf{\theta })p(\mathbf{\theta })}{\int p({\mathcal{D}}^n|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }}=\frac{p({\mathbf{x}}_n|\mathbf{\theta })p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })}{\int p({\mathbf{x}}_n|\mathbf{\theta })p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }}$
　
又：　　 $p(\mathbf{\theta }|{\mathcal{D}}^{n-1})=\frac{p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })}{\int p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }}⟹p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })=p(\mathbf{\theta }|{\mathcal{D}}^{n-1})\int p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }$
　
代入 $p(\mathbf{\theta }|{\mathcal{D}}^n)$中，可得到：
　
$\begin{array}{rl}p(\mathbf{\theta }|{\mathcal{D}}^n)& =\frac{p({\mathbf{x}}_n|\mathbf{\theta })p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })}{\int p({\mathbf{x}}_n|\mathbf{\theta })p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }}\\ & =\frac{p({\mathbf{x}}_n|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^{n-1})\int p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }}{\int p({\mathbf{x}}_n|\mathbf{\theta })p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }}\\ & =\frac{p({\mathbf{x}}_n|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^{n-1})}{\int p({\mathbf{x}}_n|\mathbf{\theta })\frac{p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })}{\int p({\mathcal{D}}^{n-1}|\mathbf{\theta })p(\mathbf{\theta })d\mathbf{\theta }}d\mathbf{\theta }}\\ & =\frac{p({\mathbf{x}}_n|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^{n-1})}{\int p({\mathbf{x}}_n|\mathbf{\theta })p(\mathbf{\theta }|{\mathcal{D}}^{n-1})d\mathbf{\theta }}\end{array}$

\qquad 如此反复下去，能够产生一系列的概率密度函数 $p(\mathbf{\theta }|{\mathbf{x}}_1),p(\mathbf{\theta }|{\mathbf{x}}_1,{\mathbf{x}}_2),\cdots ,p(\mathbf{\theta }|{\mathcal{D}}^n),\cdots$，这一过程被称为“参数估计的递归的贝叶斯方法”。

\qquad
例3：线性回归的贝叶斯学习（《Pattern Recognition and Machine Learning》3.3.1节）

\qquad 假设观测数据满足 $t_i=y(x_i,\mathbf{w})+\epsilon ,\epsilon \sim \mathcal{N}(0,{\beta }^{-1})$，其中 $\beta =\frac{1}{{\sigma }^2}=25$

\qquad 线性回归模型 $y(x,\mathbf{w})=w_0+w_{1}x,\mathbf{w}=(w_0,w_1{)}^T$，假定真实参数值为 $\mathbf{w}=(-0.3,0.5{)}^T$

\qquad 对于某个点 $x_i$，关于目标值 $t_i$ 的似然函数为 $p(t_i|x_i,\mathbf{w})=\mathcal{N}(t|y(x_i,\mathbf{w}),{\beta }^{-1})$

关于所有目标值 $\mathbf{t}$ 的似然函数为 $p(\mathbf{t}|\mathbf{w},\mathbf{X})=p(t_1,\cdots ,t_N|\mathbf{X},\mathbf{w})=\prod _{i=1}^{N}p(t_i|x_i,\mathbf{w})$

\qquad 参数向量 $\mathbf{w}$ 的先验概率为 $p(\mathbf{w})=\mathcal{N}(\mathbf{w}|{\mathbf{m}}_0,{\mathbf{S}}_0)$

\qquad 经过计算可以推导出后验概率 $p(\mathbf{w}|\mathbf{t})=\mathcal{N}(\mathbf{w}|{\mathbf{m}}_N,{\mathbf{S}}_N)$

\qquad 其中，${\mathbf{m}}_N={\mathbf{S}}_N({\mathbf{S}}_0^{-1}{\mathbf{m}}_0+\beta {\Phi }^T\mathbf{t})$
\qquad 　　　${\mathbf{S}}_N^{-1}={\mathbf{S}}_0^{-1}+\beta {\Phi }^T\Phi$

本例可参考《线性模型参数求解的最大似然估计》或者《Pattern Recognition and Machine Learning》3.3.1节

第１行 —— 中图：先验概率 $p(\mathbf{w})$，真实值 $\mathbf{w}=(-0.3,0.5{)}^T$ 位于 $‘+’$ 处
　　　 —— 右图：从 $p(\mathbf{w})$ 中抽取的、若干个不同的$\mathbf{w}$值所画出不同的直线杂乱无章（由于先验信息 $p(\mathbf{w})$ 分布广泛、不确定性大）
　
第２行 —— 右图：观测到第1个数据 $(x_1,t_1)$（蓝色圆圈）
　　　 —— 左图：计算出的似然值 $p(t_1|x_1,\mathbf{w})$ 关于 $\mathbf{w}$ 的图像
　　　 —— 中图：$p(\mathbf{w})$ 作为先验概率，计算出的后验概率 $p(\mathbf{w}|t^1)$ 开始出现尖峰，关于 $\mathbf{w}$ 的不确定性变小
　　　 —— 右图：从 $p(\mathbf{w}|t^1)$ 中抽取的、若干个不同的$\mathbf{w}$值所画出不同的直线基本都靠近点 $(x_1,t_1)$ 穿过，不再完全杂乱无章
　
第３行 —— 右图：观测到第2个数据 $(x_2,t_2)$（蓝色圆圈）
　　　 —— 左图：计算出的似然值 $p(t_2|x_2,\mathbf{w})$ 关于 $\mathbf{w}$ 的图像
　　　 —— 中图：$p(\mathbf{w}|t^1)$ 作为先验概率，计算出的后验概率 $p(\mathbf{w}|t^1,t^2)$ 的尖峰更明显，关于 $\mathbf{w}$ 的不确定性变得更小
　　　 —— 右图：从 $p(\mathbf{w}|t^1,t^2)$ 中抽取的、若干个不同的$\mathbf{w}$值所画出不同的直线的方向趋于一致
　
第４行 —— 右图：观测到第n个数据 $(x_n,t_n)$（蓝色圆圈）
　　　 —— 左图：计算出的似然值 $p(t_n|x_n,\mathbf{w})$ 关于 $\mathbf{w}$ 的图像
　　　 —— 中图：后验概率 $p(\mathbf{w}|t^1,t^2,\cdots ,t^n)$ 的尖峰非常尖锐，基本消除了 $\mathbf{w}$ 的不确定性
　　　 —— 右图：从 $p(\mathbf{w}|t^1,t^2,\cdots ,t^n)$ 中抽取的、若干个不同的$\mathbf{w}$值所画出不同的直线基本上一致，说明这些不同的$\mathbf{w}$值已经基本相同