小孩都看得懂的 HMM
全文共 7043 字，75 幅图，
预计阅读时间 36 分钟。
本文是「小孩都看得懂」系列的第十九篇，本系列的特点是内容不长，碎片时间完全可以看完，但我背后付出的心血却不少。喜欢就好！
- 小孩都看得懂的神经网络
- 小孩都看得懂的推荐系统
- 小孩都看得懂的逐步提升
- 小孩都看得懂的聚类
- 小孩都看得懂的主成分分析
- 小孩都看得懂的循环神经网络
- 小孩都看得懂的 Embedding
- 小孩都看得懂的熵、交叉熵和 KL 散度
- 小孩都看得懂的 p-value
- 小孩都看得懂的假设检验
- 小孩都看得懂的基尼不纯度
- 小孩都看得懂的 ROC
- 小孩都看得懂的 SVD
- 小孩都看得懂的 SVD 2
- 小孩都看得懂的 GMM
- 小孩都看得懂的贝塔分布
- 小孩都看得懂的多臂老虎机
- 小孩都看得懂的 GAN
- 小孩都看得懂的 HMM
0
故事背景
有个男生只要晴天就高兴，他想和远在千里的女友玩个游戏，只告诉他的心情，让女友猜他当地的天气。
女友很了解他，想既然你今天高兴，那么今天是晴天。
第二天，男生说我今天不高兴，问女朋友天气如何？
女友很了解他，想既然你今天不高兴，那么今天是雨天。
女友把这个男生卡的死死的，因为
- 如果是晴天，男友高兴
- 如果是雨天，男友烦躁
因此可以完美根据他的心情来判断天气，这个规则就是：
- 如果男友高兴，晴天
- 如果男友烦躁，雨天
如下图所示：
但人的心情不会这么简单，不可能根据天气而一成不变。让我们看看一个稍微复杂的情况：
- 如果是晴天，男友 80% 的情况下高兴，20% 的情况下烦躁
- 如果是雨天，男友 60% 的情况下烦躁，40% 的情况下高兴
如下图所示：
如果男生说这三天心情是“高兴-烦躁-高兴”，那么女友可能推断出（不是 100% 确定）天气是“晴天-雨天-晴天”。
如果男生这一周心情像过山车，隔一天高兴（Happy，H）隔一天烦躁（Grumpy，G），那么这一周他心情历程是“HGHGHG”。
他女友根据对男友的了解（即如果是晴天，男友 80% 的情况下高兴，20% 的情况下烦躁；如果是雨天，男友 60% 的情况下烦躁，40% 的情况下高兴），推断出这一周天气是“SRSRSR”，其中 S 代表 Sunny 晴天，R 代表 Rain 雨天。
Something is wrong！心情可以像过山车，但是天气一晴一雨的很少见，一般今天晴天明天大概率晴天，今天雨天明天还有可能是雨天。
因此，连续两天的天气可用一个转换概率（transition probability）来描述。如果今天是晴天，那么明天：
- 继续是晴天的概率为 80%，即 0.8
- 变成是雨天的概率为 20%，即 0.2
如果今天是雨天，那么明天：
- 继续是晴天的概率为 40%，即 0.4
- 变成是雨天的概率为 60%，即 0.6
上图就是一个 HMM 的雏形。
1
HMM 是什么
HMM 的全称是 Hidden Markov Model，中文是隐马尔可夫模型。
HMM 有四个重要的概念：观测值（observation）、隐藏状态（hidden）、转换概率（transition probability）和输出概率（emission probability）。如下图所示，让我们来一一剖析。
以女友角度出发，
- 观测值 - 就是能观测到的。比如男友的心情是高兴（H）还是烦躁（G）
- 隐含状态 - 看不到的状态，只能靠观测值来推断的。比如千里之外男友所在地的天气是晴天（S）还是雨天（R）
转换概率：隐含状态转换的概率，即
- 今天 S 到明天 S 和 R 的概率，本例分别是 0.8 和 0.2
- 今天 R 到明天 S 和 R 的概率，本例分别是 0.4 和 0.6
如下图所示。
输出概率：从隐含状态到观测值的概率，即
- S 到 H 和 G 的概率，本例分别是 0.8 和 0.2
- R 到 H 和 G 的概率，本例分别是 0.4 和 0.6
如下图所示。
总结 HMM 的示意图如下：
明晰 HMM 的概念后，让我们步入正题。
2
四个问题
为了把 HMM 讲透，接下来从易到难来分析以下四个问题。
- 如何估计转换概率和输出概率？
- 在没有任何男生情绪的信息情况下，那么晴天和雨天的概率是多少？
- 如果男生今天“高兴”，那么晴天和雨天的概率是多少？
- 如果男生连续三天是“高兴-烦躁-高兴”，那么这连续三天的天气是什么？
3
解决问题一
如何估计转换概率和输出概率？
转换概率
收集历史数据做统计，假设收集了 16 天的天气数据。
统计出“晴天-晴天”和“晴天-雨天”的个数，并标准化为概率，本例中结果为 0.8 和 0.2，两者加起来等于 1。
再统计出“雨天-晴天”和“雨天-雨天”的个数，并标准化为概率，本例中结果为 0.4 和 0.6，两者加起来等于 1。
这样就可以得到转换概率了。
输出概率
同样，收集 16 天“天气-心情”的数据。
统计出“晴天-高兴”和“晴天-烦躁”的个数，并标准化为概率，本例中结果为 0.8 和 0.2，两者加起来等于 1。
再统计出“雨天-高兴”和“雨天-烦躁”的个数，并标准化为概率，本例中结果为 0.4 和 0.6，两者加起来等于 1。
综上，我们已经计算出转换概率和输出概率，总结于下图。
第一个问题回答完毕！
4
解决问题二
在没有任何男生情绪的信息情况下，那么晴天和雨天的概率是多少？
假设今天是晴天，可能是因为昨天是晴天造成的，也可能是因为昨天是雨天造成的。如下图所示得到第一个方程
S = 0.8S + 0.4R
其中 0.8 是 “昨天 S - 今天 S” 的转换概率，0.4 是 “昨天 R - 今天 S” 的转换概率。
假设今天是雨天，可能是因为昨天是晴天造成的，也可能是因为昨天是雨天造成的。如下图所示得到第两个方程
R = 0.2S + 0.6R
其中 0.2 是 “昨天 S - 今天 R” 的转换概率，0.6 是 “昨天 R - 今天 R” 的转换概率。
仔细分析这两个方程是等价的，这样两个未知量一个方程解不出来，还需要一个，而 S + R = 1 就是我们所需的。
这样很容易求解出 S = 2/3，R = 1/3。
最终我们得到在不知道男生心情时，晴天和雨天的概率为 2/3 和 1/3。
当女友不知道男友心情时，她对天气的最优推断是 2/3 可能性是晴天，1/3 可能性是雨天。
第二个问题回答完毕！
5
解决问题三
如果男生今天“高兴”，那么晴天和雨天的概率是多少？
以上节最后结果为基准（ 2/3 可能性是晴天），先看两种情况。
情况一：当女友知道男友高兴时，，最优推断是大于 2/3 可能性（比如 4/5）是晴天。
情况二：当女友知道男友烦躁时，最优推断是小于 2/3 可能性（比如 2/5）是晴天。
那么这个具体概率值是多少呢？需要用贝叶斯定理（Bayes Theorem）来计算。这里只需要计算一天的概率，因此不需要转换概率，只需要输出概率。
因为晴天和雨天的概率为 2/3 和 1/3，不严谨地将其整数化用 2 个晴天和 1 天雨天代表概率。
根据“晴天-高兴”和“晴天-烦躁”的输出概率 0.8 和 0.2，再不严谨地将其整数化用 8 个高兴 (8/10) 和 2 天烦躁 (2/10) 代表输出概率。
根据“雨天-高兴”和“雨天-烦躁”的输出概率 0.4 和 0.6，再不严谨地将其整数化用 2 个高兴 (2/5) 和 3 天烦躁 (3/5) 代表输出概率。
回到该问题，如果男生高兴，那么晴天的概率是多少？简单，首先统计出所有男生高兴的次数，10 次。
10 次高兴中有 8 次发生在晴天，因此“当男生高兴时晴天”的条件概率为 8/10。
10 次高兴中有 2 次发生在雨天，因此“当男生高兴时雨天”的条件概率为 2/10。
同理，
- 5 次烦躁中有 2 次发生在晴天，因此“当男生烦躁时晴天”的条件概率为 2/5。
- 5 次烦躁中有 3 次发生在雨天，因此“当男生烦躁时雨天”的条件概率为 3/5。
第三个问题回答完毕！
6
解决问题四
如果男生连续三天是“高兴-烦躁-高兴”，那么这连续三天的天气是什么？
有可能是“晴天-雨天-晴天”，但怎么得到这个结果的呢？一个简单粗暴的方法就是穷举法，3 天每天 2 种天气一共 2^3 = 8 种组合，根据转换概率和输出概率计算下列 8 组概率，找一个最大值即可：
- 晴天-晴天-晴天
- 晴天-晴天-雨天
- 晴天-雨天-晴天
- 晴天-雨天-雨天
- 雨天-晴天-晴天
- 雨天-晴天-雨天
- 雨天-雨天-晴天
- 雨天-雨天-雨天
该方法简单粗暴，但是效率很低，本帖最终会介绍一个高效算法，维特比算法。首先还是看看简单粗暴法是怎么算的吧。
7
两天的情况
为了便于讲说，回顾两个观测值和两个隐含状态的缩写：
- 观测值：高兴（H）、烦躁（G）
- 隐含状态：晴天（S）、雨天（R）
先从两天开始，男生的心情是 HG，让女友推断这两天的天气是什么。
女友穷举出以下四种情况，并计算每种情况发生的概率。
以上图第二种情况 SR 为例，之前已经得到
- S 的概率是 2/3
- 输出概率 S-H 是 0.8
- 转移概率 S-R 是 0.2
- 输出概率 R-G 是 0.6
发生 SR 的总概率为
2/3 * 0.8 * 0.2 * 0.6 = 0.064
将穷举的四种组合发生概率全部计算出来，如下图所示，得到 SS 的概率最大，0.085。该方法也叫做最大似然法（maximum likelihood）。
因此当男生说两天心情为“高兴-烦躁”时，女友可推断出最有可能发生的天气是“晴天-晴天”。
8
三天的情况
再看三天，男生的心情是 HGH，让女友推断这三天的天气是什么。
女友穷举出以下八种情况，并计算每种情况发生的概率。
以上图第三种情况 SRS 为例，之前已经得到
- S 的概率是 2/3
- 输出概率 S-H 是 0.8
- 转移概率 S-R 是 0.2
- 输出概率 R-G 是 0.6
- 转移概率 R-S 是 0.2
- 转移概率 S-H 是 0.8
发生 SRS 的总概率为
2/3 * 0.8 * 0.2 * 0.6
* 0.4 * 0.8 = 0.02048
将穷举的八种组合发生概率全部计算出来，如下图所示，得到 SSS 的概率最大。
虽然可以得到正确答案，但随着天数的增多，用穷举法的情况指数性增加。
9
四天的情况
最后看四天，男生的心情是 HGHH，让女友推断这四天的天气是什么。
女友穷举出以下十六种情况，并计算每种情况发生的概率。
当天数为 N 时，需要穷举的情况为 2^N，当 N 很大时，这种暴力穷举法根本行不通，我们需要更高效的方法，一种算法就是下节要讲的维特比（Viterbi）算法。
10
Viterbi 算法
根据男生的心情链 HHGGGH，女友推断最有可能发生的天气链。
解决此问题的核心思路是迭代法，就是把第 k 步的结果和第 k-1 步的结果连立起来。
先看最后一天，星期六的天气可能是 S 和 R。
注意力先聚焦到星期六天气为 S。从星期一到星期五有很多有路径可以到 S（如下图三条），那么总有一种是最有可能发生的。
假设我们找到中间一天是“星期一到星期五”最有可能发生的，那么连上星期六的 S，就是“星期一到星期六但最后一天是 S”最有可能发生的。
这样迭代关系就建立了，当星期六是 S，
5 天最有可能天气链 + S
= 6 天最有可能天气链
当星期六是 R，
5 天最有可能天气链 + R
= 6 天最有可能天气链
两者比较找到最大值，就锁定了“6 天最有可能发生的天气链”。
接下来我们过一遍 Viterbi 算法。
11
Viterbi 算法
1. 从起点星期一开始，天气为 S 的概率为 0.67，天气为 R 的概率是 0.33（答案由问题二解答）。
2. 接下来计算当天气为 S 和 R 而导致心情 H 的概率（用输出概率）：
- P(H|S) = 0.67*0.8 = 0.533
- P(H|R) = 0.33*0.4 = 0.133
因为 0.533 > 0.133，发现 S 比 R 导致心情 H 的可能性更大。
3. 到星期二了，有两种情况，S 和 R。
先看星期二的 S，可从星期一的 S 和 R 而来（用输出概率和转换概率）：
- P(H|SS) = 0.533*0.8*0.8 = 0.341
- P(H|RS) = 0.133*0.4*0.8 = 0.043
因为 0.341> 0.043，发现 SS 比 RS 导致心情 H 的可能性更大。
再看星期二的 R，可从星期一的 S 和 R 而来（用输出概率和转换概率）：
- P(H|SR) = 0.533*0.2*0.4 = 0.043
- P(H|RR) = 0.133*0.6*0.4 = 0.032
因为 0.043 > 0.032，发现 SR 比 RS 导致心情 H 的可能性更大。
4. 接着从星期三到星期六，一直可按照上述算法，得到截止到某一天是 S 或 R 而导致当天心情是 H 或 G 的最大概率值。具体过程如下面动图所示。
5. 下图展示最终结果，这些数字的含义是：
- P(H|星期一 = S) = 0.533
- P(H|星期一 = R) = 0.133
- max( P(HH|星期二 = S) )= 0.341
- max( P(HH|星期二 = R) ) = 0.043
- max( P(HHG|星期三 = S) )= 0.0546
- max( P(HHG|星期三 = R) ) = 0.041
- max( P(HHGG|星期四 = S) ) = 0.0087
- max( P(HHGG|星期四 = R) ) = 0.0147
- max( P(HHGGG|星期五 = S) ) = 0.0014
- max( P(HHGGG|星期五 = R) ) = 0.0053
- max( P(HHGGGH|星期六 = S) )= 0.0017
- max( P(HHGGGH|星期六 = R) ) = 0.0013
6. 沿着日期，找出每个日期下的最大值连成一条线，就是要找到的天气链，SSSRRS。
7. 如果男生说这一周他心情是“高兴-高兴-烦躁-烦躁-烦躁-高兴”，他女友用 Viterbi 算法得到的一周天气是“晴天-晴天-晴天-雨天-雨天-晴天”。
12
Python 实现
回顾 HMM 示意图。
设定初始晴天和雨天概率、转换概率和输出概率，并初始化男生一周心情链。
# Initial Probabilities
p_s, p_r = 2/3, 1/3
# Transition Probabilities
p_ss, p_sr, p_rs, p_rr, = 0.8, 0.2, 0.4, 0.6
# Emission Probabilities
p_sh, p_sg, p_rh, p_rg = 0.8, 0.2, 0.4, 0.6
mood = ['H', 'H', 'G', 'G', 'G', 'H']
proba = []
weather = []在起点星期一时，计算从“晴天 S 和雨天 R” 到“高兴 H 和烦躁 G” 的概率。
if mood[0] == 'H':
proba.append((p_s*p_sh, p_r*p_rh))
else:
proba.append((p_s*p_sg, p_r*p_rg))
proba[(0.5333333333333333, 0.13333333333333333)]和上节手算结果比对没问题。
接下来用 Viterbi 算法来迭代计算从“前一天 S 和 R” 到“当天 S 和 R” 导致 H 和 G 的概率，取最大值并更新。
for i in range(1,len(mood)):
sunny_1, rainy_1 = proba[-1]
if mood[i] == 'H':
sunny0 = max(sunny_1*p_ss*p_sh, rainy_1*p_rs*p_sh)
rainy0 = max(sunny_1*p_sr*p_rh, rainy_1*p_rr*p_rh)
proba.append((sunny0, rainy0))
else:
sunny0 = max(sunny_1*p_ss*p_sg, rainy_1*p_rs*p_sg)
rainy0 = max(sunny_1*p_sr*p_rg, rainy_1*p_rr*p_rg)
proba.append((sunny0, rainy0))
proba[(0.5333333333333333, 0.13333333333333333),
(0.3413333333333334, 0.04266666666666667),
(0.05461333333333335, 0.04096000000000001),
(0.008738133333333337, 0.014745600000000001),
(0.0013981013333333341, 0.005308416),
(0.00169869312, 0.00127401984)]和上节手算结果比对没问题。
变量 proba 是一个长度为 T 的列表，第 t 个元素代表第 t 天上天气在各个转态下的概率，假设有 N 个状态，那么 proba 里每个元素是一个长度为 N 的元组，第 n 个元素代表第 n 个天气状态。在本例中 T = 6，N = 2。
根据 proba 晴天概率 p[0] 和雨天概率 p[1] 的大小，得出当天是晴天还是雨天。
for p in proba:
if p[0] > p[1]:
weather.append('S')
else:
weather.append('R')
weather['S', 'S', 'S', 'R', 'R', 'S']和上节手算结果比对没问题。
13
Python 进一步实现
上节已经用 Python 代码实现出来本帖的例子了，但是不够通用。比如所有转换概率和输出概率都是用标量表示的，一旦状态和观测值多的话，代码会非常难看。本节用 numpy array 来实现 Viterbi 算法，专门用一个函数来实现它。
代码看起来很多，但其实就是不同矩阵或向量之间相乘，只要把形状匹配就没问题了。
将本帖例子用列表或者 numpy 数组来表示：
带入 viterbi() 函数中计算结果：
probability, state_sequence = viterbi( states, observations,
prior_probability,
trans_probability,
emission_probability,
observation_sequence )
probability
state_sequence array([[0.53333333, 0.13333333],
[0.34133333, 0.04266667],
[0.05461333, 0.04096 ],
[0.00873813, 0.0147456 ],
[0.0013981 , 0.00530842],
[0.00169869, 0.00127402]])
['sunny', 'sunny', 'sunny', 'rainy', 'rainy', 'sunny']结果和上节的完全吻合！
朋友们，你们弄懂了 HMM 了吗？
腾讯云开发者
