非齊次隱馬爾可夫模型及其參數(shù)估計(jì)

徐妍

摘要：隱馬爾可夫模型是現(xiàn)今被廣泛使用的統(tǒng)計(jì)模型之一。本文在現(xiàn)有的對隱馬爾可夫模型研究的基礎(chǔ)上，通過在轉(zhuǎn)移模型的馬爾可夫鏈中增加協(xié)變量，構(gòu)建了非齊次隱馬爾可夫模型。在對模型進(jìn)行參數(shù)估計(jì)時(shí)，首先用k-means聚類分析確定了觀測數(shù)據(jù)的隱狀態(tài)，然后用極大似然估計(jì)對模型的轉(zhuǎn)移參數(shù)進(jìn)行估計(jì)。在數(shù)值模擬時(shí)，以非齊次隱馬爾可夫多元正態(tài)分布和非齊次隱馬爾可夫多元線性回歸為例，利用文章所介紹的方法對模型的參數(shù)進(jìn)行估計(jì)，驗(yàn)證了估計(jì)的可靠性。

關(guān)鍵詞：非齊次隱馬爾可夫模型 k-means聚類 極大似然估計(jì)

一、引言

隱馬爾可夫模型最早是由Leonard E.Baum等在1966年和1970年的兩篇論文中提出的[1，2]。隱馬爾可夫模型由兩部分組成，一部分是轉(zhuǎn)移模型，用來描述模型隱狀態(tài)之間的轉(zhuǎn)移關(guān)系，在實(shí)際應(yīng)用中隱狀態(tài)是不能被直接觀測到的;另一部分是發(fā)射模型，即和隱狀態(tài)一一對應(yīng)的觀測變量，每個(gè)發(fā)射模型來自一個(gè)特定的分布。20世紀(jì)70年代，隱馬爾可夫模型首先被應(yīng)用在語音識別領(lǐng)域[3]，之后，又成為分析生物遺傳信息的有用工具[4]。隨著隱馬爾可夫模型的不斷發(fā)展，對隱馬爾可夫模型的研究也越來越多。Dempster等提出了EM算法，用來對包含隱狀態(tài)的模型進(jìn)行極大似然估計(jì)[5]。Jelinek等將Baum-Welch算法與EM算法的關(guān)系進(jìn)行了完整的描述，說明了Baum-Welch算法是EM算法在隱馬爾可夫模型中的具體應(yīng)用[6]。國內(nèi)也有很多關(guān)于隱馬爾可夫模型的研究，夏葉茂等研究了隱馬爾可夫因子分析模型的半?yún)?shù)貝葉斯分析 [7]。劉鶴飛等研究了隱狀態(tài)個(gè)數(shù)未知的隱馬爾可夫多元正態(tài)分布的貝葉斯推斷，利用可逆跳躍MCMC算法對隱狀態(tài)個(gè)數(shù)進(jìn)行模型選擇，再對模型參數(shù)進(jìn)行估計(jì)[8]。王坤等將隱馬爾可夫模型與結(jié)構(gòu)方程模型相結(jié)合，并通過貝葉斯方法對模型的未知參數(shù)進(jìn)行估計(jì)[9]。

隨著大數(shù)據(jù)的發(fā)展，隱馬爾可夫模型在生活中也得到了更廣泛的應(yīng)用。金融領(lǐng)域，隱馬爾可夫模型可以用來對股市收益率波動和狀態(tài)轉(zhuǎn)換進(jìn)行建模和預(yù)測[10];互聯(lián)網(wǎng)領(lǐng)域，隱馬爾可夫模型可以用來進(jìn)行網(wǎng)絡(luò)環(huán)境監(jiān)管;安全領(lǐng)域，隱馬爾可夫模型可以和貝葉斯神經(jīng)網(wǎng)絡(luò)結(jié)合，通過分析一些先前發(fā)生的事件來預(yù)測恐怖分子在未來一段時(shí)間內(nèi)可能發(fā)動恐怖活動的概率，用來預(yù)防可能發(fā)生的恐怖活動[11]。

目前國內(nèi)外期刊發(fā)表的論文主要集中在齊次隱馬爾可夫模型的研究，即在轉(zhuǎn)移模型中只考慮了隱狀態(tài)從時(shí)刻的狀態(tài)轉(zhuǎn)移到時(shí)刻的狀態(tài)的概率，沒有考慮觀測變量自身的某些因素對狀態(tài)之間轉(zhuǎn)移概率的影響。而我們發(fā)現(xiàn)，在許多實(shí)際應(yīng)用中，模型在時(shí)刻的狀態(tài)往往不僅僅受時(shí)刻狀態(tài)的影響，還可能會受時(shí)刻自身某些協(xié)變量的影響。為此，本文提出了一種新的轉(zhuǎn)移模型用來描述隱狀態(tài)個(gè)數(shù)已知的情況下，隱狀態(tài)之間的非齊次概率轉(zhuǎn)移關(guān)系，稱為非齊次隱馬爾可夫模型，這是本文的第一個(gè)創(chuàng)新點(diǎn)。本文的第二個(gè)創(chuàng)新點(diǎn)是提出了利用k-means聚類方法確定非齊次隱馬爾可夫模型的隱狀態(tài)，在隱狀態(tài)確定的條件下，再分別估計(jì)轉(zhuǎn)移模型和發(fā)射模型中的未知參數(shù)。

二、模型描述

以上描述的就是非齊次隱馬爾可夫模型的狀態(tài)轉(zhuǎn)移過程，我們將這種形式的轉(zhuǎn)移模型稱作“連續(xù)比對數(shù)轉(zhuǎn)移模型”。

在隱狀態(tài)確定的條件下，觀測向量為服從特定分布的向量，其中p是觀測向量的維度。即時(shí)，對非齊次隱馬爾可夫模型來說，向量為模型的觀測變量，為影響轉(zhuǎn)移效果的固定協(xié)變量，q是固定協(xié)變量的維度。本文要研究的就是如何利用可觀測到的信息和去估計(jì)非齊次隱馬爾可夫模型中的轉(zhuǎn)移模型和每個(gè)隱狀態(tài)下發(fā)射模型中的未知參數(shù)。

三、估計(jì)原理

（一）隱狀態(tài)向量

隱馬爾可夫模型的觀測變量來自不同的隱狀態(tài)，然而在實(shí)際應(yīng)用中，隱狀態(tài)無法直接觀測，這也是隱馬爾可夫模型推斷的困難之處。已有的許多研究，都試圖利用觀測變量中隱含的信息首先對觀測變量的隱狀態(tài)進(jìn)行判定，在確定隱狀態(tài)向量之后，再對隱馬爾可夫模型中的轉(zhuǎn)移參數(shù)和每個(gè)隱狀態(tài)下特定分布的未知參數(shù)進(jìn)行估計(jì)。其中，最有代表性的就是向前向后遞歸算法。然而向前向后遞歸算法的理論性很強(qiáng)，計(jì)算方法非常復(fù)雜，要求使用者具有一定的統(tǒng)計(jì)學(xué)理論背景和計(jì)算機(jī)編程能力，這極大地限制了隱馬爾可夫模型在實(shí)際生產(chǎn)生活中的推廣和應(yīng)用。

本文利用k-means聚類方法，來確定隱馬爾可夫模型觀測變量的隱狀態(tài)，在隱狀態(tài)確定的條件下對模型中的未知參數(shù)進(jìn)行估計(jì)。

（二）標(biāo)簽交換

在對非齊次隱馬爾可夫模型進(jìn)行研究時(shí)，首先要解決標(biāo)簽交換問題。如果沒有解決這一問題，模型的參數(shù)估計(jì)結(jié)果就不具有可解釋性，甚至?xí)l(fā)生混淆。Richardson和 Green 在文章中通過比較均值來解決標(biāo)簽交換問題[12]。本文借鑒這一經(jīng)驗(yàn)，在發(fā)射模型為多元正態(tài)分布時(shí)，通過比較多元正態(tài)分布中均值向量第一分量的方法來解決這一問題;在發(fā)射模型為多元線性回歸時(shí)，通過比較回歸系數(shù)的方法來解決這一問題。即在參數(shù)估計(jì)之后，重新確定觀測變量隱狀態(tài)的標(biāo)簽，然后再根據(jù)隱狀態(tài)標(biāo)簽結(jié)果，重新確定每個(gè)隱狀態(tài)下參數(shù)的估計(jì)結(jié)果。

（三）轉(zhuǎn)移模型中未知參數(shù)的估計(jì)

利用極大似然估計(jì)對非齊次隱馬爾可夫模型中轉(zhuǎn)移模型的未知參數(shù)進(jìn)行估計(jì)。未知參數(shù)的似然函數(shù)如下：

最大似然估計(jì)就是要找到使得似然函數(shù)取最大值時(shí)未知參數(shù)和的值。本文中似然函數(shù)含有個(gè)未知參數(shù)。在實(shí)際應(yīng)用中，可以用統(tǒng)計(jì)軟件中求極值的函數(shù)得到未知參數(shù)的數(shù)值解。比如，本文就是通過R語言BB程序包中的fun函數(shù)，近似求解似然函數(shù)的極大值點(diǎn)。

（四）發(fā)射模型中未知參數(shù)的估計(jì)

隱馬爾可夫模型的發(fā)射模型可以來自各種不同的特定分布。本文選擇比較經(jīng)典和常用的多元正態(tài)分布和多元線性回歸作為兩個(gè)模擬實(shí)驗(yàn)的發(fā)射模型。

1.多元正態(tài)分布。假設(shè)聚類分析確定隱狀態(tài)向量后，第k個(gè)隱狀態(tài)下的觀測集合為，。則觀測模型中的待估參數(shù)為每個(gè)隱狀態(tài)下的均值向量和協(xié)方差矩陣。

用極大似然估計(jì)對發(fā)射模型中的待估參數(shù)進(jìn)行估計(jì)，結(jié)果如下：

2.多元線性回歸。假設(shè)聚類分析確定隱狀態(tài)向量后，第k個(gè)隱狀態(tài)下的觀測向量為。則觀測模型中的待估參數(shù)為每個(gè)隱狀態(tài)下自變量的系數(shù)向量。

用最小二乘估計(jì)對觀測模型中的待估參數(shù)進(jìn)行估計(jì)，結(jié)果如下：

四、實(shí)證分析

本實(shí)例的數(shù)據(jù)來源于曲靖師范學(xué)院數(shù)學(xué)與統(tǒng)計(jì)學(xué)院2014級數(shù)學(xué)與應(yīng)用數(shù)學(xué)專業(yè)兩個(gè)班91位同學(xué)八個(gè)學(xué)期的綜合測評成績。每名學(xué)生每學(xué)期的綜合測評成績分為文化知識成績和創(chuàng)新發(fā)展成績兩類，文化知識成績和創(chuàng)新發(fā)展成績均采用百分制計(jì)算。我們將每名學(xué)生每學(xué)期的文化知識成績和創(chuàng)新發(fā)展成績作為發(fā)射模型的觀測數(shù)據(jù)，將它們近似的看作一個(gè)二維的正態(tài)分布。將每名學(xué)生的性別作為固定協(xié)變量，男生設(shè)置為1，女生設(shè)置為0;將學(xué)生是否獲得獎學(xué)金作為模型的隱狀態(tài)，綜合測評分?jǐn)?shù)高的同學(xué)獲得獎學(xué)金。用k-means聚類分析確定每個(gè)觀測變量的隱狀態(tài)，再利用文章中介紹的標(biāo)簽交換的方法，聚類后，表示學(xué)生獲得獎學(xué)金，表示學(xué)生沒有獲得獎學(xué)金。

根據(jù)文中介紹的參數(shù)估計(jì)方法，得到各參數(shù)的估計(jì)值：

獲得獎學(xué)金的同學(xué)，綜合測評成績服從：

沒有獲得獎學(xué)金的同學(xué)，綜合測評成績服從：

根據(jù)發(fā)射模型參數(shù)估計(jì)結(jié)果可知：獲得獎學(xué)金的同學(xué)綜合測評成績明顯高于沒有獲得獎學(xué)金的同學(xué)。具體來說，獲得獎學(xué)金的同學(xué)，文化知識成績平均分為86.24，創(chuàng)新發(fā)展成績平均分為86.40;沒有獲得獎學(xué)金的同學(xué)，文化知識成績平均分為78.48，創(chuàng)新發(fā)展成績平均分為77.36。

根據(jù)隱狀態(tài)轉(zhuǎn)移概率公式可以知道，當(dāng)學(xué)生性別為女時(shí)，如果前一個(gè)學(xué)期該學(xué)生獲得了獎學(xué)金，那么后一個(gè)學(xué)期該學(xué)生獲得獎學(xué)金的概率為0.80，不能獲得獎學(xué)金的概率為0.20;如果前一個(gè)學(xué)期該學(xué)生沒有獲得獎學(xué)金，那么后一個(gè)學(xué)期該學(xué)生可以獲得獎學(xué)金的概率為0.45，不能獲得獎學(xué)金的概率為0.55。當(dāng)學(xué)生性別為男生時(shí)，如果前一個(gè)學(xué)期該學(xué)生獲得獎學(xué)金，那么后一個(gè)學(xué)期該學(xué)生可以獲得獎學(xué)金的概率為0.51，不能獲得獎學(xué)金的概率為0.49;如果前一個(gè)學(xué)期該學(xué)生沒有獲得獎學(xué)金，那么后一個(gè)學(xué)期該學(xué)生可以獲得獎學(xué)金的概率為0.18，不能獲得獎學(xué)金的概率為0.82。

五、結(jié)論

本文在現(xiàn)有的對齊次隱馬爾可夫模型研究的基礎(chǔ)上，通過在轉(zhuǎn)移模型的馬爾可夫鏈中增加協(xié)變量，提出了非齊次隱馬爾可夫模型，解決了當(dāng)模型的觀測變量存在固定協(xié)變量時(shí)，對隱狀態(tài)轉(zhuǎn)移關(guān)系進(jìn)行建模的問題，這是本文的第一個(gè)創(chuàng)新點(diǎn)。本文采用了“連續(xù)比對數(shù)轉(zhuǎn)移模型”來描述隱狀態(tài)之間的轉(zhuǎn)移關(guān)系，因?yàn)槟Ｐ椭形粗獏?shù)較多，且不同情況下未知參數(shù)的個(gè)數(shù)也不同，所以今后可以對非齊次轉(zhuǎn)移模型進(jìn)行相關(guān)的改進(jìn)研究，使得轉(zhuǎn)移模型形式更簡單，使用更方便。

本文的第二個(gè)創(chuàng)新點(diǎn)是利用較簡單的k-means聚類分析確定模型觀測變量的隱狀態(tài)，在隱狀態(tài)確定的情況下對非齊次隱馬爾可夫模型中的未知參數(shù)進(jìn)行估計(jì)。在今后的研究中，可以嘗試?yán)孟到y(tǒng)聚類方法對隱狀態(tài)個(gè)數(shù)未知的隱馬爾可夫模型進(jìn)行隱狀態(tài)個(gè)數(shù)的模型選擇，或者嘗試?yán)蒙窠?jīng)網(wǎng)絡(luò)、決策樹等聚類方法來確定觀測變量的隱狀態(tài)。

參考文獻(xiàn)：

[1]Baum L E，Petrie T.Statistical inference for probabilistic functions of finite state Markov chains[J].Annals of Mathematical Statistics，1966，37（1）：1554-1563.

[2]Baum L E，Petrie T，Soules G，Weiss N.A maximization technique occurring in the statistical analysis of probabilistic functions of Markov chains[J].Annals of Mathematical Statistical，1970，41（1）：164-171.

[3]Jelinek，F(xiàn).，Bahl，L.，Mercer，R.Design of a linguistic statistical decoder for the recognition of continuous speech[J].IEEE Transactions on Information Theory，1975 2（3）：250-256.

[4]Bishop，Martin J，Thompson，Elizabeth A.Maximum Likelihood alignment of DNA sequences[J].Journal of Molecular Biology，1986，190（2）：159-165.

[5]Dempster，A.P.Maximum likelihood from incomplete data via the EM algorithm[J].Journal of Royal Statistical Society B，1977，39（1）：1-38.

[6]Jelinek，F(xiàn)rederick.Statistical methods for speech recognition[M].1997.

[7]夏業(yè)茂，勾建偉，劉應(yīng)安.隱馬爾可夫因子分析模型的半?yún)?shù)貝葉斯分析[J].高校應(yīng)用數(shù)學(xué)學(xué)報(bào)，2015，30（1）：17-30.

[8]劉鶴飛，王坤，蔣成飛.隱狀態(tài)個(gè)數(shù)未知的隱馬爾可夫多元正態(tài)分布的貝葉斯推斷[J].統(tǒng)計(jì)研究，2017，34（12）：119-125.

[9]王坤，劉鶴飛，蔣成飛.隱馬爾可夫結(jié)構(gòu)方程及其貝葉斯估計(jì)[J].數(shù)理統(tǒng)計(jì)與管理，2018，37（2）：273-279.

[10]劉金全，李楠，鄭挺國.隨機(jī)波動模型的馬爾可夫鏈—蒙特卡羅模擬方法—在滬市收益率序列上的應(yīng)用[J].數(shù)理統(tǒng)計(jì)與管理，2010，29（6）：1026-1035.

[11]戰(zhàn)兵，韓銳.基于隱馬爾可夫的恐怖事件預(yù)測模型[J].解放軍理工大學(xué)學(xué)報(bào)，2015，16（4）：386-393.

[12]Richardson S，Green P J.On Bayesian analysis of mixtures with an unknown number of components[J].Journal of the Royal Statistical Society（Series B），1997，59（4）：731—758.

（作者單位：首都經(jīng)濟(jì)貿(mào)易大學(xué)統(tǒng)計(jì)學(xué)院）