【机器学习基础】一文读懂隐马尔可夫模型 HMM

闲言碎语

以前一直在看和复现 cv 方向的论文，最近有需求所以去接触了一些 DL-EEG（深度学习-脑电图）的论文和模型，由于这一方向基本没有开源代码，所以从基础学起方便以后自己搭建模型;

而 GMM-HMM 是一个语音识别系统常用的的模型，对于 EEG 这样的时序数据也同样适用，所以今天就看一下 GMM-HMM 模型；

课件地址：

链接: https://pan.baidu.com/s/1T-oApPtrDADrdRr6-BsXtw

提取码: 5q8y

什么是马尔可夫性？

在马尔可夫性的定义中，在已知“现在”的条件下，“将来”与“过去”无关，这种特性就叫强马尔可夫性，具有这种性质的马尔可夫过程叫强马尔可夫过程；

在数学上定义为：

设 $X(t),t\in T$ 是一个随机过程，$E$ 为其状态空间，若对任意 $t_1<t_2<\dots <t_n<t$，有任意的 $x_1,x_2,\dots ,x_n,x\in E$，随机变量 $X(t)$ 在已知变量 $X(t_1)=x_1,\dots ,X(t_n)=x_n$ 之下的条件分布函数就至于 $X(t_n)=x_n$ 有关，即条件分布式满足等式：

$$F(x,t|x_n,x_{n-1},\dots ,x_2,x_1,t_n,t_{n-1},\dots ,t_2,t_1)=F(x,t|x_n,t_n)$$

即

P { X ( t ) ≤ x ∣ X ( t 1 ) = x 1 , … , X ( t n ) = x n } = P { X ( t ) ≤ x ∣ X ( t n ) = x n } P\lbrace X(t)\le x|X(t_1)=x_1,…,X(t_n)=x_n\rbrace=P\lbrace X(t)\le x|X(t_n)=x_n\rbrace P{X(t)≤x∣X(t1​)=x1​,…,X(tn​)=xn​}=P{X(t)≤x∣X(tn​)=xn​}

一个简单的例子

举一个简单的股票市场的例子（来自徐亦达老师，附上链接 https://github.com/roboticcam）：

我们称在 $t$ 时刻观测到的值/事件记为 y t ∈ { Y 1 , Y 2 , … , Y L } y_t\in\lbrace Y_1,Y_2,…,Y_L\rbrace yt​∈{Y1​,Y2​,…,YL​} ，这 $L$ 个值分别代表了 $L$ 种可能事件（图中为 up,down,unchanged 三种可能事件，即 $L=3$）；

同时这个过程中也有一些隐藏的离散的状态量 q t ∈ { Q 1 , Q 2 , … , Q K } q_t\in\lbrace Q_1,Q_2,…,Q_K\rbrace qt​∈{Q1​,Q2​,…,QK​}，被称为隐状态，这 $K$ 个值分别代表了 $K$ 种可能隐状态（图中为Bull,Bear,Even 三种无法直接观测的隐状态，即 $K=3$）；

以下图的时序数据为例：

根据马尔可夫性则有：

$$P(q_t|q_1,...,q_{t-1},y_1,...,y_{t-1})=P(q_t|q_{t-1})$$

$$P(y_t|q_1,...,q_{t-1},q_t,y_1,...,y_{t-1})=P(y_t|q_t)$$

其中：

1. 称 $P(q_t|q_{t-1})$ 为转换概率（transition probability），即从隐状态 $q_{t-1}$ 到 $q_t$ 的概率；
2. 称 $P(y_t|q_t)$ 为输出概率（emission probability），即在隐状态 $q_t$ 下观测到 $y_t$ 的概率；

为了方便计算，我们将这些概率记为矩阵的形式：

$$A=\left\{\begin{array}{ccc}{a}_{0,0}& a_{0,1}& a_{0,2}\\ a_{1,0}& a_{1,1}& a_{1,2}\\ a_{2,0}& a_{2,1}& a_{2,2}\end{array}\right\}$$

$$B=\left\{\begin{array}{ccc}{b}_{0,0}& b_{0,1}& b_{0,2}\\ b_{1,0}& b_{1,1}& b_{1,2}\\ b_{2,0}& b_{2,1}& b_{2,2}\end{array}\right\}$$

其中：

1. $A$ 的大小为 $K\times K$，$B$ 的大小为 $K\times L$；
2. $a_{i,j}$ 表示从 $Q_i$ 转移到 $Q_j$ 的概率，$b_{i,j}$ 表示在 $Q_i$ 下观测到 $Y_j$ 的概率；

由第一张图我们可以分别写出这些概率，例如：

1. 从隐状态 EVEN 到 BULL 的概率为 0.4，可以写作 $P(q_t=BULL|q_{t-1}=EVEN)=0.4$，即 $a_{3,1}=0.4$；
2. 在隐状态 BEAR 下观测到 DOWN 的概率为 0.6，可以写作 $P(y_t=DOWN|q_t=BEAR)=0.6$，即 $b_{2,2}=0.6$；

为了方便书写，我们令：
q t ∈ { B U L L , B E A R , E V E N } = { 1 , 2 , 3 } q_t\in\lbrace BULL,BEAR,EVEN\rbrace=\lbrace 1,2,3\rbrace qt​∈{BULL,BEAR,EVEN}={1,2,3}
y t ∈ { U P , D O W N , U N C H A N G E D } = { 1 , 2 , 3 } y_t\in\lbrace UP,DOWN,UNCHANGED\rbrace=\lbrace 1,2,3\rbrace yt​∈{UP,DOWN,UNCHANGED}={1,2,3}

同时我们可以分别写出 $A$ 和 $B$：

$$A=\left\{\begin{array}{ccc}0.6& 0.2& 0.2\\ 0.5& 0.3& 0.2\\ 0.4& 0.5& 0.1\end{array}\right\}$$

$$B=\left\{\begin{array}{ccc}0.7& 0.1& 0.2\\ 0.1& 0.6& 0.3\\ 0.3& 0.3& 0.4\end{array}\right\}$$

并且我们记 λ = { A , B } \lambda=\lbrace A,B\rbrace λ={A,B}，称之为隐马尔可夫模型的参数；

引入隐状态的公式推理

直观上可以看到，原来相互依赖的观测值都变成相互独立了，相互依赖的量变成了 $q_{t-1}\to q_t\to q_{t+1}$，我们就称这样一条链为马尔可夫链；

那么就有：
$$\begin{array}{rl}P(y_{t+1},y_t,...,y_1)& =\sum _{q_{t+1}=1}^K\sum _{q_t=1}^K...\sum _{q_1=1}^{K}P(y_{t+1},y_t,...,y_1,q_{t+1},q_t,...,q_1)\\ & =\sum _{q_{t+1}=1}^K\sum _{q_t=1}^K...\sum _{q_1=1}^{K}P(y_{t+1}|y_t,y_{t-1},q_{t+1},q_t,q_{t-1})\times P(y_t,y_{t-1},q_{t+1},q_t,q_{t-1})\end{array}$$
由马尔可夫性可以化简为：
$$\begin{array}{rl}P(y_{t+1},y_t,...,y_1)& =\sum _{q_{t+1}=1}^K\sum _{q_t=1}^K...\sum _{q_1=1}^{K}P(y_{t+1}|q_{t+1})\times P(y_t,y_{t-1},q_{t+1},q_t,q_{t-1})\end{array}$$

同样可以继续化简为：

$$\begin{array}{rl}P(y_{t+1},y_t,...,y_1)& =\sum _{q_{t+1}=1}^K\sum _{q_t=1}^K...\sum _{q_1=1}^{K}P(y_{t+1}|q_{t+1})\times P(q_{t+1}|y_t,y_{t-1},q_t,q_{t-1})\times P(y_t,y_{t-1},q_t,q_{t-1})\\ & =\sum _{q_{t+1}=1}^K\sum _{q_t=1}^K...\sum _{q_1=1}^{K}P(y_{t+1}|q_{t+1})\times P(q_{t+1}|q_t)\times P(y_t,y_{t-1},q_t,q_{t-1})\\ & =...（注意到这是递归式）\\ & =\sum _{q_{t+1}=1}^K\sum _{q_t=1}^K...\sum _{q_1=1}^{K}P(y_{t+1}|q_{t+1})\times P(q_{t+1}|q_t)\times ...\times P(y_1|q_1)\times P(q_1)\end{array}$$

当然我们还可以注意到，式子中有一项 $P(q_1)$ 是未知的，那我们就给它一个初始化值，称之为 $\pi$，这样我们的隐马尔可夫模型完整的参数就写作：

λ = { A , B , π } \lambda=\lbrace A,B,\pi\rbrace λ={A,B,π}

那么在参数 $\lambda$ 下完整的隐马尔可夫概率模型就可以表示为：

$$\begin{array}{rl}P(y_{t+1},y_t,...,y_1|\lambda )& =\sum _{q_{t+1}=1}^K\sum _{q_t=1}^K...\sum _{q_1=1}^{K}P(y_{t+1}|q_{t+1})\times P(q_{t+1}|q_t)\times ...\times P(y_1|q_1)\times P(q_1)\end{array}$$

更具体一点，对于一组观测值 Y T = { y 1 , y 2 , . . . , y T } Y_T=\lbrace y_1,y_2,...,y_T\rbrace YT​={y1​,y2​,...,yT​}，隐马尔可夫模型对它的概率预测为：

$$\begin{array}{rl}P(y_T,y_{T-1},...,y_1|\lambda )& =\sum _{q_T=1}^K\sum _{q_{T-1}=1}^K...\sum _{q_1=1}^{K}P(y_T|q_T)\times P(q_T|q_{T-1})\times ...\times P(y_1|q_1)\times P(q_1)\\ & =\sum _{q_T=1}^K\sum _{q_{T-1}=1}^K...\sum _{q_1=1}^K\pi (q_1)\prod _{t=2}^T{a}_{q_{t-1},q_t}{b}_{q_t}(y_t)\end{array}$$

简化推理式

虽然我们得到了上述式子，但是存在的问题就是它的计算复杂度较高，对于长度为 $T$ 的序列，需要进行 $K^T$ 次计算；

为了简化推理式，我们首先定义如下概念：

定义一：${\alpha }_i(t)=P(y_1,y_2,...,y_t,q_t=i|\lambda )$，为 $q_t=i$ 和它的所有观测值 $y_1,y_2,...,y_t$ 的联合概率；

定义二：${\beta }_i(t)=P(y_{t+1},...,y_T|q_t=i,\lambda )$，为 $q_t=i$ 条件下其他的观测值 $y_{t+1},...,t_T$ 的联合概率；

那么我们可以写出：

(1)
$$\begin{array}{rl}{\alpha }_i(1)& =P(y_1,q_1=i)\\ & =P(y_1|q_1=i)\times P(q_1=i)\\ & =b_i(y_1)\times \pi \end{array}$$

(2)

$$\begin{array}{rl}{\alpha }_j(2)& =P(y_1,y_2,q_2=j)\\ & =\sum _{i=1}^{K}P(y_1,y_2,q_1=i,q_2=j)\\ & =\sum _{i=1}^{K}P(y_2|q_2=j)\times P(q2=j|q_1=j)\times P(y_1,q_1=i)\\ & =P(y_2|q_2=j)\sum _{i=1}^{K}P(q_2=j|q_1=i)\times {\alpha }_i(1)\\ & =b_j(y_2)\sum _{i=1}^K{a}_{i,j}{\alpha }_i(1)\end{array}$$

$$......$$

以此类推，我们就可以得到：

$$\begin{array}{r}{\alpha }_j(T)=b_j(y_T)\sum _{i=1}^K{a}_{i,j}{\alpha }_i(T-1)\end{array}$$

这样的递推关系中，由于每个 $\alpha (t)$ 只计算一次 $sum$，故计算复杂度为 $K\times T$；

此时我们再看一下定义一：

$${\alpha }_i(t)=P(y_1,y_2,...,y_t,q_t=i|\lambda )$$

我们可以写出：

$$P(y_T,y_{T-1},...,y_1|\lambda )=\sum _{j=1}^K{\alpha }_j(T)$$

这种算法就叫做 Forward Algorithm；

还有一种对应的算法 Backward Algorithm，可以在已知全部序列时算出所有事件的概率，即：

$$\begin{array}{rl}P(Y,q_t=i)& =P(Y|q_t=i)\times P(q_t=i)\\ & =P(y_1,...,y_t|q_t=i)\times P(y_{t+1},...,y_T|q_t=i)\times P(q_t=i)\\ & ={\alpha }_i(t)\times {\beta }_i(t)\end{array}$$

如何获取参数值？

参数值可以通过最大似然估计获取，即对于给定的 $n$ 个已知观测值的数据 Y N = { Y 1 , Y 2 , . . . , Y N } Y_N=\lbrace Y^{1},Y^{2},...,Y^{N}\rbrace YN​={Y1,Y2,...,YN}，我们可以通过以下方程求解参数 $\lambda$ 的估计值：

$$\hat{\lambda }=\underset{\lambda }{\mathrm{argmax}}logP(Y^1,Y^2,...,Y^n|\lambda )$$

通常这个估计方程的解会非常复杂，因此我们采用 EM 算法去估计参数值，EM 算法的详细解释和收敛性证明可以看我的这一篇：

【机器学习基础】EM算法详解及其收敛性证明

这里直接给出迭代步骤：

我们先假设隐藏参数 ${\epsilon }_t(i,j)=P(q_t=i,q_{t+1}|Y,\lambda )$，表示对于给定的参数 $\lambda$ 和 观测序列 $Y\in Y_N$， Y = { y 1 , . . . , y n } Y=\lbrace y_1,...,y_n\rbrace Y={y1​,...,yn​}，在 $t$ 时刻为状态 $i$ 并且在 $t+1$ 时刻为状态 $j$ 的概率；

展开为：

$$\begin{array}{rl}{\epsilon }_t(i,j)& =\frac{P(q_t=i,q_{t+1}=j,Y|\lambda )}{P(Y|\lambda )}\\ & =\frac{{\alpha }_t(i)a_{i,j}{b}_j(y_{t+1}){\beta }_{t+1}(j)}{P(Y_N|\lambda )}\\ & =\frac{{\alpha }_t(i)a_{i,j}{b}_j(y_{t+1}){\beta }_{t+1}(j)}{{\sum }_{i=1}^K{\sum }_{j=1}^K{\alpha }_t(i)a_{i,j}{b}_j(y_{t+1}){\beta }_{t+1}(j)}\end{array}$$

并定义参数 ${\gamma }_t(i)={\sum }_{i=1}^n{\epsilon }_t(i,j),1\le t\le T$，表示对于给定的 $\lambda$ 和 $Y$ 下 $t$ 时为状态 $i$ 的概率；

那么将 $t$ 带入上式，就可以得到表示为状态 $i$ 转移出去的次数的期望值 ${\epsilon }_t(i,j)$ ，和表示为从状态 $i$ 到状态 $j$ 的次数的期望值 ${\gamma }_t(i)$；

因此就可以得到参数 $\lambda$ 的估计值：

$$\widehat{{\pi }_i}={\gamma }_1(i)$$
$${\hat{a}}_{i,j}=\frac{{\sum }_{t=1}^{T-1}{\epsilon }_t(i,j)}{{\sum }_{t=1}^{T-1}{\gamma }_t(i)}$$
$${\hat{b}}_j(k)=\frac{{\sum }_{t=1}^T{\gamma }_t(j)\delta (y_t,k)}{{\sum }_{t=1}^T{\gamma }_t(j)}$$

其中：

$$\delta (x,k)=\{\begin{array}{ll}1& \text{x=k}\\ 0& \text{x!=k}\end{array}$$

迭代过程为：

1. 初始化参数 λ = { A , B , π } \lambda=\lbrace A,B,\pi\rbrace λ={A,B,π}；
2. E-Step：由参数估计值 ${\lambda }_i$ 计算 ${\epsilon }_t(i,j)$ 和 ${\gamma }_t(i)$；
3. M-Step：用期望 ${\epsilon }_t(i,j)$ 和 ${\gamma }_t(i)$ 得到参数估计值 ${\lambda }_{i+1}$；

其中:

E-step 和 M-Step 循环迭代，直到参数 $\lambda$ 收敛；

关注我的公众号：

感兴趣的同学关注我的公众号——可达鸭的深度学习教程：