统计学习方法 第十章 隐马尔可夫模型

隐马尔可夫模型（hidden Markov model,HMM）是可用于标注问题的统计学习模型，描述隐藏的马尔科夫链随机生成观测序列的过程，属于生成模型。由初始状态概率向量π、状态转移概率矩阵A和观测概率矩阵B决定的隐马尔可夫模型λ=（A,B,π）的3个基本问题，分别为：
1.概率计算问题，在给定模型和观测序列，求该模型下序列出现的概率；
2.学习问题，给定观测序列去求模型λ=（A,B,π）的参数，使得该模型下给定观测序列概率最大；
3.预测问题，已知模型λ=（A,B,π）和观测序列，求最有可能对应的状态序列。
对于问题1，常见的有前向算法及后向算法，问题2的常见方法有监督学习方法和Baum-Welch算法（由EM算法推导得出），问题3的常见算法有近似算法和维特比算法。
这里通过代码大致解释一下思路：
前向算法：

import numpy as np

#前向算法
def forwardProbCal(dataTransProb,dataObsProb,dataInitProb,obList):
    matA = np.mat(dataTransProb)
    matB = np.mat(dataObsProb)
    vecPi = np.mat(dataInitProb).T
    m,n = np.shape(matB)
    t = len(obList)
    #用于计算每个状态下的前向概率
    matAccum = np.zeros((m,t))
    #初始化
    matAccum[:,0] = np.multiply(vecPi,matB[:,obList[0]].A).T
    #递推计算
    for i in range(1,t):
        for j in range(0,m):
            prob = matAccum[:,i - 1] * matA[:,j]
            matAccum[j,i] = prob * matB[j,obList[i]]
    
    print(matAccum)
    resProb = matAccum[:,-1].sum()
    return resProb,matAccum

第一项为所得的观测概率，第二项为时序的前向概率矩阵。
以书上例10.2的数据为例，代入：

dataA = [[0.5,0.2,0.3],[0.3,0.5,0.2],[0.2,0.3,0.5]]
dataB = [[0.5,0.5],[0.4,0.6],[0.7,0.3]]
vecPi = [0.2,0.4,0.4]
obList = [0,1,0]
resProb = forwardProbCal(dataA,dataB,vecPi,obList)
print(resProb[0])

可以获得与书上一致的结果。
后向算法：

#后向算法
def backwardProbCal(dataTransProb,dataObsProb,dataInitProb,obList):
    matA = np.mat(dataTransProb)
    matB = np.mat(dataObsProb)
    vecPi = np.mat(dataInitProb).T
    m,n = np.shape(matB)
    t = len(obList)
    #用于计算每个状态下的后向概率
    matAccum = np.zeros((m,t))
    #初始化
    matAccum[:,-1] = 1
    for i in range(t-2,-1,-1):
        for j in range(0,m):
            matAccum[j,i] = np.multiply(np.multiply(matA[j,:].T, matB[:,obList[i + 1]]),np.mat(matAccum[:,i + 1]).T).sum()
    
    print(matAccum)
    resProb = np.multiply(np.multiply(vecPi,matB[:,obList[0]]),np.mat(matAccum[:,0]).T).sum()
    return resProb,matAccum

第一项为所得的观测概率，第二项为时序的后向概率矩阵。
用上面的数据测试，可以获得与前向算法一致的结果。
对于书上习题10.1，可以使用上面的后向算法得到结果（当然前向算法结果也一样），对于习题10.2，可以使用上面的两个算法各自得出前向概率矩阵及后向概率矩阵后使用np.multiply得到前向后向概率矩阵，再以第3个状态之和除以概率总和值(使用.sum(axis = 1))即可。
下面关于维特比算法的简单实现：

#维特比算法
def viterbi(dataTransProb,dataObsProb,dataInitProb,obList):
    matA = np.mat(dataTransProb)
    matB = np.mat(dataObsProb)
    vecPi = np.mat(dataInitProb).T
    m,n = np.shape(matB)
    t = len(obList)
    matProbCal = np.zeros((m,t))
    matRoute = np.zeros((m,t))
    #初始化
    matProbCal[:,0] = np.multiply(vecPi,matB[:,obList[0]]).T
    matRoute[:,0] = np.mat(list(range(1,m + 1)))
    for i in range(1,t):
        for j in range(0,m):
            #计算概率
            probCal = np.multiply(np.mat(matProbCal[:,i - 1]).T,matA[:,j])
            matProbCal[j,i] = np.max(probCal) * matB[j,obList[i]]
            matRoute[j,i] = np.argmax(probCal) + 1
    
    return matProbCal,matRoute

使用书上的例10.3的数据，可以获得对应的概率矩阵和路径矩阵
[[0.1 0.028 0.00756]
[0.16 0.0504 0.01008]
[0.28 0.042 0.0147 ]]
[[1. 3. 2.]
[2. 3. 2.]
[3. 3. 3.]]