隐马尔科夫模型数学理论推导（HMM）

一、隐马尔科夫模型基本概念

隐马尔可夫模型由初始状态概率向量 $π$ 、状态转移概率矩阵 A 和观测概率矩阵 B 决定。A，B， $π$ 称为隐马尔可夫模型的三要素。

λ = (A, B, π)

Note: Picture source

上图中是一个简单的描述股票市场的马尔科夫模型：其中隐藏状态为 ${B u l l, B e a r, E v e n}$ ，观测状态为 ${u p, d o w n, u n c h a n g e d}$

根据图模型，我们可以轻易地写出状态转移概率矩阵 A (Let Bull = 1, Bear = 2, Even = 3)：

A = [\begin{matrix} 0.6 & 0.2 & 0.2 \\ 0.5 & 0.3 & 0.2 \\ 0.4 & 0.1 & 0.5 \end{matrix}]

观测概率矩阵 B (Let up = 1, down = 2, unchanged= 3):

B = [\begin{matrix} 0.7 & 0.1 & 0.2 \\ 0.1 & 0.6 & 0.3 \\ 0.3 & 0.3 & 04 \end{matrix}]

Note: Picture source

隐马尔可夫模型作了两个 基本假设（ $y_{t}$ 为t时刻的观测， $q_{t}$ 为t时刻的状态）：

$p (q_{t} | q_{1}, \dots, q_{t - 1}, y_{1}, \dots, y_{t - 1}) = p (q_{t} | q_{t - 1})$
$p (y_{t} | q_{1}, \dots, q_{t - 1}, q_{t}, y_{1}, \dots, y_{t - 1}) = p (y_{t} | q_{t})$

对第一条假设的解释：隐藏的马尔可夫链在任意时刻t的状态只依赖于其 前一时刻 的状态，与其他时刻的状态及观测无关，也与时刻t无关；

对第二条假设的解释：任意时刻的观测只依赖于 该时刻 的马尔可夫链的状态，与其他观测及状态无关。

下面提出 HMM 的三个基本问题：

（记： $λ = (A, B, π)$ $Q = {q_{1}, \dots, q_{T}}$ $Y = {y_{1}, \dots, y_{T}}$ ）

概率计算问题（计算其产生观测序列的概率） $\Rightarrow$ 计算 $p (Y | λ)$
参数学习问题（训练模型使选取的参数能最好的描述观测数据） $\Rightarrow λ_{M L E} = \underset{λ}{a r g m a x} p (Y | λ)$
解码(decoding)问题（找到与此观测序列最匹配的隐状态序列） $\Rightarrow \underset{Q}{a r g m a x} p (Y | Q, λ)$

下面，我们主要讨论如何解决最常用到的前两个基本问题。

二、概率计算算法

1、直接计算法

\begin{aligned} (1) & p (Y | λ) & = \sum_{Q} p (Y, Q | λ) \\ (2) & = \sum_{Q} p (Y | Q, λ) p (Q | λ) \\ (3) & = \sum_{q_{1}, \dots, q_{T}} p (y_{1}, \dots, y_{T} | q_{1}, \dots, q_{T}; λ) p (q_{1}, \dots, q_{T} | λ) \\ (4) & = \sum_{q_{1}, \dots, q_{T}} & b_{q_{1}} (y_{1}) \dots b_{q_{T}} (y_{T}) \cdot p (q_{1} | λ) p (q_{2} | q_{1}, λ) p (q_{3} | q_{1}, q_{2}, λ) \\ (5) & \dots p (q_{T} | q_{1}, \dots, q_{T - 1}, λ) \\ (6) & = \sum_{q_{1}, \dots, q_{T}} ( & b_{q_{1}} (y_{1}) \dots b_{q_{T}} (y_{T})) \cdot (p (q_{1} | λ) p (q_{2} | q_{1}, λ) p (q_{3} | q_{2}, λ) \\ (7) & \dots p (q_{T} | q_{T - 1}, λ)) \\ (8) & = \sum_{q_{1}, \dots, q_{T}} (b_{q_{1}} (y_{1}) \dots b_{q_{T}} (y_{T})) \cdot (π_{q_{1}} a_{q 1, q 2} \dots a_{q_{T - 1, T}}) \\ (9) & = \sum_{q_{1}, \dots, q_{T}} π_{q_{1}} \prod_{t = 2}^{T} a_{q_{t - 1, t}} b_{q_{t}} (y_{t}) \end{aligned}

但是这种方法的计算量很大，是 $O (T N^{T})$ ，因此这种算法在实际不可行。

2、前向和后向算法

隐马尔科夫模型数学理论推导（HMM）
Note: Picture source

由上图定义了前向概率（左）和后向概率（右），前向概率描述了 $y_{1}$ 到 $y_{t}$ 和 t 时刻为第 $i$ 个状态时的联合分布，后向概率在已知t 时刻为第 $i$ 个状态时描述了 $y_{t + 1}$ 到 $y_{T}$ 的联合分布；

前向算法：

t=1 时，
$α_{i} (1) = p (y_{1}, q_{1} = i | λ) = p (q_{1} = i | λ) p (y_{1} | q_{1} = i, λ) = π_{i} b_{i} (y_{1})$

t=2 时，

\begin{aligned} (33) & α_{j} (1) & = p (y_{1}, y_{2}, q_{2} = j | λ) \\ (34) & = \sum_{i = 1}^{k} p (y_{1}, y_{2}, q_{1} = i, q_{2} = j | λ) \\ (35) & = \sum_{i = 1}^{k} p (q_{1} = i) p (y_{1} | q_{1} = i) p (q_{2} = j | q_{1} = i) \\ (36) & \cdot p (q_{2} = j | q_{1} = i) p (y_{2} | q_{2} = j) \\ (37) & = [\sum_{i = 1}^{k} α_{i} (1) a_{i, j}] b_{j} (y_{2}) \end{aligned}

……

因此，t $\geq$ 2 时， $α_{j} (t + 1) = [\sum_{i = 1}^{k} α_{i} (t) a_{i, j}] b_{j} (y_{t + 1})$

又显然有 $p (Y | λ) = \sum_{i = 1}^{k} α_{i} (T)$ ……（1）

下图可以直观地理解这个过程：

前向算法：

输入：模型参数 $λ$ ，观测 $Y$

输出： $p (Y | λ)$

(1) 初值

$α_{i} (1) = π_{i} b_{i} (y_{1})$

(2) 递推对 t=1,……,T-1

$α_{j} (t + 1) = [\sum_{i = 1}^{k} α_{i} (t) a_{i, j}] b_{j} (y_{t + 1})$

(3) 终止

$p (Y | λ) = \sum_{i = 1}^{k} α_{i} (T)$

后向算法：

t=T 时， $β_{i} (T) = 1$

t=T-1 时，

\begin{aligned} (15) & β_{i} (T - 1) & = p (y_{T} | q_{T - 1} = i) \\ (16) & = \sum_{j = 1}^{k} p (y_{T}, q_{T} = j | q_{T - 1} = i) \\ (17) & = \sum_{j = 1}^{k} p (q_{T} = j | q_{T - 1} = i) p (y_{T} | q_{T} = j, q_{T - 1} = i) \\ (18) & = \sum_{j = 1}^{k} p (q_{T} = j | q_{T - 1} = i) p (y_{T} | q_{T} = j) \\ (19) & = \sum_{j = 1}^{k} a_{i, j} b_{j} (y_{T}) \end{aligned}

t=T-2 时，

\begin{aligned} (20) & β_{i} (T - 2) & = p (y_{T}, y_{T - 1} | q_{T - 2} = i) \\ (21) & = \sum_{j = 1}^{k} p (y_{T}, y_{T - 1}, q_{T - 1} = j, q_{T} = l | q_{T - 2} = i) \\ (22) & = \sum_{j = 1}^{k} \sum_{l = 1}^{k} p (q_{T} = j | q_{T - 1} = i) p (y_{T} | q_{T} = j) \cdot \\ (23) & p (q_{T - 1} = j | q_{T - 2} = i) p (y_{T - 1} | q_{T - 1} = j) \\ (24) & = \sum_{j = 1}^{k} p (q_{T} = j | q_{T - 1} = i) p (y_{T} | q_{T} = j) \\ (25) & = \sum_{j = 1}^{k} a_{i, j} b_{j} (y_{T - 1}) β_{j} (T - 1) \end{aligned}

……

因此， $t \leq T - 1$ 时：

β_{i} (t) = \sum_{j = 1}^{k} a_{i, j} b_{j} (y_{t + 1}) β_{j} (t + 1)

又显然有 $p (Y | λ) = \sum_{i = 1}^{k} π_{i} b_{i} (y_{1}) β_{i} (1)$ ……（2）

后向算法：

输入：模型参数 $λ$ ，观测 $Y$

输出： $p (Y | λ)$

(1) 初值

$β_{i} (T) = 1$

(2) 递推对 t=T-1,……,1

$β_{i} (t) = \sum_{j = 1}^{k} a_{i, j} b_{j} (y_{t + 1}) β_{j} (t + 1)$

(3)终止

$p (Y | λ) = \sum_{i = 1}^{k} π_{i} b_{i} (y_{1}) β_{i} (1)$

前向算法和后向算法的统一

利用前面的定义可以将观测序列概率 $p (Y | λ)$ 统一：

p (Y | λ) = \sum_{i = 1}^{k} \sum_{j = 1}^{k} α_{i} (t) a_{i, j} b_{j} (y_{t + 1}) β_{j} (t + 1) ， t = 1, \dots, T - 1

此式当 t=1 和 t=T-1 时分别为式（1）和（2）。

特定状态的概率

给定模型 $λ$ 和观测 $Y$ ，在时刻 t 处于状态 $i$ 的概率如下；

$p (q_{t} = i | Y, λ) = \frac{p (Y, q_{t} = i | λ)}{p (Y | λ)} = \frac{p (Y, q_{t} = i | λ)}{\sum_{j = 1}^{k} p (Y, q_{t} = j | λ)} = \frac{α_{i} (t) β_{i} (t)}{\sum_{j = 1}^{k} α_{i} (t) β_{i} (t)}$

其中，

\begin{aligned} (26) & p (Y, q_{t} = i | λ) & = p (Y | q_{t} = i, λ) p (q_{t} = i | λ) \\ (27) & = p (y_{1}, \dots, y_{t} | q_{t} = i) p (y_{t + 1}, \dots, y_{T} | q_{t} = i) p (q_{t} = i | λ) \\ (28) & = p (y_{1}, \dots, y_{t}, q_{t} = i | λ) p (y_{t + 1}, \dots, y_{T} | q_{t} = i) \\ (29) & = α_{i} (t) β_{i} (t) \end{aligned}

三、参数学习算法

根据训练数据是包括观测序列和对应的状态序列还是只有观测序列，可以分别由监督学习与非监督学习实现。但由于监督学习需要使用训练数据，而人工标注数据代价也往往很高，因此我们会利用非监督的学习方法来学习参数。

将观测序列数据看作观测数据 Y，状态序列数据看作不可观测的隐数据 Q，那么隐马尔可夫模型事实上是一个含有隐变量的概率模型 :

p (Y | λ) = \sum_{Q} p (Y | Q, λ) p (Q | λ)

我们先回顾一下EM算法：EM算法推导（收敛性证明和在GMM中的应用）

在 HMM 中，我们可以写成如下：

\begin{aligned} (30) & λ^{(g + 1)} & = \underset{λ}{a r g m a x} (\int_{q} l n (p (Y, q | λ)) p (q | Y, λ^{(g)})) \\ (31) & = \underset{λ}{a r g m a x} (\int_{q} l n (p (Y, q | λ)) \frac{p (q, Y | λ^{(g)})}{p (Y | λ^{(g)})}) \\ (32) & = \underset{λ}{a r g m a x} (\int_{q} l n (p (Y, q | λ)) p (q, Y | λ^{(g)})) \end{aligned}

（ $p (Y | λ^{(g)})$ 为常数，因此可以省略）

E step:

求 Q 函数，见如下：

Q (λ, λ^{(g)}) = \int_{q} l n (p (Y, q | λ)) p (q, Y | λ^{(g)}) = \sum_{q_{1}, \dots, q_{T}} (l n π_{q_{1}}) p (q, Y | λ^{(g)}) + \sum_{q_{1}, \dots, q_{T}} (\sum_{t = 1}^{T} l n a_{q_{t - 1}, q_{t}}) p (q, Y | λ^{(g)}) + \sum_{q_{1}, \dots, q_{T}} (\sum_{t = 1}^{T} l n b_{q_{t}} (y_{t})) p (q, Y | λ^{(g)})

M step:

极大化 Q 函数，求模型参数 A,B, $π$ ；

观察上述 Q 函数，要极大化的参数分别单独地出现在3个项中，所以只需对各项分别极大化；

（1）第一项可以写出：

\sum_{q_{1}, \dots, q_{T}} (l n π_{q_{1}}) p (q, Y | λ^{(g)}) = \sum_{i = 1}^{k} (l n π_{i}) p (q_{1} = i, Y | λ^{(g)})

s . t . \sum_{i = 1}^{k} π_{i} = 1

利用拉格朗日乘子法，即可求解；

π_{i} = \frac{p (q_{1} = i, Y | λ^{(g)})}{\sum_{i = 1}^{k} p (q_{1} = i, Y | λ^{(g)})}

（1）第二项可以写出：

\sum_{q_{1}, \dots, q_{T}} (\sum_{t = 1}^{T} l n a_{q_{t - 1}, q_{t}}) p (q, Y | λ^{(g)}) = \sum_{i = 1}^{k} \sum_{j = 1}^{k} \sum_{t = 1}^{T} l n a_{i, j} p (q_{t - 1} = i, q_{t} = j, Y | λ^{(g)})

s . t . \sum_{i = 1}^{k} a_{i j} = 1

同样利用拉格朗日乘子法，即可求解：

a_{i, j} = \frac{\sum_{t = 1}^{T} p (q_{t - 1} = i, q_{t} = j, Y | λ^{(g)})}{\sum_{t = 1}^{T} p (q_{t - 1} = i, Y | λ^{(g)})}

（3）第二项可以写出：

\sum_{q_{1}, \dots, q_{T}} (\sum_{t = 1}^{T} l n b_{q_{t}} (y_{t})) p (q, Y | λ^{(g)}) = \sum_{j = 1}^{k} \sum_{t = 1}^{T} (l n b_{j} (y_{t})) p (q_{t} = j, Y | λ^{(g)})

s . t . \sum_{i = 1}^{k} b_{j} (y_{t}) = 1

同样利用拉格朗日乘子法，即可求解（注意：只有在 $y_{t} = v_{l}$ 时，偏导数才不为0）：

b_{j} (y_{t} = v_{l}) = \frac{\sum_{t = 1}^{T} p (q_{t} = j, Y = v_{l} | λ^{(g)})}{\sum_{t = 1}^{T} p (q_{t} = j, Y | λ^{(g)})}

上面的 EM 算法又称为 Baum-Welch 算法。

Baum-Welch 算法：

输入：观测序列Y

(1) 初始化 $λ_{0} = (A^{(0)}, B^{(0)}, π^{(0)})$

(2) 递推直至EM算法收敛

$π_{i}^{(g + 1)} = \frac{p (q_{1} = i, Y | λ^{(g)})}{\sum_{i = 1}^{k} p (q_{1} = i, Y | λ^{(g)})}$

$a_{i, j}^{(g + 1)} = \frac{\sum_{t = 1}^{T} p (q_{t - 1} = i, q_{t} = j, Y | λ^{(g)})}{\sum_{t = 1}^{T} p (q_{t - 1} = i, Y | λ^{(g)})}$

$b_{j} (l)^{(g + 1)} = \frac{\sum_{t = 1}^{T} p (q_{t} = j, Y = v_{l} | λ^{(g)})}{\sum_{t = 1}^{T} p (q_{t} = j, Y | λ^{(g)})}$

(3) 终止。得到参数 $λ^{(n + 1)} = (A^{(n + 1)}, B^{(n + 1)}, π^{(n + 1)})$

四、参考资料

[1] 李航《统计学习方法》
[2] 徐亦达教授的自视频
[3] machine-learning-notes.Professor Richard Xu .