模式识别-Ch3-参数估计
Ch3 参数估计-最大似然和贝叶斯参数估计
不考EM、HMM
基本概念
贝叶斯分类器: 已知类先验概率\(P(w_j)\)和类条件概率密度\(p(\mathbf{x}\vert w_j)\)，按某决策规则确定判别函数和决策面。
但类先验概率和类条件概率密度在实际中往往是未知的。
因此，我们要换一种处理问题的方式：“从样本出发来设计分类器”。根据设计方法，可以将分类器分为两类：
- 估计类先验概率和类条件概率密度函数（产生式方法）
- 直接估计类后验概率或判别函数（判别式方法）
| 参数估计 | 非参数估计 |
|---|---|
| 样本所属的类条件概率密度函数的形式已知，而概率密度函数的参数是未知的。 | 样本所属的类条件概率密度函数的形式和参数都是未知的。 |
| 目标是由已知类别的样本集估计概率密度函数的参数。 | 目标是由已知类别的样本集估计类条件概率密度函数本身。 |
| 例如，知道样本所属总体为正态分布，而正态分布的参数未知\(p(\mathbf{x}\vert\mu,\sigma)=\frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{1}{2}\left(\frac{\mathbf{x}-\mu}{\sigma}\right)^2\right)\) | --- |
| 基本概念 | 说明 | 例子 |
|---|---|---|
| 统计量 | 样本中包含总体的信息，我们希望通过样本集将有关信息估计出来。根据不同要求构造出有关样本的某种函数，在统计学中称为统计量\(d(\mathbf{x}_1,\mathbf{x}_2,\ldots,\mathbf{x}_n)\)。 | \(\mu=\frac{1}{n}\sum_{i = 1}^{n}\mathbf{x}_i\) |
| 参数空间 | 将未知待估计参数记为\(\theta\)，参数\(\theta\)的全部允许取值集合构成参数空间，记为\(\Theta\)。 | --- |
| 点估计 | 点估计问题就是构造一个统计量\(d(\mathbf{x}_1,\mathbf{x}_2,\ldots,\mathbf{x}_n)\)作为参数\(\theta\)的估计\(\hat{\theta}\)。 | 常用的均值估计：\(\hat{\mu}=\frac{1}{n}\sum_{i = 1}^{n} \mathbf{x}_i\) |
| 区间估计 | 与点估计不同，区间估计要求采用\((d_1,d_2)\)作为参数\(\theta\)可能取值范围的一种估计。这个区间称为置信区间。这类估计问题称为区间估计。 | --- |
最大似然估计
基本假设
- 独立同分布假设：每类样本均是从类条件概率密度\(p(x\vert w_j)\)中独立抽取出来的。
- \(p(x\vert
w_j)\)具有确定的函数形式，只是其中的参数\(\theta\)未知：
- 比如，当\(\mathbf x\)服从一维正态分布\(N(\mu,\sigma^2)\)，未知的参数为\(\theta = [\mu,\sigma]^T\)，为一个二维向量。
- 各类样本只包含本类的分布信息：即不同类别的参数是独立的。可以分别处理\(c\)个独立问题。
基本原理
已知随机抽取的\(n\)个样本(观测值)，最合理的参数估计应该是使得从该模型中能抽取这\(n\)个样本的概率最大。
直观想法：一个随机试验如有若干个可能的结果：A，B，C，…。若仅作一次试验，结果A出现，则认为试验条件(模型参数)对A出现有利，也即A出现的概率很大。
一般地，事件A发生的概率与参数\(\theta\)相关，A发生的概率记为\(P(A\vert \theta)\)，则\(\theta\)的估计应该使上述概率达到最大，这样的\(\theta\)估计意义称为极大似然估计。
设样本集包含\(n\)个样本\(D = \{\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n\}\)，这些样本是从概率密度函数\(p(x\vert \theta)\)中独立抽取的，则获得\(n\)个样本的联合概率为： \[ l(\theta)=P(D|\theta)=P(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n|\theta)=\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta) \] \(l(\theta)\)是\(\theta\)的函数，描述了在不同参数取值下取得当前样本集的可能性。
\(l(\theta)\)被称为参数\(\theta\)相对于样本集\(D\)的似然函数: 似然函数给出了从总体中抽出\(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n\)这\(n\)个样本的概率。
方法描述
令\(l(\theta)\)为样本集\(D\)的似然函数，\(D = {\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n}\)。
如果\(\theta\)是参数空间\(\Theta\)中能使\(l(\theta)\)极大化的\(\theta\)值，那么\(\theta\)就是\(\theta\)的最大似然估计量，即\(\hat{\theta}=\arg\max_{\theta\in\Theta}l(\theta)\)
为计算方便，通常采用对数似然函数： \[ H(\theta)=\ln(l(\theta))=\ln(\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta))=\sum_{i = 1}^{n}\ln(p(\mathbf{x}_i|\theta))\\ \arg \max_{\theta\in\Theta}l(\theta)=\arg \max_{\theta\in\Theta}H(\theta) \]
问题求解
\[ H(\theta)=\ln(l(\theta))=\ln(\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta))=\sum_{i = 1}^{n}\ln(p(\mathbf{x}_i|\theta))\\ \hat{\theta}=\text{arg max}_{\theta\in\Theta}l(\theta) \]
当\(l(\theta)\)可微时：\(\frac{\partial l(\theta)}{\partial\theta}=0,\text{ or }\frac{\partial H(\theta)}{\partial\theta}=0\)对于多维情形\(\theta = [\theta_1,\theta_2,\cdots,\theta_m]^T\)，梯度向量为零： \[ \nabla_{\theta}(l(\theta))=\frac{\partial l(\theta)}{\partial\theta}=\left[\frac{\partial l(\theta)}{\partial\theta_1},\frac{\partial l(\theta)}{\partial\theta_2},\cdots,\frac{\partial l(\theta)}{\partial\theta_m}\right]^T = 0 \] 用梯度上升法求解\(\theta^{t + 1}=\theta^{t}+\eta\frac{\partial l(\theta)}{\partial\theta}\)问题求解 : \[ H(\theta)=\ln(l(\theta))=\ln(\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta))=\sum_{i = 1}^{n}\ln(p(\mathbf{x}_i|\theta))\\ \hat{\theta}=\arg\max_{\theta\in\Theta}l(\theta) \]
| 当\(l(\theta)\)是可微凹函数时 | 当\(l(\theta)\)是一般可微函数时 |
|---|---|
| \(\arg\max_{\theta\in\Theta}l(\theta)\Leftrightarrow\frac{\partial l(\theta)}{\partial\theta}=0\) | \(\arg\max_{\theta\in\Theta}l(\theta)\Rightarrow\frac{\partial l(\theta)}{\partial\theta}=0\) |
| 梯度等于0是最优解的充要条件 | 梯度等于0是最优解的必要条件 |
在高维空间中，寻找全局最优解是极困难的，通常满足于局部最优解。
例1：黑白球
一个袋子里装有白球与黑球，但是不知道它们之间的比例。现有放回地抽取10次，结果获得8次黑球2次白球，估计袋子中的黑球的比例。
最大似然估计法：设抽到黑球的概率为p，取到8次黑球2次白球的概率为： \[ l(p)=\begin{pmatrix}10\\8\end{pmatrix}p^8(1-p)^2\\ \frac{\partial l(p)}{\partial p}=\begin{pmatrix}10\\8\end{pmatrix}\frac{\partial p^8(1-p)^2}{\partial p}=\begin{pmatrix}10\\8\end{pmatrix}(10p^9-18p^8+8p^7)=0\Rightarrow\hat p=0.8\\ \frac{\partial \ln l(p)}{\partial p}=\frac{\partial\ln\begin{pmatrix}10\\8\end{pmatrix}+8\ln p+2\ln(1-p)}{\partial p}=\frac 8p-\frac2{1-p}=0\Rightarrow\hat p=0.8 \]
例2：高斯分布下的最大似然估计
总的来说，这边就是一个很简单的多元正态分布求解MLE，非常简单。在多元统计分析里面也算过了。
预备公式： \[ \begin{align}&tr(A)=\sum_{i = 1}^{d}A_{ii},\text{ where }A=(A_{ij})\in R^{d\times d}\\ &s = tr(s),\ s\text{是标量}\Rightarrow x^TAx = tr(x^TAx),\ x\in R^d,A\in R^{d\times d}\\ &\frac{\partial\left|\Sigma\right|}{\partial\Sigma}=\left|\Sigma\right|(\Sigma^{- 1}),\text{ if }\Sigma\text{对称}\\ &\frac{\partial tr(A\Sigma^{-1}B)}{\partial\Sigma}=(-\Sigma^{-1}BA\Sigma^{-1})^T\Rightarrow\frac{\partial(x^T\Sigma^{-1}x)}{\partial\Sigma}=-\Sigma^{-1}xx^T\Sigma^{-1},\text{ if }\Sigma\text{对称} \end{align} \]
已知\(x = [\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_d]^T\in R^d\)，\(\mu = [\mu_1,\mu_2,\cdots,\mu_d]^T\in R^d\)，\(\Sigma\in R^{d\times d}\)
概率密度函数为： \[ p(x|\mu,\Sigma)=\frac{1}{(2\pi)^{d/2}|\Sigma|^{1/2}}\exp\left(-\frac{1}{2}(x - \mu)^T\Sigma^{-1}(x - \mu)\right) \] 对数似然函数为： \[ \ln p(\mathbf{x}_i|\mu,\Sigma)=-\frac{1}{2}\ln[(2\pi)^d|\Sigma|]-\frac{1}{2}(\mathbf{x}_i - \mu)^T\Sigma^{-1}(\mathbf{x}_i - \mu),\text{ for }i = 1,2,\cdots,n \] 似然函数求导： \[ \nabla_{\theta}(H(\theta))=\sum_{i = 1}^{n}\nabla_{\theta}\ln p(\mathbf{x}_i|\theta)=0,\text{ (where }\theta=\mu,\text{ and (or) }\Sigma) \]
估计\(\mu\)
\[ \frac{\partial(x^TAx)}{\partial x}=(A+A^T)x\\ \frac{\partial(x^TA)}{\partial x}=A\\ \frac{\partial(Ax)}{\partial x}=A^T \]
\[ \ln p(\mathbf{x}_i|\mu,\Sigma)=-\frac{1}{2}\ln[(2\pi)^d|\Sigma|]-\frac{1}{2}(\mathbf{x}_i - \mu)^T\Sigma^{-1}(\mathbf{x}_i - \mu),\text{ for }i = 1,2,\cdots,n\\ \frac{\partial H(\mu)}{\partial\mu}=\sum_{i = 1}^{n}\frac{\partial\ln p(\mathbf{x}_i|\mu,\Sigma)}{\partial\mu}=\sum_{i = 1}^{n}\Sigma^{-1}(\mathbf{x}_i - \mu)\\ \sum_{i = 1}^{n}\Sigma^{-1}(\mathbf{x}_i - \mu)=0\Rightarrow\sum_{i = 1}^{n}\mathbf{x}_i - n\mu=0\Rightarrow\hat{\mu}=\frac{1}{n}\sum_{i = 1}^{n}\mathbf{x}_i \]
估计\(\Sigma\)
\[ \ln p(\mathbf{x}_i|\mu,\Sigma)=-\frac{1}{2}\ln[(2\pi)^d|\Sigma|]-\frac{1}{2}(\mathbf{x}_i-\mu)^T\Sigma^{-1}(\mathbf{x}_i-\mu),\text{ for }i = 1,2,\cdots,n\\ \frac{\partial H(\mu,\Sigma)}{\partial\Sigma}=\sum_{i = 1}^{n}\frac{\partial\ln p(\mathbf{x}_i|\mu,\Sigma)}{\partial\Sigma}\\ =\sum_{i = 1}^{n}\left(-\frac{n}{2}\frac{1}{|\Sigma|}\frac{\partial|\Sigma|}{\partial\Sigma}-\frac{1}{2}\frac{\partial(\mathbf{x}_i-\mu)^T\Sigma^{-1}(\mathbf{x}_i-\mu)}{\partial\Sigma}\right)\\ =-\frac{n}{2}\Sigma^{-1}+\frac{1}{2}\Sigma^{-1}\sum_{i = 1}^{n}(\mathbf{x}_i-\mu)(\mathbf{x}_i-\mu)^T\Sigma^{-1}=0 \]
\[ \begin{align} \Sigma^{-1}\left(\sum_{i = 1}^{n}(\mathbf{x}_i-\mu)(\mathbf{x}_i-\mu)^T\Sigma^{-1}-n\right)&=0\\ \sum_{i = 1}^{n}(\mathbf{x}_i-\mu)(\mathbf{x}_i-\mu)^T\Sigma^{-1}&=n\\ \sum_{i = 1}^{n}(\mathbf{x}_i-\mu)(\mathbf{x}_i-\mu)^T&=n\Sigma\\ \frac{1}{n}\sum_{i = 1}^{n}(\mathbf{x}_i-\mu)(\mathbf{x}_i-\mu)^T&=\hat{\Sigma} \end{align} \]
\(\mu,\Sigma\)均未知的情况下，可以使用极大似然估计得到\(\mu,\Sigma\)的估计量： \[ \hat{\mu}=\frac{1}{n}\sum_{i = 1}^{n}\mathbf{x}_i,\hat{\Sigma}=\frac{1}{n}\sum_{i = 1}^{n}(\mathbf{x}_i - \mu)(\mathbf{x}_i - \mu)^T \] 一维情况： \[ \text{max}\sum_{i = 1}^{n}\ln p(\mathbf{x}_i|\mu,\sigma)\\ p(x|\mu,\sigma)=\frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{1}{2}\frac{(x - \mu)^2}{\sigma^2}\right)\\ \hat{\mu}=\frac{1}{n}\sum_{i = 1}^{n}\mathbf{x}_i,\hat{\sigma}^2=\frac{1}{n}\sum_{i = 1}^{n}(\mathbf{x}_i - \hat{\mu})^2 \]
贝叶斯估计
贝叶斯估计是概率密度估计中另一类主要的参数估计方法。其结果在很多情况下与最大似然法十分相似，但是，两种方法对问题的处理视角是不一样的。
| 贝叶斯估计 | 最大似然估计 |
|---|---|
| 将待估计的参数视为一个随机变量，其中的一个核心任务是根据观测数据对参数的分布进行估计。 | 将待估计的参数当作未知但固定的变量，其任务是根据观测数据估计其在参数空间中的取值。 |
| \(p(x\vert D)\sim N(\mu_{n},\sigma^{2}+\sigma_{n}^{2})\\\mu_{n}=\frac{n\sigma_{0}^{2}}{\sigma_{0}^{2}+\sigma^{2}}\hat{\mu}_{n}+\frac{\sigma^{2}}{\sigma_{0}^{2}+\sigma^{2}}\mu_{0}\\\sigma_{n}^{2}=\frac{\sigma_{0}^{2}\sigma^{2}}{n\sigma_{0}^{2}+\sigma^{2}}\) | \(p(x\vert D)\sim N(\hat{\mu}_{n},\sigma^{2})\\\hat{\mu}_{n}=\frac{1}{n}\sum_{i = 1}^{n}\mathbf{x}_{i}\\\) |
上面公式给出的是一维下估计。
基本方法
参数先验分布\(p(\theta)\)：是在没有任何数据时，有关参数\(\theta\)的分布情况（根据领域知识或经验）
给定样本集\(D = \{\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n\}\)，数据独立采样，且服从数据分布：(数据是互相独立的) \[ p(D|\theta)=p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n|\theta)=\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta) \] 利用贝叶斯公式计算参数的后验分布\(p(\theta\vert D)\)：\(p(\theta\vert D)\)中融合了先验知识和数据信息。 \[ p(\theta|D)=\frac{p(D|\theta)p(\theta)}{p(D)} \] \(p(D)\)是与参数无关的归一化因子，根据全概率公式（连续）： \[ p(D)=\sum_{\theta}p(D|\theta)p(\theta)\\ p(D)=\int_{\theta}p(D|\theta)p(\theta)d\theta\\ p(D|\theta)\Rightarrow p(x|\mu,\sigma)=\frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{1}{2}\frac{(x - \mu)^2}{\sigma^2}\right) \] 对于一般情况，计算\(p(D)\)十分困难
可得贝叶斯参数估计中的后验概率密度函数： \[ p(\theta|D)=\frac{p(D|\theta)p(\theta)}{\int_{\theta}p(D|\theta)p(\theta)d\theta}=\frac{\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta)p(\theta)}{\int_{\theta}\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta)p(\theta)d\theta}=\alpha\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta)p(\theta)\\ \alpha=\frac 1{\int_{\theta}\prod_{i = 1}^{n}p(\mathbf{x}_i|\theta)p(\theta)d\theta} \]
Q: 如何使用\(p(\theta\vert D)\)获得关于数据的分布？
得到\(p(\theta\vert D)\)只是获得了关于参数\(\theta\)的后验分布，并没有像最大似然估计那样获得参数\(\theta\)的具体取值。
方法1 方法2 方法3 对\(p(\theta\vert D)\)采样，计算平均值 最大后验估计(Maximum A Posteriori estimation, MAP) 后验数据分布(完整的贝叶斯方法) \(\hat{\theta}=\frac{1}{M}\sum_{i = 1}^{M}\theta_i,\theta_i\sim p(\theta\vert D),i = 1,\cdots,M\) \(\begin{align}&\hat{\theta}=\arg\max p(\theta\vert D)\\\Leftrightarrow&\hat{\theta}=\arg\max p(D\vert \theta)p(\theta)\\\Leftrightarrow&\hat{\theta}=\arg\max\ln p(D\vert \theta)+\ln p(\theta)\end{align}\) \(p(x\vert \mu,\Sigma)=\frac{1}{(2\pi)^{d/2}\vert \Sigma\vert ^{1/2}}\exp\left(-\frac{1}{2}(x-\mu)^T\Sigma^{-1}(x-\mu)\right)\) PR/ML方法中普遍使用的L2正则，等价于假设参数服从\(N(0,I)\)
后验数据分布
最终目的：根据\(D\)中的样本来估计概率密度函数\(p(x\vert D)\)。
比如，假定观测样本服从正态分布\(p(x\vert \mu,\Sigma)\)，给定\(D\)，可以估计得到具体的\(\mu\)和\(\Sigma\)的取值，代入如下公式可得关于样本的密度分布函数： \[ p(x\vert \mu,\Sigma)=\frac{1}{(2\pi)^{d/2}\vert \Sigma\vert ^{1/2}}\exp\left(-\frac{1}{2}(x - \mu)^T\Sigma^{-1}(x - \mu)\right) \] 但现在获得了有关\(\theta\)的后验估计\(p(\theta\vert D)\)，如何估计\(p(x\vert D)\)？
考虑全概率公式和边际分布： \[ \begin{align}p(x\vert D)&=\int_{\theta}p(x,\theta\vert D)d\theta\\ &=\int_{\theta}p(x\vert \theta)p(\theta\vert D)d\theta \end{align} \] - \(p(x\vert \theta)=p(x\vert \theta,D)\): 在给定参数\(\theta\)时，样本分布与训练集\(D\)无关 - \(\int_{\theta}p(x\vert \theta)p(\theta\vert D)d\theta\): 不同参数的密度函数的加权平均
积分通常很难计算，使用蒙特卡洛近似方法： 是\(M\)个不同参数的密度函数的平均。 \[ \hat{p}(x\vert D)=\frac{1}{M}\sum_{i = 1}^{M}p(x\vert \theta_i),\theta_i\sim p(\theta\vert D),i = 1,\cdots,M \]
一维情形：假定\(X\sim N(\mu,\sigma^2)\)且仅\(\mu\)未知
假定参数\(\mu\)的先验概率也服从正态分布：\(\mu\sim N(\mu_0,\sigma_0^2)\) \[ p(x\vert \mu)=N(\mu,\sigma^2),\ p(\mu)=N(\mu_0,\sigma_0^2) \] 第一个任务：给定样本集\(D\)，在上述条件下，估计关于参数的后验分布\(p(\mu\vert D)\)。
回顾我们前面得到的公式： \[ p(\theta\vert D)=\frac{\prod_{i = 1}^{n}p(\mathbf{x}_i\vert \theta)p(\theta)}{\int_{\theta}\prod_{i = 1}^{n}p(\mathbf{x}_i\vert \theta)p(\theta)d\theta}=\alpha\prod_{i = 1}^{n}p(\mathbf{x}_i\vert \theta)p(\theta)\\ \] （应用后验估计） \[ \begin{align} p(\mu\vert D)&=\alpha\prod_{i = 1}^{n}p(\mathbf{x}_i\vert \mu)p(\mu)\\ &=\alpha\prod_{i = 1}^{n}\frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{1}{2}\frac{(\mathbf{x}_i - \mu)^2}{\sigma^2}\right)\frac{1}{\sqrt{2\pi}\sigma_0}\exp\left(-\frac{1}{2}\frac{(\mu - \mu_0)^2}{\sigma_0^2}\right)\\ &=\alpha'\prod_{i = 1}^{n}\exp\left\{-\frac 1 2 \sum^n_{i=1}\frac{(\mathbf{x}_i-\mu)^2}{\sigma^2}-\frac n2\frac{(\mu-\mu_0)^2}{\sigma_0^2}\right\}\\ &=\alpha'\exp\left\{-\frac{1}{2}\left[\left(\frac{1}{\sigma^2}+\frac{1}{\sigma_0^2}\right)\mu^2 - 2\left(\frac{1}{\sigma^2}\sum_{i = 1}^{n}\mathbf{x}_i+\frac{\mu_0}{\sigma_0^2}\right)\mu\right]\right\}\\ &=\alpha''\exp\left\{-\frac{1}{2}\left[\left(\frac{\sigma_0^2+\sigma^2}{\sigma^2\sigma_0^2}\right)\mu^2 - 2\left(\frac{1}{\sigma^2}\sum_{i = 1}^{n}\mathbf{x}_i+\frac{\mu_0}{\sigma_0^2}\right)\mu\right]\right\}\end{align} \]
一维后验分布的性质
- \(p(\mu\vert D)\)是关于\(\mu\)的二次函数的\(\text{exp}\)函数，因此，也是一个正态分布密度函数。
- \(p(\mu\vert D)\)被称为再生密度（reproducing density），因为对于任意数量的训练样本，当样本数量\(n\)增加时，\(p(\mu\vert D)\)仍然保持正态分布。
由于\(p(\mu\vert D)\)是一个正态密度函数，我们可以将其改写为如下形式： \[ p(\mu\vert D)\sim N(\mu_{n},\sigma_{n}^{2})=\frac{1}{\sqrt{2\pi\sigma_{n}^{2}}}\exp\left(-\frac{1}{2}\frac{(\mu - \mu_{n})^{2}}{\sigma_{n}^{2}}\right) \] 同时，我们也得到其公式为 \[ p(\mu\vert D)=\alpha^{\prime}\exp\left\{-\frac{1}{2}\left[\left(\frac{n}{\sigma^{2}}+\frac{1}{\sigma_{0}^{2}}\right)\mu^{2}-2\left(\frac{1}{\sigma^{2}}\sum_{i = 1}^{n}\mathbf{x}_{i}+\frac{\mu_{0}}{\sigma_{0}^{2}}\right)\mu\right]\right\}\\ \frac{1}{\sigma_{n}^{2}}=\frac{n}{\sigma^{2}}+\frac{1}{\sigma_{0}^{2}},\quad\frac{\mu_{n}}{\sigma^2_n}=\frac{n}{\sigma^{2}}\bar{\mu}_{n}+\frac{\mu_{0}}{\sigma_{0}^{2}},\quad \bar{\mu}_{n}=\frac{1}{n}\sum_{i = 1}^{n}\mathbf{x}_{i} \] 进一步可解得： \[ \mu_{n}=\frac{n\sigma_{0}^{2}}{n\sigma_{0}^{2}+\sigma^{2}}\bar{\mu}_{n}+\frac{\sigma^{2}}{n\sigma_{0}^{2}+\sigma^{2}}\mu_{0},\quad\sigma_{n}^{2}=\frac{\sigma^{2}\sigma_{0}^{2}}{n\sigma_{0}^{2}+\sigma^{2}} \] 这些方程展示了先验信息如何与样本中的经验信息相结合以获得后验密度\(p(\mu\vert D)\)。
- \(\mu_{n}\)：代表在获得\(n\)个样本后对\(\mu\)的最佳猜测。
- \(\sigma_{n}^{2}\)：衡量对\(\mu\)猜测的不确定性。
- 因为\(\sigma_{n}^{2}\)随\(n\)单调递减，每增加一个观测值都将有助于减少我们对\(\mu\)真实值的不确定性。（这种先验起到了平滑的效果，导致了更加鲁棒的估计）
后验分布的变化趋势：因为\((\sigma_{n})^{2}\)随\(n\)单调递减，每增加一个观测值都将有助于减少我们对\(\mu\)真实值的不确定性。随着\(n\)的增加，\(p(\mu\vert D)\)变得越来越尖锐，当\(n\)趋于无穷大时，趋近于狄拉克δ函数（Dirac delta function）。
现在，我们希望获得后验数据分布 ：
可以将\(p(\mathbf{x}\vert D)\)视为服从正态分布\(N(\mu_n,\sigma^2+\sigma^2_n)\)
多元情形：高维
已知条件是： \[ \begin{align}p(\mathbf x\vert \mathbf \mu)&\sim N(\mathbf \mu,\Sigma),p(\mu)\sim N(\mu_{0},\Sigma_{0})\\ p(\theta\vert D)&=\alpha\prod_{i = 1}^{n}p(\mathbf x_{i}\vert \theta)p(\theta)\\ &=\alpha^{\prime}\exp\left\{-\frac{1}{2}\mu^{T}(n\Sigma^{- 1}+\Sigma_{0}^{-1})\mu - 2\mu^{T}(\Sigma^{-1}\sum_{i = 1}^{n}\mathbf x_{i}+\Sigma_{0}^{-1}\mu_{0})\right\}\\ &=\alpha^{\prime\prime}\exp\left\{-\frac{1}{2}(\mu - \mu_{n})^{T}\Sigma_{n}^{-1}(\mu - \mu_{n})\right\}\end{align} \] 参照上面一维的情况，可以推出： \[ \begin{align} p(\theta\vert D)&=\alpha^{\prime\prime}\exp\left\{-\frac{1}{2}(\mu - \mu_{n})^{T}\Sigma_{n}^{-1}(\mu - \mu_{n})\right\}\Rightarrow p(\theta\vert D)\sim N(\mu_{n},\Sigma_{n})\\ \Rightarrow\Sigma_{n}^{-1}&=n\Sigma^{-1}+\Sigma_{0}^{-1},\quad \Sigma_{n}^{-1}\mu_{n}=n\Sigma^{-1}\hat{\mu}_{n}+\Sigma_{0}^{-1}\mu_{0},\quad \hat{\mu}_{n}=\frac{1}{n}\sum_{i = 1}^{n}\mathbf x_{i}\\ \mu_{n}&=\Sigma_{0}(\Sigma_{0}+n^{-1}\Sigma)^{-1}\hat{\mu}_{n}+(\Sigma_{0}+n^{-1}\Sigma)^{-1}\Sigma_{0}\mu_{0}\\ \Sigma_{n}&=\Sigma_{0}(\Sigma_{0}+n^{-1}\Sigma)^{-1}\frac{1}{n}\Sigma \end{align} \] > \((A^{-1}+B^{-1})^{-1}=A(A+B)^{-1}B=B(A+B)^{-1}A\)
数据后验分布服从正态分布： \[ p(\mathbf x\vert D)=\int_{\mu}p(\mathbf x\vert \mu)p(\mu\vert D)d\mu\sim N(\mu_{n},\Sigma+\Sigma_{n}) \]
贝叶斯学习
一般情形下的贝叶斯估计（总结）
基本假设：
密度\(p(\mathbf{x}\vert \theta)\)的形式已知，但参数向量的值未知。
关于\(\theta\)的初始知识包含在已知的先验密度\(p(\theta)\)中。
关于\(\theta\)的其余知识包含在根据未知概率密度\(p(\mathbf{x})\)独立抽取的\(n\)个样本\(x_{1},x_{2},\cdots,x_{n}\)的集合\(D\)中。
基本问题：计算关于参数\(\theta\)的后验密度\(p(\theta\vert D)\)和关于数据的后验密度\(p(\mathbf{x}\vert D)\)。 \[ p(\theta\vert D)=\frac{p(D\vert \theta)p(\theta)}{\int p(D\vert \theta)p(\theta)d\theta},\quad p(D\vert \theta)=P(x_{1},x_{2},\cdots,x_{n}\vert \theta)=\prod_{i = 1}^{n}p(x_{i}\vert \theta)\\ p(\mathbf{x}\vert D)=\int_{\theta}p(\mathbf{x}\vert \theta)p(\theta\vert D)d\theta \]
遇到的困难:
除了一些特殊的分布（共轭分布）之外，对于一般情形，积分很难计算： \[ p(\theta\vert D)=\frac{p(D\vert \theta)p(\theta)}{\int p(D\vert \theta)p(\theta)d\theta},\quad p(\mathbf{x}\vert \theta)=\int_{\theta}p(\mathbf{x}\vert \theta)p(\theta\vert D)d\theta \]
参数先验\(p(\theta)\)怎么选取？对结果有何影响？
p(θ) 的选择对结果有直接影响。先验分布过于强烈可能会导致数据驱动的结果被先验主导，而过于弱的先验分布可能导致计算结果不稳定。
给定\(D\)，我们真的能通过\(p(\mathbf{x}\vert D)\)将\(p(\mathbf{x})\)估计得很好吗？或者说，随着\(D\)中样本的增多，\(p(\mathbf{x}\vert D)\)收敛于\(p(\mathbf{x})\)吗？
根据贝叶斯学习的性质，当数据量\(n \to \infty\)时，后验分布\(p(\theta\vert D)\)会集中在最大似然估计值附近，即：\(p(\theta\vert D) \to \delta(\theta-\theta_{\text{MLE}})\)这意味着后验分布的方差会逐渐缩小，预测分布 p(∣D)p(D) 也会趋近于真实分布。
贝叶斯学习的迭代计算公式:
记\(D^{n}=\{x_{1},x_{2},\cdots,x_{n}\}\)，由于样本是独立选样，则： \[ p(D^{n}\vert \theta)=p(x_{n}\vert \theta)p(D^{n - 1}\vert \theta)=p(x_{n}\vert \theta)p(x_{n - 1}\vert \theta)p(D^{n - 2}\vert \theta)=\cdots \]
于是有如下迭代公式： \[ \begin{align} p(\theta|D^{n})&=\frac{p(D^{n}|\theta)p(\theta)}{\int p(D^{n}|\theta)p(\theta)d\theta}=\frac{p(x_{n}|\theta)p(D^{n - 1}|\theta)p(\theta)}{\int p(x_{n}|\theta)p(D^{n - 1}|\theta)p(\theta)d\theta} \\ &=\frac{p(x_{n}|\theta)}{\int p(x_{n}|\theta)\frac{p(D^{n - 1}|\theta)p(\theta)}{\int p(D^{n - 1}|\theta)p(\theta)d\theta}d\theta}=\frac{p(x_{n}|\theta)}{\int p(x_{n}|\theta)p(\theta|D^{n - 1})d\theta} \\ &=\frac{p(x_{n}|\theta)p(\theta|D^{n - 1})}{\int p(x_{n}|\theta)p(\theta|D^{n - 1})d\theta} \\ p(\theta|D^{n - 1})&=\frac{p(D^{n - 1}|\theta)p(\theta)}{\int p(D^{n - 1}|\theta)p(\theta)d\theta} \end{align} \]
为统一表示，记参数先验分布\(p(\theta)\)为\(p(\theta\vert D^{0})\)，表示没有样本情形下的参数概率密度估计。
记\(D^{n}=\{x_{1},x_{2},\cdots,x_{n}\}\)，随着样本的增加，可以得到一系列对参数概率密度函数的估计： \[ p(\theta),p(\theta\vert x_{1}),p(\theta\vert x_{1},x_{2}),\cdots,p(\theta\vert x_{1},x_{2},\cdots,x_{n}),\cdots \] 一般来说，随着样本数目的增加，上述序列函数逐渐尖锐，逐步趋向于以\(\theta\)的真实值为中心的一个尖峰。当样本无穷多时，此时将收敛于一个脉冲函数（参数真值）.
例：贝叶斯估计
假设一维随机变量\(X\)服从\([0,\theta]\)上的均匀分布： \[ p(\mathbf{x}\vert \theta)=U(0,\theta)=\begin{cases} \frac 1 \theta, & 0\leq \mathbf{x}\leq\theta\\ 0, & \text{otherwise} \end{cases} \] 基于先验知识，我们知道\(0 < \theta < 10\)，并希望利用迭代的贝叶斯方法从样本\(\{4,7,2,8\}\)中估计参数\(\theta\)。
迭代过程
在任何数据到达之前，我们有\(p(\theta\vert D^{0}) = p(\theta)=U(0,10)\)。
当第一个数据点\(x_{1}=4\)到达时，则： \[ p(\theta\vert D^{1})=\frac{p(x_{1}\vert \theta)p(\theta\vert D^{0})}{\int p(x_{1}\vert \theta)p(\theta\vert D^{0})d\theta}=\alpha p(x_{1}\vert \theta)p(\theta\vert D^{0})=\alpha\frac{1}{\theta}\frac{1}{10}\\ p(\theta\vert D^{1})\propto\begin{cases} 1/\theta, & 4\leq\theta\leq10\\ 0, & \text{otherwise} \end{cases} \] 其中忽略了归一化。因为\(\theta\)一定要大于等于观测值\(\mathbf x\)。
当第二个数据点\(x_{2}=7\)到达时，我们有： \[ p(\theta\vert D^{2})\propto p(x_{2}\vert \theta)p(\theta\vert D^{1})=\frac{1}{\theta^{2}},\quad p(\theta\vert D^{2})\propto\begin{cases} 1/\theta^{2}, & 7\leq\theta\leq10\\ 0, & \text{otherwise} \end{cases} \] 当第三个数据点\(x_{3}=2\)到达时，我们有： \[ p(\theta\vert D^{3})\propto p(x_{3}\vert \theta)p(\theta\vert D^{2})=\frac{1}{\theta^{3}},\quad p(\theta\vert D^{3})\propto\begin{cases} 1/\theta^{3}, & 7\leq\theta\leq10\\ 0, & \text{otherwise} \end{cases} \] 当第四个数据点\(x_{4}=8\)到达时，我们有： \[ p(\theta\vert D^{4})\propto p(x_{4}\vert \theta)p(\theta\vert D^{3})=\frac{1}{\theta^{4}},\quad p(\theta\vert D^{4})\propto\begin{cases} 1/\theta^{4}, & 8\leq\theta\leq10\\ 0, & \text{otherwise} \end{cases} \] 当数据点\(x_{n}\)到达时，我们有： \[ p(\theta\vert D^{n})\propto p(x_{n}\vert \theta)p(\theta\vert D^{n - 1})=\frac{1}{\theta^{n}},\quad p(\theta\vert D^{n})\propto\begin{cases} 1/\theta^{n}, & \max\{D^{n}\}\leq\theta\leq10\\ 0, & \text{otherwise} \end{cases} \]
关于参数\(\theta\)的分布的调整过程：
参数\(\theta\)的最后估计结果：
\[
p(\theta\vert D^{4})=\begin{cases} 3147.5/\theta^{4}, &
8\leq\theta\leq10\\ 0, & \text{otherwise} \end{cases}
\] 最后的分布： \[
p(\mathbf{x}\vert D)=\int_{\theta}p(\mathbf{x}\vert \theta)p(\theta\vert
D)d\theta\\
p(\mathbf{x}\vert D)=\begin{cases} 0.1134, & 0\leq \mathbf{x}\leq8\\
786.875\left(\frac{1}{\mathbf{x}^{4}}-\frac{1}{10^{4}}\right), & 8
< \mathbf{x}\leq10\\ 0, & \text{otherwise} \end{cases}
\]
最大似然估计做法
对于数据，其似然函数为： \[ l(\theta)=p(x_{1},x_{2},x_{3},x_{4}\vert \theta)=\frac{1}{\theta^{4}} \] 显然，\(l(\theta)\)单调递减，\(\theta\)越小，\(l(\theta)\)越大。但同时，\(\theta\)一定要大于等于最大观测数据。在现有样本\(\{4,7,2,8\}\)中，使似然函数\(l(\theta)\)取值最大的\(\theta\)只能等于8。所以由于是均匀分布，所以\(\theta\)的最大似然估计值为8。
图中展示了最大似然估计（ML）和贝叶斯估计（Bayes）在样本后验分布上的区别。文中提到最大似然方法估计的是\(\theta\)空间中的一个点，而贝叶斯方法估计的是一个分布。
特征维数问题:分类错误率与特征的关系
模式分类与特征的关系
- 贝叶斯决策(0-1损失)：\(w^{*}=\arg\max_{j}p(w_{j}\vert \mathbf{x})\)
- 特征空间给定时，贝叶斯分类错误率就确定了，即分类性能的理论上限就确定了（与分类器、学习算法无关）
增加特征有什么好处：判别性：类别间有差异的特征有助于分类
带来什么问题：泛化性能、Overfitting（过拟合）、计算、存储
高斯分布（两类问题）
\(p(\mathbf{x}\vert w_{j})\sim N(\mu_{j},\Sigma),j = 1,2\)，等协方差矩阵
Bayes error rate（贝叶斯错误率）: \[ \begin{align} P(\text{error})&=P(\mathbf{x}\in R_{2},w_{1})+P(\mathbf{x}\in R_{1},w_{2})\\ &=\int_{R_{2}}p(\mathbf{x}\vert w_{1})P(w_{1})dx+\int_{R_{1}}p(\mathbf{x}\vert w_{2})P(w_{2})dx\\ &=P(\mathbf{x}\in R_{2}\vert w_{1})P(w_{1})+P(\mathbf{x}\in R_{1}\vert w_{2})P(w_{2})\\ P(\text{error})&=\frac{1}{\sqrt{2\pi}}\int_{r/2}^{+\infty}e^{-u^{2}/2}du,\quad r^{2}=(\mu_{1}-\mu_{2})^{T}\Sigma^{-1}(\mu_{1}-\mu_{2}) \end{align} \]
高斯分布（两类问题）
Conditionally independent case（条件独立情况）\(\Sigma=\text{diag}(\sigma_{1}^{2},\sigma_{1}^{2},\cdots,\sigma_{d}^{2})\)每
一维的二类均值之间距离反映区分度，从而决定错误率。
特征增加有助于减小错误率（因为\(r^{2}\)增大） \[ r^{2}=\sum_{i = 1}^{d}\frac{(\mu_{1i}-\mu_{2i})^{2}}{\sigma_{i}^{2}},\quad r^{2}=(\mu_{1}-\mu_{2})^{T}\Sigma^{-1}(\mu_{1}-\mu_{2}) \]
特征维数决定可分性的例子 :
- 3D空间完全可分
- 2D和1D投影空间有重叠
然而，增加特征也可能导致分类性能更差，因为有模型估计误差（wrong model）
过拟合(Overfitting)
特征维数高、训练样本少导致模型参数估计不准确
比如协方差矩阵需要样本数在\(d\)以上.
图中展示了一个10阶多项式拟合的例子，函数为\(f(\mathbf{x})=ax^{2}+bx + c+\epsilon\)，其中\(p(\epsilon)=N(0,\sigma^{2})\)，虽然完美拟合训练数据，但测试误差很大.
克服办法
特征降维：特征提取(变换)、特征选择
参数共享/平滑
方法一：共享协方差矩阵\(\Sigma_{0}\)
方法二：Shrinkage（a.k.a. Regularized Discriminant Analysis)第i类协方差矩阵给出需要分别使用第i类数据和所有数据计算\(\Sigma_i,\Sigma\)、属于启发式方法。 \[ \Sigma_{i}(\alpha)=\frac{(1-\alpha)n_{i}\Sigma_{i}+\alpha n\Sigma}{(1-\alpha)n_{i}+\alpha n},\quad \Sigma(\beta)=(1-\beta)\Sigma+\beta I \]
扩展：开放集分类的特征维数问题
开放集分类问题：
- 已知类别：\(w_{i},i = 1,\cdots,c\)
- 后验概率\(\sum_{i = 1}^{c + 1}P(w_{i}\vert \mathbf{x})=1\)
- \(w_{c + 1}\)无训练样本，测试样本作为outlier（异常值）拒识
特征维数问题：
- 区分\(c + 1\)个类别比区分\(c\)个类别需要更多的特征
- 如果分类器训练时瞄准区分\(c\)个已知类别，测试时易造成outlier与已知类别样本的混淆
- 因此，在\(c\)类样本上训练分类器时，要使特征表达具有区分更多类别的能力
- 比如，训练神经网络时加入数据重构损失(类似auto-encoder)作为正则项，生成一些假想类样本(通过组合已知类别样本)
期望最大化 (Expectation-Maximization, EM)
引入
数据完整情况下例子
数据不完整情况下例子
概率模型：
EM 算法求解
- 第四步：根据第二次掷硬币\(C\)后进行的10次试验，来猜测\(C\)的取值。
- 如果\(C\)为正面，那么掷10次硬币\(A\)，其联合概率为\(0.3^{2}×(1-0.3)^{8}=0.0052\)。
- 如果\(C\)为反面，那么掷10次硬币\(B\)，其联合概率为\(0.61^{3}×(1-0.61)^{7}=0.0003\)。
- 所以猜测\(C\)的取值为正面。
- 第五步：根据猜测的\(C\)值来更新\(p_{A}\)和\(p_{B}\)的值：\(p_{A}=5/20 = 0.25,\ p_{B}=0.61\)
- 第六步：根据第三次掷硬币\(C\)后进行的10次试验，来猜测\(C\)的取值。
- 如果\(C\)为正面，那么掷10次硬币\(A\)，其联合概率为\(0.25^{6}×(1-0.25)^{4}=7.7248×10^{-5}\)。
- 如果\(C\)为反面，那么掷10次硬币\(B\)，其联合概率为\(0.61^{6}×(1-0.61)^{4}=0.0012\)。
- 所以猜测\(C\)的取值为反面。
- 如果\(C\)为正面，那么掷10次硬币\(A\)，其联合概率为\(0.25^{6}×(1-0.25)^{4}=7.7248×10^{-5}\)。
- 第七步：根据猜测的\(C\)值来更新\(p_{A}\)和\(p_{B}\)的值： -\(p_{A}=0.25,\ p_{B}=6/10 = 0.6\)
- 第八步：根据第四次掷硬币\(C\)后进行的10次试验，来猜测\(C\)的取值。
- 如果\(C\)为正面，那么掷10次硬币\(A\)，其联合概率为\(0.25^{3}×(1-0.25)^{7}=0.0021\)。
- 如果\(C\)为反面，那么掷10次硬币\(B\)，其联合概率为\(0.6^{3}×(1-0.6)^{7}=3.5389×10^{-4}\)。
- 所以猜测\(C\)的取值为正面。
- 第九步：\(p_{A}=8/30 = 0.267,\ p_{B}=0.6\)
- 第十步：根据第五次掷硬币\(C\)后进行的10次试验，来猜测\(C\)的取值。
- 如果\(C\)为正面，那么掷10次硬币\(A\)，其联合概率为\((8/30)^{8}×(1-8/30)^{2}=1.3752×10^{-5}\)。
- 如果\(C\)为反面，那么掷10次硬币\(B\)，其联合概率为\(0.6^{8}×(1-0.6)^{2}=0.0027\)。 -所以猜测\(C\)的取值为反面。
- 如果\(C\)为正面，那么掷10次硬币\(A\)，其联合概率为\((8/30)^{8}×(1-8/30)^{2}=1.3752×10^{-5}\)。
- 第十一步：\(p_{A}=0.267,\ p_{B}=14/20 = 0.7\)
一般情况
EM是一类通过迭代实现参数估计的优化算法: 作为最大似然法的替代，用于对包含隐变量（latent variable）或缺失数据（incomplete-data）的概率模型进行参数估计。
EM 算法的核心思想是：
- E 步（Expectation）：计算隐变量的期望。
- M 步（Maximization）：最大化对数似然，更新参数。
EM算法解决的问题：包含隐变量的概率密度参数估计
- 观测变量：\(\mathbf x\)；隐含变量：\(z\)
- 任务：给定数据集\(X = \{x_1,x_2,\cdots,x_n\}\)，估计观测数据概率密度的参数
| EM算法基本要素 | 公式 |
|---|---|
| 观测数据 | \(X=\{x_1,x_2,\cdots,x_n\}\)（不完全数据） |
| 隐含数据 | \(Z = \{z_1,z_2,\cdots,z_n\}\) |
| 观测数据的概率密度函数 | \(p(\mathbf{x}\vert\theta)\) |
| 完全数据的联合概率密度函数 | \(p(\mathbf{x},z\vert \theta)\) |
| 观测数据的对数似然函数 | \(\ln\prod_{i = 1}^{n}p(\mathbf x_i\vert \theta)=\sum_{i = 1}^{n}\ln p(\mathbf x_i\vert \theta)\) |
| 完全数据的对数似然函数 | \(\ln\prod_{i = 1}^{n}p(\mathbf x_i,z_i\vert \theta)=\sum_{i = 1}^{n}\ln p(\mathbf x_i,z_i\vert \theta)\) |
算法步骤
| 算法步骤 | 对应公式 |
|---|---|
| 1. 初始化 | \(\theta^{old}\) |
| 2. (迭代)E步 | 基于当前\(\theta^{old}\)和样本，估计隐变量的后验分布\(p(z_i\vert \mathbf x_i,\theta^{old})\) |
| 3.(迭代)M步 | 基于当前所估计的\(p(z\vert \mathbf{x},\theta^{old})\)更新参数\(\theta\)（如果考虑1维离散隐变量）：\(\theta^{new}=\arg\max_{\theta}Q(\theta,\theta^{old})=\sum_i E_{p(\mathbf x_i\vert \mathbf x_i,\theta^{old})}\left[\ln(p(\mathbf x_i,_i\vert\theta))\right]\\=\arg\max_{\theta}\sum_{i}\sum_{z_i}p(z_i\vert \mathbf x_i,\theta^{old})\ln(p(\mathbf x_i,z_i\vert \theta))\) |
| 4.(迭代) t=t+1 | 回到E步，开始迭代 |
另一种等价表述
初始化：\(\theta^{old}\)
重复以下步骤：
E步（E step）：基于当前\(\theta^{old}\)和样本，估计隐变量的后验分布\(p(z_i\vert \mathbf x_i,\theta^{old})\)，并计算\(Q(\theta,\theta^{old})\)： \[ Q(\theta,\theta^{old})=\sum_{i}\sum_{z_i}p(z_i\vert x_i,\theta^{old})\ln(p(x_i,z_i\vert \theta)) \]
M步（M step）：更新参数\(\theta\)： \[ \theta^{new}=\text{arg max}_{\theta}Q(\theta,\theta^{old}) \]
更新时间步：\(t = t + 1\)
EM for Gaussian Mixture
混合密度模型：由\(K\)个不同成分组成。
| 混合密度模型 | 说明 |
|---|---|
| 权重条件 | 每个成分的权重为\(\pi_k\)，\(k = 1,2,\cdots,K\)，且满足：\(\sum_{k = 1}^{K}\pi_k = 1\)，\(\forall k:\pi_k\geq0\) |
| 概率密度函数 | 每个成分的概率密度函数为\(p(\mathbf{x}\vert \theta_k)\) |
| 混合密度模型公式 | \(p(\mathbf x\vert \pi,\theta)=\sum_{k = 1}^{K}\pi_k p(\mathbf x\vert \theta_k)\)其中，\(\theta = \{\theta_1,\theta_2,\cdots,\theta_K\}\)，\(\pi=\{\pi_1,\pi_2,\cdots,\pi_K\}\)是混合密度模型的参数 |
| 混合密度模型的参数估计 | 已知样本集\(D = \{x_1,x_2,\cdots,x_n\}\)，且样本是从以上混合密度函数中独立抽取的，通过\(D\)估计\((\pi,\theta)\) |
高斯混合模型（Gaussian Mixture Model，GMM）
GMM: 成分密度为高斯密度的混合模型=高斯混合模型。 \[ p(\mathbf{x}\vert \theta)=\sum_{k = 1}^{K}\pi_k p(\mathbf{x}\vert \theta_k)=\sum_{k = 1}^{K}\pi_k \mathcal{N}(\mathbf{x}\vert \mu_k,\Sigma_k) \] 其中\(\mathcal{N}(\mathbf x\vert \mu_k,\Sigma_k)\)是高斯密度函数。
- 权重参数：\(\pi_k\)
- 成分参数：\(\mu_k\)，\(\Sigma_k\)（\(k = 1,2,\cdots,K\)）
参数估计：最大似然（Maximum Likelihood，ML） \[ \max LL=\ln\prod_{i = 1}^{n}p(x_i)=\sum_{i = 1}^{n}\ln\sum_{k = 1}^{K}\pi_k\mathcal{N}(x_i\vert \mu_k,\Sigma_k) \\\nabla_{\pi_k}LL = 0,\quad\nabla_{\mu_k}LL = 0,\quad\nabla_{\Sigma_k}LL = 0 \] （可通过梯度下降迭代求解，但不能解析求解）
Q: 解析解？
A: 解析求解是指通过数学推导和公式运算，直接得到问题的精确解。在参数估计中，如果能够解析求解，就意味着可以通过一系列的代数运算和推导，得到模型参数的精确表达式。因此只能使用梯度下降这样的迭代算法来逐步逼近最优解。
引入隐变量分析混合密度模型
引入隐变量\(z\)用于指示不同的成分密度，\(z\in\{1,\cdots,K\}\)。
| 引入隐变量 | 公式 |
|---|---|
| 假设条件 | \(P(z = k)=\pi_k\)，\(p(\mathbf x\vert z = k)=\mathcal{N}(\mathbf x\vert \mu_k,\Sigma_k)\) |
| 观测变量\(\mathbf x\)和隐变量\(z\)联合分布 | \(p(\mathbf x,z = k)=p(\mathbf x\vert z = k)P(z = k)=\pi_k\mathcal{N}(\mathbf x\vert \mu_k,\Sigma_k)\) |
| 观测变量\(\mathbf x\)的边缘分布 | \(p(\mathbf x)=\sum_{z = 1}^{K}p(\mathbf x,z)=\sum_{k = 1}^{K}\pi_k\mathcal{N}(\mathbf x\vert \mu_k,\Sigma_k)\) |
因此，可以将混合密度模型看成包含隐变量\(z\)的概率密度函数，从而用EM算法估计参数。
EM算法估计参数
不完全数据\(X\)，完全数据\(\{X,Z\}\)：Missing the latent value\(z\)$ for each sample\(z\in\{1,2,\cdots,K\}\)。
完全数据对数似然期望：\(Q(\theta,\theta^{old})=\sum_{i}\sum_{z}p(Z\vert X,\theta^{old})\ln(p(X,Z\vert \theta))\)
- 选择初始参数集\(\theta^{old}\)
- 执行以下操作
- E-step：计算\(p(Z\vert X,\theta^{old})\)。
- M-step：更新参数：\(\theta^{new}=\arg\max_{\theta}Q(\theta,\theta^{old})\)。
- 如果收敛条件不满足：\(\theta^{old}\leftarrow\theta^{new}\)
- 结束
E步（E-Step）：固定当前估计的参数\(\{\pi_k,\mu_k,\Sigma_k\}\)，对每个样本求\(P(z_i\vert \mathbf x_i,\theta^{old})\)，\(i = 1,2,\cdots,n\)。（\(z_i = 1,2,\cdots,K\)） \[ P(z_i\vert \mathbf x_i,\theta^{old})=\frac{p(\mathbf x_i,Z_i\vert \theta^{old})}{p(\mathbf x_i\vert \theta^{old})}=\frac{\pi_{Z_i}\mathcal{N}(\mathbf x_i\vert \mu_{Z_i},\Sigma_{Z_i})}{\sum_{k = 1}^{K}\pi_k\mathcal{N}(\mathbf x_i\vert \mu_k,\Sigma_k)} \]
M步（M-Step）：固定\(\{P(z_i\vert \mathbf x_i,\theta^{old})\}\)，通过\(\max Q(\theta,\theta^{old})\)更新参数\(\{\pi_k,\mu_k,\Sigma_k\}\)。 \[ \begin{align}Q(\theta,\theta^{old})&=\sum_{i}\sum_{z_i = 1:K}p(z_i\vert \mathbf x_i,\theta^{old})\ln(p(\mathbf x_i,z_i\vert \theta)) \\ &=\sum_{i}\sum_{z_i = 1:K}p(z_i\vert \mathbf x_i,\theta^{old})\ln(\pi_{z_i}\mathcal{N}(\mathbf x_i\vert \mu_{z_i},\Sigma_{z_i}))\\ &=\sum_{i}\sum_{z_i = 1:K}p(z_i\vert \mathbf x_i,\theta^{old})(\ln\pi_{z_i}+\ln\mathcal{N}(\mathbf x_i\vert \mu_{z_i},\Sigma_{z_i}))\\ &=\sum_{i}\sum_{z_i = 1:K}p(z_i\vert \mathbf x_i,\theta^{old})\ln\pi_{z_i}+\sum_{i}\sum_{z_i = 1:K}p(z_i\vert \mathbf x_i,\theta^{old})\ln\mathcal{N}(\mathbf x_i\vert \mu_{z_i},\Sigma_{z_i}) \\ &=\sum_{i}\sum_{z_i = 1:K}p(z_i\vert \mathbf x_i,\theta^{old})\ln\pi_{z_i}+\sum_{k = 1:K}\sum_{i}p(z_i = k\vert \mathbf x_i,\theta^{old})\ln\mathcal{N}(\mathbf x_i\vert \mu_k,\Sigma_k) \\ &\frac{\partial Q(\theta,\theta^{old})}{\partial\mu_k}=0,\quad\frac{\partial Q(\theta,\theta^{old})}{\partial\Sigma_k}=0,\quad k = 1,2,\cdots,K\\ &\text{arg max}_{\pi}\sum_{i}\sum_{z_i = 1:K}p(z_i\vert \mathbf x_i,\theta^{old})\ln(\pi_{z_i})\quad\text{s.t.}\quad\sum_{k = 1}^{K}\pi_k = 1 \end{align} \]
成分 对应公式 成分权重\(\hat{\pi}_k\) \(\hat{\pi}_k = \frac{1}{n} \sum_{i = 1}^{n} P(z_i = k\) 成分均值\(\hat{\mu}_k\) \(\hat{\mu}_k = \frac{\sum_{i = 1}^{n} P(z_i = k \vert x_i, \theta^{old}) x_i}{\sum_{i = 1}^{n} P(z_i = k \vert x_i, \theta^{old})}\) 成分协方差矩阵\(\hat{\Sigma}_k\) \(\hat{\Sigma}_k = \frac{\sum_{i = 1}^{n} P(z_i = k \vert x_i, \theta^{old}) (x_i-\hat{\mu}_k)(x_i-\hat{\mu}_k)^T}{\sum_{i = 1}^{n} P(z_i = k \vert x_i, \theta^{old})}\) \(\hat\theta\)记录所有的未知参数, 所以要迭代: \[ P(Z_i = k \vert x_i, \theta^{old}) = \frac{p(x_i, Z_i = k \vert \theta^{old})}{p(x_i \vert \theta)} = \frac{p(x_i \vert \hat{\mu}_k, \hat{\Sigma}_k) \hat{\pi}_k}{\sum_{j = 1}^{K} p(x_i \vert \hat{\mu}_j, \hat{\Sigma}_j) \hat{\pi}_j} \] 在极端情况下，即当样本\(x_i\)来自于一个成分时，其后验概率\(\hat{P}(\omega_k \vert x_i, \hat\theta)=1\)，否则就为零，此时有： \[ P(z_i \vert x_i, \hat{\mu}) = \begin{cases} 1, & x_i \in \omega_k \\ 0, & x_i \notin \omega_k \end{cases}\\ \hat{\pi}_k = \frac{n_k}{n}, \quad \hat{\mu}_k = \frac{1}{n_k} \sum_{i = 1}^{n_i} x_i^{(k)}, \quad \hat{\Sigma}_k = \frac{1}{n_k} \sum_{i = 1}^{n_k} (x_i^{(k)} - \hat{\mu}_k)(x_i^{(k)} - \hat{\mu}_k)^T \] 上标\((k)\)表示属于第\(k\)个成分的样本，\(n_k\)表示属于第\(k\)个成分样本总数。
EM：数据缺失情况下的参数估计
\(\mathbf{x}=\{\mathbf{x_{ig}},\mathbf{x_{ib}}\}\):\(\mathbf{x}\)包含\(\mathbf{x_{ig}}\) (good features,没有缺失的数据部分)、以及\(\mathbf{x_{ib}}\)(bad features, 缺失的、无法观测到的特征)
对缺失数据求期望： \[ Q(\theta,\theta^{old})=\sum_i E_{p(\mathbf{x}_{ib}|\mathbf{x}_{ig},\theta^{old})}\left[\ln(p(\mathbf{x}_{ig},\mathbf{x}_{ib}\vert \theta))\right] \]
- 初始化\(\theta^{old},T,t=0\)
- E 步（期望步骤）：\(Q(\theta,\theta^{old})=\sum_i
E_{p(\mathbf{x}_{ib}\vert
\mathbf{x}_{ig},\theta^{old})}\left[\ln(p(\mathbf{x}_{ig},\mathbf{x}_{ib}\vert
\theta))\right]\)
- 用已有的参数\(\theta^{\text{old}}\) 预测缺失数据可能的取值(隐变量（缺失数据）的期望)。
- M 步（最大化步骤）：\(\theta^{new}\leftarrow\arg\max_\theta
Q(\theta,\theta^{old})\)
- 利用完整的数据（观测数据和填补的缺失数据）重新估计参数。
- 交替迭代 E 步和 M 步，直到收敛（即参数的变化很小或者对数似然函数几乎不再增加）。
例：2-D高斯EM
数据集\(D = \{x_1, x_2, x_3, x_4\}\) 是 2D 高斯分布的一部分，其中有一个样本 x4x_4 的第一个维度\(x_{41}\) 缺失（用 * 表示）。 \[ D = \{x_1, x_2, x_3, x_4\}=\left\{\begin{pmatrix}0\\ 2\end{pmatrix}, \begin{pmatrix}1\\0\end{pmatrix}, \begin{pmatrix}2\\2\end{pmatrix},\begin{pmatrix}*\\4\end{pmatrix}\right\} \]
高斯分布的参数为\(\theta = (\mu_1, \mu_2, \sigma_1, \sigma_2)\)，即两维数据的均值和方差。 \[ p(x_i | \theta) = \frac{1}{2\pi\sigma_1\sigma_2} \exp\left(-\frac{1}{2}\left(\frac{(x_{i1} - \mu_1)^2}{\sigma_1^2} + \frac{(x_{i2} - \mu_2)^2}{\sigma_2^2}\right)\right) \]
我们的目标是：在部分数据缺失的情况下，通过迭代的方法估计参数\(\theta\)。
挑战：缺失数据\(x_{41}\) 使得直接用最大似然估计参数变得困难，需要引入 EM 算法
EM算法步骤
初始化参数为\(\theta^{\text{old}} = (0, 0, 1, 1)^T\)。
E 步：计算缺失数据的期望：第四个样本\(x_4 = (*, 4)^T\)，假设其分布为：\(p(x_{41} \vert x_{42} = 4, \theta^{\text{old}})\) 根据条件高斯分布公式，计算\(x_{41}\) 的条件分布，然后求期望值和相关的加权概率。 \[ \begin{align} Q(\theta,\theta^{old})&=\sum^3_{i=1}\ln p(\mathbf{x_i}\vert\theta)+E_{p(\mathbf{x}_{41}|\mathbf{x}_{42},\theta^{old})}[\ln p(\mathbf{x}_4\vert\theta)]\\ &=\sum_{i = 1}^{3}\ln p(\mathbf{x}_i|\theta)+\int_{-\infty}^{+\infty}\ln p(\mathbf{x}_{4}|\theta)p(\mathbf{x}_{41}|\mathbf{x}_{42}=4,\theta^{old})d\mathbf{x}_{41}\\ &=\sum_{i = 1}^{3}\ln p(\mathbf{x}_i|\theta)+\int_{-\infty}^{+\infty}\ln p\left(\left(\begin{array}{c}\mathbf x_{41}\\4\end{array}\right)|\theta\right)\frac{p\left(\left(\begin{array}{c}x_{41}\\4\end{array}\right)|\theta^{old}\right)}{\int^{+\infty}_{-\infty}{p\left(\left(\begin{array}{c}\mathbf x'_{41}\\4\end{array}\right)|\theta^{old}\right)}d\mathbf{x'}_{41}}d\mathbf{x}_{41}\\ &=\sum_{i = 1}^{3}\ln p(\mathbf{x}_i|\theta)+\frac{1}{\alpha}\int_{-\infty}^{+\infty}\ln p\left(\left(\begin{array}{c}\mathbf x_{41}\\4\end{array}\right)|\theta\right)p\left(\left(\begin{array}{c}\mathbf{x}_{41}\\4\end{array}\right)|\theta^{old}\right)d\mathbf{x}_{41}\\ &=\sum_{i = 1}^{3}\ln p(\mathbf{x}_i|\theta)+\frac{1}{\alpha}\int_{-\infty}^{+\infty}\ln p\left(\left(\begin{array}{c}\mathbf x_{41}\\4\end{array}\right)|\theta\right)\frac{1}{2\pi\begin{vmatrix}1&0\\0&1\end{vmatrix}^{1/2}}\exp\left(-\frac{1}{2}(\mathbf{x}_{41}^2+4^2)\right)d\mathbf{x}_{41}\\ &=\sum_{i = 1}^{3}\ln p(\mathbf{x}_i|\theta)-\frac{1 + \mu_1^2}{2\sigma_1^2}-\frac{(4-\mu_2)^2}{2\sigma_2^2}-\ln(2\pi\sigma_1\sigma_2) \end{align} \]
M 步：最大化对数似然\(\theta=(0.75,2.0,0.938,2.0)^T\)
重复 E 步和 M 步: 迭代 3 次后，收敛到：\(\mu = \begin{pmatrix} 1.0 \\2.0 \end{pmatrix}, \quad \Sigma = \begin{pmatrix} 0.667 & 0 \\ 0 & 2.0 \end{pmatrix}\)
EM算法的理论解释
对任意分布\(q(z)\)： \[ \ln p(\mathbf{x}|\theta)=E_{q(z)}\ln\frac{p(\mathbf{x},z|\theta)}{q(z)} + KL(q(z)||p(z|\mathbf{x},\theta)) \]
KL距离（KL散度）：也称KL散度，衡量相同事件空间中两个概率分布之间的差异。KL距离恒大于等于零；当且仅当\(\forall z, p(z) = q(z)\)时，有\(KL(p(z)\vert \vert q(z)) = 0\)。
因为KL距离恒大于等于零，因此\(E_{q(z)}\ln\frac{p(\mathbf{x},z\vert \theta)}{q(z)}\)是\(\ln p(x\vert \theta)\)的下界，当且仅当\(q(z)=p(z\vert \mathbf{x},\theta)\)时，有\(E_{q(z)}\ln\frac{p(\mathbf{x},z\vert \theta)}{q(z)}=\ln p(\mathbf{x}\vert \theta)\).
令\(L(q,\theta)\equiv E_{q(z)}\ln\frac{p(\mathbf{x},z\vert \theta)}{q(z)}\)
E步（E Step）：固定\(\theta\)，最大化\(\max L(q,\theta)\) \[ \because L(q,\theta)=\ln p(\mathbf{x}|\theta)-KL(q(z)||p(z|\mathbf{x},\theta))\\ \therefore \max_{q}L(q,\theta)\Leftrightarrow\min_{q}KL(q(z)||p(z|\mathbf{x},\theta))\\ \Rightarrow q(z)=p(z|\mathbf{x},\theta) \]
M步（M Step）：固定\(q(z)\)，最大化\(\max_{\theta}L(q,\theta)\) \[ \max_{\theta}L(q,\theta)\Leftrightarrow\max_{\theta}E_{q(z)}\ln p(\mathbf{x},z|\theta) \]
因此，EM是在通过坐标轮替法最大化\(L(q,\theta)\)。但因为\(L(q,\theta)\)是对\(\ln p(\mathbf{x}\vert \theta)\)的近似（下界），也可以粗略地说EM在做极大似然估计。
隐马尔可夫模型(HMM)
时间序列模型:\(X=\{\mathbf{x}_1,\mathbf{x}_2.\dots,\mathbf{x}_n\}\)
- \(n\)是序列长度
- \(\mathbf{x}_t\in\mathbb{R}^d\)是\(X\)在t时刻的观察数据
- 不满足独立假设、观测数据间具有很强的相关性。
Q: 如何对序列数据表示、学习和推理?
首先需要引入关于数据分布和时间轴依赖关系的概率模型，即如何表示\(p(\mathbf{x}_1,\mathbf{x}_2,\dots,\mathbf{x}_n)\)
对\(P(X)\)的假定
- 方法1具有极强的灵活性、通用性，但参数量大、计算复杂度高
- 方法2具有极差的灵活性、通用性，但参数量小、计算复杂度低
- 如何平衡灵活性和复杂度?
| 方法 | 联合分布 |
|---|---|
| 方法1：不对数据做任何独立性假设，直接对条件分布建模(\(\mathbf{x}_t\)和它的全部历史相关) | \(p(\mathbf{x}_1,\dots,\mathbf{x}_n)=p(\mathbf{x_1})\prod^n_{t=2}p(\mathbf{x}_t\vert \mathbf{x}_1,\mathbf{x}_2.\dots,\mathbf{x}_{t-1} )\) |
| 方法2：假设\(\{\mathbf{x}_1,\mathbf{x}_2.\dots,\mathbf{x}_n\}\)独立，只对边缘分布建模 | \(p(\mathbf{x}_1,\dots,\mathbf{x}_n)=\prod^n_{i=1}p(\mathbf{x}_t)\) |
| 方法3：(Markov性)\(\mathbf{x}_{t}\)只与\(\mathbf{x}_{t-1}\)有关\(p(\mathbf{x}_t\vert \mathbf{x}_1,\dots,\mathbf{x}_{t-1})=p(\mathbf{x}_t\vert \mathbf{x}_{t-1})\) | \(p(\mathbf{x}_1,\dots,\mathbf{x}_n)=p(\mathbf{x}_1)\prod^n_{t=2}p(\mathbf{x}_t\vert \mathbf{x}_{t-1})\) |
HMM的表示
Markov链
静态、离散、一阶Markov链的联合分布： \[ p(\mathbf{x}_1,\dots,\mathbf{x}_n)=p(\mathbf{x}_1)\prod^n_{t=2}p(\mathbf{x}_t\vert \mathbf{x}_{t-1}) \]
| 参数 | 说明 |
|---|---|
| 离散马氏链 | \(\mathbf{x}_t\in\{1,2,\dots,K\},K\)为状态数 |
| 静态马氏链 | 转移概率\(p(\mathbf{x_t}\vert \mathbf{x_{t-1}})\)只与状态有关，与时间\(t\)无关 |
| 初始状态分布 | \(p(\mathbf{x}_1)=\pi\in\mathbb{R}^K\) |
| 状态转移概率 | \(p(\mathbf{x_t}\vert \mathbf{x_{t-1}})=A\in\mathbb{R}^{K\times K}\) |
例子
- \(\mathbf{x}_t \in \{雨天,晴天,阴天\}, K = 3\)
- 初始状态分布：\(\pi=[0.1, 0.6, 0.3]\)
- 状态转移概率：
\[ A=\begin{bmatrix} 0.1 & 0.4 & 0.5 \\ 0.1 & 0.6 & 0.3 \\ 0.2 & 0.4 & 0.4 \end{bmatrix} \]
已知第\(t\)天是雨天，第\(t + 2\)天是晴天的概率？ \[ p(\mathbf{x}_t)=[1,0,0]^T\\ p(\mathbf{x}_{t + 1})=A^T p(\mathbf{x}_t)=[0.1, 0.4, 0.5]^T\\ p(\mathbf{x}_{t + 2})=A^T p(\mathbf{x}_{t + 1})=[0.15, 0.48, 0.37]^T \]
HMM简介
- HMM的基本思想
- 观测序列由一个不可见的马尔可夫链生成。
- HMM的随机变量可分为两组：
- 状态变量\(\{z_1,z_2,\cdots,z_n\}\)：构成一阶、离散、静态马尔可夫链。用于描述系统内部的状态变化，通常是隐藏的，不可被观测的。其中，\(z_t\)表示第\(t\)时刻系统的状态。
- 观测变量\(\{\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n\}\)：其中，\(\mathbf{x}_t\)表示第\(t\)时刻的观测变量，通过条件概率\(p(\mathbf{x}_t\vert z_t)\)由状态变量\(z_t\)生成；根据具体问题，\(\mathbf{x}_t\)可以是离散或连续，一维或多维。
- 主要用于时序数据建模，在CV、NLP、语音识别中有诸多应用。
HMM的图结构
- 在图中，箭头表示依赖关系。
- \(t\)时刻的观测变量\(\mathbf{x}_t\)的取值仅依赖于状态变量\(z_t\)。当\(z_t\)已知，\(\mathbf{x}_t\)与其它状态独立。
- \(t\)时刻的状态变量\(z_t\)的取值仅依赖于\(t-1\)时刻的状态变量\(z_{t-1}\)。当\(z_{t-1}\)已知，\(z_t\)与其余\(t-2\)个状态独立。即\(\{z_t\}\)构成马尔可夫链，系统下一刻的状态仅由当前状态决定，不依赖于以往任何状态。
HMM中的条件独立性: \[ p(\mathbf{x}_1,\cdots,\mathbf{x}_n|z_t)=p(\mathbf{x}_1,\cdots,\mathbf{x}_t|z_t)p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t)\\ p(\mathbf{x}_1,\cdots,\mathbf{x}_n|z_{t - 1},z_t)=p(\mathbf{x}_1,\cdots,\mathbf{x}_{t - 2}|z_{t - 1})p(\mathbf{x}_{t - 1}|z_{t - 1})p(\mathbf{x}_t|z_t)p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t)\\\\ p(\mathbf{x}_1,\cdots,\mathbf{x}_{t - 1}|\mathbf{x}_t,z_t)=p(\mathbf{x}_1,\cdots,\mathbf{x}_{t - 1}|z_t)\\ p(\mathbf{x}_1,\cdots,\mathbf{x}_{t - 1}|z_{t - 1},z_t)=p(\mathbf{x}_1,\cdots,\mathbf{x}_{t - 1}|z_{t - 1})\\ p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|\mathbf{x}_t,z_t)=p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t)\\ p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t,z_{t + 1})=p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_{t + 1}) \]
HMM联合概率分布: \[ p(\mathbf{x}_1,\cdots,\mathbf{x}_n,z_1,\cdots,z_n)=p(z_1)\prod_{t = 2}^n p(z_t|z_{t - 1})\prod_{t = 1}^n p(\mathbf{x}_t|z_t) \]
| 基本要素 | 对应公式 |
|---|---|
| 初始状态概率向量\(\pi \in R^K\) | \(\pi_k = P(z_1 = k),\quad 1\leq k\leq K\) |
| 状态转移概率矩阵\(A\in R^{K\times K}\) | \(A_{i,j}=P(z_t = j\vert z_{t-1}=i),\quad 1\leq i,j\leq K\) |
| 发射概率矩阵\(B\in R^{K\times M}\) | 离散：\(B_{i,j}=P(\mathbf{x}_t = j\vert z_t = i),\quad 1\leq i\leq K,1\leq j\leq M\) |
例子
HMM的学习
三个基本问题
| 三个基本问题 | 简述 | 说明 |
|---|---|---|
| 给定模型\([A,B,\pi]\)，如何有效地计算其产生观测序列\(\mathbf{x}=\{\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n\}\)的概率\(P(\mathbf{x}\vert A,B,\pi)\)？ | 评估模型与观测数据的匹配程度。 | 许多任务需要根据以往的观测序列来预测当前时刻最有可能的观测值。 |
| 给定模型\([A,B,\pi]\)和观测序列\(\mathbf{x}=\{\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n\}\)，如何找到与此观测序列相匹配的状态序列\(z=\{z_1,z_2,\cdots,z_n\}\)？ | 根据观测序列推断出隐藏的模型状态。(解码问题) | 在语言识别中，观测值为语音信号，隐藏状态为文字，目标就是观测信号推断最有可能的状态。 |
| 给定观测序列\(\mathbf{x}=\{\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n\}\)，如何调整模型参数\([A,B,\pi]\)使该序列出现的概率\(P(\mathbf{x}\vert A,B,\pi)\)最大？ | 如何模型使其能够最好地描述观测数据。(参数估计-学习问题) | 在大多数实际应用中，人工指定参数已变得不可行，需要根据训练样本学习最优模型。 |
参数学习的基本任务
通过拟合观测序列，确定HMM中的参数:\(\theta=(\pi,A,B)\)
EM算法步骤:
- E Step: 对给定的\(\theta\)，估计： \[ q(z_1,\cdots,z_n) = p_{\theta}(z_1,\cdots,z_n|\mathbf{x}_1,\cdots,\mathbf{x}_n) \]
- M Step: 用估计出的\(q(z_1,\cdots,z_n)\)，更新\(\theta\)： \[ \theta = \arg \max_{\theta} \sum_{z} q(z_1,\cdots,z_n) \ln p_{\theta}(\mathbf{x}_1,\cdots,\mathbf{x}_n,z_1,\cdots,z_n) \]
- E步和M步迭代运行，直至收敛
M step: 更新\(\theta\)
\(p_{\theta}(\mathbf{x}_1,\cdots,\mathbf{x}_n,z_1,\cdots,z_n)= p_{\theta}(z_1)\prod_{t = 2}^{n} p_{\theta}(z_t\vert z_{t-1})\prod_{t = 1}^{n} p_{\theta}(\mathbf{x}_t\vert z_t)\)
\[ \begin{align*} Q(\theta,\theta^{old})&=\sum_{z} q(z_1,\cdots,z_n) \ln p_{\theta}(\mathbf{x}_1,\cdots,\mathbf{x}_n,z_1,\cdots,z_n)\\ &=\sum_{z} q(z_1,\cdots,z_n) (\ln p_{\theta}(z_1)+\sum_{t = 2}^{n} \ln p_{\theta}(z_t|z_{t - 1})+\sum_{t = 1}^{n} \ln p_{\theta}(\mathbf{x}_t|z_t))\\ &=\sum_{z_1=1}^K q(z_1) \ln p_{\theta}(z_1)+\sum_{t = 2}^{n} \sum_{z_{t - 1},z_t = 1}^{K} q(z_{t - 1},z_t) \ln p_{\theta}(z_t|z_{t - 1}) +\sum_{t = 1}^{n} \sum_{z_t = 1}^{K} q(z_t) \ln p_{\theta}(\mathbf{x}_t|z_t)\\ &=\sum_{z_1=1}^K q(z_1)\ln p_\theta(\pi_{z_1})+\sum_{t = 2}^{n} \sum_{z_{t - 1},z_t = 1}^{K} q(z_{t - 1},z_t) \ln A_{z_{t-1},z_t}+\sum_{t = 1}^{n} \sum_{z_t = 1}^{K} q(z_t) \ln B_{z_t,x_t} \end{align*} \]
用拉格朗日乘子法优化以下问题，可得：
| 参数 | 计算 | 约束 | 结果 |
|---|---|---|---|
| \(\pi\) | \(\arg\max_{\pi} \sum_{z_1 = 1}^{K} q(z_1) \ln \pi_{z_1}\) | \(\sum_{k = 1}^{K} \pi_k = 1\) | \(\pi_k = q(z_1 = k)\) |
| \(A\) | \(\arg\max_{A} \sum_{t = 2}^{n} \sum_{z_{t-1},z_t = 1}^{K} q(z_{t-1},z_t) \ln A_{z_{t-1},z_t}\) | \(\forall i \sum_{j = 1}^{K} A_{i,j} = 1\) | \(A_{i,j} = \frac{\sum_{t = 2}^{n} q(z_{t-1}=i,z_t = j)}{\sum_{t = 2}^{n} \sum_{k = 1}^{K} q(z_{t-1}=i,z_t = k)}\) |
| \(B\) | \(\arg\max_{B} \sum_{t = 1}^{n} \sum_{z_t = 1}^{K} q(z_t) \ln B_{z_t,\mathbf{x}_t}\) | \(\forall i \sum_{j = 1}^{M} B_{i,j} = 1\) | \(B_{i,j} = \frac{\sum_{t = 1}^{n} \mathbf{1}\{\mathbf{x}_t == \hat j\}q(z_t = i)}{\sum_{t = 1}^{n} q(z_t = i)}\) |
E Step: 对给定的\(\theta\)，估计\(q(z_1,\cdots,z_n)=p_{\theta}(z_1,\cdots,z_n\vert \mathbf{x}_1,\cdots,\mathbf{x}_n)\)
只需估计： \[ q(z_t) = p_{\theta}(z_t|\mathbf{x}_1,\cdots,\mathbf{x}_n)\\ q(z_{t - 1},z_t) = p_{\theta}(z_{t - 1},z_t|\mathbf{x}_1,\cdots,\mathbf{x}_n) \] 以下省略\(\theta\)： \[ \begin{align} q(z_t) &= \frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n,z_t)}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)}=\frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n|z_t)p(z_t)}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)}\\ &= \frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_t|z_t)p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t)p(z_t)}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)}\\ &=\frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_t,z_t)p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t)}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)}\\\\ q(z_{t - 1},z_t)&=\frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n,z_{t - 1},z_t)}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)}\\ &=\frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n,z_t|z_{t - 1})p(z_{t - 1})}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)}\\ &=\frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_{t - 1}|z_{t - 1})p(\mathbf{x}_t,\cdots,\mathbf{x}_n,z_t|z_{t - 1})p(z_{t - 1})}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)}\\ & =\frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_{t - 1},z_{t - 1})p(\mathbf{x}_t,\cdots,\mathbf{x}_n,z_t|z_{t - 1})}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)} \end{align} \]
前向-后向算法(forward-backward algorithm)
\[ q(z_t) = \frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_t,z_t)p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t)}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)}\\ q(z_{t - 1},z_t) = \frac{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_{t - 1},z_{t - 1})p(\mathbf{x}_t,\cdots,\mathbf{x}_n,z_t|z_{t - 1})p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t)}{p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)} \]
- 前向概率：\(\alpha_t(z_t) = p(\mathbf{x}_1,\cdots,\mathbf{x}_t,z_t)\)
- 后向概率：\(\beta_t(z_t) = p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n\vert z_t)\)
- 观测概率：\(p(\mathbf{x}_1,\mathbf{x}_2,\cdots,\mathbf{x}_n)=\sum_{k = 1}^{K} \alpha_n(z_n)\) \[ q(z_t) = \frac{\alpha_t(z_t)\beta_t(z_t)}{\sum_{z_n = 1}^{K} \alpha_n(z_n)},\quad q(z_{t - 1},z_t) = \frac{\alpha_{t - 1}(z_{t - 1})p(z_t|z_{t - 1})p(\mathbf{x}_t|z_t)\beta_t(z_t)}{\sum_{z_n = 1}^{K} \alpha_n(z_n)} \]
从前向后计算，\(t = 1,\cdots,n\)
\[ \begin{align*} \alpha_t(z_t)&=p(\mathbf{x}_1,\cdots,\mathbf{x}_t,z_t)=\sum_{z_{t - 1}} p(\mathbf{x}_1,\cdots,\mathbf{x}_t,z_{t - 1},z_t)\\ &=\sum_{z_{t - 1}} p(\mathbf{x}_1,\cdots,\mathbf{x}_t,z_t|z_{t - 1})p(z_{t - 1})\\ &=\sum_{z_{t - 1}} p(\mathbf{x}_1,\cdots,\mathbf{x}_{t - 1}|z_{t - 1})p(\mathbf{x}_t,z_t|z_{t - 1})p(z_{t - 1})\\ &=\sum_{z_{t - 1}} p(\mathbf{x}_1,\cdots,\mathbf{x}_{t - 1},z_{t - 1})p(\mathbf{x}_t,z_t|z_{t - 1})\\ &=\sum_{z_{t - 1}} \alpha_{t - 1}(z_{t - 1})p(z_t|z_{t - 1})p(\mathbf{x}_t|z_t)\\ &=p(\mathbf{x}_t|z_t)\sum_{z_{t - 1}} \alpha_{t - 1}(z_{t - 1})p(z_t|z_{t - 1}) \end{align*} \]
从后向前计算，\(t = n,n-1,\cdots,1\)
\[ \begin{align*} \beta_t(z_t)&=p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t)=\sum_{z_{t + 1}} p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n,z_{t + 1}|z_t)\\ &=\sum_{z_{t + 1}} p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_t,z_{t + 1})p(z_{t + 1}|z_t)\\ &=\sum_{z_{t + 1}} p(\mathbf{x}_{t + 1},\cdots,\mathbf{x}_n|z_{t + 1})p(z_{t + 1}|z_t)\\ &=\sum_{z_{t + 1}} p(\mathbf{x}_{t + 2},\cdots,\mathbf{x}_n|z_{t + 1})p(\mathbf{x}_{t + 1}|z_{t + 1})p(z_{t + 1}|z_t)\\ &=\sum_{z_{t + 1}} \beta_{t + 1}(z_{t + 1})p(\mathbf{x}_{t + 1}|z_{t + 1})p(z_{t + 1}|z_t) \end{align*} \]
HMM的解码
- 在实际问题中，状态变量通常有明确的含义。如语音识别中，\(z_t\)表示语音信号\(\mathbf{x}_t\)对应的文本。因此，经常需要根据观测序列推断状态序列。
- 对给定的HMM模型\(\theta = (\pi, A, B)\)和观测序列\(\{\mathbf{x}_1, \mathbf{x}_2, \cdots, \mathbf{x}_n\}\)，求解： \[ z^* = \arg\max_{z} p_{\theta}(z_1, \cdots, z_n | \mathbf{x}_1, \cdots, \mathbf{x}_n) \] \(z^*\)是最大后验概率对应的状态序列，也称为最优状态路径。
- 这对应分类问题中的最大后验概率决策，\(z_t\)对应\(\mathbf{x}_t\)的类别。
- 与分类中对\(\mathbf{x}_t\)独立解码不同，HMM需要联合解码。
状态路径：\(z_1, \cdots, z_n\)
- 对于给定的HMM模型和观测序列\(\{\mathbf{x}_1, \mathbf{x}_2, \cdots, \mathbf{x}_n\}\)，不同状态路径对应不同的后验概率：\(p_{\theta}(z_1, \cdots, _n \vert \mathbf{x}_1, \cdots, \mathbf{x}_n)\)。
- 共有\(K^n\)条可能的状态路径，对应\(K^n\)个概率值。
- 直接计算这些概率，然后选出\(z^*\)的复杂度为\(O(K^n)\)。
HMM的解码算法：维特比算法（Viterbi, 1967）
最优子问题：寻找以状态\(z_t\)结束的前\(t\)步最优状态路径 \[ w_t(z_t) \equiv \max\ln p_{\theta}(\mathbf{x}_1, \cdots, \mathbf{x}_t, z_1, \cdots, z_{t - 1}, z_t) \in R^K\\ z^* = \arg\max_{z} p_{\theta}(z_1, \cdots, z_n | \mathbf{x}_1, \cdots, \mathbf{x}_n) \Leftrightarrow \arg\max_{z} p_{\theta}(z_1, \cdots, z_n, \mathbf{x}_1, \cdots, \mathbf{x}_n) \]
动态规划算法
- For\(z_1 = 1, \cdots, K\)：\(w_1(z_1) = \ln p(z_1) + \ln p(\mathbf{x}_1 \vert z_1)\)
- For\(t = 2, \cdots, n\)：
- For\(z_t = 1, \cdots, K\)： \[ w_t(z_t) = \ln p(\mathbf{x}_t | z_t) + \max_{z_{t - 1} \in \{1, \cdots, K\}} \{w_{t - 1}(z_{t - 1}) + \ln p(z_t | z_{t - 1})\} \]
计算复杂度：\(O(nK^2)\)