HMM算法学习笔记

本文将讲解HMM的各种算法运用

背景：有３个色子ABC，编号是1到m（此例中取3），每个色子可掷1到6，但分布律各不相等，现在我掷n次，每次随机取一个色子随机掷得一个点数并记为X[i]，但是我不知道每次掷的色子是哪个色子。

概念引入：

隐藏状态数：m，现在是3种色子，所以隐藏状态数是3种

输出状态数：s，可掷到的数是1到6，所以输出状态数是6种

隐藏状态链：取得的色子编号序列

可见状态链：X[i]序列（i从1到n）

初始向量：向量π中元素π[i]表示首次取色子时，取到色子编号为i的概率

输出矩阵：矩阵B中元素B[j][k]表示这次抽到的色子是M[j]，掷得点数是k的概率

转移矩阵：矩阵A中元素A[i][j]表示上次抽到的色子是M[i]，这次抽到的色子是M[j]的概率

 

问题一：最简单的动态规划

已知：隐藏状态数，输出状态数，可见状态链，初始向量，转移矩阵，输出矩阵

求解：问得到这串可见状态链的概率

方法一：前向算法

设置变量：DPA[i]表示第i次抽到的色子是A且能掷出可见状态链前i个点数的概率，DPB与DPC同理

初始状态：DPA[1]=P(A)*P(A->X[1])，DPB[1]=P(B)*P(B->X[1])，DPC[1]=P(C)*P(C->X[1])

动态规划：for(int i=2;i<=n;i++)循环执行下面三行

DPA[i]={DPA[i-1]*P(A->A)+DPB[i-1]*P(B->A)+DPC[i-1]*P(C->A)}*P(A->X[i])

DPB[i]={DPA[i-1]*P(A->B)+DPB[i-1]*P(B->B)+DPC[i-1]*P(C->B)}*P(B->X[i])

DPC[i]={DPA[i-1]*P(A->C)+DPB[i-1]*P(B->C)+DPC[i-1]*P(C->C)}*P(C->X[i])

最终结果：max(DPA[n]，DPB[n]，DPC[n])

补充说明：

P(A)= π[1]，P(A->X[i])=B[1][i]，P(A->A)=A[1][1]，P(A->B)=A[1][2]，P(A->C)=A[1][3]

P(B)= π[2]，P(B->X[i])=B[2][i]，P(B->A)=A[2][1]，P(B->B)=A[2][2]，P(B->C)=A[2][3]

P(C)= π[3]，P(C->X[i])=B[3][i]，P(C->A)=A[3][1]，P(C->B)=A[3][2]，P(C->C)=A[3][3]

方法二：后向算法

设置变量：DPA[i]指第i次抽到的色子是A且能掷出可见状态链第i到n个点数的概率，DPB与DPC同理

初始状态：DPA[n]= P(A->X[n])，DPB[n]= P(B->X[n])，DPC[n]= P(C->X[n])

动态规划：for(int i=n;i>=2;i++)循环执行下面三行

DPA[i-1]={DPA[i]*P(A->A)+DPB[i]*P(A->B)+DPC[i]*P(A->C)}*P(A->X[i])

DPB[i-1]={DPA[i]*P(B->A)+DPB[i]*P(B->B)+DPC[i]*P(B->C)}*P(B->X[i])

DPC[i-1]={DPA[i]*P(C->A)+DPB[i]*P(C->B)+DPC[i]*P(C->C)}*P(C->X[i])

最终结果：max(P(A)*DPA[1]，P(B)*DPB[1]，P(C)*DPC[1])

补充说明：

P(A)= π[1]，P(A->X[i])=B[1][i]，P(A->A)=A[1][1]，P(A->B)=A[1][2]，P(A->C)=A[1][3]

P(B)= π[2]，P(B->X[i])=B[2][i]，P(B->A)=A[2][1]，P(B->B)=A[2][2]，P(B->C)=A[2][3]

P(C)= π[3]，P(C->X[i])=B[3][i]，P(C->A)=A[3][1]，P(C->B)=A[3][2]，P(C->C)=A[3][3]

 

 

 

问题二：对比前文把求和变求最大值，带记忆数组回溯的动态规划

已知：隐藏状态数，输出状态数，可见状态链，初始向量，转移矩阵，输出矩阵

求解：求最大概率对应的隐藏状态链

方法：动态规划－维特比算法

设置变量：

DP[i][j]表示第i次抽到的色子编号是j，能掷出可见状态链前i个点数的概率

res[i][j]表示第i次抽到的色子编号是j，是从从i-1状态中哪个色子转移过来

初始状态：

for(int j=1;j<=m;j++)DP[1][j]=π[i]*B[j][X[i]];

for(int i=2;i<=n;i++)for(int j=1;j<=m;j++)DP[i][j]=0;

动态规划：

for(int i=2;i<=n;i++)//第i次抽

  for(int j=1;j<=m;j++)//这一次抽到j号色子

    for(int k=1;k<=m;k++)//上一次抽的是k号色子

      If(DP[i][j]<DP[i-1][k]*A[k][j]*B[i][X[i]]){

              DP[i][j]=DP[i-1][k]*A[k][j]*B[i][X[i]];

              res[i][j]=k;

      }

最终结果：

Answer[n]=for(int i=1;i<=m;i++)max(Answer[n],DP[n][i]);

for(int i=n-1;i>=1;i--)Answer[i]=res[i][Answer[i+1]];

则Answer数组序列就是第1到n次取到的色子编号

 

 

问题三：

(引用自https://www.cnblogs.com/pinard/p/6972299.html)

上面是已知部分隐藏链的直接最大似然，下面是不知隐藏链的用鲍姆-韦尔奇算法

 

已知：隐藏状态数（事先人为预判），输出状态数（从可见状态链中数出来），可见状态链（相当于唯一已知条件）

求解：求初始向量，转移矩阵，输出矩阵

方法：鲍姆-韦尔奇算法

补充说明：

π[i]是初始向量，aij即A[i][j]转移矩阵，bjk即B[j][k]输出矩阵，在第一步随机赋值

一共有D组样本(用d遍历），每组样本各抽掷色子T次(用t遍历）

表示第d组样本第t次掷色子掷到点数是i概率，其实就是B[j][k]的预测值

表示第d组样本第t次取色子取到j号色子概率，已知上次取到i号色子，其实就是A[i][j]预测值

表示第d组样本第t次掷色子是否掷到vk

此处的O[t]相当于前文设的X[i]可见状态链，vk就是前面bjk中规定输出点数是k

为什么要写成vk，那是因为我前文规定点数取值范围是1开始的自然数而已，但实际上色子可掷得的点数取值范围可以是234567，这样就vk=k+1，k仍是取1到6，写成vk可以更加抽象，更加一般性

 

 

补充：鲍姆-韦尔奇算法抽象起来叫EM算法思想，在聚类时就是这种思想

不断迭代，先E步，再M，下面列举一个简单的EM算法例子演算过程

例子引用自：https://zhuanlan.zhihu.com/p/40991784