使用R语言进行机制检测的隐马尔可夫模型HMM

在本文中,将对“牛市”和“熊市”两个独立机制下的市场收益进行模拟。隐马尔可夫模型识别处于特定状态的概率。

由Kaizong Ye,Sherry Deng撰写

在概述了模拟数据的过程之后,将隐马尔可夫模型应用于股票数据,以确定基本机制。

市场体制

将隐马尔可夫模型应用于状态检测是棘手的,因为该问题实际上是无监督学习的一种形式。也就是说,没有“基础事实”或标记数据用来“训练”模型。是否有两个,三个,四个或更多个“真正的”隐藏市场机制?

×

1 隐马尔可夫模型HMM

隐马尔科夫模型(Hidden Markov Model,以下简称HMM)是比较经典的机器学习模型了,它在语言识别,自然语言处理,模式识别等领域得到广泛的应用。

当然,随着目前深度学习的崛起,尤其是RNNLSTM等神经网络序列模型的火热,HMM的地位有所下降。

但是作为一个经典的模型,学习HMM的模型和对应算法,对我们解决问题建模的能力提高以及算法思路的拓展还是很好的。

1.1 什么样的问题需要HMM模型

首先我们来看看什么样的问题解决可以用HMM模型。

使用HMM模型时我们的问题一般有这两个特征:

1)我们的问题是基于序列的,比如时间序列,或者状态序列。

2)我们的问题中有两类数据,一类序列数据是可以观测到的,即观测序列;而另一类数据是不能观察到的,即隐藏状态序列,简称状态序列。

有了这两个特征,那么这个问题一般可以用HMM模型来尝试解决。这样的问题在实际生活中是很多的。比如:我现在在打字写博客,我在键盘上敲出来的一系列字符就是观测序列,而我实际想写的一段话就是隐藏序列,输入法的任务就是从敲入的一系列字符尽可能的猜测我要写的一段话,并把最可能的词语放在最前面让我选择,这就可以看做一个HMM模型了。再举一个,我在和你说话,我发出的一串连续的声音就是观测序列,而我实际要表达的一段话就是状态序列,你大脑的任务,就是从这一串连续的声音中判断出我最可能要表达的话的内容。

从这些例子中,我们可以发现,HMM模型可以无处不在。但是上面的描述还不精确,下面我们用精确的数学符号来表述我们的HMM模型。

1.2 HMM模型的定义

对于HMM模型,首先我们假设Q是所有可能的隐藏状态的集合,V是所有可能的观测状态的集合,即:

其中,N是可能的隐藏状态数,M是所有的可能的观察状态数。

对于一个长度为 [公式] 的序列, [公式] 对应的状态序列, [公式] 是对应的观察序列,即:

HMM模型做了两个很重要的假设如下:

1) 齐次马尔科夫链假设。即任意时刻的隐藏状态只依赖于它前一个隐藏状态。当然这样假设有点极端,因为很多时候我们的某一个隐藏状态不仅仅只依赖于前一个隐藏状态,可能是前两个或者是前三个。但是这样假设的好处就是模型简单,便于求解。如果在时刻 [公式] 的隐藏状态是 [公式] ,在时刻 [公式] 的隐藏状态是 [公式] ,则从时刻 [公式] 到时刻 [公式] 的HMM状态转移概率 [公式] 可以表示为:

这样 [公式] 可以组成马尔科夫链的状态转移矩阵 [公式] :

2) 观测独立性假设。即任意时刻的观察状态只仅仅依赖于当前时刻的隐藏状态,这也是一个为了简化模型的假设。如果在时刻 [公式] 的隐藏状态是 [公式] ,而对应的观察状态为 [公式] ,则该时刻观察状态 [公式] 在隐藏状态 [公式] 下生成的概率 [公式] 满足:

这样 [公式] 可以组成观测状态生成的概率矩阵 [公式] :

除此之外,我们需要一组在时刻 [公式] 的隐藏状态概率分布 [公式] :

一个HMM模型,可以由隐藏状态初始概率分布 [公式] ,状态转移概率矩阵 [公式] 和观测状态概率矩阵 [公式] 决定。 [公式] 决定状态序列, [公式] 决定观测序列。因此,HMM模型可以由一个三元组 [公式] 表示如下: [公式]

1.3 一个HMM模型实例

下面我们用一个简单的实例来描述上面抽象出的HMM模型。这是一个盒子与球的模型,例子来源于李航的《统计学习方法》。

假设我们有3个盒子,每个盒子里都有红色和白色两种球,这三个盒子里球的数量分别是:

按照下面的方法从盒子里抽球,开始的时候,从第一个盒子抽球的概率是0.2,从第二个盒子抽球的概率是0.4,从第三个盒子抽球的概率是0.4。以这个概率抽一次球后,将球放回。然后从当前盒子转移到下一个盒子进行抽球。

规则是:如果当前抽球的盒子是第一个盒子,则以0.5的概率仍然留在第一个盒子继续抽球,以0.2的概率去第二个盒子抽球,以0.3的概率去第三个盒子抽球。如果当前抽球的盒子是第二个盒子,则以0.5的概率仍然留在第二个盒子继续抽球,以0.3的概率去第一个盒子抽球,以0.2的概率去第三个盒子抽球。如果当前抽球的盒子是第三个盒子,则以0.5的概率仍然留在第三个盒子继续抽球,以0.2的概率去第一个盒子抽球,以0.3的概率去第二个盒子抽球。如此下去,直到重复三次,得到一个球的颜色的观测序列:

O={红,白,红}

注意在这个过程中,观察者只能看到球的颜色序列,却不能看到球是从哪个盒子里取出的。

那么按照我们上一节HMM模型的定义,我们的观察集合是:

V={红,白},M=2

我们的状态集合是:

Q={盒子1,盒子2,盒子3},N=3

而观察序列和状态序列的长度为3.

初始状态分布为:

[公式]

状态转移概率分布矩阵为:

观测状态概率矩阵为:

1.4 HMM观测序列的生成

从上一节的例子,我们也可以抽象出HMM观测序列生成的过程。

1.5 HMM模型的三个基本问题

HMM模型一共有三个经典的问题需要解决:

1) 评估观察序列概率。即给定模型 [公式] 和观察序列 [公式] ,计算在模型 [公式] 下观测序列 [公式] 出现的概率 [公式] 。这个问题的求解需要用到前向后向算法,这个问题是HMM模型三个问题中最简单的。

2)模型参数学习问题。即给定观测序列 [公式] ,估计模型 [公式] 的参数,使该模型下观测序列的条件概率 [公式] 最大。这个问题的求解需要用到基于EM算法的鲍姆-韦尔奇算法, 这个问题是HMM模型三个问题中最复杂的。

3)预测问题,也称为解码问题。即给定模型 [公式] 和观测序列 [公式] ,求给定观测序列条件下,最可能出现的对应的状态序列,这个问题的求解需要用到基于动态规划的维特比算法,这个问题是HMM模型三个问题中复杂度居中的算法。

2 前向后向算法评估观察序列概率

2.1 回顾HMM问题一:求观测序列的概率

首先我们回顾下HMM模型的问题一。这个问题是这样的。我们已知HMM模型的参数 [公式] 。其中A是隐藏状态转移概率的矩阵,B是观测状态生成概率的矩阵, [公式] 是隐藏状态的初始概率分布。同时我们也已经得到了观测序列 [公式] ,现在我们要求观测序列O在模型λ下出现的条件概率 [公式] 。

乍一看,这个问题很简单。因为我们知道所有的隐藏状态之间的转移概率和所有从隐藏状态到观测状态生成概率,那么我们是可以暴力求解的。

虽然上述方法有效,但是如果我们的隐藏状态数 [公式] 非常多的那就麻烦了,此时我们预测状态有 [公式] 种组合,算法的时间复杂度是 [公式] 阶的。因此对于一些隐藏状态数极少的模型,我们可以用暴力求解法来得到观测序列出现的概率,但是如果隐藏状态多,则上述算法太耗时,我们需要寻找其他简洁的算法。

前向后向算法就是来帮助我们在较低的时间复杂度情况下求解这个问题的。

2.2 用前向算法求HMM观测序列的概率

前向后向算法是前向算法和后向算法的统称,这两个算法都可以用来求HMM观测序列的概率。我们先来看看前向算法是如何求解这个问题的。

前向算法本质上属于动态规划的算法,也就是我们要通过找到局部状态递推的公式,这样一步步的从子问题的最优解拓展到整个问题的最优解。

在前向算法中,通过定义“前向概率”来定义动态规划的这个局部状态。什么是前向概率呢, 其实定义很简单:定义时刻 [公式] 时隐藏状态为 [公式] ,观测状态的序列为 [公式] 的概率为前向概率。记为:

既然是动态规划,我们就要递推了,假设已经找到了在时刻t时各个隐藏状态的前向概率,现在需要递推出时刻t+1时各个隐藏状态的前向概率。

我们的动态规划从时刻1开始,到时刻 [公式] 结束,由于 [公式] 表示在时刻 [公式] 观测序列为 [公式] ,并且时刻 [公式] 隐藏状态 [公式] 的概率,只要将所有隐藏状态对应的概率相加,即 [公式] 就得到了在时刻 [公式] 观测序列为 [公式] 的概率。

下面总结下前向算法。

从递推公式可以看出,我们的算法时间复杂度是 [公式] ,比暴力解法的时间复杂度 [公式] 少了几个数量级。

2.3 HMM前向算法求解实例

这里用盒子与球的例子来显示前向概率的计算。

首先计算时刻1三个状态的前向概率:

时刻1是红色球,隐藏状态是盒子1的概率为:

隐藏状态是盒子2的概率为:

隐藏状态是盒子3的概率为:

现在可以开始递推了,首先递推时刻2三个状态的前向概率:

时刻2是白色球,隐藏状态是盒子1的概率为:

隐藏状态是盒子2的概率为:

隐藏状态是盒子3的概率为:

继续递推,现在我们递推时刻3三个状态的前向概率:

时刻3是红色球,隐藏状态是盒子1的概率为:

隐藏状态是盒子2的概率为:

隐藏状态是盒子3的概率为:

2.4 用后向算法求HMM观测序列的概率

熟悉了用前向算法求HMM观测序列的概率,现在我们再来看看怎么用后向算法求HMM观测序列的概率。

后向算法和前向算法非常类似,都是用的动态规划,唯一的区别是选择的局部状态不同,后向算法用的是“后向概率”,那么后向概率是如何定义的呢?

这样我们得到了后向概率的递推关系式如下:

现在我们总结下后向算法的流程,注意下和前向算法的相同点和不同点:

此时我们的算法时间复杂度仍然是 [公式] 。

2.5 HMM常用概率的计算

利用前向概率和后向概率,可以计算出HMM中单个状态和两个状态的概率公式。

上面这些常用的概率值在求解HMM问题二,即求解HMM模型参数的时候需要用到。

3 鲍姆-韦尔奇算法求解HMM参数

HMM模型参数求解的问题,这个问题在HMM三个问题里算是最复杂的。

3.1 HMM模型参数求解概述

HMM模型参数求解根据已知的条件可以分为两种情况。

可见第一种情况下求解模型还是很简单的。但是在很多时候,我们无法得到HMM样本观察序列对应的隐藏序列,只有 [公式] 个长度为 [公式] 的观测序列,即 [公式] 是已知的,此时我们能不能求出合适的HMM模型参数呢?这就是我们的第二种情况,也是我们本文要讨论的重点。

它的解法最常用的是鲍姆-韦尔奇算法,其实就是基于EM算法的求解,只不过鲍姆-韦尔奇算法出现的时代,EM算法还没有被抽象出来,所以我们本文还是说鲍姆-韦尔奇算法。

3.2 鲍姆-韦尔奇算法原理

鲍姆-韦尔奇算法原理既然使用的就是EM算法的原理,那么我们需要在E步求出联合分布 [公式] 基于条件概率 [公式] 的期望,其中 [公式] 为当前的模型参数,然后再M步最大化这个期望,得到更新的模型参数 [公式] 。接着不停的进行EM迭代,直到模型参数的值收敛为止。

3.3 鲍姆-韦尔奇算法的推导

利用前向概率的定义可得:

现在我们来看看A的迭代公式求法。方法和 [公式] 的类似。由于A只在最大化函数式中括号里的第二部分出现,而这部分式子可以整理为:

现在我们来看看B的迭代公式求法。方法和 [公式] 的类似。由于B只在最大化函数式中括号里的第三部分出现,而这部分式子可以整理为:

3.4 鲍姆-韦尔奇算法流程总结

这里我们概括总结下鲍姆-韦尔奇算法的流程。

以上就是鲍姆-韦尔奇算法的整个过程。

4 维特比算法解码隐藏状态序列

HMM模型最后一个问题的求解,即给定模型和观测序列,求给定观测序列条件下,最可能出现的对应的隐藏状态序列。

HMM模型的解码问题最常用的算法是维特比算法,当然也有其他的算法可以求解这个问题。同时维特比算法是一个通用的求序列最短路径的动态规划算法,也可以用于很多其他问题,比如用维特比算法来做分词。

本文关注于用维特比算法来解码HMM的的最可能隐藏状态序列。

4.1 HMM最可能隐藏状态序列求解概述

一个可能的近似解法是求出观测序列 [公式] 在每个时刻 [公式] 最可能的隐藏状态 [公式] ,然后得到一个近似的隐藏状态序列 [公式] 。在给定模型 [公式] 和观测序列 [公式] 时,在时刻 [公式] 处于状态 [公式] 的概率 [公式] ,这个概率可以通过HMM的前向算法与后向算法计算。这样我们有:

近似算法很简单,但是却不能保证预测的状态序列的整体是最可能的状态序列,因为预测的状态序列中某些相邻的隐藏状态可能存在转移概率为0的情况。

而维特比算法可以将HMM的状态序列作为一个整体来考虑,避免近似算法的问题,下面我们来看看维特比算法进行HMM解码的方法。

4.2 维特比算法概述

维特比算法是一个通用的解码算法,是基于动态规划的求序列最短路径的方法。

既然是动态规划算法,那么就需要找到合适的局部状态,以及局部状态的递推公式。在HMM中,维特比算法定义了两个局部状态用于递推。

4.3 维特比算法流程总结

现在我们来总结下维特比算法的流程:

4.4 HMM维特比算法求解实例

下面用盒子与球的例子来看看HMM维特比算法求解。

按照我们上一节的维特比算法,首先需要得到三个隐藏状态在时刻1时对应的各自两个局部状态,此时观测状态为1:

现在开始递推三个隐藏状态在时刻2时对应的各自两个局部状态,此时观测状态为2:

4.5 HMM模型维特比算法总结

维特比算法的核心是定义动态规划的局部状态与局部递推公式,这一点在中文分词维特比算法和HMM的维特比算法是相同的,也是维特比算法的精华所在。

维特比算法也是寻找序列最短路径的一个通用方法,和dijkstra算法有些类似,但是dijkstra算法并没有使用动态规划,而是贪心算法。同时维特比算法仅仅局限于求序列最短路径,而dijkstra算法是通用的求最短路径的方法。

5 用 [公式] 学习隐马尔可夫模型 [公式]

5.1 hmmlearn概述

hmmlearn安装很简单,"pip install hmmlearn"即可完成。

hmmlearn实现了三种HMM模型类,按照观测状态是连续状态还是离散状态,可以分为两类。GaussianHMM和GMMHMM是连续观测状态的HMM模型,而MultinomialHMM是离散观测状态的模型,也是我们在HMM原理系列篇里面使用的模型。

对于MultinomialHMM的模型,使用比较简单,"startprob_"参数对应我们的隐藏状态初始分布

[公式] , "transmat_"对应我们的状态转移矩阵A, "emissionprob_"对应我们的观测状态概率矩阵B。

对于连续观测状态的HMM模型,GaussianHMM类假设观测状态符合高斯分布,而GMMHMM类则假设观测状态符合混合高斯分布。一般情况下我们使用GaussianHMM即高斯分布的观测状态即可。以下对于连续观测状态的HMM模型,我们只讨论GaussianHMM类。

在GaussianHMM类中,"startprob_"参数对应我们的隐藏状态初始分布 [公式] , "transmat_"对应我们的状态转移矩阵A, 比较特殊的是观测状态概率的表示方法,此时由于观测状态是连续值,我们无法像MultinomialHMM一样直接给出矩阵B。而是采用给出各个隐藏状态对应的观测状态高斯分布的概率密度函数的参数。

如果观测序列是一维的,则观测状态的概率密度函数是一维的普通高斯分布。如果观测序列是

N维的,则隐藏状态对应的观测状态的概率密度函数是N维高斯分布。高斯分布的概率密度函数参数可以用μ表示高斯分布的期望向量,Σ表示高斯分布的协方差矩阵。在GaussianHMM类中,“means”用来表示各个隐藏状态对应的高斯分布期望向量μ形成的矩阵,而“covars”用来表示各个隐藏状态对应的高斯分布协方差矩阵Σ形成的三维张量。


这些问题的答案在很大程度上取决于要建模的资产类别,时间范围的选择以及所使用数据的性质。 

模拟数据

在本节中,从独立的高斯分布中生成模拟的收益率数据,每个分布都代表“看涨”或“看涨”的市场机制。看涨收益来自均值正且方差低的高斯分布,而看跌收益来自均值略为负但方差较高的高斯分布。

第一个任务是安装depmixS4和quantmod库,然后将它们导入R。

install.packages('depmixS4')
install.packages('quantmod')
library('depmixS4')
library('quantmod')
set.seed(1)


课程

R语言数据分析挖掘必知必会

从数据获取和清理开始,有目的的进行探索性分析与可视化。让数据从生涩的资料,摇身成为有温度的故事。

立即参加

牛市分布N(0.1,0.1)而空头市场分布为N(−0.05,0.2)。通过以下代码设置参数:

#创建牛市和
熊市收益分布
Nk_lower <- 50
Nk_upper <- 150
bull_mean <- 0.1
bull_var <- 0.1
bear_mean <- -0.05
bear_var <- 0.2

 所述Nk值是随机选择的:

# 为每个方案创建时间列表(以天为单位)
days <- replicate(5, sample(Nk_lower:Nk_upper, 1))

第k个周期的收益是随机抽取的:

# 创建各种牛市和熊市收益
market_bull_1 <- rnorm( days[1], bull_mean, bull_var ) 
market_bear_2 <- rnorm( days[2], bear_mean, bear_var ) 
market_bull_3 <- rnorm( days[3], bull_mean, bull_var ) 
market_bear_4 <- rnorm( days[4], bear_mean, bear_var ) 
market_bull_5 <- rnorm( days[5], bull_mean, bull_var )

创建真实状态

# 创建真实的机制状态列表和完整的收益列表
true_regimes <- c( rep(1,days[1]), rep(2,days[2]), rep(1,days[3]), rep(2,days[4]), rep(1,days[5]))
returns <- c( market_bull_1, market_bear_2, market_bull_3, market_bear_4, market_bull_5)

绘制收益图显示机制切换之间均值和方差的明显变化:

plot(returns, type="l", xlab='', ylab="Returns") 

在此阶段,可以使用Expectation Maximization算法指定隐马尔可夫模型并进行拟合:

在模型拟合之后,可以绘制处于特定状态的后验概率。post_probs包含后验概率。

财务数据

在本节中,将执行两个单独的建模任务。第一种将使HMM具有两个机制状态以拟合S&P500收益率,而第二个将利用三个状态。比较两个模型之间的结果。

使用quantmod库下载:绘制时间序列:

plot(gspcRets)

使用EM算法拟合隐马尔可夫模型。每种方案的收益率和后验概率作图:

请注意,在2004年和2007年期间,市场较为平稳,因此在此期间,隐马尔可夫模型第二种机制的可能性较高。然而,在2007年至2009年之间发生次贷危机

市场在2010年变得较为平静,但在2011年又出现了更多动荡,这导致HMM再次给第一类机制带来了较高的后验概率。2011年之后,市场再次趋于平静,HMM始终给第二种机制以高概率。2015年,市场再次变得更加混乱,这反映在HMM机制之间的切换增加。

由于该模型考虑三个单独的机制,因此在2004-2007年的平静时期导致了机制2和机制3之间的转换。但是,在2008、2010和2011年的动荡时期,机制1主导着后验概率,表明高度波动状态。在2011年之后,模型恢复为在机制2和机制3之间切换。


可下载资源

关于作者

Kaizong Ye拓端研究室(TRL)的研究员。在此对他对本文所作的贡献表示诚挚感谢,他在上海财经大学完成了统计学专业的硕士学位,专注人工智能领域。擅长Python.Matlab仿真、视觉处理、神经网络、数据分析。

本文借鉴了作者最近为《R语言数据分析挖掘必知必会 》课堂做的准备。

​非常感谢您阅读本文,如需帮助请联系我们!

 
QQ在线咨询
售前咨询热线
15121130882
售后咨询热线
0571-63341498

关注有关新文章的微信公众号


永远不要错过任何见解。当新文章发表时,我们会通过微信公众号向您推送。

技术干货

最新洞察

This will close in 0 seconds