生物序列比對中的算法課件

生物序列比對中的算法中科院計算所生物信息學實驗室華大—曙光聯合實驗室張法提綱背景知識序列相似性的比較兩條序列的比對問題多序列的比對問題一些啟發(fā)式的算法生物序列比對中的并行算法DNA(1)脫氧核糖核酸DNA的分子組成核甘(nucleotides)磷酸鹽(phosphate)糖(sugar)一種堿基腺嘌呤(Adenine)鳥嘌呤(Guanine)胞嘧啶(Cytosine)胸腺嘧啶(Thymine)DNA(2)堿基的配對原則A(腺嘌呤)—T(胸腺嘧啶)C(鳥嘌呤)—G(胞嘧啶)一個嘌呤基與一個嘧啶基通過氫鍵聯結成一個堿基對。DNA分子的方向性5'→3'DNA(3)DNA的雙螺旋結構

堿基對之間的互補能力DNA(4)DNA的復制在DNA解旋酶的作用下兩條鏈分離開，分別作為一個模板，在聚合酶的作用下合成一條新鏈。RNA、轉錄和翻譯RNA(核糖核酸)：單鏈結構、尿嘧啶U代替胸腺嘧啶T、位于細胞核和細胞質中。轉錄：DNA鏈→RNA鏈信使RNA(mRNA)，啟動子。翻譯：mRNA上攜帶遺傳信息在核糖體中合成蛋白質的過程。變異進化過程中由于不正確的復制，使DNA內容發(fā)生局部的改變。變異的種類主要有以下三種：替代(substitution)插入或刪除(insertionordeletion)

indel重排(rearrangement)蛋白質由氨基酸依次鏈接形成在生物體中總共有20種氨基酸。蛋白有十分復雜的三維結構。其三維機構決定了蛋白質的功能?；蚴裁词腔?？DNA上具有特定功能的一個片斷，負責一種特定性狀的表達。一般來講，一個基因只編碼一個蛋白質?；蚪M任何一條染色體上都帶有許多基因，一條高等生物的染色體上可能帶有成千上萬個基因，一個細胞中的全部基因序列及其間隔序列統(tǒng)稱為genomes（基因組）。DNA上的基因

基因基因的編碼基因編碼是一個邏輯的映射，表明存儲在DNA和mRNA中的基因信息決定什么樣的蛋白質序列。每個堿基三元組稱為一個密碼子(codon)堿基組成的三元組的排列共有43＝64種，而氨基酸共有20種類型，所以不同的密碼子可能表示同一種氨基酸。

帶來的問題序列排列問題基因組的重排問題蛋白質結構和功能的預測基因(外顯子、內含子)查找問題序列裝配(SequenceAssembly)問題動機在生物學的研究中，將未知序列同已知序列進行比較分析已經成為一種強有力的研究手段，生物序列相似性比較中絕大部分的問題在計算機科學領域中主要體現為字符串的匹配和查找。圖譜分析生物序列片斷的拼接基因區(qū)域的預測生物序列同源的比較蛋白質結構和功能預測兩條序列比對問題的分類全局比對(GlobalAlignment)局部比對(LocalAlignment)空位罰分(GapPenalty)全局比對(1)－定義定義1：兩個任意的字符x和y，(x,y)表示表x和y比較時的分值。

(x,x)=2,(x,y)=(x,-)=(-,y)=-1定義2：S=s1…sn和T=t1…tm，其全局比對A可以用序列S′和T′來表示，其中：(1)|S′

|=|T′|；

(2)將S′和T′中的空字符除去后所得到的序列分別為S和T；比對A的分值Score為：全局比對(2)－原始算法輸入：序列S和T，其中|S|=|T|=n

輸出：S和T的最優(yōu)比對fori=0tondofor(S的所有的子序列A，其中|A|=i)dofor(T的所有的子序列B，其中|B|=i)do

……全局比對(3)動態(tài)規(guī)劃(DynamicProgramming)Bellman在50年代提出的理論基礎－最優(yōu)化原理弱點：thecurseofdimensionalitySmith-WatermanAlgorithm全局比對(4)Smith-Waterman算法計算出兩個序列的相似分值，存于一個矩陣中。(editmatrix、DP矩陣)根據此矩陣，按照回溯的方法尋找最優(yōu)的比對序列。全局比對(5)前提條件遞歸關系計算editmatrix：fori=0tondoforj=0tomdoCalculateV(i,j)usingV(i-1,j-1),V(i,j-1),V(i-1,j)end全局比對(6)計算完editmatrix后，得到序列比對的最優(yōu)分值，通過回溯（Traceback）的方法可獲得序列的最優(yōu)比對序列。回溯的算法：for(i=|S|，j=|T|;i>0&&j>0;){if(V[i,j]=V[i–1,j–1]+σ(S[i],T[j])){i––;j––;}elseif(V[i,j]=V[i–1,j]+σ(S[i],–)){i––;insert(‘–’,T,j);}elseif(V[i,j]=V[i,j–1]+σ(–,T[i])){j––;insert(‘–’,S,i);}}例：S=“acgctg”和T=“catgt”

(x,x)=2,(x,y)=(x,-)=(-,y)=-1

ji01c2a3t4g5t00-1-2-3-4-51a-1-110-1-22c-2100-1-23g-300-1214c-4-1-1-1115t-5-2-21036g-6-3-3032三種可能的最優(yōu)比對序列：S:acgctg-T:-c–atgtS:acgctg-T:-ca–tgtS:-acgctgT:catg-t-

局部比對(1)兩條序列在一些局部的區(qū)域內具有很高的相似度。在生物學中局部比對比全局比對更具有實際的意義。局部比對(2)前提條件：

V(i,0)=0;

V(0,j)=0;遞歸關系：

找出i*和j*，使得：

局部比對(3)對全局比對策略稍作修改可得到局部最優(yōu)比對算法。比對的路徑不需要到達搜索圖的盡頭,如果某種比對的分值不會因為增加比對的數量而增加時，這種比對就是最佳的。依賴于記分系統(tǒng)的性質：因為某種路徑的記分會在不匹配的序列段減少,當分值降為零時，路徑的延展將會終止，一個新的路徑就會產生。局部比對(4)S=“abcxdex”，T=“xxxcde”記分函數(x,y)：

(x,x)=2,(x,y)=(x,-)=(-,y)=-1

ji01x2x3x4c5d6e000000001a00000002b00000003c00002104x02221105d01111326e00000257x0222114S=“abcxdex”，T=“xxxcde”

局部最優(yōu)比對是：

cxde

c-de或

x-de

xcde

空位罰分(1)空位:序列中任意連續(xù)的盡可能長的空格?？瘴坏囊胧菫榱搜a償那些插入或缺失，但是在序列的比對中引入的空位不能太多，否則會使序列的排列變得面目全非。每引入一個空位，比對的分值都會有所扣除，常見的罰分規(guī)則主要有兩種：權值恒定模型仿射罰分模型空位罰分(2)權值恒定模型：在每個空位中的空格的分值為零，即：(x,-)=(-,y)=0

比對的分值為：其中，Wg為開放一個空位的罰分空位罰分(3)仿射罰分模型(AffineGapModel)

用一個附加的罰分比例去乘空位的長度Wg+q×Ws

比對的分值為：空位罰分(4)初始條件：

V(0,0)=0;

V(i,0)=E(i,0)=Wg

+iWs;

V(0,j)=F(0,j)=Wg

+jWs;遞歸條件：

V(i,j)=max{G(i,j),E(i,j),F(i,j)};

G(i,j)=V(i-1,j-1)+σ(S[i],T[j]);

E(i,j)=max{E(i,j-1)+Ws,V(i,j-1)+Wg+Ws}

F(i,j)=max{F(i-1,j)+Ws,V(i-1,j)+Wg+Ws}.利用動態(tài)規(guī)劃計算序列最優(yōu)比對的算法的復雜度分析：時間復雜度；O(mn)

多序列比對(1)兩條序列比對問題的一般化推廣動機：攜帶更多的消息。DNA或蛋白質數據庫容量的指數級的增長，即使使用很簡單的記分函數也很難找到一種在有效時間內的解決方法，而幾乎所有的近似算法和許多的啟發(fā)式算法，實際上都是在算法的計算速度和獲得最佳比對結果的敏感性之間尋找一種權衡策略。多序列比對(2)rigorous算法兩條序列多條序列，NPhardtree_based算法和iterative算法首先從序列中找出兩條相似度最大的比對，然后按照相似度遞減的順序連續(xù)添加一些序列到最優(yōu)的比對序列中。序列非常接近原始算法基礎上的近似算法，但是它們也是非常耗時的。多序列比對(3)記分方法：用函數δ(x,y)來計算字符x和y之間的距離，兩個序列的比對的距離分值我們用V來表示：

k條序列比對的分值：為k條序列中任意兩條序列(共有Ck2條)的分值(距離)V之和，用SP(Sum_of_Pairs)來表示:動態(tài)規(guī)劃的方法在多序列比對中的時間復雜度為O(2n)k中心星算法輸入：一組含有k條序列的集合?

找出序列St，St∈?，使得∑i≠tD(

Si,St)的值最小，令A={St}逐次地添加Si∈?-{St}到A中，并使Si與St的比對的值最??；假設S1，S2，…Si-1已添加到A中，由于在分別和St進行比對的過程中需要加入一些空格，故此時A為：A={S1’，S2’，…Si-1’，S’t}。添加Si到A中，按照兩條序列比對的動態(tài)規(guī)劃算法比較S’t和Si，分別產生新的序列St’’和Si’，再按照St’’中添加空格的位置調節(jié)序列S1’，S2’，…Si+1’為S1’’，S2’’，…Si-1’’，并用St’’替換St。啟發(fā)式算法-FASTA基本思想是：一個能夠揭示出真實的序列關系的比對至少包含一個兩個序列都擁有的字(片斷)，把查詢序列中的所用字編成索引，然后在數據庫搜索時查詢這些索引，以檢索出可能的匹配，這樣那些命中的字很快被鑒定出來。

確定參數ktup，在兩個序列中查找長度為ktup的、相匹配的片段(增強點)。為了提高速度，可以通過查詢表格或hash表來完成，然后在表格中搜索與另一條序列相匹配的、長度為ktup的片段。在同一條對角線中臨近的增強點成為一個增強段。每一個增強點都賦予一個正的分值，一個增強段中相鄰的兩個增強點之間的不匹配區(qū)域賦予一定的負值。一個增強段對應于一段相匹配的子序列，分值最高的段被標記為init1。

引入indel。把那些沒有重疊(non-overlap)的增強段拼接起來(增強段的分值之和減去空位處罰)。分值最高的區(qū)域記為initn。對最有可能的匹配序列進一步評分：以增強段init1所在的對角線為中心，劃分出一個較狹窄的對角線帶，利用S-W算法，來獲得分值最高的局部比對，記作opt。決定采用initn或opt的分值，前者敏感度低但速度快。FASTA對每一個檢索到的比對都提供一個統(tǒng)計學顯著性的評估，以判斷該比對的意義。FASTA對DNA序列搜索的結果要比對蛋白質序列搜索的結果更敏感。它對數據庫的每一次搜索都只有一個最佳的比對，一些有意義的比對可能被錯過。

啟發(fā)式算法-BLAST基本思想是：通過產生數量更少的但質量更好的增強點來提高速度。BALST算法是建立在嚴格的統(tǒng)計學的基礎之上的。它集中于發(fā)現具有較高的相似性的局部比對，且局部比對中不能含有空位。由于局部比對的限制條件，在大多數情況下比對會被分解為若干個明顯的HSP(High-scoreSequencePairs)。首先確定一個終止值S、步長參數w和一個閾值t，S值通常是基于統(tǒng)計學的原理指明一個預期的終止E值，然后軟件會在考慮搜索背景性質的基礎上計算出合適的S值。使要比對的序列中包含一個分值不小于S的HSP。引入鄰近字串的思想：不需要字串確切地匹配，當有一個字串的分值高于t時，BALST就宣稱找到了一個選中的字串。為了提高速度，允許較長的字串長度W。W值很少變化，這樣，t值就成為權衡速度和敏感度的參數。一個字串選中后，程序會進行沒有空位的局部尋優(yōu)，比對的最低分值是S，當比對延伸時會遇到一些負的分值，使得比對的分值下降，當下降的分值小于S時，命中的延伸就會終止。這樣系統(tǒng)會減少消耗于毫無指望的選中延伸的時間，使系統(tǒng)的性能得以改進。BLAST的改進在1997年提出了對BLAST程序的改進算法，提高了搜索速度、敏感度和實用性?？商幚黹g隔(gap)的gappedBLAST算法PSI-BLAST算法對一個選中字串長度標準的延伸利用profile(表頭文件)的數據結構來進行搜索算法算法擴大步長，以步長為2w來搜索。允許位于不同的對角線的兩個片段拼接在一起。位于不同對角線的兩個片段拼接在一起的前提條件是：拼接后片段的分值不小于某一個終止值。執(zhí)行通常的BLAST算法，使用一種不同的記分方式，根據高度顯著比對(HSPs)的最高分值建立一個最初的profile。根據該profile反復利用BLAST算法對數據庫進行搜索，這一步實際上是根據表頭文件的統(tǒng)計結果擴展局部比對。這一過程是反復進行的，直到再沒有發(fā)現新的有意義的匹配為止。由于在每一輪都會有新的片段加入，因此在操作過程中profile需要在每一個循環(huán)結束之后更新。為適應同源比較的不通情況，BLAST得到了各種各樣的補充和發(fā)展，形成了一個程序家族。程序名檢索序列數據庫序列BLASTP蛋白質蛋白質BLASTN核苷酸核苷酸BLASTX核苷酸蛋白質TBLASTN蛋白質核苷酸

SubstitutionMatrixPAM250BLOSUM62在SUN10000上進行BLAST的測試：8CPU8GmemoryVERSION:BLASTN2.2.1[Jul-12-2001]blastall-pblastn-i$query.seq–d/mnt/database/public/NCBI/ftp.ncbi.nlm.nih.gov/BLAST/nt/2001-12-02/nt-e1e-10-o$result.out-a4>rzy_63001.y1.abdCHROMAT_FILE:rzy_63001.y1.abdPHD_FILE:rzy_63001.y1.abd.phd.1CHEM:termDYE:ETTIME:ThuJun2116:03:572001TAAAGGAGAATGACTCTTTTAAGTCATAGCTTGCTGNCATTAGTCGTTTGAAAGGTCGAGTGCTCAATATCTTAACCCTTTTACCTTTTTTGAAAAACACTAGGCTAGATAGATGCCATATCGGATTGAGGGGTTATTATTGGAGATTTGGAGGGGCATTTTTTTTTAGTTATTTTTTTTCTAAATTTTTTGGATGGTGAGGCTATATGAGGTGCTCTTGAAGATGCTCTCCCTTTATGTCCGAGTGATTGAATGGTATAATTTGCAACACACAAGCTTGCTTTTACAGAACTCATGGACTTCTTTACTTATATCGTATAATTCGCAGGAATATTTCGCGAAATGGGGAGAGCAAGGCACGGTGGATCTGAAACGCGAGCTCGACCTCCTCATCCTGACGATCGCGAGCCGCGTCCTCCTCGGCAAGGAGGTTAGGGAGACCATGTTCGCCGACGTCGTGGCATCCTTCCACGAGCTCATGGACAACAGCATGCACCTCATCAGCCTCTGCTTCCCCAACCTCCCGATCCCGCGGCACCGCCGGCGCGACACGGCGTCCGCTAGGCTCAAGGAGCTCTTCTCCCGTGCCATCCAGCTTCGCCGCGGCTCCGGCCGCGCCGAGGACGACGTGCTCCAGCGGTTCTTGGATCCAGGTACAGGGACGGCntdatabase:

mnt/database/public/NCBI/ftp.ncbi.nlm.nih.gov/BLAST/nt/2001-12-02/ntSize:4.8G1′030′000個reads

借助于一些比最初的動態(tài)規(guī)劃算法更快的啟發(fā)式算法。借助硬件來完成動態(tài)規(guī)劃的算法采用高性能計算系統(tǒng)，把問題有效地分配到多個處理器上運行，最后再把各處理器的運算結果收集起來生物序列比對中的

并行算法兩條序列比對的并行算法根據序列的相似性比較，找出兩者的最佳匹配找出從一條序列轉化到另一條序列的最佳路徑核心：動態(tài)規(guī)劃動態(tài)規(guī)劃中的數據相關Lander的處理方法如果已知第i-1和i-2行反對角線上的各項值，那么第i行反對角線上各項的值可以并行計算?！鳌睢鳌睢鳌睢鳌睢??????D12D21D31D11D22D13Dm,n●●●處理器01m時間步12m+n-2m+n-1Smith-Waterman的并行處理Row-Wavefront算法和DiagonalWavefront算法

Row-Wavefront算法是讓每個處理器順序處理動態(tài)規(guī)劃矩陣中的每一行，當一個處理器檢測到它所需要的上一行矩陣的元素還沒有計算出來時，該處理器就阻塞自己，直到所需要的元素計算出來。

處理器1處理器2處理器3序列S序列TSmith-Waterman的并行處理Row-Wavefront算法和DiagonalWavefront算法

DiagonalWavefront算法中所有的處理器同時沿著矩陣的一條反對角線進行計算，當一條反對角線的元素全部計算完，才轉到下一條反對角線計算。當一個處理器空閑的時候，它從當前的反對角線上尋找一個還沒有計算的元素進行計算。這種算法沒有嚴格的處理器計算順序序列T序列SHuang的處理方法采用了分而治之的策略對回溯算法進行了改進:按照中間的一條反對角線來分割動態(tài)規(guī)劃矩陣。(0,0)r(