羅志兵

（西南林業(yè)大學大數(shù)據(jù)與智能工程學院，昆明 650224）

基于動態(tài)規(guī)劃的基因雙序列比對研究

羅志兵

（西南林業(yè)大學大數(shù)據(jù)與智能工程學院，昆明 650224）

在生物信息學的研究中，最基礎的問題是對生物序列進行分析，序列比對是最基本的操作。通過比對找出序列間的相似性，根據(jù)相似性進行同源性分析，推導出演化過程。因此，序列比對常用于研究由共同祖先進化而來的序列，特別是如蛋白質序列或DNA序列等生物序列。它對于發(fā)現(xiàn)生物序列中的功能、結構和進化信息具有非常重要的意義。如何開發(fā)出高效、準確的序列比對算法是目前序列比對的一個重點。經過對雙序列比對算法的研究，探索動態(tài)規(guī)劃算法在全局比對和局部比對中的應用以及算法的實現(xiàn)過程。

生物信息學；序列比對；動態(tài)規(guī)劃

1 序列比對的概念及其意義

序列比對的定義是：根據(jù)特定的計分規(guī)則，兩個或多個符號序列按位置比較后排列，盡可能反映序列間的相似性，這一過程稱為序列比對。

生物學上一個基本問題是，一個基因或者蛋白是否和別的基因或者蛋白存在聯(lián)系。蛋白質在序列層次的關聯(lián)暗示了其同源性，同時暗示了它們可能有類似的功能。通過分析多個DNA或者蛋白質序列，我們可能會發(fā)現(xiàn)一些保守序列和區(qū)域。這種基因序列或者蛋白質序列之間關聯(lián)性的分析就是通過序列比對來完成的[2]。如今，人們完成了許多不同生物體的基因組分析，獲得了大量的基因和蛋白質序列數(shù)據(jù)。序列比對能夠向我們揭示某一生物體中蛋白質之間的關聯(lián)以及蛋白質在不同生物體中的關聯(lián)，從而進一步理解生命和進化。

序列比對的主要用途：

（1）用于搜索相似序列,當你獲得一段DNA序列或氨基酸序列后，發(fā)現(xiàn)對它一無所知時，可以在已知的核酸序列數(shù)據(jù)庫中搜索關于這一序列的信息，一個有效的方法是采用比對算法在數(shù)據(jù)庫中找到一系列與該序列有相似性的序列，并按相似程度由高到低排列。

（2）在基因組測序中主要用于，通過比對算法找到有重疊區(qū)的另外的序列片段，把它們合并起來還原成原來的長核酸序列，得到長核酸序列的堿基排列順序。

（3）尋找序列中的特定位點。當一個基因的某一位點發(fā)生突變時，它與原基因進行比對時就能發(fā)現(xiàn)這個位點，這在尋找致病基因時尤為重要。

（4）預測蛋白質結構，通過把待測的序列與已知結構的序列進行比對，根據(jù)相似性去預測待測序列局部或全部的結構。

2 序列比對與記分規(guī)則

2.1 全局比對與局部比對

序列比對根據(jù)對比對后要排列的片斷范圍可將比對分為全局比對與局部比對。全局比對是全部待研究的全部序列的全部符號參加比較，最后也是全部序列的全部符號進行排列與計分，比對的結果中各序列長度相同；而局部比對是全部序列的全部符號參加比較，最后只將各序列中得分高的片斷中的符號進行排列與計分，即只排列局部的序列片斷[5]。

相比于全局比對，局部比對更能反映序列間的相似性。因為，蛋白質功能位點往往是由較短的序列片段組成的，盡管在序列的其它部位可能有插入、刪除等突變，但這些關鍵的功能部位的序列往往具有相當大的保守性。而局部比對往往比整體比對對這些功能區(qū)段具有更高的靈敏度，因此其結果更具生物學意義。

2.2 最優(yōu)匹配與記分規(guī)則

序列比對的結果有很多種可能，而目標是要找出最優(yōu)比對結果，即使匹配項盡可能地多一點，而引入的空格盡可能地少一些。為此，我們引入了記分規(guī)則，用以給每一個匹配結果打分。而分數(shù)越高說明把序列轉換為對應序列的代價越小，相似性也就越大，所以最后得分最高的匹配即為最優(yōu)匹配。

計分規(guī)則是比對的重要條件，計分方法的生物學意義常常就決定了比對所反映的生物學特征。在使用差異較大的不同計分方法時將會產生不同的比對結果。根據(jù)所代表的生物學意義可以粗略地將計分方法分為三類：匹配計分、結構與性質計分、可觀察變換計分。匹配計分的規(guī)則是字符進行比較時只有3至4個分值：兩個字母相同（匹配）獎勵一個分值、兩個字母不同（錯配）罰一個分值、字母對空格罰1至2個分值。例如常用的生物信息學軟件BLAST中的核酸比對計分就是采用匹配計分。由于這種方法簡單，較容易用它說明比對的一般原理，所以本文的核酸序列比對都采用這種方法。

3 比對算法

有不少的序列比對算法已出現(xiàn)在文獻及應用軟件中，其中一些得到廣泛的應用，如動態(tài)規(guī)劃法、累進方法等。兩序列比對與多序列比對的算法有差異，兩序列比對的經典方法是動態(tài)規(guī)劃法。全局比對動態(tài)規(guī)劃法是Needleman與Wunsch在1970年提出[1]，一直沿用至今，這個算法是生物信息學的基礎算法之一。局部比對的動態(tài)規(guī)劃法是Temple F.Smith和Michael S.Waterman于1981年提出的[3]。Smith-Waterman算法實際上是Needle-Wunsch算法的一個變體。二者都是動態(tài)規(guī)劃算法。這一算法的優(yōu)勢在于可以在給定的打分方法下找出兩個序列的最優(yōu)的局部比對。

動態(tài)規(guī)劃算法的思想大致上是，若要解一個給定問題，我們需要解其不同部分（即子問題），再合并子問題的解以得出原問題的解。通常許多子問題非常相似，為此動態(tài)規(guī)劃法試圖僅僅解決每個子問題一次，從而減少計算量：一旦某個給定子問題的解已經算出，則將其記憶化存儲，以便下次需要同一個子問題解之時直接查表[4]。這種做法在重復子問題的數(shù)目關于輸入的規(guī)模呈指數(shù)增長時特別有用。動態(tài)規(guī)劃高效地解決了如何將一個龐大的數(shù)學問題分解為一系列小問題，并且從一系列小問題的解決方法重建大問題的解決方法的過程。

基于動態(tài)規(guī)劃的這一基本思想，序列比對算法的主要過程大致是：

首先，一個DNA螺旋由一串被稱為基的分子組成，可能的基包括腺嘌呤（Adenie）、鳥嘌呤（Guanine）、胞嘧啶（Cytosine）和胸腺嘧啶（Thymine）。分別以它們的首字母來代表這些基，于是，一個DNA螺旋就可以表示為在有窮集合A,G,C,T上的一個串。

以全局比對為例進行說明，對序列X=ACT與Y=CGG進行比對，比對的結果有很多種可能，其中一部分的可能比對結果如下：（規(guī)定空格對空格無效?。?/p>

將比對全過程分為若干步，每一步增加一個位置。因為空格對空格無效，所以增加一個位置時有三種情況：第一個序列增加一個字母而第二個序列增加一個空格；第一個序列增加一個空格而第二個序列增加一個字母；兩個序列都增加一個字母。這樣要進行n步的話就可能有3n種可能。動態(tài)規(guī)劃算法的巧妙之處是把第一序列已比對字母且第二序列已比對字母都相同的各種比對結果放在一起進行判斷，只留最優(yōu)結果。用這種篩選方面一直進行下去，直到所有的字母都進行過比對為止。最后所得的最優(yōu)解就是動態(tài)規(guī)劃算法的最后結果。以序列GC和序列AT的比對中間過程為例：

以上三種中間比對都是：第一序列的G已比對且第二序列的A已比對。而顯然這三種可能中，第三種的比對最優(yōu)，因為它的空格引入最少。于是舍棄了前兩個可能，只保存第三中可能的結果。

因此，用動態(tài)規(guī)劃算法進行兩序列比對的過程可用矩陣顯示，矩陣中的每一元素可表示第一序列已比對字母且第二序列已比對字母相同的各種比對結果的最優(yōu)者，最后的一格（即右下格）的最優(yōu)結果就是整個比對的最優(yōu)結果。在具體算的過程中，每一格只用最優(yōu)比對的得分來表示。矩陣的計算過程可表示如下：

給出的長度為m的序列x和長度為n的序列y，建立一個（n+2）×（m+2）的矩陣，函 數(shù) F（i,j）遞歸求序列 x[1，…，i]和 y[1，…，j]的最佳比配的得分。

其中：

而g=-2（與空格對齊的話罰兩分）。

4 比對算法的實現(xiàn)

4.1 動態(tài)規(guī)劃的全局比對的過程

算法的整體功能是：

·輸入：將要比對的兩條序列，x（A,G,C）,y（A,A,A,C）,以及打分函數(shù)

·輸出：最大比對得分所對應的比對結果。主要步驟如下：

（1）建立得分矩陣

建立一個（n+2）×（m+2）的矩陣，其中，m，n分別為X,Y的長度。將序列x（欲比較的序列）從第一行的第三列開始填入，而序列y（參照序列）從第一列的第三行開始填入，如表1所示。

表1 得分矩陣

（2）選擇得分體系

字符的比對有三種情況：匹配，錯配，插入或刪除。因為插入和刪除都是引入空位的操作（在不同序列上）。序列1上某一位置的插入相當于在序列2上對應位置的刪除，所以，在實際計算中，只需要考慮何時引入空位即可。

這三種得分情況有很多打分標準，不同的打分標準會得出不同的比對結果。這里將選用的是最簡單的打分體系，即，匹配得一分，錯配罰一分，空格罰兩分。

（3）填充得分矩陣

開始于第二行中的第二列，初始值為0。通過矩陣一行一行,從左向右地移動，計算每個位置的分數(shù)。每一單元格的最佳得分（即最高）是從已有的左側，頂部或左上方（對角線）最佳得分計算而來，最后取這三者中的最大值與該單元格的得分值相加，并用一個箭頭指向得分的來源處（在某些位置，可能存在著兩個或三個相同的最大得分）。當一個得分從頂部計算，或從左邊這代表在對準時發(fā)生的插入缺失。當我們從對角線計算分數(shù)這表示兩個字母所在位置匹配的對準。因為不存在“向上”或“左上”的位置對第二行，我們只能從左側現(xiàn)有單元求和。因此，我們每一步分數(shù)都-2，代表著相對于上一個的比分發(fā)生了插入或缺失。這導致在第一行是0，-2，-4，-6。這同樣適用于第二列，因為對于第二列而言，只存在向上的分數(shù)。因此，我們有圖1。

圖1 矩陣填充圖

（4）回溯得出最優(yōu)匹配

從矩陣的最右下角的單元格開始，根據(jù)指向上一個位置的箭頭來回溯取得最佳匹配的路徑。然后根據(jù)從起始位置指向最終位置的路徑來完成匹配的建立。如果最佳匹配路徑有多個的話，在回溯開始前，先決定好回溯時優(yōu)先選擇追蹤箭頭的方向，如優(yōu)先選指向高處的箭頭或指向低處的箭頭。如圖2所示：

圖2 回溯圖

如果選擇第一種回溯方式，表明在回溯時，優(yōu)先選擇指向高處的箭頭（即首先選1箭頭，如果沒有就選擇2箭頭，如果沒有則選3箭頭）。選擇第二種回溯方式則反之。選擇不同的回溯方式，最后得出的最佳匹配結果也不一樣。下面給出以上兩種回溯方式下得出的匹配結果，如圖3所示。

圖3 結果圖

若以上面第一種從高到低的偏好來選擇回溯方式的話，得出的匹配結果為（圖3綠線途徑）：

最佳匹配得分為-1。而以第二種從低到高的順序選擇回溯方式的話，得出的匹配結果則為（圖3黃線路徑）：

最佳匹配得分為-1。

（5）具體代碼實現(xiàn)

首先填充矩陣，即計算每一個單元格的得分值。在具體代碼實現(xiàn)的時候，實際上數(shù)組并沒有保存指向上一個最優(yōu)值的箭頭，而只是記錄了分值而已。由于矩陣的第一行和第一列分別用來存儲序列x,y，所以計算從第二行第二列開始。m ←length（x）

n ← length（y）

for j=1 to（m+1）

F（1,j）← g*（j-1）

for i=1 to（n+1）

F（i,1）← g*（i-1）

for i=2 to n+1

for j=2 to m+1

{

Match ← F（i-1,j-1）+S（x_i,y_j）

Delete←F（i-1,j）+g /*代表序列 x在此處發(fā)生了deleted*/

Insert←F（i,j-1）+g /*代表序列 x在此處發(fā)生了inserted*/

F（i,j）← max（Match,Insert,Delete）

}

以下代碼將實現(xiàn)如何回溯得出最佳匹配

AlignmentX←""

AlignmentY←""

i←n+1

j←m+1

while（i＞ 1 or j＞ 1）

{

if（i＞ 1 and F（i,j）==F（i-1,j）+g）

{

AlignmentX←"-"+AlignmentX

AlignmentY←y_i+AlignmentY

i←i-1

}

else if（j＞ 1 and F（i,j）==F（i,j-1）+g）

{

AlignmentX←x_j+AlignmentX

AlignmentY←"-"+AlignmentY

j←j-1

}

else

{

AlignmentX←x_i+AlignmentX

AlignmentY←y_j+AlignmentY

i←i-1

j←j-1

}

}

4.2 動態(tài)規(guī)劃的局部比對的實現(xiàn)

基于動態(tài)規(guī)劃的局部比對算法最著名的算法是Smith-Waterman算法。它是一種進行局部序列比對（相對于全局比對）的算法，用于找出兩個核苷酸序列或蛋白質序列之間的相似區(qū)域。該算法的目的不是進行全序列的比對，而是找出兩個序列中具有高相似度的片段。與全局比對的主要區(qū)別在于該算法不存在負分（負分被替換為零），它使得對序列局部的比對成為可能。當任何位置分數(shù)低于0，即表示此前的序列不具有相似性；而此時將此元素分數(shù)設為0，即達到了“重置”的效果，使得從此位置開始的比對不受之前位置的影響。初始化階段的區(qū)別使得該算法可從一序列任意位置開始匹配另一序列不受罰分。最后，回溯是從分數(shù)最高的矩陣元素開始，直到遇到分數(shù)為零的元素停止。分數(shù)最高的局部比對結果在此過程中產生。

以序列 x（A,A,G,A）和序列 y（T,T,A,A,G）為例進行局部比對，直觀效果圖如下：從圖4可看出，局部比對的結果為：

圖4 局部比對圖

得分為3。而如果進行全局比對的話，最大得分則為-3。

（1）局部比對分值計算過程

局部比對的遞歸函數(shù)與全局比對的稍有不同，因為不存在負分，因此將會引入零進行比較。如果從已有的相鄰三個分值計算來的得分均小于零的話，表示此前的序列不具有相似性。于是重置為零，接下來將重新計分。

（2）局部比對的回溯過程

首先需要找出得分矩陣里的最大分值，然后從最大分值處回溯直到遇上零停止回溯，即得出一個局部最佳匹配。

max←0

i_{max}←0

j_{max}←0

for（i=2;i＜=n+1;i++）

{

for（j=2;j＜=m+1;j++）

{

if（F[i,j]＞ max）

{

max=F[i,j];i_{max}=i;j_{max}=j;

}

}

}

AlignmentX←""

AlignmentY←""

cur_i←i_{max}

cur_j←j_{max}

while（F[cur_i,cur_j]!=0）

{

if（ cur_i＞1 and F[cur_i,cur_j]==F[cur_i-1,cur_j]+g）

{

AlignmentX←"-"+AlignmentX

AlignmentY←y_i+AlignmentY

cur_i←cur_i-1

}

else if（ cur_j＞1 and F[cur_i,cur_j]==F[cur_i,cur_j-1]+g）

{

AlignmentX←x_j+AlignmentX

AlignmentY←"-"+AlignmentY

cur_j←cur_j-1

}

else

{

AlignmentX←x_i+AlignmentX

AlignmentY←y_j+AlignmentY

cur_i←cur_i-1

cur_j←cur_j-1

}

}

5 結語

全局比對與局部比對都各有優(yōu)勢，因此，在現(xiàn)實生活中，如果進行匹配的兩條序列的長度相差很大的話，當進行全局性配對時，匹配程度將會很低。因為該配對方法要求將兩個序列的所有堿基都進行配對，所以當兩序列長度相差很大時，將會造成配對結果含有很長的間隔。盡管在某處，配對序列和參考序列的匹配程度相當高，在全局性配對時，將不會把這兩個片段直接進行配對，而是盡量嘗試讓參考序列內的所有堿基都參與匹配，造成結果與局部性配對方法有極大差別。因此，在選擇配對算法的時候，必須參考兩序列的具體長度和配對的實際意義和需要，選擇適當?shù)姆椒?，才能得到有生物學意義的結果。一般來說，全局性配對主要應用于比較同源性的基因的相似性而局部性配對主要應用于在非同源性的基因序列中尋找具有相似性的同源的基因域。

當堿基對或殘基對之間匹配會導致更低分數(shù)時，需要空位的引入，即讓堿基或殘基與空位匹配?？瘴涣P分決定了插入或者刪除的分值。無論是全局比對或局部比對，本文在選擇打分機制時都是采用的是最簡單的線性打分機制，即：

w（k）=k × g

其中，k為空位的個數(shù)，g為每一個空位的罰分值。這種線性空位罰分方式為每次插入或者刪除的得分相同。然而在生物學上，兩序列之間一般存在著具有不同相似度的區(qū)域。另外，單個基因突變事件可能導致一長串空位的插入。因此，一個連續(xù)的較長的空位優(yōu)于多個分散的小的空位。雖然多個分散的小的空位可以產生更多匹配，但一個連續(xù)的較長的空位代表這個區(qū)域只在一個序列中出現(xiàn)，使用后者可以避免為了得到高分而過度匹配這段序列。

[1]Saul B.;Needleman and Christian D Wunsch.A General Mmethod Applicable to the Search for Similarities in the Amino Acid Sequence of Two Proteins.Journal of Molecular Biology,1970.

[2]Durbin R,Eddy S,Krogh A et al.生物序列分析,蛋白質和核酸的概率論模型[M].北京清華大學出版社,2002.

[3]Michael S.Waterman Temple F.Smith.Identi fi cation of Common Molecular Subsequences.Journal of Molecular Biology,1981.

[4]Charles E.Leiserson Thomas H.Cormen.Introduction to algorithms.機械工業(yè)出版社出版,2002.

[5]David W.Mount.Bioinformatics Sequence and Genome Analysis.科學出版社,2002.

Dynamic Programming in Pair Sequences Alignment

LUO Zhi-bing
（College of Big Data and Intelligence Engineering,Southwest Forestry University,Kunming 650224）

DNA or Protein sequences analysis is the most important job in bioinformatics and sequences alignment is the most fundamental operation.Similarities between the sequences are found by sequences alignment.Those similarities can be used to analyze homology in DNA or protein sequences and derivate evolution information.Therefore,sequences alignment is really important in derivation of information about protein function and structures.Nowadays,improving correction and reducing time consumption in searching algorithms has become the main problem in this field.By doing research in algorithms about solving pair sequences alignment,puts forward dynamic programming algorithm in solving both global alignment and local alignment,and also gives an implementation of dynamic programming.

Bioinformatics;Sequences Alignment;Dynamic Programming

1007-1423（2017）32-0028-07

10.3969/j.issn.1007-1423.2017.32.007

羅志兵（1989-），女，廣西壯族自治區(qū)貴港人，在讀碩士研究生，研究方向為操作系統(tǒng)研發(fā)、算法分析

2017-08-29

2017-10-25