數(shù)據(jù)壓縮在序列比對(duì)中的應(yīng)用

1、數(shù)據(jù)壓縮在序列比對(duì)中的應(yīng)用    數(shù)據(jù)壓縮在序列比對(duì)中的應(yīng)用呼廣躍(南開(kāi)大學(xué)數(shù)學(xué)科學(xué)學(xué)院,LMPC,天津，300071)  摘要同源或非同源長(zhǎng)基因組序列的分析比較需要高效率的比對(duì)算法。論文開(kāi)發(fā)出一個(gè)新的兩兩比對(duì)工具“超級(jí)壓縮比對(duì)”（簡(jiǎn)稱SCA），該系統(tǒng)是建立在Sequitur編碼理論專為長(zhǎng)基因組序列的兩兩比對(duì)設(shè)計(jì)。SCA是一個(gè)線性算法，并且能夠處理序列重排。實(shí)驗(yàn)證明SCA算法能夠準(zhǔn)確快速的完成長(zhǎng)基因組序列的比對(duì)。關(guān)鍵詞序列比對(duì),Sequitur碼，錨定文獻(xiàn)標(biāo)識(shí)碼A中圖分類號(hào)TP301.6Datacompressionsapplicat

2、ioninalignmentHuGuangyueSchoolofmathematcscienceofNankaiUniversity,LMPC,TianJin,300071AbstractTocompareandanalyzelargegenomicDNAsequencesofrelatedorganismsandofdifferentspecies,researchersneedefficientmethodstoalignlongsequences.WedevelopedanewtoolSuperCompressionAlignment(SCAforshort),anewsystemspe

3、ciallyforrapidglobalalignmentofgenomicsequences.ThenewsystemisbasedonSequiturcodingtheory.SCAisalinearalgorithmanditcandealwithrearrangements.SCAhasbeenprovedtoaligngenomicsequencesefficientlyandaccurately.keywordsequencealignment，Sequiturcode，anchor1引言序列比對(duì)問(wèn)題是生物信息學(xué)的一個(gè)基本問(wèn)題，它是近代生命科學(xué)研究的基礎(chǔ)，也是數(shù)學(xué)和計(jì)算機(jī)科學(xué)與生物

4、學(xué)的一個(gè)重要匯合點(diǎn)，是數(shù)學(xué)理論和計(jì)算機(jī)理論在生物學(xué)中得到成功應(yīng)用的典范6。近年來(lái)，隨著測(cè)序技術(shù)的發(fā)展，越來(lái)越多的完整同源生物的基因組序列出現(xiàn)了，這些序列的長(zhǎng)度規(guī)模極大，有的甚至數(shù)十兆到數(shù)百兆，因此迫切需要快速準(zhǔn)確的序列比對(duì)算法分析比較它們的變化關(guān)系。序列比對(duì)可以尋找序列共有的穩(wěn)定區(qū)，從總體上觀察序列的變化特征與變化趨勢(shì)，解決關(guān)于突變前序列結(jié)構(gòu)的類型與預(yù)測(cè)問(wèn)題。此外序列比對(duì)在蛋白質(zhì)空間結(jié)構(gòu)研究中也有重要應(yīng)用。早期的序列比對(duì)算法主要是為短序列設(shè)計(jì)的，比如Needleman-Wunsch算法1和Smith-Waterman算法2,它們能夠得到兩條序列的最優(yōu)比對(duì)，但算法的時(shí)間復(fù)雜度是,其中n是輸入序列

5、的平均長(zhǎng)度。Needleman-Wunsch算法和Smith-Waterman算法都是動(dòng)態(tài)規(guī)劃算法，但由于時(shí)間復(fù)雜度的限制，它們不適用于長(zhǎng)序列比對(duì)。超級(jí)兩兩比對(duì)（SPA）算法4綜合運(yùn)用概率統(tǒng)計(jì)和組合分析，把兩兩比對(duì)的時(shí)間和空間復(fù)雜度降為線性并得到序列的最優(yōu)或次優(yōu)比對(duì)，但必須經(jīng)過(guò)改進(jìn)才能用于長(zhǎng)基因組比對(duì)5。目前適合長(zhǎng)基因組序列兩兩比對(duì)的工具有很多，比如AVID，MUMer，LAGAN等58。這些算法大都是啟發(fā)式的，它們使用哈希表或者后綴樹(shù)對(duì)序列進(jìn)行錨定處理從而加快比對(duì)速度。但是這類算法處理同源性差的序列的時(shí)候效率往往不高，并且盡管空間復(fù)雜度為線性，但當(dāng)序列很長(zhǎng)時(shí)耗費(fèi)內(nèi)存仍然很大。本文的SCA算法

6、使用Sequitur壓縮編碼理論降低了空間復(fù)雜度，并能很好的處理同源性差的序列的比對(duì)問(wèn)題。2.SCA算法21Sequitur碼理論 Sequitur碼3,7是一種基于語(yǔ)法的無(wú)失真碼，它分析消息序列的語(yǔ)法結(jié)構(gòu)，把重復(fù)出現(xiàn)的片斷進(jìn)行歸類，最終實(shí)現(xiàn)數(shù)據(jù)壓縮的目標(biāo)。如果信源字符集合為，而要壓縮的消息序列記為。下面是一個(gè)消息序列的語(yǔ)法： 序列x的語(yǔ)法結(jié)構(gòu)如下： 上面就是的一個(gè)語(yǔ)法，它包含一個(gè)狀態(tài)集，和三個(gè)產(chǎn)生式規(guī)則（1），（2），（3）。使用算術(shù)碼可以將這個(gè)語(yǔ)法編譯成一個(gè)二進(jìn)制流。YK編碼最關(guān)鍵的步驟是根據(jù)約簡(jiǎn)規(guī)則對(duì)輸入序列進(jìn)行語(yǔ)法分析，找出所有的產(chǎn)生式規(guī)則和最后的狀態(tài)集。

7、YK碼的約簡(jiǎn)規(guī)則見(jiàn)文獻(xiàn)3。例：令A(yù)0,1,x=00110010100011111011，n=20，求x的YK語(yǔ)法的過(guò)程如下：(1)在時(shí)，無(wú)重復(fù)的編碼狀態(tài)，從m6開(kāi)始，這時(shí)S=001100其語(yǔ)法如下：（2）讀數(shù)據(jù)，得到，出現(xiàn)新的產(chǎn)生式規(guī)則，語(yǔ)法調(diào)整為   （3）讀入，語(yǔ)法調(diào)整為   （4）讀入，語(yǔ)法調(diào)整為  （5）如此繼續(xù)，最后得到語(yǔ)法為得到序列x的語(yǔ)法之后，譯碼工作非常簡(jiǎn)單，只要按照各個(gè)狀態(tài)的產(chǎn)生式規(guī)則去翻譯即可。YK碼的編、譯有唯一可譯性，并且有線性的計(jì)算復(fù)雜度7。22SCA的算法步驟設(shè)A、B為要比對(duì)的兩條基因組

8、序列，SCA的算法步驟如下：步驟-1：用Sequitur語(yǔ)法變化的方法對(duì)序列A、B做語(yǔ)法分析，求出A、B的聯(lián)合語(yǔ)法，并記下各個(gè)狀態(tài)變量出現(xiàn)的位置。步驟-2：利用步驟1得到的狀態(tài)變量集，可以直接求出兩個(gè)序列的匹配子串，并將匹配串向左右擴(kuò)展，直到匹配串的相似度降到某個(gè)值為止，將每個(gè)匹配記為向量，其中分別代表該匹配在A、B兩個(gè)序列中的坐標(biāo)，代表匹配子串的長(zhǎng)度。步驟-3：將所有匹配按照坐標(biāo)遞增排列，找出一個(gè)沒(méi)有重疊的匹配串的最大子集作為錨定。步驟-4：用SPA算法結(jié)合模結(jié)構(gòu)完成錨定間隔的比對(duì)。下面我們用兩條序列來(lái)解釋以上步驟.A=accgttatatB=acccatattt1.A、B聯(lián)合的YK語(yǔ)法分析

9、結(jié)果為：  2.根據(jù)出現(xiàn)的位置可以找到下面3個(gè)匹配(1,1,3,3),(5,8,2,2),(7,5,4,4).3.本例比較簡(jiǎn)單，錨定和最終的比對(duì)結(jié)果大致一致，根據(jù)第二步得到的匹配直接可以確定比對(duì)結(jié)果，序列A在位置4發(fā)生一個(gè)替換突變，位置5發(fā)生一個(gè)置換突變就得到序列B。長(zhǎng)序列比對(duì)的時(shí)候錨定間隔都有一定的長(zhǎng)度，使用SPA算法能很快將序列全部對(duì)齊。  3.結(jié)果31SCA算法分析 眾所周知，基因組序列一般都很長(zhǎng)，比如人類基因組總長(zhǎng)度大概3.0G，由于PC機(jī)的空間限制，使用哈希表和后綴樹(shù)往往會(huì)消耗太多主存從而降低運(yùn)算速度。對(duì)DNA序列預(yù)先使用YK碼編碼，

10、能夠大大降低內(nèi)存空間，并且得出DNA序列唯一的一個(gè)語(yǔ)法結(jié)構(gòu)，即使對(duì)于相似性較低的序列仍然能夠快速而精確的得到比對(duì)的錨定。Sequitur碼的編碼和譯碼過(guò)程的計(jì)算和空間復(fù)雜度都是線性的。SPA算法是一種綜合運(yùn)用概率統(tǒng)計(jì)和組合分析的新型算法，在保證比對(duì)精確度的損失很小的前提下，把兩兩比對(duì)的時(shí)間和空間的復(fù)雜度降為線性。SCA在完成錨定后，使用SPA算法并采用了模結(jié)構(gòu)大大簡(jiǎn)化了間隔比對(duì)過(guò)程。綜合上面的分析SCA是線性的全基因組比對(duì)算法。與以往的算法相比，SCA能夠更有效的處理同源性較差的基因組序列的比對(duì)。通常由于相似度的降低，MUMmer，AVID等算法的效率會(huì)降低，但是SCA算法不會(huì)受到影響，主要是

11、因?yàn)閷?duì)序列用YK算法對(duì)序列編碼后，相似片斷的位置就已經(jīng)標(biāo)出，不需要在花費(fèi)時(shí)間做大量搜索工作。把輸入序列中的一條倒轉(zhuǎn)以后再編碼，就能找出基因組倒轉(zhuǎn)重排的變換。圖1：兩株大腸桿菌序列的比對(duì)坐標(biāo)圖圖2：兩株大腸桿菌序列模擬重排后的比對(duì)坐標(biāo)圖32試驗(yàn)數(shù)據(jù)和結(jié)果 我們使用了大量基因組數(shù)據(jù)測(cè)試SCA算法，我們使用主頻3.0G主存1.0G的PC機(jī)測(cè)試SCA算法。SCA完成擬南芥的2號(hào)和4號(hào)染色體總長(zhǎng)度為39M的基因組比對(duì)只需要40秒。部分?jǐn)?shù)據(jù)和程序運(yùn)行時(shí)間見(jiàn)表1。表2是比對(duì)的“相似度”和“基因覆蓋率”（包括功能已知的基因和功能不確定的基因）。表1,表2中用NCBI的genbank序列編號(hào)來(lái)代替序列

12、名稱，其中NC_003098和NC_003028是兩條S.pneumoniae基因組，NC_004741和NC_004337是兩條Shigellaflexneri基因組，NC_002655和NC_002695是Escherichiacoli的兩條基因組，hum_chr和chimp_chr是人類和大猩猩的兩個(gè)染色體片斷。這些數(shù)據(jù)可以在NCBI的官方網(wǎng)站/下載，圖1，圖2是比對(duì)結(jié)果坐標(biāo)圖，其中縱坐標(biāo)軸代表U00096序列，橫坐標(biāo)代表BA000007序列。圖1兩條原始序列的比對(duì)圖，圖2是我們對(duì)序列進(jìn)行局部倒轉(zhuǎn)、易位重排后的比對(duì)圖，由圖2可以看出SC

13、A算法可以將序列重排找出。由表1，表2的數(shù)據(jù)可以看出：（1） SCA算法的速度大大高于BLASTZ算法，最高相差的達(dá)100倍之多。BLASTZ的速度受序列相似度的影響較大，而SCA算法的的時(shí)間復(fù)雜度是線性的，只跟輸入序列的長(zhǎng)度有關(guān)。（2） SCA算法得到的比對(duì)序列的相似度比BLASTZ稍低，而基因覆蓋率并沒(méi)有損失。序列1 序列2 長(zhǎng)度（M） 時(shí)間（秒）   SCA blastzHum_chr Chimp_chr 1.9，1.4 5.0 43NC_003098

14、0;NC_003028 2.0，2.1 6.0 612NC_004741 NC_004337 4.6,4,7 10 1033NC_002655 NC_002695 5.6,5.5 12 670 表1：四組基因組序列組成的測(cè)試集序列1 序列2 基因覆蓋率 相似度  SCA Blastaz SCA blastzHum_chr Chimp_chr 98 99 75

15、% 78%NC_003098 NC_003028 98 98 92% 93%NC_004741 NC_004337 90 90 80% 81%NC_002655 NC_002695 93 93 91% 91%表2：四組基因組序列比對(duì)的相似度33相關(guān)工作和展望 使用壓縮算法做生物序列的比對(duì)和分析是生物信息的一個(gè)新的研究方向，SCA算法經(jīng)過(guò)改進(jìn)擴(kuò)展可以用于多重序列比對(duì)，序列搜索和基因識(shí)別等生物信息領(lǐng)域。致謝：.參考文獻(xiàn)1 

16、;Needleman,S.B.andWunschC.D.Ageneralmethodapplicabletothesearchforsimilaritiesintheaminoacidsequenceoftwoproteins,J.Mol.Biol,1970;48(3):443453.2 Smith,T.F.andWatermanM.S.Identificationofcommonmolecularsubsequences,J.Mol.Biol,1981;147(1):195197.3 E.-H.YangandJ.C.Keiffer.Efficientuniversallo

17、sslessdatacompressionalgorithmsbasedonagreedysequentialgrammartransformPartone:Withoutcontextmodels.IEEETrans.Inform.Theory,2000;147(1):755777.4 ShenSY,YangJYaoAHwangPI.SuperPairwiseAlignment(SPA):AnEfficientApproachtoGlobalAlignmentforHomologousSequences,JComputBiol,2002;9(3):477486.5 GuangyueHuandShiYiShen.A