首页 > 其他分享 >R语言中的隐马尔可夫HMM模型实例|附代码数据

R语言中的隐马尔可夫HMM模型实例|附代码数据

时间:2023-06-10 23:57:10浏览次数:84  
标签:## 模型 iteration HMM 马尔可夫 实例 软糖

原文链接:http://tecdat.cn/?p=17592

最近我们被客户要求撰写关于隐马尔可夫HMM模型的研究报告,包括一些图形和统计输出。

最近,我们使用隐马尔可夫模型开发了一种解决方案,并被要求解释这个方案

HMM用于建模数据序列,无论是从连续概率分布还是从离散概率分布得出的。它们与状态空间和高斯混合模型相关,因为它们旨在估计引起观测的状态。状态是未知或“隐藏”的,并且HMM试图估计状态,类似于无监督聚类过程。

例子

在介绍HMM背后的基本理论之前,这里有一个示例,它将帮助您理解核心概念。有两个骰子和一罐软糖。B掷骰子,如果总数大于4,他会拿几颗软糖再掷一次。如果总数等于2,则他拿几把软糖,然后将骰子交给A。现在该轮到A掷骰子了。如果她的掷骰大于4,她会吃一些软糖,但是她不喜欢黑色的其他颜色(两极分化的看法),因此我们希望B会比A多。他们这样做直到罐子空了。

现在假设A和B在不同的房间里,我们看不到谁在掷骰子。取而代之的是,我们只知道后来吃了多少软糖。我们不知道颜色,仅是从罐子中取出的软糖的最终数量。我们怎么知道谁掷骰子?HMM。

在此示例中,状态是掷骰子的人,A或B。观察结果是该回合中吃了多少软糖。如果该值小于4,骰子的掷骰和通过骰子的条件就是转移概率。由于我们组成了这个示例,我们可以准确地计算出转移概率,即1/12。没有条件说转移概率必须相同,例如A掷骰子2时可以将骰子移交给他,例如,概率为1/36。

模拟

首先,我们将模拟该示例。B平均要吃12颗软糖,而A则需要4颗。

 
# 设置
simulate <- function(N, dice.val = 6, jbns, switch.val = 4){

  #模拟变量
    #可以只使用一个骰子样本
    #不同的机制,例如只丢1个骰子,或任何其他概率分布
  
    b<- sample(1:dice.val, N, replace = T) + sample(1:dice.val, N, replace = T)
    a <- sample(1:dice.val, N, replace = T) + sample(1:dice.val, N, replace = T)
    bob.jbns <- rpois(N, jbns[1])
    alice.jbns <- rpois(N, jbns[2])

    # 状态 
    draws <- data.frame(state = rep(NA, N), obs = rep(NA, N), 





    # 返回结果
    return(cbind(roll = 1:N, draws))


# 模拟场景


draws <- simulate(N, jbns = c(12, 4), switch.val = 4)

# 观察结果

ggplot(draws, aes(x = roll, y = obs)) + geom_line()

图片

如您所见,仅检查一系列计数来确定谁掷骰子是困难的。我们将拟合HMM。由于我们正在处理计数数据,因此观察值是从泊松分布中得出的。

 
fit.hmm <- function(draws){

  # HMM 
  mod <- fit(obs ~ 1, data = draws, nstates = 2, family = poisson()



  # 通过估计后验来预测状态
  est.states <- posterior(fit.mod)
  head(est.states)

  # 结果
 
  hmm.post.df <- melt(est.states, measure.vars = 


  # 输出表格
  print(table(draws[,c("state", "est.state.labels")]))
 
## iteration 0 logLik: -346.2084 
## iteration 5 logLik: -274.2033 
## converged at iteration 7 with logLik: -274.2033 
##        est.state.labels
## state   alice bob
##   a     49   2
##   b      3  46

模型迅速收敛。使用后验概率,我们估计过程处于哪个状态,即谁拥有骰子,A或B。要具体回答该问题,我们需要更多地了解该过程。在这种情况下,我们知道A只喜欢黑软糖。否则,我们只能说该过程处于状态1或2。下图显示了HMM很好地拟合了数据并估计了隐藏状态。

 
# 绘图输出


    g0 <- (ggplot(model.output$draws, aes(x = roll, y = obs)) + geom_line() +
        theme(axis.ticks = element_blank(), axis.title.y = element_blank())) %>% ggplotGrob
    g1 <- (ggplot(model.output$draws, aes(x = roll, y = state, fill = state, col = state)) + 
      


    g0$widths <- g1$widths
    return(grid.arrange(g0, g1


plot.hmm.output(hmm1)

图片

令人印象深刻的是,该模型拟合数据和滤除噪声以估计状态的良好程度。公平地说,可以通过忽略时间分量并使用EM算法来估计状态。但是,由于我们知道数据形成一个序列,因为观察下一次发生的概率取决于前一个即\(P(X_t | X_ {t-1})\),其中\(X_t \ )是软糖的数量。


点击标题查阅往期内容

图片

隐马尔可夫模型(HMM)识别不断变化的股市状况股票指数预测实战

图片

左右滑动查看更多

图片

01

图片

02

图片

03

图片

04

图片

考虑到我们构造的问题,这可能是一个相对简单的案例。如果转移概率大得多怎么办?

 
 simulate(100, jbns = c(12, 4), switch.val = 7)
 
## iteration 0 logLik: -354.2707 
## iteration 5 logLik: -282.4679 
## iteration 10 logLik: -282.3879 
## iteration 15 logLik: -282.3764 
## iteration 20 logLik: -282.3748 
## iteration 25 logLik: -282.3745 
## converged at iteration 30 with logLik: -282.3745 
##        est.state.labels
## state   alice bob
##   alice    54   2
##   bob       5  39
 
plot(hmm2)

图片这有很多噪音数据,但是HMM仍然做得很好。性能的提高部分归因于我们对从罐中取出的软糖数量的选择。分布越明显,模型就越容易拾取转移。公平地讲,我们可以计算中位数,并将所有低于中位数的值都归为一个状态,而将所有高于中位数的值归为另一状态,您可以从结果中看到它们做得很好。这是因为转移概率非常高,并且预计我们会从每个状态观察到相似数量的观察结果。当转移概率不同时,我们会看到HMM表现更好。

如果观察结果来自相同的分布,即A和B吃了相同数量的软糖怎么办?

 
hmm3 <- fit.hmm(draws)
plot(hmm3)

图片不太好,但这是可以预期的。如果从中得出观察结果的分布之间没有差异,则可能也只有1个状态。

实际如何估算状态?

首先,状态数量及其分布方式本质上是未知的。利用对系统建模的知识,用户可以选择合理数量的状态。在我们的示例中,我们知道有两种状态使事情变得容易。可能知道确切的状态数,但这并不常见。再次通过系统知识来假设观察结果通常是合理的,这通常是合理的。

从这里开始,使用 Baum-Welch算法 来估计参数,这是EM算法的一种变体,它利用了观测序列和Markov属性。除了估计状态的参数外,还需要估计转移概率。Baum-Welch算法首先对数据进行正向传递,然后进行反向传递。然后更新状态转移概率。然后重复此过程,直到收敛为止。

在现实世界

在现实世界中,HMM通常用于

  • 股票市场预测,无论市场处于牛市还是熊市
  • 估计NLP中的词性
  • 生物测序
  • 序列分类

仅举几例。只要有观察序列,就可以使用HMM,这对于离散情况也适用。

图片

点击文末 “阅读原文”

获取全文完整资料。

本文选自《R语言中的隐马尔可夫HMM模型实例》。

点击标题查阅往期内容

python中使用马尔可夫决策过程(MDP)动态编程来解决最短路径强化学习问题
隐马尔可夫模型(HMM)识别不断变化的股市状况股票指数预测实战
马尔可夫Markov区制转移模型分析基金利率
马尔可夫区制转移模型Markov regime switching
时变马尔可夫区制转换MRS自回归模型分析经济时间序列
马尔可夫转换模型研究交通伤亡人数事故时间序列预测
如何实现马尔可夫链蒙特卡罗MCMC模型、Metropolis算法?
Matlab用BUGS马尔可夫区制转换Markov switching随机波动率模型、序列蒙特卡罗SMC、M H采样分析时间序列
R语言BUGS序列蒙特卡罗SMC、马尔可夫转换随机波动率SV模型、粒子滤波、Metropolis Hasting采样时间序列分析
matlab用马尔可夫链蒙特卡罗 (MCMC) 的Logistic逻辑回归模型分析汽车实验数据
stata马尔可夫Markov区制转移模型分析基金利率
PYTHON用时变马尔可夫区制转换(MRS)自回归模型分析经济时间序列
R语言使用马尔可夫链对营销中的渠道归因建模
matlab实现MCMC的马尔可夫转换ARMA - GARCH模型估计
R语言隐马尔可夫模型HMM识别不断变化的股票市场条件
R语言中的隐马尔可夫HMM模型实例
用机器学习识别不断变化的股市状况—隐马尔科夫模型(HMM)
Matlab马尔可夫链蒙特卡罗法(MCMC)估计随机波动率(SV,Stochastic Volatility) 模型
MATLAB中的马尔可夫区制转移(Markov regime switching)模型
Matlab马尔可夫区制转换动态回归模型估计GDP增长率
R语言马尔可夫区制转移模型Markov regime switching
stata马尔可夫Markov区制转移模型分析基金利率
R语言如何做马尔可夫转换模型markov switching model
R语言隐马尔可夫模型HMM识别股市变化分析报告
R语言中实现马尔可夫链蒙特卡罗MCMC模型

标签:##,模型,iteration,HMM,马尔可夫,实例,软糖
From: https://www.cnblogs.com/tecdat/p/17472219.html

相关文章

  • imessages数据检测,imessages过蓝检测,用applescript检测手机号码是否注册imessage实
    一、检测iMessage发送数据的2种方式:1.人工筛选,将要验证的号码输出到文件中,以逗号分隔。再将文件中的号码粘贴到iMessage客户端的地址栏,iMessage客户端会自动逐个检验该号码是否为iMessage账号,检验速度视网速而定。红色表示不是iMessage账号,蓝色表示iMessage账号。2.编写脚本控制......
  • 数据分析实例
    1、导入用于分析和可视化作图的库importpandasaspdimportmatplotlib.pyplotaspltimportseabornassns#seaborn也很强大,可以小试一下da=pd.read_csv('D:/datasource/mycrawldata/dataanylist.csv')da.head()da.columnsda.shape#查看数据的行列数,这里是(465,7)da.desc......
  • dubbo+spring+zookeeper的集成入门实例
    一、启动zookeeper我用的kafka自带的zookeeper任务管理器输入bin\windows\zookeeper-server-start.batconfig\zookeeper.properties开启zookeeperCtrl+c输入Y关闭服务 二、安装dubbo—admin管理控制台1、打开https://github.com/apache/dubbo-admin/tree/......
  • DevExpress 动态创建实例化类 (xpo)
    使用xpo(devexpress)时动态创建一个持久化类。这样方便访问数据库。/*使用DevExpress控件xpoXPObject持久化对象数据库访问表XPObject*///z2011-07-2722:06:[email protected]转载请注明出处classProgram{staticvoidMain(string[]args){XpoD......
  • boost.array 使用实例
    #include<iostream>//z包含array相关头文件。#include<boost/array.hpp>usingnamespacestd;usingnamespaceboost;//z仿函数,输出array各元素。classPrintInt{private:intsum;intcnt;public:PrintInt(intval):sum(......
  • emoji食用实例
    \(emoji\)食用实例使用方法:1.使用系统自带的emoji,快捷键\(win+;\)或\(win+.\)当然你也可以用微软的拼音输入法点击那个笑脸,就像这样......
  • 实例讲解Flink 流处理程序编程模型
    摘要:在深入了解Flink实时数据处理程序的开发之前,先通过一个简单示例来了解使用Flink的DataStreamAPI构建有状态流应用程序的过程。本文分享自华为云社区《Flink实例:Flink流处理程序编程模型》,作者:TiAmoZhang。在深入了解Flink实时数据处理程序的开发之前,先通过一个简单......
  • Wpf(Storyboard)动画简单实例
    Wpf(Storyboard)动画简单实例动画的三种变换方式RotateTransform:旋转变换变化值:CenterX围绕转的圆心横坐标      CenterY纵坐标       Angle旋转角度(角度正负表示方向) ScaleTransform:缩放变换变化值:ScaleX横向放大倍数 ScaleY纵向(负值时翻转)  TranslateTransform......
  • 用Spring MVC实现用户登录的完整实例
    用SpringMVC实现用户登录的完整实例本例子是再Eclipse中建立一个Tomcat工程,来讲解SpringMVC的全过程,实例代码如下:<一>编写日记文件放在myMVC/WEB-INF/src下#指定日志输入文件的大小log4j.appender.stdout.MaxFileSize=500KBlog4j.appender.stdout.MaxBackupI......
  • Mybatis框架及原理实例分析
    摘要本篇文章只是个人阅读mybatis源码总结的经验或者个人理解mybatis的基本轮廓,作为抛砖引玉的功能,希望对你有帮助,如果需要深入了解细节还需亲自去阅读源码。mybatis基本架构mybatis的源码应该算是比较容易阅读的,首先mybatis核心功能就是执行Sql语句,但在其基础上又有许多增强的地方......