（计算机软件与理论专业论文）基于fft和kmer划分的多序列比对.pdf

摘要 序列比对是现代生物信息学中一个最基本的研究课题。随着生物数据库快 速持续的增长，对多序列比对算法的敏感性和运算速度提出了更高的要求，开 发具有高敏感性和高效率的算法成为当今研究的重点。本文提出一个新的算法 k m dm s a ，它能在保持原来比对精度的前提下，降低比对的时间复杂度。 本文首先介绍了序列比对涉及的基本问题：空位罚分，替换矩阵和比对结 果评价标准。接着对基于渐进方法构建的多序列比对算法c l u s 诅1 w 进行了深入 的研究。然后通过对这些算法的分析，对当前最流行的渐进比对提出了改进。 m a f f t 最早将分治思想应用到序列比对。它把快速傅立叶变换应用到序列 比对，使得两序列比对的时间复杂度从o ( l 2 ) 降低到o ( l l g l ) 。通过把k - m e r 应 用在序列比对中，并结合分治思想，本文提出了一个新的多序列比对方法： k m d 。它通过快速的识别同源区域，把比对问题分成子问题予以解决。_msa 而寻找同源片断的时间复杂度也从d l 班) 降低到0 ( ) 。为了降低时间和空间 复杂度，该方法包含两个技术：k - m e t 查找和压缩字母。为了验证比对效果，把 k m d 的性能同大多数当前流行的方法作以比较。实验结果表明_msa k m dm s a 和其他方法具有可比的精度，而时间上却有更小的开销。这也证明 该方法是有效的多序列比对方法。 关键词：序列比对m a f f t 分治快速傅立叶变换詹1 t l e rk m d s a a b s t r a c t s e q u e n c ea l i g n m e n t i st h em o s tc o l m n o nf u n d a m e n t a ls u b j e c ti nm o d e m b i o i n f o r m a t i c s a sb i o l o g i c a ld a t a b a s e sg r o w i n ge x p o n e n t i a l l y ，i ti sp o s e dh i g h e r r e q u e gf o rt h es e n s i t i v i t ya n de f f i c i e n c yo fm u l 矗p l ea l i g n m e n ta l g o r i t h m s t h e r e s e a r c ho ff a s ta n ds e n s i t i v eb i o l o g ys e q u e n c ea l i g n m e n ta l g o r i t h mi sac u r r e n th o t t o p i co fb i o i n f o r m a t i c s w ep r e s e n tan o v e la p p r o a c hc a l l e dk m d _ m s a ，w h i c h d r a s t i c a l l yr e d u c e st h et i m ec o m p l e x i t yw i t h o u ts a c r i f i c i n gt h ea c c u r a c y f i r s t l y , w ed e s c r i b et h eb a s i cp r o b l e ma b o u ts e q u e n c ea l i g n m e n t ，g a pp e n a l t y , s u b s t i t u t i o nm a t r i xa n ds t a n d a r do fa s s e s s i n ga l i g n m e n tr e s u l ta n ds t ) o n ，s e c o n d l y ， w em a i n l ys t u d ya n dd e s c d b et h ea l g o r i t h mc l u s t a l ww h i c hi sb a s e do nt h e p r o g r e s s i v ea l i g n m e n ts t r a t e g y t h e nt h r o u g ht h ea n a l y s i so ft h e s ea l g o r i t h m s ，w e i m p r o v et h em u l t i p l es e q u e n c ea l i g n m e n t sa l g o r i t h m sb a s e d o nt h e p r o g r e s s i v e a l g o r i t h ms t r a t e g y t h ed i v i d e a n d c o n q u e ra r ea p p l i e di nm a f f te a r l i e s t ，w h i c hr e d u c e dt h ec p u t i m ef r o mo ( l 2 ) t od ( l 1 9 l ) b yt h ea p p l i c a t i o no ff f ti na l i g n m e n t i nt h i sr e p o r t , c o m b i n i n gt h ed i v i d e a n d - c o n q u e ra p p r o a c h , k - m e ri sa p p l i e dt ot h es e q u e n c e s a l i g n m e n t ，an o v e lm u l t i p l es e q u e n c ea l i g n m e n tp r o g r a m ，k m d _ m s a ，h a sb e e n d e v e l o p e d h o m o l o g o u sr e g i o n s a r er a p i d l yi d e n t i f i e d ；t h eo r i g i n a la l i g n m e n t p r o b l e mi sd i v i d e di n t os o m es u b - p r o b l e m st os o l v e a n di tr e d u c e st h et i m eo f i d e n t i f i e dh o m o l o g o u ss e g m e n t sf r o mo ( l l # ) t od ) i nt h i ss t u d y , w ec o m b i n e t w ot e c h n i q u e s ，k - m e rl o o k u pa n dc o m p r e s s e da l p h a b e t s ，t oa c h i e v er e d u c t i o ni nt h e t i m ea n ds p a c ec o m p l e x i t yo fp a l r w i s ea l i g n m e n t t h ep e r f o r m a n c eo fk m d m s a i sc o m p a r e dw i t ht h em o s tp o p u l a rm e t h o d sb yb e n c h m a r kt e s t s t h et e s ts h o w st h a t t h et i m eo fk m d _ m s ai sd r a s t i c a l l yr e d u c e da sc o m p a r e dw i t ho t h e rm e t h o d sw i t h c o m p a r a b l ea c c u r a c y t h e s ea l g o r i t h m sa r ea l s op r o v e df e a s i b l ea n de f f i c i e n ti n b i o l o g ys e n s i t i v i t ya n dc o m p u t i n ge f f i c i e n c y k e y w o r d ：s e q u e n c ea f i g n m e n tm a f f td i v i d e - a n d c o n q u e r f f tk - m e r k m dm s a 创新性声明 本人声明所呈交的论文是我个人在导师指导下进行的研究工作及取得的研 究成果。尽我所知，除了文中特别加以标注和致谢中所罗列的内容以外，论文 中不包含其他人已经发表或撰写过的研究成果；也不包含为获得西安电子科技 大学或其它教育机构的学位或证书而使用过的材料。与我一同工作的同志对本 研究所做的任何贡献均已在论文中做了明确的说明并表示了谢意。 申请学位论文与资料若有不实之处，本人承担一切相关责任。 差 本人签名：u 势a 日期呻7 彳 关于论文使用授权的说明 本人完全了解西安电子科技大学有关保留和使用学位论文的规定，即：研 究生在校攻读学位期间论文工作的知识产权单位属西安电子科技大学。本人保 证毕业离校后，发表论文或使用论文工作成果时署名单位仍然为西安电子科技 大学。学校有权保留送交论文的复印件，允许查阅和借阅论文；学校可以公布 论文的全部或部分内容，可以允许采用影印、缩印或其它复制手段保存论文。 ( 保密的论文在解密后遵守此规定) 本学位论文属于保密，在年解密后适用本授权书。 立职知 本人签名：。 导师签名： 率_ l 日期矽7 形 日期纽q z ：碰 ， 第一章绪论 第一章绪论 i 1 引言 生物技术是二十世纪最伟大的科学成就之一，正在引发生命科学和技术产 业革命，并改变着人类的生活方式。一方面，生物科学技术迅猛发展，无论从 数量上还是从质量上都极大地丰富了生物科学的数据资源。数据资源的急剧膨 胀迫使人们寻求一种强有力的工具去组织这些数据，以利于储存、加工和进一 步利用，协助人脑完成对这些数据的分析工作。另一方面，以数据分析、处理 为本质的计算机科学技术和网络技术迅猛发展并日益渗透到生物科学的各个领 域。于是，二者结合，一门崭新的、拥有巨大发展潜力的新学科生物信息 学兴起。生物信息学是一门交叉科学，它包含了生物信息的获取、加工、存储、 分配、分析、解释等在内的所有方面，它综合运用数学、计算机科学和生物学 的各种工具，来阐明和理解大量数据所包含的生物学意义【1 1 。 序列比较是生物信息学中最基本、最重要的操作，通过序列比较可以发现 生物序列中的功能、结构和进化的信息。序列比较的根本任务是：通过比较生 物分子序列，发现它们的相似性，找出序列之间共同的区域，同时辨别序列之 间的差异。在分子生物学中，d n a 或蛋白质的相似性是多方面的，可能是核酸或 氨基酸序列的相似，或者是结构的相似，也可能是功能的相似。一个普遍的规 律是序列决定结构，结构决定功能。研究序列相似性的目的之一是，通过相似 的序列得到相似的结构或相似的功能；另一个目的是通过序列的相似性，判别 序列之间的同源性，推测序列之间的进化关系。 在生物信息学研究中，最常用的研究手段是比对( a l i g n ) 。序列的比对是 一种关于序列相似性的定性描述，它反映在什么部位两条序列相似，在什么部 位两条序列存在差别。它是一种基本的信息处理方法，对于发现生物序列中的 功能、结构和进化的信息具有非常重要的意义。序列比对就是运用某种特定的 数学模型或算法，找出两个或多个序列之间的最大匹配的碱基或残基数，比对 结果体现了算法在多大程度上反映了序列之间的相似性关系以及它们的生物学 特征。这个过程往往需要通过数据库搜索，找出具有相似性的同源序列。对于 d n a 序列，可推测该未知序列可能属于哪个基因家族，具有哪些生物学功能：而 对蛋白质序列来说，有可能找到已知三维结构的同源蛋白质，而推测其可能的 空间结构与功能。 在序列的比对中，对各种大分子序列进行分析是一项非常基本的工作。序 2 基于f f t 和k - m e r 划分的多序列比对 列分析的主要任务就是从大量的序列信息中获取基因结构、功能和进化等知识。 例如，在发现一个基因或蛋白质序列后，为了推测这个序列的性质、结构、功 能，同源性、变异性等等，首先必须应用序列数据库进行相似性检索，得到与 未知序列相似的目的序列，从而为后续实验确定初步的研究方向。很多有价值 的发现，都有赖于利用计算机进行的序列同源性比较与分析。双序列比对是序 列分析的基础。而对于构成基因家族的成组序列，只有建立多个序列之间的关 系，才能揭示整个基因家族的特征，这需要用到多序列比对。因此，多序列比 对在阐明一组相关序列的重要生物学模式方面起着相当重要的作用。 i 2 生物信息学序列比对现状 双序列比对( p a i r w i s es e q u e n c ea l i g n m e n t ) 是指通过特定的算法对两个d n a 或蛋白质序列进行比较，找出两者之间最大相似性匹配。双序列比对是序列分 析常用方法之一，也是多序列比对和数据库搜索的基础。 双序列比对又分为全局比对和局部比对算法，典型的全局比对算法是 n e e d l e m a n w u n s c h 算法【2 】；局部比对算法的基础是s m i t h w a t e r m a n 算法【3 】双 序列比对的算法很多，但目前大都是基于运筹学中动态规划的算法思想。 多序列比对是现代生物信息学中一个最有挑战性的课题，也是本文研究的 主要问题。现有的多序列比对方法大体可以分为：精确比对方法，渐进比对方 法，迭代比对方法以及基于图论的比对方法。 精确比对方法完全基于动态规划算法，最为经典的是多维 n e e d l m a n - w u n s c h 算法，但其可行的计算维数为3 1 4 】。这是因为多维的动态规 划算法的时间复杂度的数量级是以序列长度为底，维数为指数的幂。 c a r r i l l o - l i p m a n 算法p j 通过减小计算空间，将计算维数提高到1 0 。s t o y e 描述 了一种新的分治算法d c a t 6 j 将序列在不影响其特征表现的前提下，切割成完全 满足c a r r i l l o - l i p m a n 算法的长度要求的子序列集，利用c a r r i l l o l i p m a n 算法进 行序列比对。2 0 0 0 年又提出了一种减少了存储空间的对d c a 算法的递归调用 算法0 m a 算法【7 】 当前，多序列比对大都是采用启发式算法。这类算法包含很多种类，如： 渐进多序列比对、迭代多序列比对、基于图论的比对等方法。 渐进比对方法首先由h o g e w e g 和h e s p e # 8 1 提出，f e n g g 茅1 t a y l o r 1 0 】又加以 完善。其基本思想是重复地对两条序列使用动态规划算法，先由两条序列的比 对开始，逐渐添加新序列，直到所有序列都加入为止。基于这种方法的软件很 多，主要代表是c l u s t a l w 【l j 】。 迭代比对方法是使用比对计分函数反复添加一条附加的序列到已比对的序 第一章绪论 3 列中。这类算法根据改善比对的策略可以分为确定型和随机型迭代比对方法。 最简单的迭代比对类型是确定型。随机迭代方法包括p r r p l l 2 】，隐马尔可夫 ( h m m ) 【l 列模型，模拟退火，遗传算法等其他方法。主要优点在于将优化过程 和目标函数在概念上进行了分离。虽然这类算法并不能提供获得优化比对结果 的保证，但却具有对比对序列个数不敏感的特征。 很多方法都是渐进方法和迭代方法的混合。p r r p 是一个完善的基于动态规 划的迭代比对算法；m u l t a l i g n 方法是一个典型的基于渐进与迭代混合的多序列 比对方法。 基于图论的比对方法的研究起步比较晚，其主要代表就是偏序比对( p o a ) 0 4 ，一种以有向无环图( d a g ) 的表示方式取代行列表示的新多序列比对方法。 p o a 方法利用图的结构，在多序列比对领域打开了一个全新的视角。它允许结 构域重组，使之成为多结构域蛋白质和表达序列标签比对的有效工具。 自从启发式的渐进比对提出以来，很少有关于速度上的改进。k a t o h 等人 提出的m a f f t i ”1 算法在原来的比对基础上，把分治的思想应用到序列比对当 中，是当前最快的方法之一。r o b e r tc e d g a r 将k - m e r 的概念和压缩氨基酸字母 表豹方法【l6 】应用到序列比对中计算距离矩阵，并提出了线性时间内计算两条序 列间的距离的算法。它也可用在所序列比对中搜索同源片断，例如基于字符串 匹配算法的f a s t a l l 7 1 和b l a s t t l 8 】比对。 专一般来说，评价生物序列比对算法的标准有两个：一为算法的运算速度， 二为获得最佳比对结果的敏感性( s e n s i t i v e ) 或准确性( a c c u r a c y ) ，简称为精度。 当前人们虽然已经提出了众多的多序列比对算法，但由于比对问题的计算复杂 性，目前还没有快速且十分有效的方法。所以，序列比对问题尤其是多序列比 对问题已经成为生物信息学中一个非常重要且具有挑战性的研究课题。 1 3 本文所做的工作 如何提高比对的速度和改进比对的精度已经成为序列比对改进的两个重要 方向。本人在阅读了大量的生物信息学资料、期刊、文献之后，针对当前国际 上流行的多序列比对算法进行了细致分析和深入研究，提出了一种改进的多序 列比对的方法k m dm s a 。该方法结合m a f f t 分治的思想，m a f f t 最早将两 条序列比对的时间复杂度从o ( l 2 ) 降低到o ( l i g l ) ”】，而k m dm s a 又将时间复 杂度降低到d ( ) ，并且比对的精度保持不变。 在提出了改进两条序列比对方法后，本文又成功将该方法应用到两组序列 间的比对。因为在渐进比对的方法中，主要是以两两比对构造距离矩阵和后来 的组闻比对为主要计算量，所以，本文结合渐进比对算法，将分治的思想可以 4 基于f f t 和k - m e r 划分的多序列比对 完全应用到渐进比对中的每一步，从而将渐进比对的时间复杂度降低了一个数 量级。 接下来，针对上述算法，设计了一个实用的计算机系统。该系统包括 m a f f t 、k m d _ m s a 两种方法，通过输入标准的序列组文件( m s f 格式文件) ， 输出比对结果。最后以标准的b a l i b a s e 蛋白质库1 9 1 为例，给出比对效果的测 量值：s p s 值和c s 僮【2 0 】。通过实验和现有的多序列比对方法进行了比较，验证 了本方法的有效性和可行性。 最后，在上述工作的基础上，把这些工作整理成论文。 论文内容具体安排如下： 第一章主要介绍序列比对问题的背景知识，序列比对的研究意义和现状， 最后对本人所做工作进行了简要的介绍。 第二章主要介绍序列比对中涉及到的一些基本概念和问题，包括：序列比 对问题的描述，序列比对中涉及的空位罚分，替换矩阵，评价标准等问题。然 后介绍了序列比对的过程，并结合当前多序列比对中的主流算法进行详细介绍。 第三章主要介绍了多序列比对方法m a f f t ，它利用f f t 把序列比对划分 为更小的子问题，即划分为若干段，分段比对，最后把比对的结果合并得到原 问题的结果。 第四章把k - m e r 应用在多序列比对当中，提出了k m dm s a 比对方法。 文中阐述了k - m e r 在序列比对中的两个主要作用：求距离矩阵和划分比对序列 分段比对，并利用压缩字母表改进距离矩阵和划分的性能。这是在第三章介绍 的m a f f t 的划分基础进一步改进。最后，用b a l i b a s e 库中的蛋白质数据对 该算法进行了测试，并且和现有的方法进行了比较分析。 第五章对整篇文章进行了总结，分析方法中的不足，并指出了今后的研究 方向。 第二章生物序列比对基础 5 第二章生物序列比对基础 2 1 序列比对概述 序列比对的理论基础是进化学说。许多生物学的事实表明：不同的核酸或 蛋白质序列可能源于同一原始序列，经过序列内残基的取代，残基或序列片断 的缺失，以及序列重组等遗传变异过程分别演化而来。这些信息对于揭示序列 的结构和功能是至关重要的。因此，序列比对可用于蛋白质的功能域识别，二 级结构预测，基因识别，以及分子系统发育分析等方面的研究。 2 1 i 序列比对的基本概念 在生物分子信息处理过程中，将生物分子序列抽象为字符串，其中的字符 取自特定的字母表。字母表是一组符号或字符，序列由字母表中的字符组成。 如d n a 序列由四种核苷酸组成，用“a ，t ，c ，“g ”代表四种碱基， 复杂度为4 ，“c c a t c _ , - c t a g a t ”可代表一个简单的d n a 序列。d n a 序列中字 母的意义如表2 1 所示。 表2 1 d n a 序列各字母意义 蛋白质序列则是由2 0 种氨基酸组成，由f a ，c ，d ，e ，f ，g h ，i ，k ，l ，m ，n ，p ， q ，r ，s ，t ，v ，w ，y ) 代表不同的残基。另外，用“x ”表示某个不确定的残基，“b ” 表示天冬胺或天冬胺酸，“z ”表示谷氨酰胺或谷氨酸，则完整的蛋白质字母表为： a ，b ，c ，d ，e ，f ，q h ，i ，k l ，m ，n ，p ，q ，rs ，t ，v w ，x ，y z ) ，其复杂度为2 3 ， “b e g s s l 烈m a b n n m a ”可表示一个简单的蛋白质序列。因此生物序列比对可 以看作字符串的比对。这里字符不区分大小写。蛋白质中具体的字母代表的意 义如表2 2 所示。 在理想情况下，同源基因或蛋白质序列在相互比较时，残基之间相互对应， 从而使替代的情况很明显地表现出来。在某些位置，一个序列中拥有某些残基 6基于h 可和k o m e r 划分的多序列比对 而另一个序列中缺少这些残基，表明这些残基是插入到前者或是从后者中丢失 的。这些空位在序列比对时用连续的短线填补。 表2 2 蛋白质序列各字母意义 因此，和上面生物的特性对应，对字符串的比对操作有以下三种： 1 ) 插入：在序列中插入一个或多个字符。 2 ) 删除：在序列中删除一个或多个字符。 3 ) 替换：用另一个字符替代某个字符。 插入和删除互为逆过程。已知两个序列，如果在一个序列中插入一个( 或 多个) 字符得到另一个，那么同样从后一个序列中删除这个字符也会得到前一 个。如果用“一”来表示插入一个字符，则序列比对问题就可以描述为：对两个 或多个字符串序列通过匹配相对应的字符或插入“”来显示插入或删除而得到 序列之间的最大相似性排列，佼得两个序列之间所对应的相同字符最多。 如果给每一种操作( 插入、删除、替换) 赋予一定的分值，则序列比对问 题就是要找出两条序列间各字符对齐得分最优的一种匹配模式。最优的匹配可 以不唯一。 假设比对的两序列分别为s l 和，小写字母a 和b 代表生物分子字母表中 任意字符。序列匹配之间常用的表示方法如下： 第二章生物序列比对基础 7 表示完全匹配。 ，一) 表示从西中删除字符口，或是在中插入空位。 ( a ，b ) 表示用岛中的字符b 替代& 中的a ( a 6 ) 。 ( - ：6 ) 表示在蜀中插入空位，或是从& 中删除字符b 。 s i 和中的问题是对称的，因此在两中的删除也可以看作是在岛中的插 入，反之亦然。 例如下面就是一个简单的双序列比对的例子。 输入序列：s l = a t c g a g c t g g t ，& = a t c g a g c g g t 。 插入空位前两序列最优匹配如下： 两：at cg a gctg gt 兕：at co a gcggt 插入空位后两序列最优匹配如下： s l ：at co a gc tg gt & ：at cg a gc ggt 插入空格前各种对齐方式最多有8 个相同残基，插入空位后有l o 个相同残 基，在插入空位后残基的匹配数比原来没有空位的时候增加了，这就说明空位 插入的必要性，也反应了序列进化过程中的插入及删除过程。但是如何进行空 位的插入就是一个很困难的问题，也是本文要研究的重点。 2 1 2 打分系统 相对于上面所述的每一个操作赋予一个分值，通常插入和删除的分值是相 同的，这样，可以利用比对的算法求出最优得分值。如果每一个操作用得分来 衡量，则求出最大值，反之，则求出最小值。这个分值通常以矩阵的形式存储。 如果在序列比对过程中，只考虑残基的同一性，即两个序列之间完全相同的匹 配残基数目，则这种单一的计分矩阵具有很大的局限性。对于不同的替换，应 该具有不同的分值；同样，对于不同特征的序列比对，即使是相同的操作这个 得分值也应该是灵活变化的。为了得到那些潜在的具有生物学意义的最佳匹配， 就需要改进计分矩阵的性能，因此提出了构建相似性替换矩阵。 在进行序列两两比对对，有两方面的问题直接影响比对得分值：替换矩阵 和空位罚分。粗糙的比对方法仅仅用相同或不同来描述两个残基的关系，显然 这种方法无法描述残基取代对结构和功能的不同影响效果，缬氨酸对异亮氨酸 的取代与谷氨酸对异亮氨酸的取代应该给予不同的打分。因此如果用一个替换 矩阵来描述氨基酸残基两两取代的分值会大大提高比对的敏感性和生物学意 义。在进行序列比对时，需要通过插入空位( g a p ) 使其具有最佳匹配，即序列 8 基于f f t 和k - m e t 划分的多序列比对 之间所对应的相同残基最多。但是。空位的弓l 入意味着两个序列之间残基的插 入或删除，如果引入空位不加限制，所得的比对结果即使分值很高，也缺乏生 物学意义。因此，必须有一种机制来对空位的引入加以限制，常用的方法就是 空位罚分。 相似性替换矩阵 序列比对时的计分方法最常见的是用替换矩阵来表示。用于d n a 序列的替 换矩阵相对比较简单，如常用的有单位矩阵。而2 0 种氨基酸之间的替换，远比 核苷酸复杂，因为没有一个矩阵可以适合各种情况。 相似性替换矩阵的构建，是基于远距离进化过程中观察到的残基替换率， 并用不同的分值表征不同残基之间相似性程度。恰当选择相似性替换矩阵，可 以提高序列比对的敏感度，特别是在两条序列之间完全相同的残基数比较少的 情况下。 虽然针对不同的研究目标和对象应该构建适宜的替换矩阵，但国际上常用 的替换矩阵是突变数据矩阵( m u t a t i o nd a t a ，m d ) 和残基片段替换矩阵( b l o c k s s u b s t i t u t i o nm a t r i x 矩阵，b l o s u m 矩阵) 等，它们来源于不同的构建方法和不 同的参数选择。 突变数据计分方法是基于蛋白质序列中单点可接受突变( p o i n ta c c e p t e d m u t a t i o n ，p a m ) 的概念，通过对蛋白质进化模式的研究而建立的。p a m 矩阵系 列有p a m 2 0 ，p a m 6 0 ，p a m l 2 0 ，p a m 2 5 0 等。数字表示进化距离，数字越小， 距离越近，反之距离越远。比如p a m 2 5 0 相当于序列闻有2 0 相同残基，更适 用于相似度在2 0 左右的序列比对。距离直接从向导树被测量。相似度范围和 使用p a m 矩阵系列对应为：8 0 1 0 0 ：p a m 2 0 ，6 0 8 0 ：p a m 6 0 ，4 0 6 0 ： p a m l 2 0 ，0 4 0 ：p a m 3 5 0 。 残基片段替换矩阵( b l o s u m 矩阵) 是由h e n i k o f f 夫妇于1 9 9 2 年提出的， 针对序列不同的距离，也构造了一系列的b l o s u m 矩阵。如高于或等于8 0 相同残基的序列组成序列模块可用于产生b l o s u m 8 0 矩阵等。使用b l o s u m 系列的范围是：8 0 1 0 0 ：b l o s u m 8 0 ，6 0 8 0 ：b l o s u m 6 2 ，3 0 6 0 ： b l o s u m 4 5 ，o 3 0 ：b l o s u m 3 0 。 空位罚分 空位罚分是为了补偿插入和缺失对序列相似性的影响。由于没有什么合适 的理论模型能很好地描述空位问题，因此空位罚分缺乏理论依据而更多的带有 主观特色。常用的空位罚分规则有三种： 1 ) 常量空位罚分 常量空位罚分是对插入比对序列中的每个空格都赋予一个常数量的罚分 ，整个比对的空位罚分就是插入的全部空格数k 的罚分之和，即玎。七。如下 第二章生物序列比对基础 9 面比对的两条序列蜀和岛： s l ： ct te f fedt ta cc s 2 ：c 一一edfctf dtt 一一c 则& 和是序列中的空格个数分别为3 和5 ，所以为与岛序列比对的空位 罚分就为( 3 + 5 ) 。 这种罚分策略的优点是简单，不会增加额外的时间复杂度。但是其缺点也 是显而易见的。因为在实际的生物分子进化中，基因中不同位点的突变概率是 不同的，对每个空位使用相同的常量罚分显然是不能准确刻画序列比对的生物 意义。 2 ) 恒定空位罚分 这种罚分策略从整体上处理插入的空格，即序列中插入的相连空格被看作 一个整体的空位进行罚分。这样每个空位的罚分畋都是与其长度无关的。如上 面的两和序列的比对，在& 序列中有2 个整体空位，在& 中有3 个整体空 位，这样这比对的罚分值就为( 2 + 3 ) 。恒定空位罚分虽然可以避免常量空 位罚分将罚分单纯依赖于空位长度的缺点，但是却忽略了空位长度对比对结果 的影响，可能会因为插入过多空格而导致割裂整个比对片段。因此更适当的罚 分策略应当既考虑引入空位长度的影响，又要避免单纯的长度依赖。 3 ) 仿射空位罚分 仿射型的空位罚分把整个罚分分成两个罚分值，起始空位罚分( g o p ) 和 扩展空位罚分( g e p ) 。所谓起始空位，是指序列比对时，在一条序列中插入一 个空位，使两条序列之间有更好的匹配；所谓扩展空位，是指在引入一个或几 个空位后，继续引入下一个连续的空位，使两条序列之间有更好的匹配。任一 空位的出现均需给以起始空位罚分，而任一空位的延伸必须给以扩展空位罚分。 对于一个空位长度为k 的罚分阢可以用如下式表示： 巩= 口+ b k 式( 2 1 ) 其中口表示起始空位罚分。b 表示延伸空位罚分，k 为空位长度。a 和b 两 个参数设置的不同将产生不同的比对结果，其影响如表2 3 所示。 表2 3 起始空位罚分和扩展空位罚分对匹配的影响 起始空位罚分扩展空位罚分 说明 通过对上面公式的改进得到了一种更为理想的空位罚分，改进的公式如下 式所示： 1 0 基于f f t 和k - m e r 划分的多序列比对 = 仃l o g ( k + 1 )式( 2 2 ) 大多数比对算法针对特定的替换矩阵设定了空位罚分的默认值( d e f a u l t ) ，如 果使用者希望使用不同的替换矩阵，则原来的空位罚分设置不一定适合。如何 设置罚分并无明确的理论可循，一般普遍认可的思路是：起始空位罚分值设置 的较大，而扩展空位罚分值则设置的很小。 2 1 3 标准测试集与评判标准 目前还没有一个公认的方法来评价多序列比对结果的优劣。本文为了评价 多序列比对结果的好坏，结合实验所得结果和标准b a l i b a s e 库h 9 1 中的参考数 据进行比较，分别计算两个分值t h es u mo fp a i rs e o r e ( s p s ) 和t h ec o l u m n s c o r e ( c s ) t 2 0 1 。假定实验所得序列个数为，每条序列有m 列，而参考序列的列 数为j j l 磊，第i 列的残基表示为：c f l ，c 1 2 ，g i n ，则这两个值的计算方法分别描述 如下： s p s 。对一列上的每一对残基勺和c 砖，定义p 归如果勺与c i k 比对上即相 同，则p 础为1 ，反之为0 。则每一列的值昌的计算公式如式( 2 3 ) ： s = ，。：。p 拈 式( 2 3 ) 假设s 值为参考数据对应的墨值。 s p s = 盯s i = l ：s 。 z 一 j f i l ” 则s p s 值计算公式如式( 2 4 ) ： 式( 2 - 4 ) c & 如果每一列上的所有残基都檩等，则c 产l ，否则c 产吣，则c s 值计算 公式如式( 2 - 5 ) ： 四= 崩m q 肛 式( 2 5 ) 显然，s p s 是残基对准确对齐的比率，而c s 值是所有序列准确对齐的比率。 实验中将采用这两个值来评估本算法的比对结果。 2 2 序列比对算法 目前，进行序列比对的算法很多，而这些算法大多是基于运筹学中动态规 划的算法思想，只是在其基础上进行了不同程度的改进而已。根据同时进行比 对的序列条数，序列比对分为双序列蛾j ( p a i r w i s es e q u e n c ea l i g n m e n t ) 和序 列比对( m u l t i p l es e q u e n c e a l i g n m e n t ) 。双序列比对是多序列比对的基础。 第二章生物序列比对基础 2 2 1 双序列比对 双序列比对是指通过一定算法对两条d n a 或蛋白质序列进行比较，找出两 者之间的最大相似性匹配。它已经成为序列比对闯题和数据库搜索的基础。 自从f i t c h 提出基于统计的方法，利用计算机来自动的比较蛋白质序列以取 代人眼的观察比较以来，国际上对序列比对的研究已有几十年的历史。最有代 表性的双序列比对算法有点阵图法【2 1 l 和动态规划算法。 点阵图法 点阵图法是一种最简单且容易实现的双序列比对方法。该方法是通过将一 条序列排在上首，另一条序列纵列排在左端，两个序列在任何位置上若出现相 同残基，就在两个序列对应的交叉位置上标注一个点。结果排列成对角线的点 列体现出两条序列间具有相同的残基，从而形象的表明序列间的相似性。如图 2 1 所示。 图2 1 点阵图 点阵图法能较容易的显示插3 缺失( 对角线水平或垂直偏移) 以及正向和反 向重复片断的存在。这种方法的主要优点在于可以找到两条序列间的所有可能 的残基匹配，但局限是大部分的点阵计算机程序不能显示真实的比对序列。因 此，实际比对中需要用其它的比对方法来检测，如动态规划算法，这些方法是 自动的。 动态规划算法 自从n e e d l e m a n 和w u n s c h 首次提出动态规划算法以来，双序列比对算法 在其后的三十多年中得到了广泛的应用和改进，成为序列分析的一个重要理论 基础。下面就详细介绍动态规划算法。 动态规划算法非常适用于序列比对。长度为所和栉的序列s l 【1 m 】和 【l 珂】在全长范围的比对结果，包含了两和的前缀子序列研 1 司和& 1 胡 ( 1 s 赶抑，1 9 鲫) 的比对结果，这是一个递归的关系。因此可以先求解子序列的最 1 2 基于f f t 和k - m e r 划分的多序列比对 优值，根据递归关系求解更大规模的予问题( 更长的子序列) ，直到求得序列 & 【l m 】和【l 川的比对最优值。然后根据各阶段最优值信息，采用回溯方法， 构造出最优序列比对结果。 假设对序列岛【1 胁】和【1 川利用动态规划算法进行比对。则首先构建 一个大小为( m + 1 ) 仍+ 1 ) 的矩阵，矩阵中的元素m i ，刀( 0 k ，一个k 次多项式也是一个次数范围为以的多项式。 如果彳和口都是次数界为抑的多项式，则说他们的积是一个次数界为 2 阶l 的多项式c ) 。其方法是把爿( d 中的每一项与b ) 中的每一项相乘然后再 把同类项合并。 n - i b ( z ) = b j x 7 ；o 2 n 2 c ( 石) = c j x 7 = 0 c i ：皂a - 0 式( 3 2 ) 式( 3 3 ) 式( 3 - 4 ) 在某种意义上说，多项式的系数表示法与点值表示法是等价的，即用点值 形式表示的多项式都对应一个系数形式的多项式。下面将介绍这两种方法，并 基于f f t 和k - m e r 划分的多序列比对 阐述如何把这两种表示结合起来，从而使两个次数界为月的多项式乘法运算在 0 ( n l g n ) 的时间内完成。 系数表示法 一个次数界为甩的多项式a o ) 的系数表示法就是一个由系数组成的向量口= ( a o ，口1 ，口 1 ) 。 采用系数表示法对于某些多项式的运算是很方便的。例如，对两个分别用 系数向量表示的多项式进行相加所需的时间是e ( n ) 。 对于两个用系数形式表示的次数界为聆的多项式彳和鼢) 的乘法运算， 如果用式( 3 2 ) 和式( 3 3 ) 所描述的方法，完成多项式乘法所需要的时间就是 0 ( 矿) 。由式( 3 3 ) 给出的结果系数向量c 也称输入向量a 和b 的卷积，表示成 c = 口圆b 。 点值表示法 个次数界为1 1 的多项式4 0 ) 的点值表示就是栉个点值对所组成的集合： ( c 功，y o ) ，( x l ，y o ，t ，y - 0 其中所有a k 各不相同，并且对k = 0 ，l ，n - 1 有： y k = a ( x k ) 一个多项式可以有很多不同的点值表示，这是由于任意胛个相异点x o ，期， 稚l 组成的集合都可以作为这种表示法的基础。 对于一个用系数形式表示的多项式来说，计算其点值表示是简单易行的， 因为所要做的就是选取r 1 个相异点x o , x 1 x n 1 ，然后对k ；0 ，l ，舻1 ，求出 a ( x i ) 的值。使用霍纳法则，求出这阼个点的值所需要时间为0 ( 矿) 。在稍后将看 到如果巧妙地选取x k ，来加速这计算过程，从而使其时间复杂度为 ( n l g n ) 。 这一过程也叫做f f t 变换。 求值运行的逆运算( 从一个多项式的点值表示确定其系数表示中的系数) 称为插值，也叫做f f t 反变换。可以证明，对于任意n 个点值对组成的集合 ， y o ) ，o l ，y 1 ) ，( h i ，肋1 ) ) ，存在唯一的次数界为r 的多项式4 ，满足y k = a ( x k ) ， = 0 ，1 ，舻1 。 因此，r 个点的求值运算与插值运算是有完备定义的互逆运算，运用这两 种运算就可以把多项式的系数表示与其点值表示相互进行转化。这两种运算用 上述算法的时间复杂度都为o ( n 2 ) 。 点值表示法对许多关于多项式的操作是很方便的。对于加法，如果a ) = 爿0 ) + 戢) ，则对于任意点x k ，c ( x k ) = a ( x k ) + b ( x k ) 。类似地，点值表示法对于 多项式乘法也是很方便的。如果c ( x ) = a 口，则对任意点x k ，有c ( x k ) = a ( x o 曰陬) ，把4 的点值表示点乘口的点值表示，就可以获得c 的点值表示。不过， 因为c 的次数界是4 的次数界与b 的次数界的和，4 和b 的标准点值表示法是 第三章基于f f t 的多序列比对方法 由n 个点值对所组成，而c 的次数界为2 甩，要获得c 的点值表示需要2 甩个点 值对，必须对a 和b 的点值进行“扩充”，使每个多项式都包含2 力个点值对。 如果给定两个扩充点值形式的输入多项式，就会看到，使其相乘从而得到 结果的点值形式所需要的运行时间为e ( n ) ，这要比采用系数形式的两个多项式 相乘所需的时间要少得多。 能否利用关于点值形式表示的多项式来加快系数形式表示的多项式乘法运 算的速度呢? 答案依赖于能否把一个多项式从系数形式转化为点值形式( 求值) 和从点值形式转化为系数形式( 插值) 。 虽然可以用任何点作为求值点，但通过精心地挑选求值点，可以把两种表 示法之间的转化所需的时间压缩为o ( n l g n ) 。 单位元素的复根 单位元素的n 次复根是满足w n = l 的复数h ，。实际上，单位元素的n 次复根 有，1 个，它们是e 捌k n ，j = 0 ，l ，一1 。为了解释这一式子，利用复数的幂的 定义： 矿= c o s ( u ) + is i n ( u ) 。 式( 3 - 5 ) 单位元素的疗个复根均等地分布在以复平面的原点为圆心的单位半径的圆 周上。值w n ；e 2 x v n 称为单位元素的主1 次根，单位元素的所有其他刀次复根都 是w 。的幂，如图3 1 所示。 f 】啊 噬一嵋一 - 1 研j 一0增 图3 1 在复平面上w o ，以一，以的值，其中w 8 ：e 抽4 是单位元的主8 次根 d f t 回顾一下式( 3 1 ) ，在”，以，嵋。处的值( 即在单位元素的h 个，1 次复 根处) ，不失一般性，假定珂是2 的幂，因为给定的次数界是可以增大的。假定 已知彳的系数形式：口= ( a 0 ，a l ，a n 1 ) 。对后= 0 ，1 ，i i - 1 ，定义结果肌如下： 2 4 基于f f t 和k - m e r 划分的多序列比对 n - i m = 4 ( 砭) = 吩钟 j 。o 式( 3 - 6 ) 向量y = ( ，j ，i ，弘1 ) 是系数向量a = ( a o ，a l 一，a n - 1 ) 的离散傅立叶变换 ( d i s c r e t ef o u r i e rt r a n s f o r m ，简称d f t ) ，也写作，= d f t 。( a ) 。