柳 楠,朱永琦,李勝華,崔曉宇
(山東建筑大學(xué) 計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院,山東 濟(jì)南 250101)
隨著二十世紀(jì)三大科學(xué)計(jì)劃之一的人類基因組計(jì)劃的實(shí)施,大量的生物學(xué)數(shù)據(jù)有待處理[1-4],如何利用計(jì)算機(jī)建模、仿真等技術(shù)去提取其中有用的數(shù)據(jù),進(jìn)而研究其中所蘊(yùn)含的生物學(xué)意義,對計(jì)算機(jī)科學(xué)技術(shù)來說是一項(xiàng)嚴(yán)峻的挑戰(zhàn)[5-6]。因此在二十世紀(jì)提出了一門新興交叉學(xué)科—計(jì)算生物學(xué)。計(jì)算生物學(xué)運(yùn)用數(shù)學(xué)、計(jì)算機(jī)和生物學(xué)相關(guān)理論解決生物學(xué)問題,已經(jīng)成為目前最活躍的研究領(lǐng)域之一[7-9]。
基因組片段填充問題[10-12]是計(jì)算生物學(xué)極其經(jīng)典的問題之一,其中含重復(fù)基因的基因組片段填充問題已經(jīng)被證明為NP-完全問題[13-14],如何優(yōu)化基因組片段填充近似算法是近些年來的討論熱點(diǎn)。依據(jù)基因樣本序列中是否含有重復(fù)基因,將該基因組填充問題分為含重復(fù)基因的基因組片段填充問題和無重復(fù)基因的基因組片段填充問題;或依據(jù)基因樣本序列不完整數(shù)量,將該基因組填充問題分為單面基因組片段填充和雙面基因組片段填充,其中一條序列完整,另一條序列缺失,稱為單面基因組序列,兩條基因序列均為不完整的,則為雙面基因組序列[15]。
該文重點(diǎn)討論單面重復(fù)基因組片段填充問題。Munoz和D. Sankoff等人[12-13]首次提出了基于最小重組距離(DCJ距離)的單面基因組填充方法,使用斷點(diǎn)圖設(shè)計(jì)了多項(xiàng)式時間算法,并證明了基于DCJ距離的單面基因組片段填充算法是多項(xiàng)式可解的。對于單面無重復(fù)基因組片段填充問題,H. Jiang等人提出了使用DCJ距離或斷點(diǎn)距離為度量的算法,并證明了其是多項(xiàng)式可解的[14];對于含重復(fù)基因的基因組片段填充問題,H. Jiang等人證明了其是NP-完全的,并提出了4/3-近似算法[14-16]。隨后N. Liu等人采用局部優(yōu)化和貪婪算法將該類問題近似度改善到1.25[17-18];J. Ma等人采用非盲局部搜索策略將該類問題近似度進(jìn)一步改善到1.2[19-20]。
在許多應(yīng)用中,基因組序列通常被定義為一系列連續(xù)的片段重疊群(contig)[21],其中任何一個contig都不能被破壞,缺失基因的插入只能在contig的兩端執(zhí)行。在此約束下,當(dāng)不存在重復(fù)基因時,單面基因組片段填充問題是多項(xiàng)式可解的;當(dāng)存在重復(fù)基因時,H. Jiang等人通過最大化公共鄰接證明了該類問題是NP-完全的,并提出了一個近似值為2的近似算法[22-23]和一個雙參數(shù)的FPT算法[22](k,公共鄰接數(shù),d,基因最大重復(fù)數(shù));L. Bulteau等人給出了一種基于最大鄰接數(shù)和最小斷點(diǎn)距離的k-Mer參數(shù)的FPT算法[24];Q. Feng等人通過構(gòu)造輔助圖和二次尋找最大匹配給出了2.57-近似算法[25]。
該文的主要工作有以下三個方面:系統(tǒng)歸納了基于contig的單面基因組片段填充問題的現(xiàn)有算法并通過實(shí)例實(shí)現(xiàn)了算法,有助于讀者對此類問題的進(jìn)一步了解;在技術(shù)應(yīng)用和時間復(fù)雜度等方面對現(xiàn)有算法做了對比,并分析了該些算法存在的一些弊端;分析接下來研究工作中面對的挑戰(zhàn)和可能的解決方案。
該文只關(guān)注基于contig的單面基因組片段填充算法,但其結(jié)果可以推廣到多染色體或環(huán)狀基因組。
首先,給出一些必要的定義。不失一般性,假設(shè)所有的基因和基因組都由無符號的字母和整數(shù)組成,給定一個集合Σ和一個基因序列S,使用c(S)表示基因序列S中所有符號的集合。如果Σ中的符號在基因序列S中出現(xiàn)且只出現(xiàn)一次,則稱S是Σ上的一個排列,否則稱為序列。對于Σ中的任意兩個符號x,y,如果基因序列S至少包含{xy,yx}中的任意一個子集,那么則稱x,y在S中鄰接,令P(S)為S中所有鄰接的集合。設(shè)A和B是Σ中的兩個基因序列,A={a1a2…an},B={b1b2…bm}。對于P(A)中的任意一個鄰接aiai+1和P(B)中的任意一個鄰接bjbj+1,如果aiai+1=bjbj+1(或aiai+1=bj+1bj),則稱aiai+1與bjbj+1構(gòu)成了公共鄰接,a(A,B)表示A和B的公共鄰接集合,同時稱(aiai+1,bjbj+1)為一個匹配對。如果P(A)和P(B)中不存在aiai+1=bjbj+1(或aiai+1=bj+1bj),則稱aiai+1相對于bjbj+1構(gòu)成了斷點(diǎn),bp(A,B)和bp(B,A)分別表示A和B的斷點(diǎn)集合,如圖1所示。
定義一個基因組序列是由一系列contig構(gòu)成的,且contig內(nèi)部不能插入缺失基因,即S=
下面具體給出One-Sided-SF-max問題的概念:
定義1:One-Sided-SF-max問題。
輸入:一個完整的基因組序列G和一個不完整的基因組序列S,其中S=
輸出:將X=c(G)-c(S)≠?插入S得到S',使得|a(S',G)|最大。
One-Sided-SF-max問題已經(jīng)被證明為NP-完全的[22],此類問題不能在有效時間內(nèi)求出精確解,因此設(shè)計(jì)近似算法更具有實(shí)際意義。本節(jié)主要對One-Sided-SF-max問題進(jìn)行簡要介紹,概括分析了國內(nèi)外經(jīng)典的三種算法:2-近似算法、2.57-近似算法以及k-Mer算法。
該算法由H. Jiang等人提出,主要使用了貪婪和最大匹配的思想來實(shí)現(xiàn)基于contig的單面基因組片段填充。首先在該算法中給出以下定義:對于基因序列S=
定義缺失基因集合X中有一個長度為n的子串,如果插入到slot <βi,αi+1>中(1≤i≤m-1),產(chǎn)生n+1個新公共鄰接,稱子串為n-Type-1類型串。同樣的,產(chǎn)生n個新公共鄰接,稱為n-Type-2類型串;產(chǎn)生n-1個新公共鄰接,稱為n-Type-3類型串。算法大體流程如下:
(1)計(jì)算缺失基因集合X=c(G)-c(S);
(2)對于缺失基因集合,采用貪婪策略將1-Type-1類型串插入到相應(yīng)的slot中,并將該slot鎖定,不允許其他缺失基因插入;
(3)構(gòu)造二分圖并求其最大匹配,將1-Type-2類型串插入到可構(gòu)成外部鄰接的slot處,并對該slot進(jìn)行更新:如果xj插入到slot ?ai前面,那么將此slot更新為?xj,如果?xj插入到slotβi?后面,那么將此slot更新為xj?;
(4)以步驟2后的缺失基因?yàn)轫旤c(diǎn)構(gòu)造多重圖:若x∈X,y∈X且xy為G中一個內(nèi)部鄰接,那么x和y之間添加一條邊,尋找最大匹配M。對于最大匹配M中的所有匹配對xy,如果x為步驟3中插入的元素,則將y插入到相應(yīng)slot處使得xy構(gòu)成鄰接;將其余匹配對xy任意插入到未鎖定slot中,且不能破壞現(xiàn)有鄰接;
(5)在不破壞現(xiàn)有鄰接關(guān)系的前提下,將所有剩余缺失基因任意插入到S中未鎖定的slot處;
(6)得到近似解S'。
對于該算法,下面通過一個實(shí)例(如圖2所示)來說明算法的執(zhí)行過程:
(2)搜尋1-Type-2類型串,找到2,d,g為1-Type-2類型串,并以1-Type-2類型串集合和未鎖定slot集合為頂點(diǎn),構(gòu)造二分圖,建立的二分圖BG1如圖3所示。
(3)以步驟2之后的剩余缺失基因集合為頂點(diǎn)構(gòu)造多重圖Q,建立的多重圖Q并求得最大匹配,將a插入到g之前,將7插入到d之后,將52插入到a之前。
(5)算法結(jié)束,得到填充后的基因序列S'=
<52ag4k1acbd21d727>。
以上可以看出,通過此算法可以得到8個新公共鄰接,同時此例的最優(yōu)解為S*=<2ag4k1acb1d7527d2>,有13個新公共鄰接。
所以,
所以有近似解:
該算法為一個2-近似算法,同時該算法的運(yùn)行時間主要由步驟2中計(jì)算O(n)個頂點(diǎn)的二分圖中的最大匹配以及步驟3中計(jì)算O(n)個頂點(diǎn)的多重圖中的最大匹配決定,兩者都需要O(n2.5)個時間,所以2-近似算法的時間復(fù)雜度為O(n2.5)。
Q. Feng提出的2.57-近似算法繼續(xù)考慮了冗余塊對填充過程存在的影響。該算法主要使用了最大匹配算法構(gòu)造簡單路徑來具體實(shí)現(xiàn)基于contig的單面基因組片段填充問題。
首先在該算法中給出以下定義:令F(S)為contigCi的首尾元素集合,F(xiàn)(S)=(α1,β1,…,αm,βm)。如果最大匹配M中有塊xy,xy在最大匹配M中出現(xiàn)的次數(shù)稱為xy的指示數(shù)。設(shè)xy與bp(G,S)中塊ab可構(gòu)成匹配對(xy,ab),若xy在M中出現(xiàn)次數(shù)大于ab在bp(G,S)中出現(xiàn)次數(shù),則稱塊xy為一個冗余塊。定義⊕為對稱差,A⊕B=(AB)∪(BA),令K'為此算法中的兩次最大匹配M1和M2的對稱差,則K'中每個連通分量必為簡單路徑或簡單循環(huán)。算法的大體流程如下:
(1)計(jì)算缺失基因集合X=c(G)-c(S)和斷點(diǎn)集合bp(G,S);
(2)基于缺失基因集合X、斷點(diǎn)集合bp(G,S)和S中每個contig首尾元素集合F(S)構(gòu)造一般圖Γ1,尋找最大匹配M1;
(3)對于M1中的任意塊xy,有以下三種情況:如果x和y分別為同一個slot的前后兩端,則合并此相鄰的兩個contig;如果x(或y)屬于F(S),將y(或x)插入相應(yīng)的slot中使得xy鄰接;如果x和y均不屬于F(S),則將其置于圖H'中頂點(diǎn)。依據(jù)以上更新基因序列為S1;
(4)基于G、S1和H',求得斷點(diǎn)集合bp(G,S1)和F(S1);
(5)更新圖H':刪除可與bp(G,S1)中斷點(diǎn)構(gòu)成匹配對的邊;
(6)基于缺失基因集合X、斷點(diǎn)集合bp(G,S1)和集合F(S1)構(gòu)造一般圖Γ2,尋找最大匹配M2;
(7)Δ=H'⊕M2;
(8)對于圖Δ中的任意路徑k=p1p2…pt-1pt,判斷其是否為簡單路徑:若為簡單路徑,插入到相應(yīng)slot中,反之刪除路徑k中任意一條邊得到新的路徑p1p2…pt-1pt,將路徑p1p2…pt-1pt插入到基因序列的最右側(cè);更新基因序列為S2;
(9)統(tǒng)一將c(G)-c(S2)插入到序列S2的最右側(cè);
(10)得到填充完成后的基因序列S'。
下面通過上述2-近似算法的同一個實(shí)例(見圖2)來說明算法的執(zhí)行過程:
(1)計(jì)算斷點(diǎn)集合bp(G,S)=<1a,ac,1d,d7,75,52,24,4g,ga,a2,2k,k7,7d,d2>,計(jì)算S中每個contig首尾元素集合F(S)=<4,1,c,b,1,1>;
(2)構(gòu)造圖Γ1:X∪F(S)中的所有元素被視為頂點(diǎn),對于其中任意兩個元素x,y,如果有x∈X,y∈X或y∈F(S)且存在一個斷點(diǎn)β使得與xy構(gòu)成一個匹配對,則在x與y之間添加一條邊;如果有x,y∈F(S),假設(shè)x在contigC1中且為C1中最后一個元素,y在contigC2中且為C2中第一個元素,C1和C2相鄰且存在一個斷點(diǎn)β使得與xy構(gòu)成一個匹配對,則在x與y之間添加一條邊。在圖Γ1中尋找最大匹配M1,如圖4所示;
(3)刪除M1中的冗余塊:對于M1中的塊xy,xy與斷點(diǎn)ω可構(gòu)成匹配對(xy,ω),如果xy在M1中的出現(xiàn)次數(shù)大于ω在bp(G,S)中的出現(xiàn)次數(shù),則稱塊xy為冗余塊并將其刪除;
(4)令{α1,α2,…,αr}為M1中刪除冗余塊后剩余塊的集合,H'為{α1,α2,…,αr}中的塊與bp(G,S)中的斷點(diǎn)構(gòu)成匹配對的集合,則此例中H'={(ac?ac),(d1?1d),(52?52),(g4?4g),(a2?a2),(d2?d2)};
(7)計(jì)算X'=c(G)-c(S1)得到缺失基因集合X'=<7,5,2,a,2,7,d,2>,計(jì)算S1中每個contig首尾元素集合F(S1)=
(8)計(jì)算新的斷點(diǎn)集合bp(G,S1)=
(9)使用構(gòu)造圖Γ1的同樣方法構(gòu)造圖Γ2;
(10)如果圖Γ2與圖H'中存在相同邊,則在Γ2中刪除此條邊,此例中,無此類邊;
(11)在圖Γ2中求得最大匹配M2,如圖5所示。
(12)令Δ=H'⊕M2,此例中Δ=(a2,ag,7d);
(15)將剩余缺失基因插入到序列S2的最右側(cè);
(16)得到填充完成后的基因序列S'=
該算法可以得到8個公共鄰接,同時其中一個最優(yōu)解為<2ag4k1acb1d7527d2>,有13個新公共鄰接。
該算法的近似性能比為2.57。
k-Mer算法從參數(shù)化復(fù)雜性的角度研究基因組填充問題,相較于H. Jiang等人提出的2-近似算法,主要有以下三個方面的不同:
(1)不再限制插入的基因集合為c(G)-c(S),插入集合可以包含比c(G)-c(S)更多或更少的基因集合;
(2)允許將要插入的字符串?dāng)?shù)量預(yù)先指定為輸入約束,t1為要插入的字符串?dāng)?shù)量的下限,t2為要插入的字符串?dāng)?shù)量的上限(t1≤t2);
(3)作為相似性度量依據(jù),不局限于最大化公共鄰接的數(shù)量,相反,對于一個預(yù)定的參數(shù)k,最大化公共k-mers的數(shù)目,k值越高,結(jié)果越準(zhǔn)確。
L.Bulteau等人對于此類問題給出以下定義:對于兩個基因序列G和S,G°S表示二者的串聯(lián)。存在一個正整數(shù)k,使得ak(G)={S[i,i+k]|i∈[n-k]}為序列G中k-mers的集合,則ak(G,S)=ak(G)∩ak(S)。設(shè)S[i]表示S中第i個元素,S[i,j]表示序列S中從位置i到j(luò)的基因元素。對于一個完整基因序列G和一個不完整基因序列S,令pk(S,G)=ak(G)ak(S)表示存在于G中但不存在于S中的k-mers集合,并將此類k-mers稱為潛在的公共k-mers。
定義2:k-Mer Scaffold Filling(k-Mer-SF)。
輸入:一個完整的基因組序列G和一個不完整的基因組序列S,其中S=
輸出:找到T'?T,t1≤|T'|且填充后的S'∈S+T',使得|ak(S',G)|最大。
該文給出一個k-Mer-SF實(shí)例,在此只舉例說明了k=2和k=3的填充情況(如圖6所示)。
H. Jiang等人提出的2-近似算法和Q. Feng等人提出的2.57-近似算法均為k-Mer-SF中t1=t2=|T|且k=2時的特殊情況,并可在多項(xiàng)式時間O(3·(+m)2)內(nèi)計(jì)算。k-Mer-SF相較于以上兩個近似算法,不再僅僅參考公共鄰接數(shù)目的多少,主要使用以下評估參數(shù):k,k-mers的長度;:=ak(S*,G)-ak(S,G),匹配后帶來的額外公共k-mers數(shù)目;d,一個基因在G中出現(xiàn)的最大次數(shù);m,S中重疊群的個數(shù);t2,要插入的字符串?dāng)?shù)量的上限;λ,T中字符串長度的上界。
k-Mer算法主要解決了基于動態(tài)規(guī)劃如何在2O()·nO()時間內(nèi)求得近似解并給出其參數(shù)為k+的FPT算法。首先使用著色法對k-mers進(jìn)行分類,β:T→[t2]表示潛在的公共k-mers,β:T→[t2]表示可能插入字符。圖著色后,使用動態(tài)規(guī)劃算法重建序列S,使得S中有個潛在的公共k-mers轉(zhuǎn)為已實(shí)現(xiàn)的k-mers,在動態(tài)規(guī)劃過程中逐步找到大小遞增的局部最優(yōu)解,從左到右依次將缺失基因插入到序列S中,并使用局部優(yōu)化策略避免一些缺失基因重復(fù)插入。
對于k-Mer-SF問題,首先對序列G,S及插入基因集合T做約簡操作。刪除T中多余基因元素,如果T中存在一個基因元素出現(xiàn)次數(shù)大于t2,將刪除一個元素;其次,對G中的鄰接關(guān)系分類,并只保留潛在公共鄰接,假設(shè)x為不出現(xiàn)在S和T中的基因,y為不出現(xiàn)在G中的基因,令P1°x°P2°…°Pq-1°x°Pq替換G且用Ci[1]°y°Ci[|Ci|]表示S。如果有一個潛在鄰接在G中出現(xiàn)次,則在G中刪除一個該鄰接,刪除后若鄰接滿足以下條件之一,則稱為可實(shí)現(xiàn)鄰接:
(1)b∈T且c∈T;
(2)存在一個contigCi使得b∈Ci[|Ci|]且c∈T;
(3)存在一個contigCi使得b∈Ci[|Ci|]和c∈Ci[1];
(4)存在一個contigCi使得b∈T且c∈Ci[1]。
建立鄰接圖H=(V,E):令T,G和S中的基因元素作為H中的頂點(diǎn),如果鄰接bc或鄰接cb均為可實(shí)現(xiàn)鄰接,則令兩個頂點(diǎn)b和c相鄰,求得最大匹配M。
令V(M)表示匹配的端點(diǎn),建立兩個二分圖H1和H2且頂點(diǎn)分別為B:=V(M)和C:=(VV(M))。在H1中,當(dāng)bc是一個可實(shí)現(xiàn)鄰接時,在b∈B和c∈C之間添加一條邊。在H2中,當(dāng)cb是一個可實(shí)現(xiàn)鄰接時,在b∈B和c∈C之間添加一條邊。如果H1中存在頂點(diǎn)b∈B且度數(shù)至少為2+m+1,則將鄰接bc從G中移除;如果H2中存在頂點(diǎn)b∈B且度數(shù)至少為2+m+1,則將鄰接bc從G中移除,其中c是b的任意鄰接。
完成約簡操作后,這些頂點(diǎn)的數(shù)量最多為|V(M)|·2·(2+m+1)。這給出了G中頂點(diǎn)數(shù)量的界限,從而給出了實(shí)例大小的界限且所有約簡規(guī)則可以在多項(xiàng)式時間內(nèi)執(zhí)行。
從更廣泛的角度來看,k-Mer-SF考慮了字符串的基本插入問題。事實(shí)上,可以將該算法擴(kuò)展到更一般的情況,即給定一個字符串G和一個部分字符串S,完成部分字符串S的插入得到新的字符串S',使得G和S'的相似度最優(yōu),因此其包含了H. Jiang等人的問題作為特例。
該文發(fā)現(xiàn)2-近似算法沒有考慮斷點(diǎn)對填充過程的影響和長度大于1的缺失串的插入情況,2.57-近似算法則沒有考慮n-Type-3類型串的插入情況,k-Mer算法也僅僅介紹了固定基因子串長度的一般處理情況,然而不同長度以及不同類型串的插入都會對公共鄰接數(shù)造成影響,從而影響到算法近似比。
序列中存在連續(xù)長缺失基因串,在完成1-Type-1類型串的插入后,假設(shè)存在可構(gòu)成n-Type-1類型串的缺失基因串a(chǎn)1a2…an,如若按照2-近似算法中分別處理,則該缺失基因串中1-Type-1類型串被破壞,有以下2種情況:
(1)最多有n個1-Type-2類型串,產(chǎn)生n個鄰接;
(2)最少有2個1-Type-2類型串,剩余為n-Type-3類型串,產(chǎn)生2個鄰接。
產(chǎn)生鄰接數(shù)k為2≤k≤n,若將其合并后插入基因序列中,則會產(chǎn)生n+1個鄰接。以上證明合并插入n-Type-1類型串具有更優(yōu)效果。
在處理n-Type-2類型串時,假設(shè)存在可構(gòu)成n-Type-2類型串的缺失基因串b1b2…bn,如若按照2-近似算法分別處理,則該缺失基因串中n-Type-2類型串被破壞,有以下2種情況:
(1)最多有2個1-Type-2類型串,同時bi(2≤i≤n-1)均可與其相鄰缺失基因bi-1或bi+1(2≤i≤n-1)構(gòu)成內(nèi)部鄰接,產(chǎn)生n個鄰接;
(2)最少有1個1-Type-2類型串,同時剩余n-1個缺失基因均為n-Type-3類型串,產(chǎn)生1個鄰接;產(chǎn)生鄰接數(shù)k為1≤k≤n,如若將其合并后插入基因序列中,會產(chǎn)生n個鄰接。以上證明合并插入n-Type-2類型串具有更優(yōu)效果。
對于2.57-近似算法,雖然n-Type-3類型串會產(chǎn)生n-1個鄰接,但是如果不對其進(jìn)行處理,任其插入未鎖定的slot,存在破壞已有鄰接的可能,其中最多會破壞2個鄰接,n-Type-3類型串產(chǎn)生鄰接數(shù)k為n-3≤k≤n-1,明顯降低最后的填充效率,所以對其進(jìn)行插入處理十分必要。
k-Mer算法固定了基因子串的長度,2-近似算法就是一個特例,其固定長度為2的公共鄰接數(shù)目作為算法性能參考依據(jù),以上也證明了此時并不是最優(yōu)算法,進(jìn)而說明了基因填充過程中限制基因子串長度存在影響近似性能比的可能。
重點(diǎn)介紹了基于contig的單面重復(fù)基因組片段填充問題的研究現(xiàn)狀。對該些算法在近似性能比等方面做出了詳細(xì)的對比(見表1),直觀地看出各個算法現(xiàn)存在的一定不足,說明此類算法仍有改進(jìn)空間,同時提出了該類問題的改進(jìn)思路。
表1 三種算法分析比較
(1)現(xiàn)有算法依賴于公共鄰接數(shù)目,鄰接的定義沒有考慮基因序列不存在逆序關(guān)系的情況,因此存在兩個基因序列的鄰接均為公共鄰接但二者完全不相似的情況,不利于后續(xù)對算法近似比的研究。
(2)現(xiàn)有算法對度量依據(jù)的選擇較少,并過度依賴于最大匹配算法,因此對此類問題的研究較為片面。
針對目前單面重復(fù)基因組片段填充問題的研究工作,發(fā)現(xiàn)基因組填充問題有以下發(fā)展前景。
(1)目前One-Sided-SF-max問題的最佳性能近似比為2,日后還需要進(jìn)一步優(yōu)化。
(2)現(xiàn)有算法均是在最大化鄰接基礎(chǔ)上考慮其近似比,基于最小斷點(diǎn)數(shù)層面有待研究。
(3)雙面基因組填充問題也被證為NP-完全的,但沒有提出近似比算法,此類算法可推廣到雙面基因組填充問題,有利于雙面此問題的近似比優(yōu)化。