Minimum Edit Distance
Minimum Edit Distance (最小編輯距離),是用來計算兩個字串的相似程度。
可應用於拼字校正、或是計算兩個DNA序列的相似程度。
例如,假設有三個DNA序列:
1. AGCCT
2. AACCT
3. ATCT
直覺上,會認為,DNA序列 2. 和 1. 的相似程度比 3. 和 1. 的相似程度還要高。
但為何是 2. 和 1. 較相似?相似程度又是多少?
若要把兩序列的相似程度量化成數值,就要計算「最小編輯距離」,最小編輯距離越小,就表示兩序列越相似。
Defining of Minimum Edit Distance
我們採用 Levenshtein 的最小編輯距離的定義,這是一種最簡單的定義方式。依此定義,若一個字串,編輯成另外一個字串,可以進行下面三種動作。
其中,定義 a. 和 b. 的 cost 為 , c. 的 cost 為 (因為「替換」相當於「刪除」後再「插入」)。
例如, 1. 和 2. 的最小編輯距離,可以這樣計算,如下:
上圖中,將兩字串對齊,發現要把 1. 編輯成 2. ,需要對第二個字 G 進行替換( s ),把 G 替換成 A ,而替換的 cost 為 ,故 1. 和 2. 的最小編輯距離為 。
再來,把 1. 編輯成 3. 的最小編輯距離,如下:
上圖中, 表示該處沒字元, 而要把 1. 編輯成 3. ,需要把第二個字 G 刪除( d ),並替換( s )第三個字 C ,這兩個動作的 cost 分別為 和 ,故最小編輯距離為
Computing Minimum Edit Distance
接著來看,要用什麼演算法,來計算最小編輯距離?可以用 Dynamic Programming 的方式來計算。
所謂的 Dynamic Programming 簡而言之,就是把一個問題,一層一層分解下去,分成許多子問題,再算出這些子問題,用這些子問題的解答,求出原本問題的答案。
首先,定義字串長度為 的兩字串,最小編輯距離為
在分解過程中,需要計算長度為 的兩個子字串,最小編輯距離為 ,其中 。
接著,定義起始狀態,就是把兩個給定的字串,各自編輯成空字串的距離,分別為 和 。
由於,把字串長度為 的字串,編輯成空字串長度為 ,即 , 同理,
演算法如下:從起始狀態開始,由 及 ,依序去計算長度為 字串的最小編輯距離。
其中, 為第一個字串中第 個字元, 為第二個字串中第 個字元。
用以上演算法,舉個例子來計算,若兩字串分別為:
- AGCCT
- ATCT
輸入字串長度分別為 和 ,故 。
演算法開始的第一步,先建立以下表格,右下方空格部份即為 的值。
接著,填入初始狀態的值,即
再來,計算 的值,用演算法的公式,如下:
其中,可從表格中得知 ,又因為 ,故 的最小編輯距離為 。將算出來的 結果,填入表格,如下:
用此公式,依序算出表格中其他空格的值。
全部算完以後,表格中最右下方的數字,即為 AGCCT 和 ATCT 兩序列的最小編輯距離,其值為3。
Implementation
再來是實作的部份,建立一個python程式檔 minEditDist.py ,並將以下程式複製貼上。
def minEditDist(sm,sn):
m,n = len(sm),len(sn)
D = map(lambda y: map(lambda x,y : y if x==0 else x if y==0 else 0,
range(n+1),[y]*(n+1)), range(m+1))
for i in range(1,m+1):
for j in range(1,n+1):
D[i][j] = min( D[i-1][j]+1, D[i][j-1]+1,
D[i-1][j-1] + apply(lambda: 0 if sm[i-1] == sn[j-1] else 2))
for i in range(0,m+1):
print D[i]
return D[m][n]
其中, sm
和 sn
分別為輸入的字串,而 D[i][j]
為 最小編輯距離 。
到 python interactive mode 載入模組
>>> from minEditDist import minEditDist
接著,用此程式計算 AGCCT 和 ATCT 的最小編輯距離,方法如下:
>>> minEditDist("AGCCT","ATCT")
[0, 1, 2, 3, 4]
[1, 0, 1, 2, 3]
[2, 1, 2, 3, 4]
[3, 2, 3, 2, 3]
[4, 3, 4, 3, 4]
[5, 4, 3, 4, 3]
3
程式印出 表格中的數值,並回傳最後結果,結果為 。
Reference
本文參考至coursera線上課程