【Hidden Markov Model】

隐马尔科夫模型(Hidden Markov Model)

与马尔科夫模型相比，隐马尔科夫模型有 隐藏状态 和 观测状态。两者并非对应关系，即隐藏状态数$n$和观测状态数$m$无关。

• 当前时刻隐藏状态$x_t$只与上一时刻隐藏状态$x_{t-1}$有关：
  状态转移概率$a_{ij}=P(x_t=i|x_{t-1}=j)$
• 某个观测状态$o_t$只和生成它的隐藏状态$x_t$有关：
  观测状态概率$b_j(k)=P(o_t=k|x_t=j)$
• 初始隐藏状态：
  ${\pi }_i=P(x_1=i)$

总而言之，隐马尔科夫模型表示如下：
$$HMM=(\pi ,\mathcal{A},\mathcal{B})$$
其中$\pi$为初始状态概率，$\mathcal{A}$为隐藏状态转移概率，$\mathcal{B}$为观测状态生成概率。

隐马尔可夫模型解决的问题

求解参数

已知观测状态序列 O = { o 1 , o 2 , . . . , o T } \mathbf{O} = \left \{ o_1,o_2,...,o_T \right \} O={o1​,o2​,...,oT​}、隐藏状态序列 X = { x 1 , x 2 , . . . , x T } \mathbf{X} = \left \{ x_1,x_2,...,x_T \right \} X={x1​,x2​,...,xT​}，求解隐马尔科夫模型参数$\pi ,\mathcal{A},\mathcal{B}$。

根据$\pi ,\mathcal{A},\mathcal{B}$定义，利用大数定律估计$\pi ,\mathcal{A},\mathcal{B}$即可。

滤波(filtering)

已知模型$\lambda =(\pi ,\mathcal{A},\mathcal{B})$及观测序列 O = { o 1 , o 2 , . . . , o T } \mathbf{O} = \left \{ o_1,o_2,...,o_T \right \} O={o1​,o2​,...,oT​}，求取该观测序列出现的概率$P(\mathbf{O}|\lambda )$

直接求解

为求给定模型下观测序列出现的概率，因观测状态只与其对应的隐藏状态有关，则通过所有可能的隐藏状态序列 X = { x 1 , x 2 , . . . , x T } \mathbf{X} = \left \{ x_1,x_2,...,x_T \right \} X={x1​,x2​,...,xT​}即可求取$P(\mathbf{O}|\lambda )$：
$$P(\mathbf{O}|\lambda )=\sum _{\mathbf{X}}P(\mathbf{O},\mathbf{X}|\lambda )=\sum _{\mathbf{X}}P(\mathbf{X}|\lambda )P(\mathbf{O}|\mathbf{X},\lambda )\text{\hspace{1em}(1)}$$
利用$\pi ,\mathcal{A}$可求取：
$P(\mathbf{X}|\lambda )=P(x_1|\lambda )P(x_2,...,x_T|x_1,\lambda )={\pi }_{x_1}P(x_2|x_1,\lambda )P(x_2,...,x_T|x_1,x_2,\lambda )$
利用$\mathcal{B}$可求取：
$P(\mathbf{O}|\mathbf{X},\lambda )={\prod }_{k=1}^{T}P(o_k|x_k,\lambda )$
则式$(1)$为：
$$P(\mathbf{O}|\lambda )=\sum _{\mathbf{X}}P(\mathbf{O},\mathbf{X}|\lambda )=\sum _{x_1,x_2,...,x_T}{\pi }_{x_1}{b}_{x_1}(o_1)a_{x_1{x}_2}{b}_{x_2}(o_2)a_{x_2{x}_3}...a_{x_{T-1}{x}_T}{b}_{x_T}(o_T)\text{\hspace{1em}(2)}$$
若隐藏状态数为$n$，则式$(2)$复杂度为$O(T{n}^T)$

前向算法

假设$t$时刻隐藏状态$x_t=i$，观测序列 O = { o 1 , o 2 , . . . , o t } \mathbf{O} = \left \{ o_1,o_2,...,o_t \right \} O={o1​,o2​,...,ot​}的概率被称为前向概率：
$${\alpha }_i(t)=P(o_1,o_2,...,o_t,x_t=i|\lambda )\text{\hspace{1em}(3)}$$
当$t=T$时，根据式$(3)$有${\alpha }_i(T)=P(o_1,o_2,...,o_T,x_T=i|\lambda )$，则此时可求取$P(\mathbf{O}|\lambda )$：
$$P(\mathbf{O}|\lambda )=\sum _{i=1}^n{\alpha }_i(T)\text{\hspace{1em}(4)}$$
利用${\alpha }_i(1)=P(o_1,x_t=i|\lambda )={\pi }_i{b}_i(o_1)$、${\alpha }_i(t+1)=({\sum }_{j=1}^n{\alpha }_j(t)a_{ji})b_i(o_{t+1})$可求取式$(4)$，复杂度减少为$O(T{n}^2)$

与之类似，后向算法:
$${\beta }_i(t)=P(o_{t+1},o_{t+2},...,o_T|x_t=i,\lambda )\text{\hspace{1em}(5)}$$

训练(training)

已知观测序列 O = { o 1 , o 2 , . . . , o T } \mathbf{O} = \left \{ o_1,o_2,...,o_T \right \} O={o1​,o2​,...,oT​}，求解模型参数$(\pi ,\mathcal{A},\mathcal{B})$使得$P(\mathbf{O}|\lambda )$最大

EM算法：Baum-Welch算法

期望最大化(Exoectation maximization, EM)算法可以用于含有隐变量的统计模型的参数最大似然估计，其的核心思想是：

1. 初始时随机地给模型的参数赋值λ0，该赋值遵循模型对参数的限制；
2. 根据λ0求得模型中隐变量的期望值；
3. 由该期望值又可以重新得到模型参数的新估计值λ1；
4. 循环进行该过程，直到参数收敛于最大似然估计值。

Baum-Welch算法或前向后向算法是EM算法的一种具体实现。

给定HMM的参数$\lambda =(\pi ,\mathcal{A},\mathcal{B})$及观测状态序列 O = { o 1 , o 2 , . . . , o T } \mathbf{O} = \left \{ o_1,o_2,...,o_T \right \} O={o1​,o2​,...,oT​}

时间$t$时，隐藏状态$x_t=i$的概率为：
$${\gamma }_i(t)=P(i_t=q_i|\mathbf{O},\lambda )=\frac{P(i_t=q_i,\mathbf{O}|\lambda )}{P(\mathbf{O}|\lambda )}=\frac{{\alpha }_i(t){\beta }_i(t)}{{\sum }_{j=1}^n{\alpha }_j(t){\beta }_j(t)}\text{\hspace{1em}(6)}$$

时间$t$时隐藏状态$x_t=i$、$t+1$时隐藏状态$x_{t+1}=j$的概率为：
$$\begin{gathered} {\xi }_{ij}(t)=P(x_t=i,x_{t+1}=j|\mathbf{O},\lambda ) \\ =\frac{P(x_t=i,x_{t+1}=j,\mathbf{O}|\lambda )}{P(\mathbf{O}|\lambda )} \\ =\frac{{\alpha }_i(t)a_{ij}{b}_j(o_{t+1}){\beta }_j(t+1)}{P(\mathbf{O}|\lambda )} \\ =\frac{{\alpha }_i(t)a_{ij}{b}_j(o_{t+1}){\beta }_j(t+1)}{{\sum }_{i=1}^n{\sum }_{j=1}^n{\alpha }_i(t)a_{ij}{b}_j(o_{t+1}){\beta }_j(t+1)}\text{\hspace{1em}(7)} \end{gathered}$$
根据式$(6)、(7)$可得参数估计：
$$\{\begin{array}{c}\widehat{{\pi }_i}={\gamma }_i(1)\\ \widehat{a_{ij}}=\frac{{\sum }_{t=1}^{T-1}{\xi }_{ij}(t)}{{\sum }_{t=1}^{T-1}{\gamma }_i(t)}\\ \widehat{b_j(k)}=\frac{{\sum }_{t=1,o_t=k}^T{\gamma }_i(t)}{{\sum }_{t=1}^T{\gamma }_i(t)}\end{array}\text{\hspace{1em}(8)}$$

解码(decodin)

已知模型$\lambda =(\pi ,\mathcal{A},\mathcal{B})$及观测序列 O = { o 1 , o 2 , . . . , o T } \mathbf{O} = \left \{ o_1,o_2,...,o_T \right \} O={o1​,o2​,...,oT​}，求取最有可能的状态序列 X = { x 1 , x 2 , . . . , x T } \mathbf{X} = \left \{ x_1,x_2,...,x_T \right \} X={x1​,x2​,...,xT​}，即使$P(\mathbf{X}|\mathbf{O},\lambda )$最大

Viterbi算法

对于时刻$t$，隐藏状态$x_t=i$，所有可能的情况中概率最大路径：
$${\delta }_t(i)=\underset{x_1,x_2,...,x_{t-1}}{\max }P(x_t=i,x_{t-1},...,x_1,o_t,...,o_1|\lambda )\text{\hspace{1em}(9)}$$
由此可知：
$$\begin{gathered} {\delta }_{t+1}(i)=\underset{x_1,x_2,...,x_t}{\max }P(x_{t+1}=i,x_t,...,x_1,o_{t+1},...,o_1|\lambda ) \\ =\underset{1\le j\le N}{\max }[{\delta }_t(j)a_{ji}]b_i(o_{t+1})\text{\hspace{1em}(10)} \end{gathered}$$
概率最大的路径中$t-1$时刻的隐藏状态：
$${\Psi }_t(i)=\arg \underset{1\le j\le n}{\max }[{\delta }_{t-1}(j)a_{ji}]\text{\hspace{1em}(11)}$$

1. 根据式$(10)$计算$T$时刻所有隐藏状态对应的值，若${\delta }_T(i)$为其中最大的一项，则$x_T=i$；
2. 则$T-1$时刻隐藏状态$x_{T-1}={\Psi }_T(i)=k$，$T-2$时刻时刻隐藏状态$x_{T-2}={\Psi }_{T-1}(k)$，以此类推。

算法实现

创建一个隐马尔科夫模型的类：

class HMM:
    """
    状态取值有N个，观测序列取值有M个
    """
    def __init__(self, A: np.ndarray, B: np.ndarray, pi):
        self.A = A  # 状态转移概率矩阵 [N, N]
        self.B = B  # 输出观测概率矩阵 [N, M]
        self.pi = pi  # 初设状态概率向量 [N]

        self.N = len(pi)
        self.M = self.A.shape[1]

        self.map_i = {}  # 存储状态取值对应的编码, key为编码，value为状态值
        self.map_o = {}  # 存储观测值对应的编码， key为观测值，value为编码

前向概率计算：

def forward(self, t, i, O):
    """
        前向算法
        :param t: 时刻t，以0为起点
        :param i: 状态值，[0, 1, ...., N-1]
        :param O: 观测序列
        :return:
        """
    if t == 0:
        return self.pi[i] * self.B[i][O[t]]

    return self.B[i][O[t]] * sum([self.forward(t - 1, j, O) * self.A[j][i] for j in range(self.N)])

后向概率计算：

def backward(self, t, i, O):
    """
        后向算法
        :param t: 时刻t，以0为起点
        :param i: 状态值，[0, 1, ...., N-1]
        :param O: 观测序列
        :return:
        """
    if t == len(O) - 1:
        return 1

    return sum([self.A[i][j] * self.B[j][O[t+1]] * self.backward(t+1, j, O) for j in range(self.N)])

前向后向算法
计算期望：

def gamma(self, t, i, O):
    """
    根据已知隐马尔科夫模型的参数和观测序列O，时刻t处于状态i的概率
    :param t: 时刻t，以0为起点
    :param i: 状态值，[0, 1, ...., N-1]
    :param O: 观测序列
    :return:
    """
    molecule = self.forward(t, i, O) * self.backward(t, i, O)
    denominator = sum([self.forward(t, j, O) * self.backward(t, j, O)] for j in range(self.N))
    return molecule / denominator

def xi(self, t, i, j, O):
    """
    根据已知隐马尔科夫模型的参数和观测序列O，时刻t处于状态i且时刻t+1处于状态j的概率
    :param t: 时刻t，以0为起点
    :param i: 状态值，[0, 1, ...., N-1]
    :param O: 观测序列
    :return:
    """
    molecule = self.forward(t, i, O) * self.A[i][j] * self.B[j][O[t+1]] * self.backward(t+1, j, O)
    denominator = 0
    for n in range(self.N):
        for m in range(self.N):
            denominator += self.forward(t, n, O) * self.A[n][m] * self.B[m][O[t+1]] * self.backward(t+1, m, O)

    return molecule / denominator

参数更新：

def train(self, O):
    """
    根据隐马尔科夫模型参数的估计值，极大化模型参数
    :param O:
    :return:
    """
    pi = self.pi.copy()
    A = self.A.copy()
    B = self.B.copy()

    for i in range(self.N):
        pi[i] = self.gamma(0, i, O)

        # 根据公式计算a(i,j)
        for j in range(self.N):
            A[i][j] = sum([self.xi(t, i, j, O) for t in range(len(O)-1)]
                          ) / sum([self.gamma(t, i, O) for t in range(len(O)-1)])

        # 根据公式计算b_i(k)
        B_denominator = sum([self.gamma(t, i, O) for t in range(len(O))])  # 公式中的分母
        # 求得观测序列的所有不同的取值，和对应的所有时刻。
        # 比如观测序列(1, 2, 1)，第一个和最后一个时刻的取值为1，第二个时刻的取值为2
        Map = {}
        for t in range(len(O)):
            if O[t] in Map.keys():
                Map[O[t]].append(t)
            else:
                Map[O[t]] = [t]
        # 这里k对应观测序列中的一个不同的值
        for k, tl in Map.items():
            B[i][k] = sum([self.gamma(t, i, O) for t in tl]) / B_denominator

    # 更新模型参数
    self.pi = pi
    self.A = A
    self.B = B
   
def run(self, train_nums, O):
    """
    对训练样本进行训练
    :param train_nums: 训练轮数
    :param O: 所有训练样本，二维数组
    :return:
    """
    for i in range(train_nums):
        for o in O:
            # 将观测文本值映射为对应编码
            o = [self.map_o[m] for m in o]
            self.train(o)  # 对一个观测序列训练，进行参数更新

维特比算法

def viterbi(self, O):
    """
    维特比编码求最优路径
    :param O:
    :return:
    """
    pro = np.ones([len(O), self.N])  # 记录累计概率
    backpointers = np.ones([len(O), self.N], dtype=np.int32)  # 记录最优路径

    # 计算时刻t=1，各个状态下观测序列为o_1的概率
    pro[0] = [self.pi[i] * self.B[i][O[0]] for i in range(self.N)]
    # 时刻t=2,3,....,
    for t in range(1, len(O)):
        b = [self.B[i][O[t]] for i in range(self.N)]  # 计算各个状态下观测序列为o_t的概率
        v = np.expand_dims(pro[t - 1], 1) * self.A  # 计算前一个时刻状态为j，当前时刻状态为i的概率矩阵
        pro[t] = np.max(v, 0) * b  # 当前状态为i的最大概率 * 对应状态下观测序列为o_t的概率
        backpointers[t] = np.argmax(v, 0)  # 记录，当前状态为i的概率最大时，上个时刻的状态

    # 开始回溯
    viterbi = [np.argmax(pro[-1])]  # 最终累计概率最大时，最后一个时刻的状态
    # 根据当前状态从最优路径中挑选上个时刻的状态
    for bp in reversed(backpointers[1:]):
        viterbi.append(bp[viterbi[-1]])
    viterbi.reverse()

    # 将状态编码转化为对应的状态值
    viterbi = [self.map_i[v] for v in viterbi]

    viterbi_score = np.max(pro[-1])  # 最高的累计概率
    return viterbi, viterbi_score

应用举例

拼音输入/语音识别

• 输入法把拼音看做是观察状态，需要得到的汉字为隐藏状态，这样，输入法的整句解码就变成了维特比解码，其转移概率即是二元语言模型，其输出概率即是多音字对应不同拼音的概率。
• 将上图中的拼音换成语音，就成了语音识别问题，转移概率仍然是二元语言模型，其输出概率则是语音模型，即语音和汉字的对应模型。