數學之美中的馬爾科夫鏈

與HMM模型相關的「有用」的問題是評估（前向演算法）和解碼（維特比演算法）——它們一個被用來測量一個模型的相對適用性，另一個被用來推測模型隱藏的部分在做什麼（「到底發生了」什麼）。可以看出它們都依賴於隱馬爾科夫模型（HMM）參數這一先驗知識——狀態轉移矩陣，混淆（觀察）矩陣，以及向量（初始化概率向量）。

　　然而，在許多實際問題的情況下這些參數都不能直接計算的，而要需要進行估計——這就是隱馬爾科夫模型中的學習問題。前向-後向演算法就可以以一個觀察序列為基礎來進行這樣的估計，而這個觀察序列來自於一個給定的集合，它所代表的是一個隱馬爾科夫模型中的一個已知的隱藏集合。

　　一個例子可能是一個龐大的語音處理資料庫，其底層的語音可能由一個馬爾可夫過程基於已知的音素建模的，而其可以觀察的部分可能由可識別的狀態（可能通過一些矢量數據表示）建模的，但是沒有（直接）的方式來獲取隱馬爾科夫模型（HMM）參數。

　　前向-後向演算法並非特別難以理解，但自然地比前向演算法和維特比演算法更複雜。由於這個原因，這裡就不詳細講解前向-後向演算法了（任何有關HMM模型的參考文獻都會提供這方面的資料的）。

　　總之，前向-後向演算法首先對於隱馬爾科夫模型的參數進行一個初始的估計（這很可能是完全錯誤的），然後通過對於給定的數據評估這些參數的的價值並減少它們所引起的錯誤來重新修訂這些HMM參數。從這個意義上講，它是以一種梯度下降的形式尋找一種錯誤測度的最小值。

　　之所以稱其為前向-後向演算法，主要是因為對於網格中的每一個狀態，它既計算到達此狀態的「前向」概率（給定當前模型的近似估計），又計算生成此模型最終狀態的「後向」概率（給定當前模型的近似估計）。這些都可以通過利用遞歸進行有利地計算，就像我們已經看到的。可以通過利用近似的HMM模型參數來提高這些中間概率進行調整，而這些調整又形成了前向-後向演算法迭代的基礎。

要理解前向-後向演算法，首先需要了解兩個演算法：後向演算法和EM演算法。後向演算法是必須的，因為前向-後向演算法就是利用了前向演算法與後向演算法中的變數因子，其得名也因於此；而EM演算法不是必須的，不過由於前向-後向演算法是EM演算法的一個特例，因此了解一下EM演算法也是有好處的，說實話，對於EM演算法，我也是雲里霧裡的。好了，廢話少說，我們先談談後向演算法。

1、後向演算法（Backward algorithm）

　　其實如果理解了前向演算法，後向演算法也是比較好理解的，這裡首先重新定義一下前向演算法中的局部概率at(i)，稱其為前向變數，這也是為前向-後向演算法做點準備：

　　相似地，我們也可以定義一個後向變數Bt(i)（同樣可以理解為一個局部概率）：

　　後向變數(局部概率)表示的是已知隱馬爾科夫模型及t時刻位於隱藏狀態Si這一事實，從t+1時刻到終止時刻的局部觀察序列的概率。同樣與前向演算法相似，我們可以從後向前（故稱之為後向演算法）遞歸地計算後向變數：

　　1）初始化，令t=T時刻所有狀態的後向變數為1：

　　2）歸納，遞歸計算每個時間點，t=T-1,T-2,…,1時的後向變數：

　　這樣就可以計算每個時間點上所有的隱藏狀態所對應的後向變數，如果需要利用後向演算法計算觀察序列的概率，只需將t=1時刻的後向變數（局部概率）相加即可。下圖顯示的是t+1時刻與t時刻的後向變數之間的關係：

　　上述主要參考自HMM經典論文《A tutorial on Hidden Markov Models and selectedapplications in speechrecognition》。下面我們給出利用後向演算法計算觀察序列概率的程序示例，這個程序仍然來自於UMDHMM。

後向演算法程序示例如下（在backward.c中):

void Backward(HMM *phmm, int T, int *O, double **beta, double*pprob)

{

　　int i, j;

　　int t;

　　double sum;

　　for (i = 1; i <= phmm->N;i++)

　　　　beta[T][i] = 1.0;

　　for (t = T - 1; t >= 1; t--)

　　{

　　　　for (i = 1; i <= phmm->N;i++)

　　　　{

　　　　　　sum = 0.0;

　　　　　　for (j = 1; j <= phmm->N;j++)

　　　　　　　　sum += phmm->A[i][j] *

　　　　　　　　　　　　　　(phmm->B[j][O[t+1]])*beta[t+1][j];

　　　　　　beta[t][i] = sum;

　　　　}

　　}

　　*pprob = 0.0;

　　for (i = 1; i <= phmm->N;i++)

　　　　*pprob += beta[1][i];

}

隱馬爾科夫模型（HMM）的三個基本問題中，第三個HMM參數學習的問題是最難的，因為對於給定的觀察序列O，沒有任何一種方法可以精確地找到一組最優的隱馬爾科夫模型參數（A、B、）使P(O|)最大。因而，學者們退而求其次，不能使P(O|)全局最優，就尋求使其局部最優（最大化）的解決方法，而前向-後向演算法（又稱之為Baum-Welch演算法）就成了隱馬爾科夫模型學習問題的一種替代（近似）解決方法。

　　我們首先定義兩個變數。給定觀察序列O及隱馬爾科夫模型，定義t時刻位於隱藏狀態Si的概率變數為：

　　回顧一下第二節中關於前向變數at(i)及後向變數Bt(i)的定義，我們可以很容易地將上式用前向、後向變數表示為：

　　其中分母的作用是確保：

給定觀察序列O及隱馬爾科夫模型，定義t時刻位於隱藏狀態Si及t+1時刻位於隱藏狀態Sj的概率變數為：

　　該變數在網格中所代表的關係如下圖所示：

　　同樣，該變數也可以由前向、後向變數表示：

　　而上述定義的兩個變數間也存在著如下關係：

　　如果對於時間軸t上的所有相加，我們可以得到一個總和，它可以被解釋為從其他隱藏狀態訪問Si的期望值（網格中的所有時間的期望），或者，如果我們求和時不包括時間軸上的t=T時刻，那麼它可以被解釋為從隱藏狀態Si出發的狀態轉移期望值。相似地，如果對在時間軸t上求和（從t=1到t=T-1），那麼該和可以被解釋為從狀態Si到狀態Sj的狀態轉移期望值。即：

上一節我們定義了兩個變數及相應的期望值，本節我們利用這兩個變數及其期望值來重新估計隱馬爾科夫模型（HMM）的參數，A及B：

　　如果我們定義當前的HMM模型為，那麼可以利用該模型計算上面三個式子的右端；我們再定義重新估計的HMM模型為，那麼上面三個式子的左端就是重估的HMM模型參數。Baum及他的同事在70年代證明了因此如果我們迭代地的計算上面三個式子，由此不斷地重新估計HMM的參數，那麼在多次迭代後可以得到的HMM模型的一個最大似然估計。不過需要注意的是，前向-後向演算法所得的這個結果（最大似然估計）是一個局部最優解。

　　關於前向-後向演算法和EM演算法的具體關係的解釋，大家可以參考HMM經典論文《A tutorial on Hidden MarkovModels and selected applications in speechrecognition》，這裡就不詳述了。下面我們給出UMDHMM中的前向-後向演算法示例，這個演算法比較複雜，這裡只取主函數，其中所引用的函數大家如果感興趣的話可以自行參考UMDHMM。

前向-後向演算法程序示例如下（在baum.c中):

void BaumWelch(HMM *phmm, int T, int *O, double **alpha, double**beta, double **gamma, int *pniter, double *plogprobinit, double*plogprobfinal)

{

　　int i, j, k;

　　int t, l = 0;

　　前向-後向演算法就到此為止了。

　　double logprobf, logprobb, threshold;

　　double numeratorA, denominatorA;

　　double numeratorB, denominatorB;

　　double ***xi, *scale;

　　double delta, deltaprev, logprobprev;

　　deltaprev = 10e-70;

　　xi = AllocXi(T, phmm->N);

　　scale = dvector(1, T);

　　ForwardWithScale(phmm, T, O, alpha, scale,&logprobf);

　　*plogprobinit = logprobf;

　　BackwardWithScale(phmm, T, O, beta, scale,&logprobb);

　　ComputeGamma(phmm, T, alpha, beta, gamma);

　　ComputeXi(phmm, T, O, alpha, beta, xi);

　　logprobprev = logprobf;

　　do

　　{

　　　　for (i = 1; i <= phmm->N;i++)

　　　　　　phmm->pi[i] = .001 + .999*gamma[1][i];

　　　　for (i = 1; i <= phmm->N;i++)

　　　　{

　　　　　　denominatorA = 0.0;

　　　　　　for (t = 1; t <= T - 1; t++)

　　　　　　　　denominatorA += gamma[t][i];

　　　　　　for (j = 1; j <= phmm->N;j++)

　　　　　　{

　　　　　　　　numeratorA = 0.0;

　　　　　　　　for (t = 1; t <= T - 1; t++)

　　　　　　　　　　numeratorA += xi[t][i][j];

　　　　　　　　phmm->A[i][j] = .001 +

　　　　　　　　　　　　　　　　　.999*numeratorA/denominatorA;

　　　　　　}

　　　　　　denominatorB = denominatorA + gamma[T][i];

　　　　　　for (k = 1; k <= phmm->M;k++)

　　　　　　{

　　　　　　　　numeratorB = 0.0;

　　　　　　　　for (t = 1; t <= T; t++)

　　　　　　　　{

　　　　　　　　　　if (O[t] == k)

　　　　　　　　　　　　numeratorB += gamma[t][i];

　　　　　　　　}

　　　　　　　　phmm->B[i][k] = .001 +

　　　　　　　　　　　　　　　　　.999*numeratorB/denominatorB;

　　　　　　}

　　　　}

　　　　ForwardWithScale(phmm, T, O, alpha, scale,&logprobf);

　　　　BackwardWithScale(phmm, T, O, beta, scale,&logprobb);

　　　　ComputeGamma(phmm, T, alpha, beta, gamma);

　　　　ComputeXi(phmm, T, O, alpha, beta, xi);

　　　　delta = logprobf - logprobprev;

　　　　logprobprev = logprobf;

　　　　l++;

　　}

　　while (delta > DELTA);

　　*pniter = l;

　　*plogprobfinal = logprobf;

　　FreeXi(xi, T, phmm->N);

　　free_dvector(scale, 1, T);

}

總結(Summary)

　　通常，模式並不是單獨的出現，而是作為時間序列中的一個部分——這個過程有時候可以被輔助用來對它們進行識別。在基於時間的進程中，通常都會使用一些假設——一個最常用的假設是進程的狀態只依賴於前面N個狀態——這樣我們就有了一個N階馬爾科夫模型。最簡單的例子是N= 1。

　　存在很多例子，在這些例子中進程的狀態（模式）是不能夠被直接觀察的，但是可以非直接地，或者概率地被觀察為模式的另外一種集合——這樣我們就可以定義一類隱馬爾科夫模型——這些模型已被證明在當前許多研究領域，尤其是語音識別領域具有非常大的價值。

　　在實際的過程中這些模型提出了三個問題都可以得到立即有效的解決，分別是：

　　* 評估：對於一個給定的隱馬爾科夫模型其生成一個給定的觀察序列的概率是多少。前向演算法可以有效的解決此問題。

　　* 解碼：什麼樣的隱藏（底層）狀態序列最有可能生成一個給定的觀察序列。維特比演算法可以有效的解決此問題。

　　*學習：對於一個給定的觀察序列樣本，什麼樣的模型最可能生成該序列——也就是說，該模型的參數是什麼。這個問題可以通過使用前向-後向演算法解決。

　　隱馬爾科夫模型（HMM）在分析實際系統中已被證明有很大的價值；它們通常的缺點是過於簡化的假設，這與馬爾可夫假設相關——即一個狀態只依賴於前一個狀態，並且這種依賴關係是獨立於時間之外的（與時間無關）。

　　關於隱馬爾科夫模型的完整論述，可參閱：

L R Rabiner and B H Juang, `An introduction to HMMs』, iEEEASSP Magazine, 3, 4-16.

推薦閱讀：