R語言實作viterbi算法

資料集介紹

最近初學HMM(Hidden Markov Model)，老師讓我自己試着用Java去實作viterbi算法，結果試了一下，發現資料的輸入太麻煩，因為資料集是一些各種字元串和标點符号的矩陣吧，然後才開始學習R語言的， 首先先看看資料集。

http://www.clips.uantwerpen.be/conll2000/chunking/

資料集是從上面這個網站下載下傳的，有一個train.txt還有一個test.txt。

格式大概如下：

He        PRP  B-NP
   reckons   VBZ  B-VP
   the       DT   B-NP
   current   JJ   I-NP
   account   NN   I-NP
   deficit   NN   I-NP
   will      MD   B-VP
   narrow    VB   I-VP
   to        TO   B-PP
   only      RB   B-NP
   #         #    I-NP
   1.8       CD   I-NP
   billion   CD   I-NP
   in        IN   B-PP
   September NNP  B-NP
   .         .    O      

這是一整句話的資料。 第一列我們稱為symbol或者是word,第二列稱為tag或者是标簽(tag又稱作狀态state)，第三列在這個算法中沒用到，是以先不提及。      

然後首先要知道viterbi算法的輸入和輸出，      

算法輸入輸出

Input:

tp: transition probability matrix  (不同tag之間的轉移機率矩陣, tp[i,j]表示從第i個标簽，下一個單詞的标簽是第j個标簽的機率 )

ep: emission probability matrix (表示tag和word之間的發射機率矩陣，ep[v,k]表示第v個标簽下，這個标簽會輸出第k個單詞(word/symbol)的機率 )

startp: start state probability vector (表示某個tag作為一句話開頭的機率, starp[i]表示第i個tag作為句子sequence開頭的機率)

test: the observation sequence need to be tagged (test 就是輸入進去的一句話，是向量形式輸入)

Output:

outp: it’s also a probability matrix, each row denotes each word(i.e. symbol) of thetest sequence and the column contains all the tags.(43 tags in all)

  Each element(outp[i,j]) in this matrix denotes the probability of the ithword would be tagged in the jth tag.

(我這裡定義的輸出可能和其他人的viterbi算法輸出不太一樣，這個輸出了之後還要做一些後續的處理，這裡的outp是一個矩陣，行是代表輸入的test的每一個單詞，列是代表每一個tag，然後outp[i,j]就表示第i個單詞用第j個标簽進行tag的機率，是以最後隻需要再做一步統計每一行中最大的那一列就是這一行的那個單詞的tag)

算法輸入參數預測和計算

熟悉可以跳過下面這三個代碼段，因為下面我要首先計算出 tp,ep,startp這三個矩陣裡面的數值，首先是計算轉移機率矩陣。

# This part is used to calculate the transition probability matrix tag.tp 計算轉移機率矩陣
data <- read.table('train.txt',sep=" ",quote="")  #read the data
tags_num <- table(data[,2])                       #table函數可以統計在data中每個tag出現了多少次，相當于SQL裡面的group的作用
tag.transform<-table(data[,2],data[,2])		      #transform這裡是想弄成一個矩陣，行和列都是tag，tag.transform[i,j]表示第i行轉移到第j列的頻數
for(i in 1:length(tag.transform[1,])){            #循環初始化每一個頻數都是0
  tag.transform[i,i] <- 0
}

for(i in 1:length(data[,1])){			  #這裡通過掃描data裡面總的單詞，統計前一個tag和後一個tag的連接配接頻數。
  if(data[i,2]=='.') next
  tag.transform[data[i,2],data[i+1,2]]<-tag.transform[data[i,2],data[i+1,2]]+1
}

tag.tp<-tag.transform				  #tag.tp就是最後的轉移機率矩陣，通過頻數除以總數得到轉移機率
for(i in 1:43){
  for(j in 1:43){
    tag.tp[i,j]<-tag.tp[i,j]/tags_num[[i]]
  }
}
#######################End of this part#####################
           

然後下面這部分計算輸出機率矩陣

#This part is used to calculate the emission probability using both train data and test data
symbol_data <- read.table('train_and_test.txt',sep=" ",quote="")    #這個輸出機率矩陣的計算把訓練集和測試集的資料都拿了出來，因為assignment的文檔說有些單詞train data上沒有，是以必須拿出來
symbol<-table(symbol_data[,1])                              #統計相同單詞的數量
symbol_tags_num<-table(symbol_data[,2])                     #統計相同tag的數量
symbol_name<-names(symbol)                                  #names函數的作用比較重要，把下标提取出來，不要數字在裡面
symbol.epmatrix<-table(symbol_data[,2],symbol_data[,1])     #這是用來計算的，行是tag,列是單詞，ep[i,j]，第i個tag輸出第j個單詞的機率
symbol.ep<-symbol.epmatrix
for(i in 1:43){                                             #這裡是43的原因是我之前debug程式早就知道總共43個tag，就懶得寫length(tag)了
  for(j in 1:length(symbol)){
    symbol.ep[i,j]<-symbol.ep[i,j]/symbol_tags_num[[i]]
  }
}
####################End of This part#########################
           

最後計算起始機率向量

#calculate the start probability of state(i.e. tags) 這部分看上去比較複雜，但隻是有些地方要注意
startp<-tags_num                            #整個程式是連續的，是以tags_num是引用上面的
for(i in 1:length(startp)){                 #初始化每個tag的作為起始的機率為0
  startp[i]<-0
}
sumstart<-0                                 #sumstart作為統計train的資料集中，總的作為開始的tag數量或者說sentence的數量
for(i in 1:length(symbol_data[,2])){
  if(i==1){                                 #如果是i=1因為沒有‘.’去訓示下一個就是開頭，是以這裡先判斷i=1也是開頭
    startp[symbol_data[i,2]]<-1
    sumstart<-1
    next
  }
  if(i==length(symbol_data[,2])){
    break
  }
  if(symbol_data[i,2]=='.'){                #如果這個tag='.' 表示下一個tag就是下一個sentence的開頭
    startp[symbol_data[i+1,2]]<-startp[symbol_data[i+1,2]]+1
    sumstart<-sumstart+1
  }
}
for(i in 1:length(startp)){
  startp[i]<-startp[i]/sumstart             #最後再除以總數就得到機率
}
########################End of This Part#####################
           

Viterbi算法

然後就是重點的viterbi算法了，我用了一個function去表示這個算法。

#The start of viterbi algorithm
firsttest <- read.table('test.txt',sep=" ",quote="")        #這裡讀進來的就是測試集的資料，不大
words<-as.character(firsttest[,1])                          #把單詞轉換成字元型的
p=array(0,dim=c(length(words),43))                          #再強調一下43是tag的數量，下面就不強調了，outp是輸出的機率矩陣

hmm<-function(tp,ep,p,startp,test){                     #整個算法的複雜度和理論相同O(單詞數*tag數的平方)
  outp<-p
  for(k in 1:length(test)){
    for(v in 1:43){
      if(k==1 || test[k-1]=='.'){
        outp[k,v]<-startp[v]*ep[v,test[k]]              #這裡是計算每句話初始的機率因為當k=1時不能用k=0算，此時隻能用startp去計算
      }
      else{
        for(u in 1:43){                                 #如果并不是開始的那個單詞，就可以用公式算找出最大的那個outp，在u中周遊
          maxvalue<-outp[k-1,u]*tp[u,v]*ep[v,test[k]]
          if(maxvalue>outp[k,v]){
            outp[k,v]<-maxvalue
          }
        }
      } 
    }
  }
  return(outp)                                          #最後得到的是一個結果的機率矩陣，outp[i,j]就是第i個單詞被第j個标簽tagged的機率，是以隻需要統計出每一行的最大的那個機率的那一列，就是這一行的單詞對應的tag了
}
           

因為是老師叫我寫的自己練習練習的，當初就不知道R竟然有hmm包，然後後來用了hmm包發現結果也是差不多。 最後統計出的錯誤 tag:381個 在訓練集中總共有47366個詞， 是以 準确率應該就是 381/47366= 99.19% 吧。

結尾

因為在這裡我不太能區分準确率和召回率的概念，是以就簡單的一除了，希望有能人幫我解答一下哈！！！

by the way, 原來發現java也有那個hmm的包，才發現自己實作和debug的時候是有多累啊。。 不過。歡迎讨論哈。。

還有整段代碼都是可以運作的哈~~~