基于块排序索引的生物序列局部比对查询技术

1、基于块排序索引的生物序列局部比对查询技术*block sorting index-based techniques for local alignment searches on biological sequences李永光 王 镝 王国仁 马宜菲 (东北大学信息科学与工程学院计算机系统研究所 沈阳 110004)abstract a common query against large protein and gene sequence data sets is to locate targets that are similar to an input query sequence. t

2、he current set popular search tools, such as blast, employs heuristics to improve the speed of such searches. however, such heuristics can sometimes miss targets, which in many cases is undesirable. the alternative to blast is to use an accurate algorithm, such as smith-waterman(s-w) algorithm. howe

3、ver, these accurate algorithms are computationally very expensive. recently, a new technique, oasis, has been proposed to improve the efficiency and accuracy by employing dynamical programming during traversing suffix tree and its speed is comparable to blast. but its main drawback is too much memor

4、y consuming. we propose an efficient and accurate algorithm for locally aligning genome sequences. we construct a block sorting index structure for the large sequence. the index structure is less than the suffix tree index and can be fit for large data size. experimental results show that our algori

5、thm has better performance than oasis.keywords sequence; block sorting index; accuracy1. 引言随着人类基因组计划的不断发展，在结构基因学和功能基因学的研究过程中，生物序列的相似性分析成为一种有效的手段1。通常情况下，两条dna长序列，可能只在很小的区域内(密码区)存在关系；不同家族的蛋白质往往具有功能和结构上的相同的一些区域。因而，研究序列局部相似性比研究全序列相似性往往更有意义。序列局部比对是一种关于片段相似性的定性描述1。通常的方法是通过动态规划2进行精确查找两个序列的最优局部比对，但因其代价太大，对超

6、长生物序列直接应用这种技术是不可行的。为了提高查询速度，启发式算法当前被广泛应用，这些算法可以分为三类：基于哈希表的算法、基于频率空间过滤的算法、基于后缀树索引的算法。基于哈希表的典型算法有fasta5和blast6。fasta的核心思想是对数据库序列中的所有模式进行哈希索引。在查询时基于哈希索引可以快速检索出可能的匹配模式。然后根据这些模式的扩展得到更长的序列。blast算法是建立在严格的统计学的基础之上，它主要用于发现具有较高的相似性的局部比对。在大多数情况下，根据局部比对参数会产生若干个hsp(high-score sequence pairs)。然后把所有匹配分值大于阈值的hsp作为结

7、果。这种基于启发式算法的过滤原理尽管在一定程度上尽量减少丢失局部比对的结果，但精确率却不可能达到100%，一些符合条件的匹配序列可能不在查询结果中。mrs4是基于频率空间过滤的代表方法，它采用滑动窗口和小波变换技术，设计了一种高效的生物序列搜索方法。根据频率距离是编辑距离的上界的过滤原理，可以节省大约50%的编辑距离计算。这种技术比较适合长的dna序列，但不适合蛋白质序列的比对分析。另外一些方法是基于后缀树技术的, 该类方法一般是通过对一条序列构建后缀树，然后利用后缀树快速查找到匹配的模式。对于后缀树中的某个分支，如果完全匹配的模式的个数超过给定阈值，那么就对该分支使用blast作进一步的查询

8、。这类算法对阈值比较低的查询处理效率太低。oasis3通过在遍历后缀树的同时运用动态规划算法，提供了一种精确的序列局部比对的搜索方法。其实验结果表明，它比动态规划算法快很多，而且可以和blast相提并论，但存储后缀树需要很大的存储空间，其构造过程也是复杂的。本文首先提出了一种基于块排序索引算法来计算序列比对，我们引入基于这种索引结构的片断向量的概念，通过片断向量的扩展完成动态规划算法，从而完成序列比对查找过程。算法采用a*8搜索策略，进行过滤，当计算查询序列部分匹配时，搜索算法同时计算匹配剩余查询所得分值的上界，根据阈值，如果上界分值小于阈值，那么该片段向量就被过滤掉，从而减少搜索空间的计算。

9、在搜索的每一步，算法首先扩展可能产生最优局部比对的片断向量，这就能保证返回的结果按照比对的分值由达到小排序。实验结果表明，我们的算法在性能和过滤能力上都优于oasis。本文的主要贡献在于通过运用较小的索引结构，我们提出一种精确有效的局部比对搜索算法，因而这种算法为在大的生物数据集提供了一种可行有效的精确局部比对查找策略。本文余下部分的组织结构如下：第2部分给出了本文的背景知识；第3部分具体介绍了核心数据结构和算法；第4部分给出了测试结果和性能分析；第5部分总结全文。2. 背景知识动态规划算法能精确计算生物序列全局和局部比对，在这一部分，我们首先介绍生物序列的局部比对的概念，然后讨论动态规划算法

10、的相关工作，最后介绍块排序索引结构。2. 1 序列局部比对给定两条字符序列q=q1q2qm, t=t1t2tn,这两条序列的局部比对是指通过某种方式是“排列”q和t的子序列。图1给出一个局部比对例子，同时给出三种类型的局部比对操作。l 替换：用相同或者不同字符进行替换（如图1中的1和2标记）。l 删除：允许在目标序列中忽略一个字符（如图1中的3标记）。l 插入：允许在查询序列中忽略一个字符（如图1中的4标记）。31124target：c a b i n 1query： d r a i n 图1. 局部比对的样例序列是通过三种操作所得到的分值之和来确定比对得分。每种操作可以形象描述为替换操作“”

11、，其中插入操作表示为“ ”,删除操作表示为“ ”，操作的罚分可以标记为s（）。可以通过替换矩阵s存储这些分值，所以“”的操作代价可以简单表示为s,。表1给出一个通用的替换矩阵的样例。称为单位编辑距离矩阵（因为完全匹配分值为1，否则为-1）。acgta1-1-1-1-1c-11-1-1-1g-1-11-1-1t-1-1-11-1-1-1-1-1-表1：单位编辑距离矩阵2. 2 smith-waterman(s-w)算法smith-waterman算法通过动态规划2的方法计算两条序列局部比对的最大的可能分值，计算过程需要o(mn)空间和时间的复杂度。算法形成一个mn的矩阵g，gi,j保存查询序列和

12、目标序列分别结束在qi和tj的最大比对分值。每一个gi,j是通过下面的公式(1)计算出来的： （1）2. 3块排序索引结构对于长度为n的序列s（我们的应用中，s以非序列中出现的字符作为结束字符，我们选用结束字符$，并规定$的字典序小于序列中的所有其它字符），进行块排序时对序列s循环移动0n-1位，得到一个n行列表，将各行按升序排序，就得到了矩阵m，如图2所示，左侧是对序列agtacgcctag$的块排序结果, 这里我们不详细说明构造的过程，细节可看见7等文献。图2：序列agtacgcctag的块排序索引结构下面我们描述基于块排序的索引结构。为了描述的方便我们将矩阵m的第一列称为f列，最后一列成

13、为l列。fi和lj分别表示f列中的第i个字符和l列中的第j个字符。另外，我们引入函数count (f,i)在f列的前i个字符中，fi出现的次数。类似的，我们也定义count (l,j)。现在，我们给出在块排序结构中一个字符的后继字符的表达。对于一个字符fi，如果fi=lj并且count（f，i）=count（l，j），则fj是fi在给定序列中的后继字符。块排序结构由f列和t列构成，ti是一个整数，它表示fi的后继字符在f列中的位置，即fti是fi在原序列中的后继字符。例如，图2给出序列agtaccgcctag$的块排序结果。f3的后继字符是ft3，即f5。块排序索引就有如下的保序性：性质1对于

14、f中的任意字符fi和fj，如果fi=fj并且ij，则titj。证明：令si是矩阵m中的第i行代表的序列，则fi是si的第一个字符。对于两个序列si和sj，它们的先后关系是由两序列中的第一个不相同的字符决定的。如果fi=fj，则si和sj的顺序与sti和stj的顺序是相同的。因此，若ij，则titj。例如，在图2中 f1=f2=f3=a，则（t1=3） (t2=5) successor.maxscore ) then7: successor.maxscore mi,1;8: successor.maxf max (mi,1 + hi );9: end if10: end for11: if (s

15、uccessor.maxscore = successor.maxf) & (successor.maxscore minscore ) ) then 12: successor. tag accepted;13: return successor; 14: else if (successor.maxscore maxscore）,则要看有无必要在该路径上作进一步扩展：如果maxf小于minscore，则返回标记为“unviable”的searchvector，否则标记为“viable”。 4. 性能测试和分析为了测试基于块排序过滤器的局部比对搜索算法的性能，我们在同样的平台上实现了我们的算

16、法和基于后缀树的算法（oasis），所有的测试都是在一台intel pentium 4 *2.0ghz cpu，2g内存，1t硬盘ibm x series 345服务器上进行的。我们使用比较流行的人类染色体chr22 和chr18的片断和drosophila melanogeneses作为测试数据，使用当前比较流行的单位编辑距离替代矩阵进行比对分值计算。系统使用的操作系统为windows 2003，使用microsoft visual studio.2003编译器编译。4.1阈值minscore对算法性能的影响oasis 是基于后缀树技术进行精确局部比对查询的算法，在实验中，我们算法与oasi

17、s性能比较，我们通过minscore 参数控制查询的敏感度。公式（2）给出期望得到的比对个数e与比对得分s的关系： （2） 其中，m是查询序列的长度，n为数据库的大小，k和是通过blast计算的常数，所以算法中的minscore可以用公式（3）进行计算： （3） 图 5 给出根据minscore的影响下，两种算法的比较。对于查询序列，我们用drosophila melanogenetic chromosome nt_033778 的子序列（taaacttggcctgctt）,用drosophila melanogenetic chr4 全序列作为目标序列。查询时间在图中是取以2为底数的对数数值

18、，查询序列的长度为16，minscore的范围412，因为我们算法的启发性比较好，所以从图中不难看出，minscore对我们的算法的影响比较小。图5：阈值minscore对查询性能的影响4.2 查询长度和数据库数据大小对算法性能的影响我们根据查询序列的长度和目标序列的大小进行性能比较。图6说明了不同查询长度对查询性能的影响，图7给出不同目标序列大小对查询性能的影响。我们使用8m的目标序列，目标序列是人类的基因序列chr22片段，查询序列为随机截取chr18的子序列。因为在我们实验机器的内存限制下，oasis无图6：查询长度对查询性能的影响法为这些数据集大于10m的目标序列构建后缀树索引，而bl

19、ock sorting (bs)我们可以为17m的数据构造索引结构。所以在图6中我们的目标序列大小为8m的数据，图7中oasis算法我们最大的目标序列为9m的数据，bs算法则测试到17m的数据。扩展对于oasis算法对于返回的accepted 结果，如果为未非叶子节点，则必须通过遍历后缀树后，然后根据得到的叶子节点才能得到子序列在目标序列中的位置，而基于block sorting（bs）算法中，返回accepted结果的vec中就已经保存了子序列的位置信息。所以在一定程度上提高了速度。图7： 目标长度对查询性能的影响4.3过滤效率的性能比较我们同时计算所要求扩展的列，即目标序列中需要进行动态规

20、划计算的部分，从而反映出两种算法的过滤效率。图8给出两种算法的过滤能力比较，对于block sorting（bs）算法每次扩展一步后就会加以判断当前可能的最优扩展的片段向量，然后对此片段向量进行继续扩展，而对于oasis，根据后缀树结构，需要根据节点的对应的子序列的长度进行多步扩展后才能对当前的所有可扩展后缀树节点做出最优的判断，而算法所有扩展步数之和即为动态规划算法所需要计算的部分。因而bs算法过滤效率比oasis算法要好的多。图8：过滤效率比较5. 结论本文提出了一种基于块排序索引的序列局部比对查找快速有效的算法。该方法利用块排序结构索引超常生物序列，根据块索引性质，提出了基于块索引结构完

21、成动态规划算法，在计算过程中，对进行计算的目标子序列的进行削减，以提高计算效率。性能测试结果表明，构建和存储块排序索引结构所需要的代价比后缀树要小的多，块排序索引是一种非常适合生物序列计算的数据结构。本文主要针对dna序列局部比对查询处理。在以后的工作中，我们将利用块排序技术解决多序列比对等复杂的查询问题。参 考 文 献1 the human genome project (hgp). /hgp/.2 t.smith and waterman. identification of common molecular subsequences.journal of molecular biology 147:195-197,1981.3 colin meek .jignesh m.p, shr