用hmmlearn學習隱馬爾科夫模型HMM

劉建平Pinard發表於2017-06-13
HMM馬爾科夫模型

　　　　在之前的HMM系列中，我們對隱馬爾科夫模型HMM的原理以及三個問題的求解方法做了總結。本文我們就從實踐的角度用Python的hmmlearn庫來學習HMM的使用。關於hmmlearn的更多資料在官方文件有介紹。

1. hmmlearn概述

　　　　hmmlearn安裝很簡單，"pip install hmmlearn"即可完成。

　　　　hmmlearn實現了三種HMM模型類，按照觀測狀態是連續狀態還是離散狀態，可以分為兩類。GaussianHMM和GMMHMM是連續觀測狀態的HMM模型，而MultinomialHMM是離散觀測狀態的模型，也是我們在HMM原理系列篇裡面使用的模型。

　　　　對於MultinomialHMM的模型，使用比較簡單，"startprob_"引數對應我們的隱藏狀態初始分佈$\Pi$, "transmat_"對應我們的狀態轉移矩陣$A$, "emissionprob_"對應我們的觀測狀態概率矩陣$B$。

　　　　對於連續觀測狀態的HMM模型，GaussianHMM類假設觀測狀態符合高斯分佈，而GMMHMM類則假設觀測狀態符合混合高斯分佈。一般情況下我們使用GaussianHMM即高斯分佈的觀測狀態即可。以下對於連續觀測狀態的HMM模型，我們只討論GaussianHMM類。

　　　　在GaussianHMM類中，"startprob_"引數對應我們的隱藏狀態初始分佈$\Pi$, "transmat_"對應我們的狀態轉移矩陣$A$, 比較特殊的是觀測狀態概率的表示方法，此時由於觀測狀態是連續值，我們無法像MultinomialHMM一樣直接給出矩陣$B$。而是採用給出各個隱藏狀態對應的觀測狀態高斯分佈的概率密度函式的引數。

　　　　如果觀測序列是一維的，則觀測狀態的概率密度函式是一維的普通高斯分佈。如果觀測序列是$N$維的，則隱藏狀態對應的觀測狀態的概率密度函式是$N$維高斯分佈。高斯分佈的概率密度函式引數可以用$\mu$表示高斯分佈的期望向量，$\Sigma$表示高斯分佈的協方差矩陣。在GaussianHMM類中，“means”用來表示各個隱藏狀態對應的高斯分佈期望向量$\mu$形成的矩陣，而“covars”用來表示各個隱藏狀態對應的高斯分佈協方差矩陣$\Sigma$形成的三維張量。

2. MultinomialHMM例項

　　　　下面我們用我們在HMM系列原理篇中的例子來使用MultinomialHMM跑一遍。

　　　　完整程式碼參見我的github:https://github.com/ljpzzz/machinelearning/blob/master/natural-language-processing/hmm.ipynb

　　　　首先建立HMM的模型：

import numpy as np
from hmmlearn import hmm

states = ["box 1", "box 2", "box3"]
n_states = len(states)

observations = ["red", "white"]
n_observations = len(observations)

start_probability = np.array([0.2, 0.4, 0.4])

transition_probability = np.array([
  [0.5, 0.2, 0.3],
  [0.3, 0.5, 0.2],
  [0.2, 0.3, 0.5]
])

emission_probability = np.array([
  [0.5, 0.5],
  [0.4, 0.6],
  [0.7, 0.3]
])

model = hmm.MultinomialHMM(n_components=n_states)
model.startprob_=start_probability
model.transmat_=transition_probability
model.emissionprob_=emission_probability

　　　　現在我們來跑一跑HMM問題三維特比演算法的解碼過程，使用和原理篇一樣的觀測序列來解碼，程式碼如下：

seen = np.array([[0,1,0]]).T
logprob, box = model.decode(seen, algorithm="viterbi")
print("The ball picked:", ", ".join(map(lambda x: observations[x], seen)))
print("The hidden box", ", ".join(map(lambda x: states[x], box)))

　　　　輸出結果如下：

('The ball picked:', 'red, white, red')
('The hidden box', 'box3, box3, box3')

　　　　可以看出，結果和我們原理篇中的手動計算的結果是一樣的。

　　　　也可以使用predict函式，結果也是一樣的，程式碼如下：

box2 = model.predict(seen)
print("The ball picked:", ", ".join(map(lambda x: observations[x], seen)))
print("The hidden box", ", ".join(map(lambda x: states[x], box2)))

　　　　大家可以跑一下，看看結果是否和decode函式相同。

　　　　現在我們再來看看求HMM問題一的觀測序列的概率的問題，程式碼如下：

print model.score(seen)

　　　　輸出結果是：

-2.03854530992

　　　　要注意的是score函式返回的是以自然對數為底的對數概率值，我們在HMM問題一中手動計算的結果是未取對數的原始概率是0.13022。對比一下：$$ln0.13022 \approx -2.0385$$

　　　　現在我們再看看HMM問題二，求解模型引數的問題。由於鮑姆-韋爾奇演算法是基於EM演算法的近似演算法，所以我們需要多跑幾次，比如下面我們跑三次，選擇一個比較優的模型引數，程式碼如下：

import numpy as np
from hmmlearn import hmm

states = ["box 1", "box 2", "box3"]
n_states = len(states)

observations = ["red", "white"]
n_observations = len(observations)
model2 = hmm.MultinomialHMM(n_components=n_states, n_iter=20, tol=0.01)
X2 = np.array([[0,1,0,1],[0,0,0,1],[1,0,1,1]])
model2.fit(X2)
print model2.startprob_
print model2.transmat_
print model2.emissionprob_
print model2.score(X2)
model2.fit(X2)
print model2.startprob_
print model2.transmat_
print model2.emissionprob_
print model2.score(X2)
model2.fit(X2)
print model2.startprob_
print model2.transmat_
print model2.emissionprob_
print model2.score(X2)

　　　　結果這裡就略去了，最終我們會選擇分數最高的模型引數。

　　　　以上就是用MultinomialHMM解決HMM模型三個問題的方法。

3. GaussianHMM例項

　　　　下面我們再給一個GaussianHMM的例項，這個例項中，我們的觀測狀態是二維的，而隱藏狀態有4個。因此我們的“means”引數是$4 \times 2$的矩陣，而“covars”引數是$4 \times 2 \times 2$的張量。

　　　　建立模型如下：

startprob = np.array([0.6, 0.3, 0.1, 0.0])
# The transition matrix, note that there are no transitions possible
# between component 1 and 3
transmat = np.array([[0.7, 0.2, 0.0, 0.1],
                     [0.3, 0.5, 0.2, 0.0],
                     [0.0, 0.3, 0.5, 0.2],
                     [0.2, 0.0, 0.2, 0.6]])
# The means of each component
means = np.array([[0.0,  0.0],
                  [0.0, 11.0],
                  [9.0, 10.0],
                  [11.0, -1.0]])
# The covariance of each component
covars = .5 * np.tile(np.identity(2), (4, 1, 1))

# Build an HMM instance and set parameters
model3 = hmm.GaussianHMM(n_components=4, covariance_type="full")

# Instead of fitting it from the data, we directly set the estimated
# parameters, the means and covariance of the components
model3.startprob_ = startprob
model3.transmat_ = transmat
model3.means_ = means
model3.covars_ = covars

 　　　　注意上面有個引數covariance_type，取值為"full"意味所有的$\mu,\Sigma$都需要指定。取值為“spherical”則$\Sigma$的非對角線元素為0，對角線元素相同。取值為“diag”則$\Sigma$的非對角線元素為0，對角線元素可以不同，"tied"指所有的隱藏狀態對應的觀測狀態分佈使用相同的協方差矩陣$\Sigma$

　　　　我們現在跑一跑HMM問題一解碼的過程，由於觀測狀態是二維的，我們用的三維觀測序列， 所以這裡的 輸入是一個$3 \times 2$的矩陣，程式碼如下：

seen = np.array([[1.1,2.0],[-1,2.0],[3,7]])
logprob, state = model.decode(seen, algorithm="viterbi")
print state

 　　　　輸出結果如下：

[0 0 1]

　　　　再看看HMM問題一對數概率的計算：

print model3.score(seen)

　　　　輸出如下：

-41.1211281377

　　　　以上就是用hmmlearn學習HMM的過程。希望可以幫到大家。

 

（歡迎轉載，轉載請註明出處。歡迎溝通交流： liujianping-ok@163.com） 

相關文章