## 隐马尔可夫模型的简单应用
## 状态序列:隐藏马尔科夫链随机生成的状态序列，称为i状态序列
## 观测序列：每个状态生成一个观测，由此产生的观测是随机的序列，称为观测序列
## 隐马尔科夫模型由，初始概率分布、状态转移概率分布以及观测概率分布确定
## 问题：我们有三个盒子，每个盒子都放有红球和白球
## 问题来自《统计学习方法》
## HMM - Hidden Markov Models
library(HMM)
## 所有可能状态集合
state <- c("1","2","3")
## 所有可能观测集合
Symbols <- c("红","白")
## 初始概率分布：第一次抽取到相应盒子的概率
startprob <- c(0.2,0.4,0.4)
## 状态转移概率矩阵
transpro <- matrix(c(0.5, 0.3, 0.2, 0.2, 0.5, 0.3,0.3,0.2,0.5),nrow = 3)
transpro
## [,1] [,2] [,3]
## [1,] 0.5 0.2 0.3
## [2,] 0.3 0.5 0.2
## [3,] 0.2 0.3 0.5
## 观测概率分布：每个盒子里不同颜色球的概率
emisspr <- matrix(c(0.5,0.4,0.7,0.5,0.6,0.3),nrow = 3)
emisspr
## [,1] [,2]
## [1,] 0.5 0.5
## [2,] 0.4 0.6
## [3,] 0.7 0.3
## 初始化隐马尔科夫模型
hmm <- initHMM(States = state,Symbols = Symbols,startProbs = startprob,transProbs = transpro,
emissionProbs = emisspr)
## 的到的一个观测
ober <- c("红","白","红","红")
## 1:概率计算问题
## 使用向前算法的到相应观测的概率
## 包含向前概率的矩阵。第一个维度指的是状态和第二维度的时间。
exp(forward(hmm, ober))
## index
## states 1 2 3 4
## 1 0.10 0.0770 0.041870 0.02107790
## 2 0.16 0.1104 0.035512 0.01679232
## 3 0.28 0.0606 0.052836 0.03225698
## 状态求和，得到每次个时间为该序列的概率
colSums(exp(forward(hmm, ober)))
## 1 2 3 4
## 0.5400000 0.2480000 0.1302180 0.0701272
## 使用向后算法的到相应观测的概率
## 包含向后概率的矩阵。第一个维度指的是状态和第二维度的时间。
exp(backward(hmm, ober))
## index
## states 1 2 3 4
## 1 0.132638 0.2939 0.54 1
## 2 0.140463 0.2588 0.49 1
## 3 0.122819 0.3123 0.57 1
colSums(exp(backward(hmm, ober)))
## 1 2 3 4
## 0.39592 0.86500 1.60000 3.00000
## 预测问题：给定观测序列，求最有可能的对应的状态序列
## 维特比算法求解
## 输出：一个字符串向量，包含最可能的状态路径。
viterbi(hmm,ober)
## [1] "3" "3" "3" "3"
## 此函数计算在一个给定的序列的观测和给定的Hidden Markov模型在时间k在状态x的后验概率。
posterior(hmm,ober)
## index
## states 1 2 3 4
## 1 0.1891392 0.3227036 0.3224113 0.3005667
## 2 0.3204759 0.4074242 0.2481331 0.2394552
## 3 0.4903849 0.2698722 0.4294556 0.4599782
## 从结果可以得出，最有可能的状态序列为：3，2，3，3
colSums(posterior(hmm,ober))
## 1 2 3 4
## 1 1 1 1
## 给定一个隐马尔科夫模型的模拟
simHMM(hmm,4)
## $states
## [1] "3" "2" "2" "1"
##
## $observation
## [1] "红" "白" "白" "红"
## 学习问题：已知观测序列，学习模型的参数
hmm <- initHMM(States = state,Symbols = Symbols,startProbs = startprob,
transProbs = transpro,emissionProbs = emisspr)
hmm
## $States
## [1] "1" "2" "3"
##
## $Symbols
## [1] "红" "白"
##
## $startProbs
## 1 2 3
## 0.2 0.4 0.4
##
## $transProbs
## to
## from 1 2 3
## 1 0.5 0.2 0.3
## 2 0.3 0.5 0.2
## 3 0.2 0.3 0.5
##
## $emissionProbs
## symbols
## states 红 白
## 1 0.5 0.5
## 2 0.4 0.6
## 3 0.7 0.3
ober <- c(sample(c(rep("红",100),rep("白",200))),
sample(c(rep("白",300),rep("红",100))))
bw <- baumWelch(hmm,observation = ober)
bw$hmm
## $States
## [1] "1" "2" "3"
##
## $Symbols
## [1] "红" "白"
##
## $startProbs
## 1 2 3
## 0.2 0.4 0.4
##
## $transProbs
## to
## from 1 2 3
## 1 0.5162805 0.2293114 0.2544081
## 2 0.2987254 0.5283639 0.1729107
## 3 0.2226074 0.4078684 0.3695243
##
## $emissionProbs
## symbols
## states 红 白
## 1 0.2942211 0.7057789
## 2 0.2399061 0.7600939
## 3 0.3448797 0.6551203
## 连续迭代的差异
bw$difference
## [1] 0.711033336 0.006630626 0.006323676 0.006074176 0.005842350
## [6] 0.005626430 0.005424877 0.005236343 0.005059639 0.004893716
## [11] 0.004737642 0.004590589 0.004451817 0.004320666 0.004196539
## [16] 0.004078903 0.003967275 0.003861218 0.003760336 0.003664268
## [21] 0.003572688 0.003485297 0.003401820 0.003322007 0.003245628
## [26] 0.003172473 0.003102346 0.003035068 0.002970474 0.002908409
## [31] 0.002848733 0.002791312 0.002736026 0.002682760 0.002631408
## [36] 0.002581872 0.002534060 0.002487884 0.002443265 0.002400128
## [41] 0.002358401 0.002318018 0.002278918 0.002241040 0.002204332
## [46] 0.002168739 0.002134215 0.002100712 0.002068187 0.002036598
## [51] 0.002005908 0.001976079 0.001947077 0.001918868 0.001891421
## [56] 0.001864706 0.001838696 0.001813363 0.001788683 0.001764630
## [61] 0.001741182 0.001718317 0.001696014 0.001674253 0.001653014
## [66] 0.001632281 0.001612035 0.001592261 0.001572941 0.001554062
## [71] 0.001535608 0.001517565 0.001499921 0.001482663 0.001465778
## [76] 0.001449255 0.001433083 0.001417251 0.001401748 0.001386564
## [81] 0.001371691 0.001357119 0.001342838 0.001328841 0.001315120
## [86] 0.001301666 0.001288472 0.001275531 0.001262836 0.001250379
## [91] 0.001238155 0.001226157 0.001214379 0.001202815 0.001191460
## [96] 0.001180307 0.001169353 0.001158591 0.001148016 0.001137625