HMM 隐马尔可夫模型(下)

使用 Viterbi 算法寻找隐藏状态

Posted by Xiaosheng on March 24, 2017

在上一篇《HMM 隐马尔可夫模型(上)》中,我们通过一个掷骰子的例子,简单地阐述了一下 HMM 的基本概念,以及 HMM 关注问题的解决方法。文本将正式介绍隐马尔可夫模型的数学定义,并通过一个实例,演示具体问题的求解过程。

1. 隐马尔可夫模型

隐马尔可夫模型中的变量可分为两组。第一组是状态变量 ${y_1,y_2,\cdots,y_n}$,其中 $y_i\in Y$ 表示第 $i$ 时刻的系统状态。通常假定状态变量是隐藏的、不可被观测的,因此状态变量亦称隐变量(hidden variable)。第二组是观测变量 ${x_1,x_2,\cdots,x_n}$,其中 $x_i \in X$ 表示第 $i$ 时刻的观测值。在隐马尔可夫模型中,系统通常在多个状态 ${s_1,s_2,\cdots,s_N}$ 之间转换,因此状态变量 $y_i$ 的取值范围 $Y$ (称为状态空间)通常是有 $N$ 个可能取值的离散空间。观测变量 $x_i$ 可以是离散型也可以是连续型。为方便讨论,我们仅考虑离散型观测变量,并假定其取值范围 $X$ 为 ${o_1,o_2,\cdots,o_M}$。隐马尔可夫模型的图结构如下所示:

1

图中的箭头表示了变量间的依赖关系。在任一时刻,观测变量的取值仅依赖于状态变量,即 $x_t$ 由 $y_t$ 确定,与其他状态变量及观测变量的取值无关。同时,$t$ 时刻的状态 $y_t$ 仅依赖于 $t-1$ 时刻的状态 $y_{t-1}$,与其余 $n-2$ 个状态无关。这就是所谓的马尔可夫链(Markov chain),即:系统下一时刻的状态仅由当前状态决定,不依赖于以往的任何状态。基于这种依赖关系,所有变量的联合概率分布为:

从马尔可夫链到马尔可夫模型

隐马尔可夫模型其实并不是 19 世纪俄罗斯数学家马尔可夫(Andrey Markov)发明的,而是美国数学家鲍姆(Leonard E. Baum)等人在 20 世纪六七十年代发表的一系列论文中提出的,隐马尔可夫模型的训练方法(鲍姆-韦尔奇算法)也是以他的名字命名的。

到了 19 世纪,概率论的发展从对随机变量的研究发展到对随机过程的研究。但是,随机过程要比随机变量复杂得多。首先,在任一时刻 $t$,对应的状态 $s_t$ 都是随机的;第二,任一状态 $s_t$ 的取值都可能和周围其他状态相关。这样随机过程就有了两个维度的不确定性。马尔可夫为了简化问题,提出一种假设:随机过程中各个状态 $s_t$ 的概率分布只与它前一个状态 $s_{t-1}$ 有关,即 $P(s_t\mid s_1,s_2,\cdots,s_{t-1})=P(s_t\mid s_{t-1})$。这种假设未必适合所有的应用,但是至少对以前很多不好解决的问题给出了近似解。这个假设后来被命名为马尔可夫假设,又称马尔可夫链。

隐马尔可夫模型是马尔可夫链的一个扩展,其中隐含的状态序列是一个典型的马尔可夫链。鲍姆把这种模型称为“隐含”马尔可夫模型。

除了结构信息,欲确定一个隐马尔可夫模型还需以下三组参数:

  • 状态转移概率:模型在各个状态间转换的概率,通常记为矩阵 $\mathbf{A}=[a_{ij}]_{N\times N}$,其中:

    表示在任意时刻 $t$,若状态为 $s_i$,则在下一时刻状态为 $s_j$ 的概率。

  • 输出观测概率:模型根据当前状态获得各个观测值的概率,通常记为矩阵 $\mathbf{B}=[b_{ij}]_{N\times M}$,其中

    表示在任意时刻 $t$,若状态为 $s_i$,则观测值 $o_j$ 被获取的概率。

  • 初始状态概率:模型在初始时刻各状态出现的概率,通常记为 $\boldsymbol{\pi}=(\pi_1,\pi_2,\cdots,\pi_N)$,其中:

    表示模型的初始状态为 $s_i$ 的概率。

通过指定状态空间 $Y$、观测空间 $X$ 和上述三组参数,就能确定一个隐马尔可夫模型,通常用其参数 $\lambda=[\mathbf{A},\mathbf{B},\boldsymbol{\pi}]$ 来指代。给定隐马尔可夫模型 $\lambda$,它按如下过程产生观测序列 ${x_1,x_2,\cdots,x_n}$:

  1. 设置 $t = 1$,并根据初始状态概率 $\boldsymbol{\pi}$ 选择初始状态 $y_1$;
  2. 根据状态 $y_t$ 和输出观测概率 $\mathbf{B}$ 选择观测变量取值 $x_t$;
  3. 根据状态 $y_t$ 和状态转移矩阵 $\mathbf{A}$ 转移模型状态,即确定 $y_{t+1}$;
  4. 若 $t<n$,设置 $t = t + 1$,并转到第 2 步,否则停止。

其中 $y_t\in{s_1,s_2,\cdots,s_N}$ 和 $x_t\in{o_1,o_2,\cdots,o_M}$ 分别为第 $t$ 时刻的状态和观测值。

在实际应用中,人们常关注隐马尔可夫模型的三个基本问题:

  • 给定模型 $\lambda=[\mathbf{A},\mathbf{B},\boldsymbol{\pi}]$,如何有效计算其产生观测序列 $\mathbf{x}={x_1,x_2,\cdots,x_n}$ 的概率 $P(\mathbf{x}\mid \lambda)$?换言之,如何评估模型与观测序列之间的匹配程度?
  • 给定模型 $\lambda=[\mathbf{A},\mathbf{B},\boldsymbol{\pi}]$ 和观测序列 $\mathbf{x}={x_1,x_2,\cdots,x_n}$,如何找到与此观测序列最匹配的状态序列 $\mathbf{y}={y_1,y_2,\cdots,y_n}$?换言之,如何根据观测序列推断出隐藏的模型状态?
  • 给定观测序列 $\mathbf{x}={x_1,x_2,\cdots,x_n}$,如何调整模型参数 $\lambda=[\mathbf{A},\mathbf{B},\boldsymbol{\pi}]$ 使得该序列出现的概率 $P(\mathbf{x}\mid\lambda)$ 最大?换言之,如何训练模型使其能最好地描述观测数据?

上述问题在现实应用中非常重要。例如许多任务需根据以往的观测序列 ${x_1,x_2,\cdots,x_{n-1}}$ 来推测当前时刻最有可能的观测值 $x_n$,这可以转化为求取概率 $P(\mathbf{x}\mid\lambda)$,即上述第一个问题;在语音识别等任务中,观测值为语音信号,隐藏状态为文字,目标就是根据观测信号来推断最有可能的状态序列(即对应的文字),即上述第二个问题呢;在大多数现实应用中,人工指定模型参数已变得越来越不行,如何根据训练样本学得最优的模型参数,恰是上述第三个问题。

2. 使用 Viterbi 算法寻找隐藏状态

正如前面所说,HMM 有三个金典问题。下面我们对第二个问题(即已知模型参数,寻找最可能的能产生某一特定输出序列的隐含状态的序列)给出一个实例,并通过 Viterbi 算法解决。

维特比算法(Viterbi algorithm)是一种动态规划算法。它用于寻找最有可能产生观测事件序列的维特比路径——隐含状态序列,特别是在马尔可夫信息源上下文和隐马尔可夫模型中。维特比算法由安德鲁·维特比于 1967 年提出,用于在数字通信链路中解卷积以消除噪音。现今也被常常用于语音识别关键字识别计算语言学生物信息学中。

假设我有一个住的很远的朋友,她根据每天的天气情况(雨天、晴天)来决定当天的活动(公园散步、购物及清理房间)。我对于她所住地方的天气情况并不了解,但是我知道总的趋势。这个朋友每天会在微博上发布她做的事:“我今天出去散步了”、“我今天去超市购物了”、“今天清理了房间”等等,我希望通过她每天发的微博来推断她所在地的天气情况。

可以把天气情况看成是一个马尔可夫链,其有两个状态“雨天”和“晴天”,但是我们无法直接观察,或者说它们是隐藏的。因为这个朋友会把每天的活动发在微博上,所以这些活动就是可见的观测数据。这整个系统就是一个隐马尔可夫模型。

因为我知道这个地区的总的天气趋势,并且平时知道你朋友会做的事情,所以这个隐马尔可夫模型的参数是已知的。例如,可以用 Python 语言写下来:

states = ('Rainy', 'Sunny') #状态
observations = ('walk', 'shop', 'clean') #观测值序列
#初始状态概率
start_probability = {'Rainy': 0.6, 'Sunny': 0.4}
#状态转移概率
transition_probability = {
    'Rainy' : {'Rainy': 0.7, 'Sunny': 0.3},
    'Sunny' : {'Rainy': 0.4, 'Sunny': 0.6},
}
#输出观测概率
emission_probability = {
    'Rainy' : {'walk': 0.1, 'shop': 0.4, 'clean': 0.5},
    'Sunny' : {'walk': 0.6, 'shop': 0.3, 'clean': 0.1},
}

2.1 Viterbi 算法

假设状态变量 $y_i$ 的取值范围是一个有 $N$ 个可能取值的离散空间 ${s_1,s_2,\cdots,s_N}$,且初始状态 $y_1=s_i$ 的概率为 $\pi_i$,从状态 $s_i$ 到状态 $s_j$ 的转移概率为 $a_{ij}$。令观察到的输出序列为 $x_1,x_2,\cdots,x_T$。则产生观察结果的最有可能的状态序列 $y_1,y_2,\cdots,y_T$ 由递推关系给出:

此处 $V_{t,s_i}$ 代表的是观测结果可能对应的前 $t$ 个最终状态为 $s_i$ 的状态序列的最大概率。通过保存向后指针记住在第二个等式中选择的状态 $s_k$ 可以获得维特比路径。声明一个函数 $Ptr(s_i,t)$,当 $t>1$ 时返回计算 $V_{t,s_i}$ 时用到的上一状态 $s_j$;当 $t=1$ 时,返回 $s_i$。这样:

简单一点来说,在任意时刻 $t$,我们都要计算 $y_t$ 取每一个状态 $s_i$ 时的最大概率,并且要记录下使得 $s_i$ 概率最大的上一状态,以便反推获得维特比路径。

2.2 完整 Python 代码实现

首先,我们需要一个路径概率表来存储任意时刻 $t$,每一个状态 $s_i$ 对应的状态序列的最大概率。我们用 [{s1:prob, s2:prob, ...}, ...] 这样的结构来存储路径概率表,列表 $t$ 位置的字典存储 $t$ 时刻对应的各状态 $s_i$ 的最大概率。我们首先编写路径概率表的展示函数:

# 打印路径概率表
def print_dptable(V):
	print "    ",
	for i in range(len(V)): print "%7d" % i,
	print
 
	for state in V[0].keys():
		print "%.5s: " % state,
		for t in range(len(V)):
			print "%.7s" % ("%f" % V[t][state]),
		print

接下来我们编写核心的 Viterbi 算法函数:

def viterbi(obs, states, start_p, trans_p, emit_p):
	'''
	:param obs:观测序列
	:param states:状态
	:param start_p:初始状态概率
	:param trans_p:状态转移概率
	:param emit_p:输出观测概率
	:return:
	'''
	V = [] #路径概率列表
	path = {} #维特比路径
	#初始状态 (t = 0)
	V.append({}) #状态最大概率用字典表示
	for state in states:
		V[0][state] = start_p[state] * emit_p[state][obs[0]]
		path[state] = [state]
	#对 t > 0 跑一遍维特比算法
	for t in range(1, len(obs)):
		V.append({})
		newpath = {}
		for s in states:
			(prob, state) = max([(V[t-1][s0] * trans_p[s0][s], s0) for s0 in states])
			V[t][s] = prob * emit_p[s][obs[t]]
			newpath[s] = path[state] + [s]
		#不需要保留旧路径
		path = newpath
 
	print_dptable(V)
	prob, state = max([(V[-1][s], s) for s in states])
	return prob, path[state]

运用之前我们提供的 HMM 各项参数,我们可以测试一下上面编写的算法:

print viterbi(observations,states,start_probability,transition_probability,emission_probability)

输出:

           0       1       2
Rainy:  0.06000 0.03840 0.01344
Sunny:  0.24000 0.04320 0.00259
(0.01344, ['Sunny', 'Rainy', 'Rainy'])

参考

  1. 周志华《机器学习》
  2. 维基百科《隐马尔可夫模型》《维特比算法》
  3. 码农场《HMM与分词、词性标注、命名实体识别》
  4. 吴军《数学之美》