陳鳳珍，李　玲, 操利超, 嚴(yán)志祥

( 深圳華大基因研究院, 深圳 518083)

?

四種常用的生物序列比對(duì)軟件比較

陳鳳珍，李玲, 操利超, 嚴(yán)志祥*

( 深圳華大基因研究院, 深圳 518083)

摘要：隨著高通量測(cè)序技術(shù)的快速發(fā)展，下一代測(cè)序技術(shù)也迅速發(fā)展為生物領(lǐng)域中的主流技術(shù)，而理解下一代測(cè)序數(shù)據(jù)最重要的一步是比對(duì)。比對(duì)是進(jìn)行后續(xù)生物信息分析的基石，也因此催生了很多比對(duì)軟件。本文主要選取了四種常用的比對(duì)軟件Bowtie2、BWA、MAQ和SOAP2，對(duì)這四種軟件及算法進(jìn)行綜述，并通過(guò)實(shí)際測(cè)序數(shù)據(jù)對(duì)四種軟件進(jìn)行比較和評(píng)估，為生物學(xué)研究者選擇最佳的短序列比對(duì)軟件提供理論和實(shí)踐依據(jù)。

關(guān)鍵詞：下一代高通量測(cè)序; 比對(duì)軟件; 生物信息

1引言

隨著新一代測(cè)序(Next-generation sequencing，NGS)的蓬勃發(fā)展，核酸測(cè)序成本已大大降低，高通量測(cè)序方法已被廣泛應(yīng)用到DNA測(cè)序[1]、RNA測(cè)序[2]、表觀遺傳測(cè)序[3-4]等研究。然而，無(wú)論使用何種生物測(cè)序技術(shù)和研究方法，理解這些數(shù)據(jù)的最重要的一步是序列比對(duì)分析。序列比對(duì)是將已有基因組序列作為參考基因序列(Reference)，將短序列與參考基因序列進(jìn)行序列比對(duì), 并在參考基因序列上進(jìn)行精確定位。通過(guò)序列比對(duì)可以發(fā)現(xiàn)生物序列中的功能、結(jié)構(gòu)和進(jìn)化的信息。目前已有上百種序列比對(duì)工具，面對(duì)如此多的比對(duì)工具，很多生物信息分析人員通常自由的選擇比對(duì)工具，而沒(méi)有考慮到比對(duì)工具的特點(diǎn)，準(zhǔn)確性等。然而，不同的比對(duì)軟件，對(duì)同一個(gè)數(shù)據(jù)集都有可能得出大相徑庭的結(jié)果[5]；同一算法設(shè)置不同的參數(shù)，其結(jié)果也相差很巨大。如果選擇了一個(gè)不合適的工具，將導(dǎo)致結(jié)果偏差甚至是錯(cuò)誤，可能得到錯(cuò)誤的研究結(jié)論。因而選擇合適的比對(duì)工具，對(duì)于生物研究而言顯得特別重要。

在Nuno A. Fonseca等人[6]的對(duì)60多種比對(duì)軟件比較統(tǒng)計(jì)分析中，發(fā)現(xiàn)Bowtie2[7]、BWA[8]、MAQ[9]和SOAP2[10]被引用的次數(shù)相對(duì)其他幾十種軟件較多，其中Bowtie2引用率(Citations/Years)為363.42，BWA為224.20，MAQ為251.66，而SOAP2為99.38，SOAP2的前版本SOAP為104.41。因而在本研究中，主要選取了這四種常見(jiàn)的比對(duì)工具進(jìn)行評(píng)估比較。根據(jù)比較結(jié)果分析，Bowtie2、BWA和SOAP2處理高通量短序列數(shù)據(jù)比對(duì)問(wèn)題時(shí)，計(jì)算速度快，內(nèi)存使用量低，具有高效的實(shí)用性；在同等條件下，MAQ的運(yùn)行速度較慢。 Bowtie2、BWA的比對(duì)率相比于SOAP2和MAQ高。BWA軟件與Bowtie2軟件比對(duì)的重復(fù)率較高，MAQ較低。

2四種比對(duì)軟件及算法

2.1四種比對(duì)軟件介紹

Bowtie2是一個(gè)超高速的，節(jié)約內(nèi)存且靈活與成熟的短序列比對(duì)軟件，比較適合下一代測(cè)序技術(shù)。通常使用全文分索引(FM-index)以及Burrows-Wheeler 變換(BWT)索引基因組使得比對(duì)非?？焖偾覂?nèi)存高效，但是這種方法不適合于找到較長(zhǎng)的、帶缺口的序列比對(duì)。

BWA主要應(yīng)用二代測(cè)序后的大量短小片段與參考基因組之間的定位比對(duì)。需要先對(duì)參考序列建建立索引，BWA也是基于 BWT和 FM-Index 理論來(lái)對(duì)參考基因組做索引。根據(jù)測(cè)序方法的不同，有單末端序列(Single-end，SE)比對(duì)和雙末端序列(Pair-end，PE)比對(duì)。

MAQ是使用質(zhì)量分?jǐn)?shù)推導(dǎo)序列和比對(duì)序列的一致性的短序列比對(duì)工具，并且MAQ充分利用配對(duì)信息，估計(jì)每個(gè)比對(duì)read的錯(cuò)誤的概率，同時(shí)也使用貝葉斯統(tǒng)計(jì)模型來(lái)評(píng)估最后的基因型錯(cuò)誤概率。

SOAP2是短寡核苷酸比對(duì)程序(Short Oligonucleotide Alignment Program)的一個(gè)顯著改進(jìn)版本，它減少了計(jì)算機(jī)內(nèi)存使用，并極大地提高了比對(duì)速度。SOAP2使用一個(gè)Burrows Wheeler Transformation(BWT)壓縮索引替代種子策略在主存儲(chǔ)器中索引參考序列。SOAP2適合于單末端片段和雙末端片段。此外，該工具也支持多種文本和壓縮文件格式。

2.2四種比對(duì)軟件算法

對(duì)于成千上萬(wàn)條的短序列的比對(duì)分析，目前，大多數(shù)算法是通過(guò)建立索引來(lái)加快比對(duì)的速度。常用的數(shù)據(jù)結(jié)構(gòu)有哈希表法和基于BWT (Burrows-wheeler transform) 的后綴樹(shù)兩種。

哈希表法的算法核心思想是采用種子序列定位及延伸算法(Seed-and-extend algorithm)[11]，通過(guò)掃描參考基因組序列，對(duì)參考基因組序列建立哈希表，將序列分成一定長(zhǎng)度的小片段，這種小片段也被稱之為種子。然后，在目標(biāo)序列中查找和種子序列相同的片段并標(biāo)記，以這些標(biāo)記點(diǎn)為錨點(diǎn)向左右按一定規(guī)律延伸比對(duì)，將不合條件的舍棄，符合條件的結(jié)果將輸出保存。采用基于哈希表數(shù)據(jù)結(jié)構(gòu)的比對(duì)算法的軟件包括MAQ。

羅樓組以泥巖為主，兼有火山碎屑巖，可分下段和上段。下段主要由泥巖組成。底部夾少量泥質(zhì)灰?guī)r，下部夾有一層厚約1 m的玻屑凝灰?guī)r，中、上部夾凝灰質(zhì)泥巖，泥巖呈灰-深灰色，薄層狀，主要由水云母(95%～98%)組成，含少量石英(1%～5%)和微量鐵質(zhì)和有機(jī)質(zhì)；下部泥巖時(shí)見(jiàn)生物潛穴遺跡。

后綴樹(shù)法是一種n叉樹(shù)，n為字母表大小。每個(gè)節(jié)點(diǎn)表示從根節(jié)點(diǎn)到此節(jié)點(diǎn)所經(jīng)過(guò)的所有字符組成的字符串，它的根結(jié)點(diǎn)不包含任何信息，是一種以犧牲存儲(chǔ)空間來(lái)降低序列查詢時(shí)間的字符串預(yù)處理方式。為了提高空間利用率，F(xiàn)erragina 和Manzini 提出了FM (Full-text minute-space)-index 算法，F(xiàn)M是一種基于BWT (Burrows-wheeler transform)的全文本壓縮索引結(jié)構(gòu)， BWT算法是通過(guò)統(tǒng)計(jì)基因組序列中各個(gè)堿基出現(xiàn)次數(shù)，將相同堿基盡量排列在一起，壓縮基因組序列的索引數(shù)據(jù)結(jié)構(gòu)，將基因組序列的索引數(shù)據(jù)結(jié)構(gòu)重排列，實(shí)現(xiàn)短序列在基因組中候選位點(diǎn)的快速搜索，減少內(nèi)存占用率。例如人類基因組約3GB,若不使用FM-index將要用12GB內(nèi)存存儲(chǔ)，超過(guò)了計(jì)算機(jī)內(nèi)存使用限度，而如果使用FM-index， 每隔數(shù)行建立一個(gè)索引，人類基因組占用的內(nèi)存可縮小到約1.3GB，這樣普通的計(jì)算機(jī)就可以進(jìn)行分析。采用BWT轉(zhuǎn)換的軟件有Bowtie2和SOAP2，BWA。

雖然Bowtie2、SOAP2和BWA都采用了BWT算法，然而三種軟件還有差別。其中Bowtie2采用Ferragina 和Manzini 提出的FM (Full-text minute-space)-index算法，為基因組序列創(chuàng)建具有后綴矩陣特性的 FM 索引數(shù)據(jù)結(jié)構(gòu)，實(shí)現(xiàn)短序列的快速搜索；SOAP2則采用的是 BWT 算法壓縮基因組序列哈希表索引數(shù)據(jù)結(jié)構(gòu)進(jìn)行精確匹配，采用“分割短序列策略”(Split-read strategy)進(jìn)行不精確匹配，比對(duì)速度顯著提高且內(nèi)存使用量顯著地降低。最后，BWA 軟件是采用 BWT 算法壓縮來(lái)構(gòu)建基因組序列前綴樹(shù)(Prefix tree)數(shù)據(jù)結(jié)構(gòu)，通過(guò)對(duì)壓縮數(shù)據(jù)結(jié)構(gòu)自頂向下遍歷進(jìn)行反向搜索，其比對(duì)計(jì)算過(guò)程中內(nèi)存覆蓋區(qū)域相對(duì)較小，計(jì)算時(shí)間并不隨基因組的大小而變化。

基于哈希表法和基于BWT的后綴數(shù)法數(shù)據(jù)結(jié)構(gòu)的算法都有利于提高比對(duì)效率，區(qū)別在于哈希表法占用的內(nèi)存空間大，產(chǎn)生的種子匹配多，然而哈希表法具有較高的匹配敏感性和準(zhǔn)確性。有利于發(fā)現(xiàn)SNPs和突變?？捎糜诰植科ヅ浠驈拇罅繑?shù)據(jù)中搜索匹配點(diǎn)以及跨物種序列間的比對(duì)。而后綴樹(shù)法可以有效減少不精確匹配，并可避免比對(duì)過(guò)程中做無(wú)用功，這個(gè)特點(diǎn)適用于相同物種之間相似性高的序列比對(duì)和尋找保守區(qū)。

2.3四種比對(duì)軟件比較

選擇合適的軟件要根據(jù)軟件適用的數(shù)據(jù)類型，適宜測(cè)序平臺(tái)，數(shù)據(jù)格式，適宜的reads長(zhǎng)度等進(jìn)行全面考慮，做出選擇。表1中對(duì)四種比對(duì)軟件分析的序列類型，可用于分析的測(cè)序平臺(tái)，輸入和輸出數(shù)據(jù)格式，最小和最大reads長(zhǎng)度及軟件是否開(kāi)源進(jìn)行了詳細(xì)的分析和比較。從表中可以看出在適宜測(cè)序平臺(tái)方面，SOAP2就受到限制，只適用于Illumina平臺(tái)，BWA適用的平臺(tái)最廣。在適宜的reads長(zhǎng)度方面，BWA、MAQ適用的范圍較窄。最后，根據(jù)軟件的輸入輸出格式，MAQ的適用范圍更廣。

表1　四種比對(duì)軟件比較

3軟件評(píng)估實(shí)驗(yàn)

3.1實(shí)驗(yàn)數(shù)據(jù)

3.2軟件運(yùn)行環(huán)境

32G內(nèi)存，16核處理器，linux操作系統(tǒng)服務(wù)器。

3.3結(jié)果評(píng)估

四種軟件的比對(duì)率和時(shí)間消耗如表3。從表3可以看出BWA和Bowtie2的比對(duì)率較高，而SOAP2的時(shí)間更高效，MAQ相對(duì)來(lái)說(shuō)較慢。

表2　實(shí)驗(yàn)數(shù)據(jù)

表3　四種比對(duì)軟件比對(duì)率和比對(duì)時(shí)間

從四種軟件比對(duì)的reads重復(fù)數(shù)兩兩比較可以看出，Bowtie2和BWA比對(duì)上的reads重復(fù)數(shù)較高，Maq和其他三種軟件比對(duì)上的reads重復(fù)數(shù)較低，如圖1。將四種軟件同時(shí)比較時(shí)，發(fā)現(xiàn)BWA比對(duì)軟件和其他三種軟件不重復(fù)的reads數(shù)最少，只有62 134條，Bowtie2和其他三種軟件不重復(fù)的reads數(shù)最多，為466 792條，如圖2。

從實(shí)驗(yàn)結(jié)果看出Bowtie2和BWA的比對(duì)率相比于SOAP2和MAQ高。BWA軟件與Bowtie2軟件比對(duì)的重復(fù)率較高，MAQ較低，可能與選取的實(shí)驗(yàn)數(shù)據(jù)相關(guān)，本實(shí)驗(yàn)選取的是高覆蓋度的HPV全基因組測(cè)序數(shù)據(jù)，BWA比對(duì)工具比較適合全基因組測(cè)序數(shù)據(jù)的比對(duì)分析。

4討論

通過(guò)比較和實(shí)驗(yàn)研究發(fā)現(xiàn)，Bowtie2、BWA、MAQ和SOAP2四種軟件在處理高通量短序列數(shù)據(jù)比對(duì)問(wèn)題時(shí)，計(jì)算速度較快，內(nèi)存使用量較低，具有高效的實(shí)用性。 但是，這四種常用的分析軟件都只對(duì)短序列分析較為適合，然而，第三代測(cè)序技術(shù)正在快速的發(fā)展，必將成為未來(lái)的主流技術(shù)。第三代測(cè)序技術(shù)相比于第二代測(cè)序技術(shù)特點(diǎn)之一是讀長(zhǎng)長(zhǎng)。因而開(kāi)發(fā)高準(zhǔn)確性的適合第三代測(cè)序數(shù)據(jù)的長(zhǎng)序列比對(duì)工具是未來(lái)研究的主題。

對(duì)于比對(duì)分析一個(gè)常見(jiàn)的問(wèn)題是，哪一個(gè)分析工具是本研究最適合的。一個(gè)最好最適合的比對(duì)工具不光要考慮數(shù)據(jù)的類型，一個(gè)重要的方面包含比對(duì)工具是否和比對(duì)下游的分析和分析工具結(jié)合緊密，更包含比對(duì)的工具的速度和準(zhǔn)確性。但是目前，評(píng)估一個(gè)比對(duì)工具的準(zhǔn)確性和速度仍然很難，主要的困難是缺乏不同測(cè)序技術(shù)和研究方法的金標(biāo)準(zhǔn)數(shù)據(jù)集，因?yàn)椴煌谋葘?duì)軟件，不同的數(shù)據(jù)集，數(shù)據(jù)類型，數(shù)據(jù)大小等都有可能導(dǎo)致比對(duì)準(zhǔn)確度和時(shí)間偏差。因而創(chuàng)建適合的金標(biāo)準(zhǔn)數(shù)據(jù)集對(duì)于比對(duì)工具的評(píng)估和研究特別重要。

圖1　四種比對(duì)軟件比對(duì)上的讀長(zhǎng)兩兩比較Fig. 1　The multiple compration of mapped reads using the selected software

圖2　四種軟件比對(duì)上的讀長(zhǎng)比較Fig. 2　The compration of mapped reads using the selected software

5結(jié)論

對(duì)二代測(cè)序的四種常用比對(duì)軟件的算法進(jìn)行了總結(jié)，并對(duì)四種軟件的適用性和性能等方面進(jìn)行了對(duì)比，同時(shí)利用實(shí)際的基因組數(shù)據(jù)進(jìn)行測(cè)試分析，歸納總結(jié)，給出軟件選擇的參考建議，為研究人員選擇適合的比對(duì)分析工具提供參考。

參考文獻(xiàn)

[1]MARDIS E R. Next-generation DNA sequencing methods[J]. Annual Review of Genomics and Human Genetics,2008, 9: 387-402.

[2]WANG ZHONG, GERSTEIN M, SNYDER M. RNA-Seq: a revolutionary tool for transcriptomics[J]. Nature Reviews Genetics, 2009, 10:57-63.

[3]PARK P J. ChIP-seq: advantages and challenges of a maturing technology[J]. Nature Reviews Genetics,2009, 10(10): 669-680.

[4]MEISSNER A, MIKKELSEN T S, GU H, et al. Genome-scale DNA methylation maps of pluripotent and differentiated cells[J]. Nature,2008, 454(7205):766-770.

[5]NEKRUTENKO A, TAYLOR J. Next-generation sequencing data interpretation: enhancing reproducibility and accessibility[J]. Nature Reviews Genetics, 2012, 13(9):667-672.

[6]FONSECA N A, RUNG J, BRAZMA A, et al. Tools for mapping high-throughput sequencing data[J]. Bioinformatics,2012, 28(24):3169-3177.

[7]LANGMEAD B, SALZBERG S L. Fast gapped-read alignment with Bowtie 2[J]. Nature Methods,2012, 9(4):357-359.

[8]LI HENG, DURBIN R. Fast and accurate short read alignment with Burrows-Wheeler transform[J]. Bioinformatics,2009, 25(14):1754-1760.

[9]LI HENG, RUAN JUE, DURBIN R. Mapping short DNA sequencing reads and calling variants using mapping quality scores[J]. Genome Research,2008, 18(11):1851-1858.

[10]LI Ruiqiang, YU Chang, LI Yingrui, et al. SOAP2: an improved ultrafast tool for short read alignment[J]. Bioinformatics,2009, 25(15):1966-1967.

[11]LI Heng, HOMER N. A survey of sequence alignment algorithms for next-generation sequencing[J]. Briefings in Bioinformatics,2010, 11(5):473-483.

Comparison of four common biological sequence alignment tools

CHEN Fengzhen, LI Ling, CAO Lichao,YAN Zhixiang*

(BGI-Shenzhen,Shenzhen518083,China)

Abstract：With the rapid development of high-throughput sequencing technology,Next-generation sequencing technology has rapidly developed into a mainstream technology in the biological field. Alignment is the key step in understanding the sequence data and also it is the cornerstone for bioinformatics analysis. And thus gave birth to a lot of alignment tools. In this paper, four common biological sequence alignment tools Bowtie2, BWA, MAQ and SOAP2 were selected to evaluate and compare using the whole genome sequencing data of HPV. And a comparison of four tools from many perspectives such as algorithm and suitable sequencing platforms was given. Hopefully the research can provide theoretical and practical basis for researchers to select the best biological sequence alignment tools.

Keywords：Next generation sequencing; Alignment tools; Bioinformatics

中圖分類號(hào)：Q-31

文獻(xiàn)標(biāo)志碼：A

文章編號(hào)：1672-5565(2016)01-056-05

doi:10.3969/j.issn.1672-5565.2016.01.10

作者簡(jiǎn)介：陳鳳珍，女，生物信息工程師；E-mail:chenfengzhen@genomics.cn.*通信作者：嚴(yán)志祥，男，博士，生物信息高級(jí)工程師；E-mail:yanzhixiang@genomics.cn.

基金項(xiàng)目：國(guó)家自然科學(xué)基金資助項(xiàng)目(U1301252)。

收稿日期：2016-01-19；修回日期：2016-03-08.