生物信息學(xué)序列拼接課件

1、生物信息學(xué)序列拼接基因組序列拼接基因組序列拼接生物信息學(xué)序列拼接序列拼接序列拼接 序列拼接任務(wù)即將測(cè)序生成的reads短片段拼接起來，恢復(fù)出原始的序列。該問題是序列分析的最基本任務(wù)，是基因組研究成功與失敗的關(guān)鍵，拼接結(jié)果直接影響到序列標(biāo)注，基因預(yù)測(cè)、基因組比較等后續(xù)任務(wù)。 基因組序列的拼接也是基因組研究必須解決的首要難題。其困難不僅來自它的海量數(shù)據(jù)(以人類基因組序列為例，從數(shù)量為10兆級(jí)的片斷恢復(fù)出長度為億級(jí)的原始序列)，而且源于它含有高度重復(fù)的序列。生物信息學(xué)序列拼接生物信息學(xué)序列拼接拼接問題的難點(diǎn)拼接問題的難點(diǎn) DNA測(cè)序數(shù)據(jù)有其固有的四個(gè)的特點(diǎn)，他們也正是解決實(shí)際的序列拼接問題的難點(diǎn)所在

2、： 1測(cè)序有誤差 2不完全覆蓋性 3序列所在鏈不確定 4重復(fù)序列的干擾生物信息學(xué)序列拼接1測(cè)序有誤差測(cè)序有誤差 由于測(cè)序技術(shù)的局限，難免會(huì)出現(xiàn)測(cè)序錯(cuò)誤，尤其是在序列的末端，一般錯(cuò)誤率可控制在1以下。所以對(duì)每個(gè)堿基一般有一個(gè)正確概率，以質(zhì)量打分的形式給出。因此每個(gè)ri都有個(gè)可信度。而read與read之間有不同程度的重疊，由此導(dǎo)致有的重疊可信度高，有的重疊可信度低。生物信息學(xué)序列拼接2不完全覆蓋性不完全覆蓋性 不是所有的堿基被測(cè)序的次數(shù)都等于平均測(cè)序覆蓋度。極端的情況，可能會(huì)出現(xiàn)源基因組序列上部分區(qū)域未被測(cè)序的情況(這段區(qū)域稱為gap)。即，測(cè)序的reads集合不是原始基因組序列一個(gè)完整覆蓋。此

3、時(shí)需要借助于各種圖譜如：基因組指紋圖譜(genome fingerprint map)， 基因組級(jí)物理圖譜(genome-wide physical map)，細(xì)胞發(fā)生圖譜(cytogenetic maps)等協(xié)助對(duì)reads進(jìn)行定位生物信息學(xué)序列拼接3序列所在鏈不確定序列所在鏈不確定 由于測(cè)序過程中無法確定特定片斷屬于DNA雙鏈中的哪一條鏈上，所以我們?cè)谄唇舆^程中并不清楚使用的是read的正義鏈，還是其互補(bǔ)鏈。4重復(fù)序列的干擾重復(fù)序列的干擾 DNA序列自身含有高度重復(fù)的子序列，它們一種表現(xiàn)為短序列的串級(jí)重復(fù)，比如：(GGAA)n?；駻mTn等。另一種表現(xiàn)為大量相似序列(其拷貝數(shù)可達(dá)幾十萬)散

4、布在基因組的各個(gè)地方。Repeat的存在，將導(dǎo)致fragments間overlap的不真實(shí)性，進(jìn)而產(chǎn)生錯(cuò)拼的結(jié)果。因此在拼接過程中耍確定這些序列的形式及大小，才能保證以高概率恢復(fù)出其在原始真實(shí)序列中的位置生物信息學(xué)序列拼接拼接算法評(píng)價(jià)拼接算法評(píng)價(jià) 以上拼接問題的四個(gè)難點(diǎn)不僅極大的增加了解決實(shí)際拼接問題的難度，而且從某種程度上說無法完整地恢復(fù)出原始DNA序列來。即實(shí)際上僅能構(gòu)建出若干個(gè)contig(重建的fragments的一種排列形式，它覆蓋基因組上一段連續(xù)區(qū)域)這些contig將指導(dǎo)測(cè)序項(xiàng)目finishing階段的實(shí)驗(yàn)方法最終構(gòu)建DNA完整序列。生物信息學(xué)序列拼接 目前，國際上對(duì)拼接軟件的公

5、認(rèn)評(píng)價(jià)標(biāo)準(zhǔn)包括兩方面，即重建出的contig的數(shù)目和準(zhǔn)確度。我們發(fā)展的基因組序列拼接新算法的目標(biāo)是在確保準(zhǔn)確性的前提下，構(gòu)建盡量少的contig，以減少測(cè)序后期大量的人力和財(cái)力的投入。生物信息學(xué)序列拼接基因組序列拼接算法研究現(xiàn)狀基因組序列拼接算法研究現(xiàn)狀 現(xiàn)在最常用的拼接程序使用的拼接算法可分成兩類，一類是將拼接問題轉(zhuǎn)化為在圖中尋找的Hamilton路徑的問題；另一類是將拼接問題在某種特殊情況下轉(zhuǎn)化成尋求圖中的Euler路徑的問題。他們均有其成功的典型算法。生物信息學(xué)序列拼接1.轉(zhuǎn)化為轉(zhuǎn)化為Hamilton Path問題問題 每個(gè)DNA片段(read)相當(dāng)于圖中一個(gè)結(jié)點(diǎn)，如果兩個(gè)片段之間存在著

6、重疊(overlap)關(guān)系，則在兩個(gè)結(jié)點(diǎn)之間定義一條邊，而沿著DNA原始序列從頭到尾，則必然經(jīng)過每個(gè)結(jié)點(diǎn)一次且僅一次，即是一條Hamilton路徑。一條contig表示圖中一條簡單路，此類算法以Phrap，TIGR Assembler，CAP3，GigAssemble等為代表。生物信息學(xué)序列拼接 他們都是遵循“overlap-layout-consensus”的框架。首先，為了構(gòu)建圖。計(jì)算任意兩個(gè)read間可能的比對(duì)情況。其次，通過去除歧義的或者不確信的邊得到較為準(zhǔn)確的圖，并在其上尋找非交叉的簡單路的集合，該集合對(duì)應(yīng)于contig的集合。最終，通過對(duì)包含在一個(gè)簡單路上的所有read進(jìn)行多序列比

7、對(duì)，為每一個(gè)contig構(gòu)建一個(gè)一致性序列(consensus sequence)。生物信息學(xué)序列拼接2.轉(zhuǎn)化為轉(zhuǎn)化為Euler Path問題問題 EULER是這類算法的代表。與傳統(tǒng)方法沿著“OverlapLayoutConsensus”路線不同，它不計(jì)算各個(gè)read之間的Overlap，即沒有Overlap步驟。 它的大致想法如下： 為了排除read中的錯(cuò)誤，獲得Error-Free的read，將所有的read切割成小片n-mers。生物信息學(xué)序列拼接 將每個(gè)read和Gk的近似進(jìn)行比對(duì)，尋求read的最小改變能夠使得read的所有n-mers包含在Gk的近似集合中。從而構(gòu)建了高質(zhì)量序列，而

8、對(duì)于Poor read，直接拋棄，對(duì)Chimeric read(兩端在n-mers中但整體不在的reads)進(jìn)行特殊處理。生物信息學(xué)序列拼接 初始的想法是要實(shí)現(xiàn)去除reads中的測(cè)序錯(cuò)誤的目的，如果知道原始序列G，那么直接使用測(cè)序獲得的read和G進(jìn)行比較即可。 但是實(shí)際上G并不可知，那么退而求其次， G的序列片斷Gk亦可，事實(shí)上Gk亦不可知。所以將所有的read切割成小片n-mers，所有Solid的n-mers形成的集合稱為Gk的近似。最后，構(gòu)造De Bruijn圖。生物信息學(xué)序列拼接現(xiàn)有算法的主要問題現(xiàn)有算法的主要問題 雖然已經(jīng)開發(fā)了以上的算法，基因組序列拼接問題尚未徹底解決，以上兩類算

9、法都存在著各自的缺陷。生物信息學(xué)序列拼接 對(duì)于第一類算法來說，實(shí)際上是在圖中尋找一條使得評(píng)價(jià)函數(shù)值最優(yōu)的Hamilton路徑，這是一個(gè)NP完全問題。 一般都采用greedy-merging的算法近似求解。由于這種step-by-step的局部貪心算法，其明顯的局部特性忽略了reads間“長距離”或者整體性的聯(lián)系，從而導(dǎo)致了拼接錯(cuò)誤，即拼接結(jié)果和真實(shí)的DNA原始序列不同。最近研究指出，在對(duì)已知序列的流行性感冒嗜血桿菌基因組的拼接過程中，無論是Phrap，TIGR Assembler，還是CAP3，都發(fā)生了拼接錯(cuò)誤的現(xiàn)象。生物信息學(xué)序列拼接 對(duì)于第二類算法來說，它只能在特殊的情況下，才能將問題簡化

10、成尋求一條Euler路徑，最終的結(jié)果是從多條候選的Euler超路中選擇出來的。EULER算法依然存在拼接錯(cuò)誤，且結(jié)果選擇的過程沒有理論依據(jù)。EULER軟件在實(shí)際數(shù)據(jù)集的運(yùn)行速度上和第一類算法相當(dāng)。更重要的是，EULER采用的算法過于獨(dú)立，很難利用其他輔助生物信息，導(dǎo)致其實(shí)用性和流行性大打折扣。生物信息學(xué)序列拼接局部搜索局部搜索(Local Search)方法方法胡杰胡杰生物信息學(xué)序列拼接 將局部搜索方法運(yùn)用于一個(gè)具體的問將局部搜索方法運(yùn)用于一個(gè)具體的問題，需要對(duì)以下四項(xiàng)進(jìn)行明確的定義。題，需要對(duì)以下四項(xiàng)進(jìn)行明確的定義。 1將原同題表示成一個(gè)優(yōu)化問題，即定義一個(gè)可行域以及在該可行域上的一個(gè)目標(biāo)函

11、數(shù)。 2定義可行域中的鄰域結(jié)構(gòu)，即說明滿足什么條件的兩個(gè)點(diǎn)相鄰。生物信息學(xué)序列拼接3確定在鄰域中的搜索方式。4局部極值點(diǎn)的處理。如果當(dāng)前解點(diǎn)鄰域中的所有點(diǎn)的目標(biāo)函數(shù)值都比當(dāng)前點(diǎn)大，那么這點(diǎn)就稱為局部極值點(diǎn)。在一些問題中，局部極值點(diǎn)就是全局最優(yōu)點(diǎn)。而對(duì)另一些問題而言，局部極值點(diǎn)已經(jīng)滿足求解實(shí)際問題的需求。生物信息學(xué)序列拼接基于局部搜索的序列拼接算法框架基于局部搜索的序列拼接算法框架 主要目標(biāo)主要目標(biāo)是在layout階段盡可能準(zhǔn)確的前提下，獲得更長的contig。具體的，使用局部搜索算法求得數(shù)據(jù)集上近似全局的最短超串；然后根據(jù)求得的最短公共超串對(duì)應(yīng)的fragment的排列關(guān)系為基礎(chǔ)獲得“conse

12、nsus segment”生物信息學(xué)序列拼接1.Overlap定義定義 如果一個(gè)串的前綴是另一個(gè)串的后綴則認(rèn)為這兩個(gè)串之間存在overlap，并根據(jù)over-lap構(gòu)建超串。對(duì)給定的串f和g，存在多個(gè)可能的overlap關(guān)系 比如說，若f=ACTGGGAGCAGC， g=AGCAGCTTTTACT， 那么他們之間至少存在兩個(gè)overlap形式。生物信息學(xué)序列拼接生物信息學(xué)序列拼接 在我們的算法中，僅考慮兩個(gè)串之間最大的overlap情況，并定義overlap(f，g)表示f和g之間存在的多個(gè)overlap關(guān)系中最長的一個(gè)overlap所包含的字符個(gè)數(shù)。 在上面的例子中overlap(f，g)=

13、6。如果f和g之間overlap區(qū)域長度小于M(M是一個(gè)足夠小的正整數(shù))，則overlap(f，g)=0。生物信息學(xué)序列拼接2.優(yōu)化目標(biāo)定義優(yōu)化目標(biāo)定義 我們對(duì)reads集合S中的每個(gè)元素按照其左端在超串t中首次出現(xiàn)的位置進(jìn)行排序，并沿超串t從左向右依次讀入，并將其對(duì)應(yīng)為序S=sl，s2，Sn)。我們用P(S)表示這個(gè)由字符串為元素構(gòu)成的序。在序在序P(S)中，對(duì)于每一個(gè)連續(xù)的字符串元素對(duì)中，對(duì)于每一個(gè)連續(xù)的字符串元素對(duì)si， si+1，都存在都存在overlap(si，Si+1)。生物信息學(xué)序列拼接 因此，字符串的一個(gè)排列等價(jià)于一個(gè)超串t及作用在其上的函數(shù)overlap，超串t的長度leng

14、th(t)=t =ssIs|- overlap(si，Si+1), 為了求超串t具有最小的長度。因?yàn)樵诮o定集合s中ssIs是確定的值，我們可以進(jìn)一步把問題轉(zhuǎn)化成尋找一個(gè)排列P使得 最大化。ni1生物信息學(xué)序列拼接3.鄰域定義鄰域定義 在使用局部搜索方法前，我們首先定義排列的鄰域這一概念。我們稱一個(gè)排列為問題的一個(gè)解。設(shè)reads集合S是一個(gè)具有n個(gè)元素的字符串集合，P(S)是reads集合S所有解的集合。定義操作定義操作rshift(i，j) ，即，即對(duì)一個(gè)解對(duì)一個(gè)解PP(S)，將，將P第第i位置的元素向右位置的元素向右移至第移至第j個(gè)位置個(gè)位置(1ijn) 。生物信息學(xué)序列拼接也就是說，如果

15、我們應(yīng)用這一操作于給定的解P=sl，s2，sn)，則結(jié)果P=rshift(i，j)(P)定義為；生物信息學(xué)序列拼接 類似的，定義操作lshift(i,j)將P第i位置的元素向左移至第j個(gè)位置(1jin)，定義P=lshift(i，j)(P)為；生物信息學(xué)序列拼接 這樣，一個(gè)解PP(s)的鄰域可以定義為N(P)=rshift(i，j)(P)| |1ijn) lshift(i，j)(P)|1|1j0時(shí)，p比p更好，則由p轉(zhuǎn)移到p生物信息學(xué)序列拼接 此外，對(duì)給定的解P=s1，s2Sn,內(nèi)部的元素均包含前趨和后繼。但是，頭元 素僅有后繼而尾元素僅有前驅(qū)。在算法中，為了消除頭尾元素差異，我們將排列的頭尾

16、 元素連結(jié)起來并使用局部搜索方法尋找最優(yōu)的“Loop”超串。接著在overlap最小的地方切分該環(huán)狀超串，最終還原成一條線性超串。生物信息學(xué)序列拼接5.局部極值點(diǎn)的處理局部極值點(diǎn)的處理 當(dāng)經(jīng)過以read為單元的搜索后，可以獲得一個(gè)當(dāng)前鄰域內(nèi)的局部最優(yōu)解1,k1,k1+1, ,k2,k2+1，k3km。它對(duì)應(yīng)集合s上的一個(gè)superstring。該解對(duì)應(yīng)的reads的排列中，任意相鄰兩個(gè)reads間的overlop關(guān)系薄弱，即overlap(i，i+1)M(其中，M是一個(gè)足夠小的正整數(shù)，1i n)。生物信息學(xué)序列拼接例如： 如果overlap(rkl，rkI+1)M， overlap(rk2，r

17、k2+1)M,， overlap(rkm，rkm+1)M1 ， overlap(r1，r3)M1，但overlap(r2，r3)M2和overlap(r3，r2)50k)Contigs平均長度平均長度錯(cuò)誤錯(cuò)誤contigs數(shù)數(shù)目目錯(cuò)誤錯(cuò)誤contigs總長總長(bp)AE0078725428697.5CAP35414652132876300PHRAP5405881932845200B-LSA5410821442864800AE0006577.5CAP3777199886144539PHRAP717215003104164B-LSA5011309192101102INRUENSS8CAP3681

18、2264085227941PHRAP6314287484188705B-LSA44124116215219AE0007827.5CAP310311215096279380PHRAP9411223814212060B-LSA751328949339270AE0078697.5CAP31261422436136274PHRAP1161524272229589B-LSA99203148500AL0091267.5CAP317110243586138549PHRAP15519256955168824B-LSA12321314849122422CAP3、PHRAP、B-LAS的性能比較的性能比較生物信息

19、學(xué)序列拼接結(jié)果一contigs的個(gè)數(shù)contig平均長度B-LSA40833.661PHRAP53325.483CAP356623.035-結(jié)果：結(jié)果：產(chǎn)生更少的產(chǎn)生更少的contig，而且，而且contig更長更長有效的減少了測(cè)序項(xiàng)目后續(xù)工作的工作量有效的減少了測(cè)序項(xiàng)目后續(xù)工作的工作量生物信息學(xué)序列拼接結(jié)果二總contigs數(shù)錯(cuò)誤數(shù)目錯(cuò)拼率B-LSA408153.67CAP3745526.97PHRAP566244.24-結(jié)果：結(jié)果：包含更少的錯(cuò)誤拼寫現(xiàn)象包含更少的錯(cuò)誤拼寫現(xiàn)象生物信息學(xué)序列拼接 基于原型算法框架實(shí)現(xiàn)的基于原型算法框架實(shí)現(xiàn)的Basic LSA系系統(tǒng)在針對(duì)實(shí)際序列模擬的統(tǒng)在針對(duì)

20、實(shí)際序列模擬的shotgun reads數(shù)數(shù)據(jù)集上的拼接結(jié)果表明：據(jù)集上的拼接結(jié)果表明： 在不引入其他輔助信息的前提下，在不引入其他輔助信息的前提下，Basic LSA比比PHRAP(1999)及及CAP3都有都有明顯的改進(jìn)這也證明了原型算法的有效明顯的改進(jìn)這也證明了原型算法的有效性和可行性性和可行性生物信息學(xué)序列拼接基于原型算法的優(yōu)化基于原型算法的優(yōu)化 原型算法的運(yùn)行速度過慢，以至成為原型算法的運(yùn)行速度過慢，以至成為整個(gè)系統(tǒng)運(yùn)行的整個(gè)系統(tǒng)運(yùn)行的瓶頸瓶頸，從而導(dǎo)致該系統(tǒng)無，從而導(dǎo)致該系統(tǒng)無法滿足實(shí)驗(yàn)和可能的進(jìn)一步應(yīng)用的需法滿足實(shí)驗(yàn)和可能的進(jìn)一步應(yīng)用的需要這也成了將原型系統(tǒng)實(shí)用化必須解決要這也成

21、了將原型系統(tǒng)實(shí)用化必須解決的問題的問題.提速優(yōu)化的目標(biāo)是在不改變運(yùn)行結(jié)提速優(yōu)化的目標(biāo)是在不改變運(yùn)行結(jié)果的前提下，提高系統(tǒng)運(yùn)行速度果的前提下，提高系統(tǒng)運(yùn)行速度生物信息學(xué)序列拼接“鄰域剪枝鄰域剪枝”-(Neighborhood Pruning) 盡可能通過排除那些目標(biāo)函數(shù)值低于盡可能通過排除那些目標(biāo)函數(shù)值低于當(dāng)前值的鄰接點(diǎn)來縮小搜索空間這里的當(dāng)前值的鄰接點(diǎn)來縮小搜索空間這里的關(guān)鍵在于找到了一個(gè)關(guān)于能夠增加目標(biāo)函關(guān)鍵在于找到了一個(gè)關(guān)于能夠增加目標(biāo)函數(shù)值的變換的必要條件數(shù)值的變換的必要條件生物信息學(xué)序列拼接After transposition“領(lǐng)域剪枝領(lǐng)域剪枝”（Neighborhood Pruni

22、ng）方）方法法Befor transposition生物信息學(xué)序列拼接 在搜索過程中，每一個(gè)在搜索過程中，每一個(gè)transposition將導(dǎo)將導(dǎo)致目標(biāo)函數(shù)改變致目標(biāo)函數(shù)改變overlap，Aoverlap=Overlap4+Overlap5+Overlap6 Overlapl-Overlap2-Overlap3 若該若該transposition操作能改進(jìn)當(dāng)前解，必定操作能改進(jìn)當(dāng)前解，必定有有Aoverlap0 因此，也必然滿足以下條件因此，也必然滿足以下條件 (Overlap4-Overlapl-Overlap2) +(Overlap5+Overlap6)Overlap30生物信息學(xué)序列拼接兩種可能性：兩種可能性：1Overlap4Overlapl+Overlap2生物信息學(xué)序列拼接 基于上述基于上述Neighborhood Pruning方法，方法，我們優(yōu)化了原型算法，由于局部搜索我們優(yōu)化了原型算法，由于局部搜索(Local Search)是一種啟發(fā)式