ｍＲＮＡ可變剪接問題的并行化研究

摘要：對mRNA可變剪接問題進行了并行化分析和研究，解決了原AltSplice程序中Nbest參數選取欠合理的問題，結合實驗實現了基于MPI平臺的并行AltSplice版本MPI_AltSplice，并在5萬億次的聯想深騰6800高性能計算機上獲得了較好的運行性能。

關鍵詞：可變剪接； AltSplice； 并行化； 表達序列標簽(EST)序列庫； 消息傳送接口

中圖分類號：TP311.52文獻標志碼：A

文章編號：1001－3695(2008)03－0705－04

面對人類基因組計劃產生的海量分子生物學數據，生物信息學的重要性越來越突出，它無疑將會為生命科學的研究帶來革命性變革[1]。隨著測序技術的發展與成熟，對海量生物信息數據進行分析和處理已經演變為當前或未來生物信息學的中心任務，從而也促進了高性能計算技術在生物信息領域的應用與發展。人類的基因數目只有3萬多個，這甚至比一些植物還少。那么，如何理解人比其他生物高得多的進化程度？一個重要原因是高等動物的很多基因(如人的基因有一半以上)都涉及mRNA可變剪接。mRNA可變剪接使一個基因能表達出多個不同的蛋白質，這使得人在蛋白質組水平上表現出極高的復雜性。mRNA可變剪接涉及生命現象的很多方面，如細胞分化、個體發育、免疫機理、某些遺傳病的發生等[2]。例如，最新研究表明人類14號染色體上一個基因是neurexin 3，該基因的跨度達到170萬個堿基對，是14號染色體上最大的基因，其大小是所有基因平均的30倍。由于可變剪接，這個基因可能存在上千個不同的微小變化，該基因與腦中的突觸連接有重要作用。利用高性能計算技術進行這方面的研究有極其重要的科學和現實意義。

1可變剪接

基因轉錄時，先是產生與基因整體序列完全相同的RNA序列，稱為前體mRNA；接著前體mRNA中一些稱為內含子(introns)的片段被剪掉，剩下的稱為外顯子(exons)的片段被順序連接在一起；再經過對頭尾兩端的一定處理，就成為成熟mRNA，這個過程叫做剪接[3]。一般地，一個基因的外顯子和內含子的數目及所在位置是固定的，但是也可能出現外顯子和內含子數目及所在位置不固定的情況。這時，一個基因可以轉錄成多個不同的成熟mRNA，稱mRNA可變剪接，其直接后果是一個基因產生多個不同功能的蛋白質[4]。由于EST序列是成熟mRNA的一些鏡像片段，如果某個基因在轉錄時存在可變剪接，而且相應的EST測序是充分的，就能通過EST與EST、EST與基因組序列間的比較找出各種mRNA可變剪接形式[5]。現有的實驗證實，在mRNA剪接過程中會有四種規則的可變剪接形式，如圖1所示。

1．1生物信息學方法研究可變剪接的可行性

用生物信息學方法研究可變剪接的可行性源于高通量的基因組測序和EST測序。EST來源于完全加工的mRNA，它們提供了一個廣泛的mRNA多樣性的樣品庫。計算機為這種多樣性分析提供了更為便捷的工具。最近兩年，多個研究小組通過不同的生物信息學方法從整個人基因組的水平進行分析，結果一致顯示，約35%~60%的人基因有可變剪接形式。對大多數基因來說，每個基因只測了很少幾個EST甚至沒有EST；EST不是全長的mRNA，多位于mRNA的5’和3’端；EST來源于有限的組織和發育階段；很有可能存在更多的可變剪接，而現在的EST庫中沒有顯示。因此，實際可變剪接的頻率可能比預測的更高。將高性能計算技術應用于充分測序的大容量EST序列庫進行可變剪接分析，將有助于得到更多的可變剪接位點，從而對所表達的蛋白質以及相應物種的結構和功能分析提供便利。

1．2可變剪接的生物信息學分析方法

在生物信息學中，分析可變剪接的常用方法有兩種：EST序列間相互比較；EST序列和基因組序列間比較。通過序列間相互比較來分析兩條序列間相互匹配的部分。根據匹配部分的長度和位置來決定是不是存在可變剪接以及可變剪接的形式[6]。

如圖2所示，通過ID為AA054824的EST序列和ID為AV730044的EST序列的匹配情況，可認為這種類型屬于圖1中的optional exon類型的可變剪接形式。說明在剪接過程中存在剪接變體。

2AltSplice簡介

mRNA可變剪接分析軟件程序AltSplice是由中國科學院生物物理研究所自主開發的。它通過EST序列之間的相互比較，找出那些涉及mRNA可變剪接的EST序列。該軟件的核心算法是利用BLAST(basic local alignment search tool，基本局部比對搜索工具)相識性比對確定EST序列間或EST序列與基因組序列間相識性比對的長度和區域。BLAST由Altschul等人[7]于1990年提出，其基本要點是序列片段比對。所謂序列片段對是指兩個給定序列中的一對子序列，它們長度相等，且可以形成無空位或有空位插入的完全匹配。BLAST算法首先找出待檢序列和目標序列間所有匹配程度超過一定閾值的序列片段對；然后對具有一定長度的片段對根據給定相似性閾值延伸，得到一定長度的相似性片段，稱之為高分值片段對（highscoring pairs，HSPs）。AltSplice分析得到的高分值片段對，根據該相似性片段對長度和區域特點決定是否屬于已經由生物學實驗所證實的四種規則可變剪接形式，從而找到可變剪接位點并根據這種可變剪接的特點進一步對這種可變剪接形式進行分類和輸出。除四種規則的可變剪接形式外，本軟件還得到一些不規則的剪接形式。該形式的剪接是否屬于可變剪接尚未得到生物學實驗的證實，其結果只能作為一種參考。

AltSplice主要是對EST或基因組序列庫進行分析，查找可變剪接位點，進而輸出可變剪接的形式。先獲得要分析的EST序列庫（必須為FASTA格式），選擇進行分析的序列所在EST序列庫和進行可變剪接位點查找的序列個數及位置。例如要進行分析的三個EST序列庫為Hum_E74_EST1.cl和Hum_E74_EST2.cln、Hum_E74_EST3.cln。首先，從序列庫Hum_E74_EST1.cln中選擇前500條序列進行可變剪接分析；然后，遍歷EST序列庫得到總的EST序列數目N_all，格式化所有EST序列庫，并對每條序列作索引index[1]~index[N_all]；其次，將所要進行可變剪接分析的序列和所有EST序列庫進行BLAST比對分析，獲得BLAST比對結果文件；再根據各個序列比對的長度和位置判斷是否存在可變剪接以及該可變剪接的形式；最后，輸出各種形式的可變剪接結果。程序簡單的串行流程如圖3所示。

3AltSplice 的并行優化

在AltSplice的串行流程中對所要分析的EST序列庫均采用了連續遍歷的方法，這些地方均具有明顯方便并行處理的特征。如果對n個EST序列庫，假設庫中的序列數目依次為{a1，a2，…，an}，單處理器運行時的負載明顯為S=a1+a2+…+an。由于各個EST庫之間以及各條EST序列間并無任何相關性，設想將此n個EST庫運用多處理器進行并行處理，以此提高程序運行效率。

3．1AltSplice的并行優化算法

根據AltSplice流程中對每個EST序列庫進行串行處理的特點，可以有以下兩種并行優化方案：

a）對所要處理的EST序列庫不作任何修改和重新部署，簡單將每個EST序列庫分給不同處理器進行處理。這種方案的優點是省去了對所有EST序列庫進行重新部署的開銷；程序運行過程中的參數不作任何修改，并行算法設計比較簡單。但該方案的缺點也顯而易見：首先，用戶提供的各個EST序列庫的大小并不能保證均等，這無疑會在并行程序運行過程中出現負載失衡的情況。假設有三個EST序列庫，其序列數目分別為100、1 000、100 000 000或更多，分配的進程分別為P1、P2和P3，那么，P3和P1的負載會嚴重失衡。其次，如果用戶只提供一個較大的EST序列庫，就只能利用一個處理器來處理，并不能體現并行的效果，即NP=EST序列庫的個數（NP為并行運行時選取的處理器個數），缺乏必要的靈活性。

在對一條序列和庫中所有序列進行BLAST比對時，有一個重要參數Nbest。Nbest參數表示從BLAST比對結果中取Nbest個最匹配的序列作為比對結果輸出。若對上述三個EST序列庫進行BLAST比對，Nbest取相同的值，那么對這三個EST序列庫就會有不同的意義。假設Nbest=300，那么對第一個序列數目為100的EST序列庫，Nbest參數取大于100的任何值都已經沒有意義；但是對于第二、三個EST序列庫，取 Nbest=300也許并不能將所有滿足條件的比對結果全部取出來。因此對于不同大小的EST序列庫采用相同Nbest參數進行BLAST比對是不恰當的。不僅在這種并行方案中存在Nbest參數選取欠合理的問題，在原AltSplice軟件中同樣存在類似問題。

b)采取重新部署EST序列庫的方法以克服負載失衡和Nbest參數的問題。將進行可變剪接分析的所有EST或基因組序列庫看做一個整體，均勻分割為NP個新EST序列庫。具體分割方法為Snum % NP==Libnum。其中：Snum為某條序列在所有庫序列數總和N_all中的序號；NP為并行運行時的處理器個數；Libnum為分割出來的新EST序列庫序號。新序列庫中序列數目基本相同，序列個數最多相差為1，避免了負載失衡。

對于EST序列庫的大小不一所產生的Nbest參數選取問題，由于對原EST序列庫進行了大小均衡的重新部署，保證了Nbest參數的一致性。但是，對原始EST序列庫進行部署同樣也會影響到Nbest參數選取的合理性。如果對一個或多個原EST序列庫，采取順序分割部署為大小均等的新EST序列庫，并不能解決Nbest參數恰當選取的問題，如圖4所示。

對原EST序列庫OL1和OL2進行順序分割部署。圖4中淺灰色代表進行BLAST比對結果中匹配結果較好的序列。為敘述方便，以下稱這種序列為HighScore序列，即匹配分數較高的序列。若進行順序分割部署，得到的new library 1 EST 庫中全部為HighScore序列，而new library 2 EST序列庫中則無HighScore序列。如果取Nbest參數為N<3 000 000，則在對new library 1進行BLAST比對分析時必然會丟掉本該進行AltSplice分析的HighScore序列，直接導致丟掉更多的可變剪接位點。合理提高NP和Nbest參數時，HighScore序列仍然會集中到通過順序部署得到的前幾個新EST序列庫，可能丟掉很有生物學意義的可變剪接位點。

若采用上述的Snum % NP==Libnum方法，暫稱為定步長部署法，對原EST序列庫進行分割部署，同樣地OL1和OL2原EST序列庫會被分割部署為圖5所示的狀態，HighScore序列就會均勻分布到各個新EST序列庫。如果Nbest參數選取適當就不會丟掉具有生物學意義的剪接位點。當NP>2時會將HighScore序列分配得更加均勻，反而能獲得更多甚至是意想不到的剪接位點。

原EST序列庫中HighScore序列的分布本身是沒有什么規律可言的，但若用同樣的Snum % NP==Libnum分割部署方法對圖6所示的EST序列庫進行分割部署會同樣產生圖5所示的順序部署產生的結果，同樣會丟掉有意義的HighScore序列。但在適當改變NP值的情況下，取NP=NP+1，那么HighScore序列仍會均勻分布于新EST序列庫中。在兩個不同NP值的情況下，HighScore序列集中于同一個新EST序列庫的概率為0。不僅這種HighScore序列的分布規律出現概率極低，定步長部署法本身可以取不同步長，同樣可以解決HighScore序列集中的問題。這種對原EST序列庫的分割部署策略會對AltSplice運行結果產生很大影響，同樣也會影響到并行優化的性能。

3．2AltSplice的并行實現

AltSplice串行流程中對所分析的原EST序列庫采取了順序處理方式。在遍歷EST序列庫，得到總的EST序列個數N_all、對每條EST序列生成index索引、formatdb EST序列庫，以及對每條EST序列均需要與所有EST序列庫進行BLAST分析，明顯具有可并行處理的特點。利用優化方案b），首先對所有原EST序列庫進行均衡的重新分割部署；然后，對得到的NP個新EST序列庫進行處理器與EST序列庫的一一對應，每個處理器處理一個EST序列庫；最后，對得到的可變剪接形式進行匯總輸出得到總的可變剪接分析結果。

MPI有兩種最基本的并行程序設計模式，即對等模式和主從模式[8~10]。本文采取對等模式與主從模式結合的方法對AltSplice軟件進行并行程序的設計。在程序的初始和結尾處有一個進程為主進程負責程序運行參數的廣播和程序運行結果的收集。其具體并行程序的流程如圖7所示。

每個進程處理一個新部署的EST序列庫，并動態生成與此進程對應的EST序列庫名，如NL0、NL1，…，NLNP。進程個數與新部署的EST序列庫個數相等。每個進程都需要完成的任務如下：

a）生成自己要處理的新EST序列庫及其庫名。

b）格式化EST序列庫，即formatdb(NL pid)。

c）EST序列中對每條序列作索引，即makeIndex。

d）取得要進行BLAST比對的序列，進行BLAST比對。

e）分析比對結果，找到alternative sites。

主進程除完成上述任務外還需要進行程序初始時參數的廣播，以及對各進程進行index(needed by output)、sites等的收集，最后進行可變剪接形式的輸出。并行AltSplice整個流程中涉及到的進程間通信只有主進程對各個進程的參數廣播以及index和sites的收集，其他各個進程間并不涉及通信。并行優化性能與主進程執行收集時所進行通信的時間開銷有直接關系。若進行index和sites收集的時間分別表示為T(index，pid )和T(sites，pid)。其中pid為進程序號。各個進程執行對等模式的任務時間表示為PT(sametask，pid)，則只有當

∑NP-1pid=0T(index，pid)+T(sites，pid)≤

∑NP-1pid=0PT(sametasks，pid)（1）

才會有明顯并行優化性能的體現。具體實驗環境的并行機系統互聯網絡對進程間通信效率起著決定性的作用。

3．3并行AltSplice測試結果

試驗環境：聯想深騰6800超級計算機系統；系統結構COW；運算速度為4.183萬億次/s；265個四路節點機；內存總容量2．6 TB，磁盤總容量80 TB；高速連接網絡QsNet，點對點通信帶寬大于300 MB，延遲時間小于7 μs。

試驗數據：Hum_E74_EST1.cln.001 Seqs為249 940；Hum_E74_EST1.cln.002 Seqs為233 960。

上述實驗數據為人類的EMBL release 74 EST 去雜質后的FASTA格式EST序列庫，可變剪接分析序列為Hum_E74_EST1.cln.001中的 1～100。并行AltSplice、MPI_AltSplice的測試結果如表1所示。原串行AltSplice程序的執行結果如表2所示（試驗環境和原EST序列庫不變）。

3．4并行性能分析

上述測試結果的加速比如圖8所示。NP=2時出現超線性加速比，且NP=1時，并行AltSplice比原串行程序運行時間長。原因是并行AltSplice要對原 EST序列庫進行重新分割部署。當NP=1時，會將所有原EST序列庫合并成一個新EST序列庫，增加了時間開銷。隨著Nbest參數和NP的變化，并行AltSplice所得到的可變剪接位點個數也隨之變化。Nbest×NP越大，得到的可變剪接位點越多，避免丟掉具有生物學意義的可變剪接位點。當NP繼續增大，對應的新EST序列庫數目也繼續增大時，通信開銷的增加大于各個處理器單獨處理每個新EST序列庫所帶來的好處，并行效率開始下降。

并行AltSplice在對原EST序列庫進行大小均衡的重新部署基礎上，使各個處理器的負載達到了只相差一條EST序列的程度。

4結束語

分子生物學研究的重點已經由如何得到大量生物學數據轉變為如何解析這些海量數據，進而獲得這些數據背后隱含的實際生物學意義，將高性能計算技術應用于生物信息學是發展的必然趨勢。本文系統地分析了AltSplice軟件的并行優化策略及對原 EST序列庫的分割部署和Nbest參數合理選取問題；實現了AltSplice軟件的并行化，為大容量的EST序列庫進行可變剪接分析提供了幫助。通過尋找和注釋更多的可變剪接位點，可以逐步建立世界上首個專門的人類和其他物種可變剪接位點數據庫[11]，從而大大加深對人類或其他物種結構和功能的了解。本文的工作具有廣闊的應用前景。

參考文獻:

[1]MOUNT D W. Sequence and genome analysis[M]. New York: Cold Spring Harbor Laboratory Press， 2002:323－350.

[2]MODREK B， LEE C. A genomic view of alternative splicing[J]. Nature Genet， 2002，30（1）:1319.

[3]ATTWOOD K，PARRYSMITH D J.生物信息學概論[M]. 羅靜初，等譯.北京:北京大學出版社，2002:11197.

[4]李衍達，孫之榮.生物信息學基因和蛋白質分析的實用指南[M]. 北京:清華大學出版社，2000:103－211.

[5]XIE Hanqing，ZHU Weiyong，WASSERMAN A， et al. Computational analysis of alternative splicing using EST tissue information[J]. Genomics， 2002，80(3):326－330.

[6]BLACK D L. Mechanisms of alternative premessenger RNA splicing[J]. Annu Rev Biochem， 2003，72（2）:291－336.

[7]ALTSCHUL S F，GISH W，MILLER W，et al. Basic local alignment search tool[J]. J Mol Biol， 1990，215(3):403－410.

[8]莫則堯，袁國興.消息傳遞并行編程環境MPI[M].北京:科學出版社，2001:111.

[9]陳國良.并行計算——結構·算法·編程[M].北京:高等教育出版社，1999:83179.

[10]都志輝.高性能計算并行編程技術——MPI并行程序設計[M].北京:清華大學出版社，2001:22－98.

[11]JI Hongkai， ZHOU Qing， WEN Fang， et al. AsMamDB: an alternative splice database of mammals[J]. Ucleic Acids Res， 2001，29(1):260－263.

“本文中所涉及到的圖表、注解、公式等內容請以PDF格式閱讀原文”