【优秀硕士博士论文】基于q-gram过滤的海量基因序列比对算法研究最终版

硕士学位论文题目基于QGRAM过滤的海量基因序列比对算法研究专业计算机科学与技术研究方向生物信息学授予学位日期2014年6月30日基于QGRAM过滤的海量基因序列比对算法研究计算机科学与技术专业生物信息学是21世纪生物学的核心，通过研究与探索生物信息学，不仅能发现生命科学所包含的本质问题，而且能对人类发展做出巨大的贡献。序列比对作为生物信息学的基础，具有重要的研究价值，通过序列比对，可以得到序列在结构上的相似性，近一步推测出其在功能上的相似性。然而，当今世界的基因数据库过于庞大并且还在以爆炸式的速度逐年增长，对生物数据库的研究需要去粗取精，去伪存真，才能对基因序列做出精准的预测。本文以QGRAM过滤算法以及NEEDLEMANWUNSCH算法为研究对象，详细分析了过滤阶段基于QGRAM的各种过滤器以及下一阶段NEEDLEMANWUNSCH算法的并行实现，在充分研究各种过滤器能力以及基因序列比对并行特点的基础上，实现了基于QGRAM过滤的海量基因序列比对算法的并行实现。本文的核心思想是首先对海量基因序列进行过滤，得到与查询序列相接近的候选序列库，然后将查询序列与已得候选库进行比对，将NEEDLEMANWUNSCH算法并行实现。相比于已有算法，本文工作的优点在于对海量基因序列作了预处理，防止进行大量无用的计算，在保持精度的前提下提高了效率。本文工作的主要内容包括1本文实现了QGRAM过滤算法，通过QGRAM与编辑距离上下界的关系，过滤掉了大量与查询序列相似度较低的序列，得出满足一定条件的候选基因序列集，为其后续算法做准备。2本文对整个基因序列库进行预处理来加快过滤速度，实现了多种过滤器，并且对其组合以及顺序做了分析，在此基础上进一步利用CPU多核心硬件的特点，利用OPENMP来加快过滤阶段所用时间。3对于得出的候选集，本文利用基于GPU的CUDA开发工具，利用其强大的并行处理的能力，将NEEDLEMANWUNSCH算法移植到GPU上，得出最终结果，为近一步对序列功能的分析提供比对结果。实验结果首先对过滤阶段的算法做了对比，结果表明基于QGRAM的过滤算法相比于传统的编辑距离算法，其提高的效率达到两个数量级，并且对各种过滤器及其组合做了分析，得出基于QGRAM的高效过滤顺序。在比对阶段，实验将基于CUDA平台的NEEDLEMANWUNSCH算法与传统的基于CPU的比对算法相比较，实验结果表明，基于GPU的CUDA平台算法实现相比CPU的OPENMP平台而言，有极大的优势。关键词QGRAM过滤，序列比对，并行THERESEARCHOFMASSIVEGENESEQUENCEALIGNMENTALGORITHMBASEDONQGRAMFILTERINGMETHODASISKNOWNTOALL,BIOINFORMATICSISTHECOREOFBIOLOGYINTHE21STCENTURYTHROUGHBIOINFORMATICSRESEARCHANDEXPLORATION,WECANNOTONLYDISCOVERTHEESSENCEOFISSUECONTAINEDINTHEBIOSCIENCE,BUTALSOMAKEAGREATCONTRIBUTIONTOHUMANDEVELOPMENTASTHEBASISFORBIOINFORMATICS,SEQUENCEALIGNMENTHASIMPORTANTRESEARCHVALUEBYSEQUENCEALIGNMENT,SEQUENCESIMILARITYINSTRUCTURECANBEOBTAINEDANDTOGOASTEPFURTHER,ITSSIMILARITYINFUNCTIONCANALSOBESPECULATEDHOWEVER,INTODAYSWORLD,THEGENETICDATABASEISTOOLARGEANDHASBEENGROWINGINANEXPLOSIVESPEEDYEARAFTERYEARRESEARCHONBIOLOGICALDATABASESNEEDTODISCARDTHEDROSSANDSELECTTHEESSENTIALTOMAKEACCURATEPREDICTIONSABOUTGENESEQUENCESTOTAKEQGRAMFILTERINGALGORITHMANDNEEDLEMANWUNSCHALGORITHMASTHEOBJECTOFSTUDY,THISPAPERMAKESADETAILEDANALYSISABOUTTHEVARIOUSFILTERSBASEDONQGRAMINTHEFILTERINGSTAGEANDTHEPARALLELREALIZATIONOFNEEDLEMANWUNSCHALGORITHMINTHENEXTSTAGEBASEDONTHEFULLRESEARCHESOFTHEABILITYOFVARIOUSFILTERSANDTHECHARACTERISTICSOFTHEPARALLELALIGNMENTOFTHEGENESEQUENCE,THEPARALLELREALIZATIONOFMASSIVEGENESEQUENCEALIGNMENTALGORITHMBASEDONQGRAMALGORITHMCOMESTRUETHECOREIDEAOFTHISPAPERISTOFILTERTHEMASSIVEGENESEQUENCESANDTOOBTAINACANDIDATESEQUENCEDATABASECLOSERTOQUERYSEQUENCE,THENCOMPARETHEQUERYSEQUENCEWITHTHECANDIDATEDATABASE,ANDFINALLYMAKENEEDLEMANWUNSCHALGORITHMFORPARALLELIMPLEMENTATIONSUCCESSCOMPAREDWITHEXISTINGMETHODS,THEADVANTAGEOFTHISWORKISTHEPRETREATMENTFORMASSIVEGENESEQUENCE,PREVENTINGALARGENUMBERSOFUSELESSCALCULATIONS,IMPROVINGEFFICIENCYWHILEMAINTAININGTHEACCURACYREQUIREMENTTHEMAINCONTENTSOFTHISWORKINCLUDE1THISPAPERACHIEVESTHEQGRAMFILTERINGALGORITHMFILTERINGOUTALOTOFDISSIMILARSEQUENCEBYTHERELATIONSHIPBETWEENQGRAMANDEDITDISTANCE,ANDGETSACANDIDATESETMEETINGCERTAINCONDITIONSINPREPARATIONFORTHESUBSEQUENTALGORITHM2ACCELERATETHEFILTERINGSPEEDTHROUGHTHEPRETREATMENTFORTHEENTIREGENESEQUENCEDATABASETHISPAPERACHIEVESAVARIETYOFFILERSITHASALSOMADEANANALYSISOFFILTERSCOMBINATIONSANDORDERSFURTHERMORE,BYUTILIZINGMULTICOREHARDWAREFEATURESOFCPU,THISPAPERUSESOPENMPTOSPEEDUPTHEFILTERINGSTAGE3FORTHERESULTINGCANDIDATESET,THISPAPERWITHCUDADEVELOPMENTTOOLSBASEDONGPU,USINGITSPOWERFULPARALLELPROCESSINGCAPABILITIESTOTRANSPLANTEDTHENEEDLEMANWUNSCHALGORITHMTOTHEGPU,ANDGETTHEFINALRESULT,PROVIDECOMPARATIVERESULTSFORCLOSERANALYSISOFTHESEQUENCEFUNCTIONFIRSTLY,THERESULTSOFRESEARCHMADECOMPARISONSBETWEENALGORITHMSINTHEFILTERINGSTAGEITISSHOWNTHATTHEEFFICIENCYOFALGORITHMBASEDONQGRAMIMPROVESBYTWOORDERSOFMAGNITUDETHANTHATOFOTHERSBASEDONCONVENTIONALMETHODSTHENTHISPAPERALSOMADEANANALYSISOFVARIOUSFILTERSANDTHEIRCOMBINATIONANDGOTAHIGHEFFICIENCYFILTERINGORDERBASEDONQGRAMALGORITHMINTHECOMPARISONPHASE,THEEXPERIMENTCOMPAREDNEEDLEMANWUNSCHALGORITHMBASEDONCUDASPLATFORMWITHTRADITIONALNEEDLEMANWUNSCHALGORITHMBASEDONCPUSPLATFORMTHERESULTSSHOWTHATALGORITHMBASEDONGPUSCUDAPLATFORMHASAGREATADVANTAGETHANTHATBASEDONCPUSOPENMPPLATFORMKEYWORDSQGRAMFILTERING,SEQUENCEALIGNMENT,PARALLEL目录目录I1绪论111课题研究背景112国内外研究现状3121序列的相似性查询研究现状3122基因序列比对研究现状713本文主要工作1014本文组织结构102序列比对相关算法概述1221序列相似性查询12211QGRAM定义12212QGRAM拆分及其与编辑距离之间的关系13213倒排索引结构1322双序列比对算法研究14221NEEDLEMANWUNSCH算法16222HIRSCHBERG算法1823并行技术介绍19231OPENMP技术介绍19232CUDA技术2124本章小结233基于QGRAM过滤的基因序列比对算法实现2431QGRAM过滤算法实现24311算法流程24312过滤算法描述26313过滤器介绍以及实现27314过滤顺序模型及实现32315基于OPENMP的实现36316过滤算法的数据结构和UML图3832基于GPU的双序列全局比对算法39321算法实现39322算法性能优化4433本章小结484实验结果及其分析4941实验环境4942实验数据4943实验结果与分析50431QGRAM过滤算法结果及分析50432基于GPU的双序列全局比对算法5644本章小结605总结与展望6151本文总结6152研究展望62参考文献63作者在读期间的科研成果69声明70致谢71图21QGRAM拆分13图22倒排列表14图23两个序列的比较15图24序列比对16图25得分矩阵计算描述图17图26FORKJOIN并行执行模型20图27CPU与GPU的晶体管比较21图28CUDA编程模型22图29KERNEL函数的调用23图31算法流程24图32索引创建实例25图33过滤验证模式26图34前缀过滤28图35基于内容误配的过滤器31图36过滤阶段UML设计38图37将矩阵Q划分为象限40图38计算矩阵的边界值41图39回溯阶段43图310GPU计算单元简图45图311CUDA存储器模型46图41数据库长度分布图50图42序列数据库不同基本过滤器组合的过滤能力比较51图43不同Q值对序列数据库A过滤能力及过滤时间的影响52图44不同Q值对序列数据库B过滤能力及过滤时间的影响53图45序列数据库建立索引的串行与并行时间对比（时间MS）55图46GPU实现的效率比对图（STRINGLENGTH5）57图47效率比对图（STRINGLENGTH5）57图48GPU实现的效率比对图（STRINGLENGTH30）58图49效率比对图（STRINGLENGTH30）59图410GPU实现的效率比对图（STRINGLENGTH75）59图411效率比对图（STRINGLENGTH75）601绪论11课题研究背景21世纪是信息的时代，同时也是生命科学的时代。作为生命科学和信息科学相结合的新兴学科，生物信息学综合利用应用数学、信息学、统计学以及计算机科学等方法来研究生物学中的问题1，其研究的对象是经过各种方式获取的生物学数据，研究的方法包括对生物学数据的收集、筛选、处理以及利用。生物信息学研究是从理论上认识生物本质的必要途径，通过生物信息学研究和探索，可以更为全面和深刻地认识生物科学中的本质问题，了解生物分子信息的组织和结构，破译基因组信息，阐明生物信息之间的关系。生物信息学研究在医学上也有重要的意义。通过生物信息学分析，可以了解基因与疾病的关系，了解疾病产生的机理，为疾病的诊断和治疗提供依据。蛋白质组学和基因组学是生物信息学的两个主要研究方面，具体说就是对核酸序列以及蛋白质序列进行分析，获取序列的结构和功能。在序列分析中，最常用的方法便是将未知序列同已知序列进行相似性比较，近而可由其结构的相似性推测出它们在功能上的相似性，由已知基因的功能预测出未知基因的功能。这种相似性比较的方法在生物信息学中被称为序列比对（SEQUENCEALIGNMENT），有非常重要的理论意义和实践意义。序列比对的理论基础是进化学说2生物序列在进化过程中，某些基因序列保持不变，而另一些则可能会发生基因序列的插入、缺失和取代。如果两条基因序列拥有共同的祖先，则称为同源序列，其可能具有相似的功能。然而，进化历程已经过去，现在序列的同源关系只能通过序列之间的相似程度来推断，在理想情况下，序列比对可以真实反映了基因组的进化历程。所以我们可以通过将各种操作（插入、删除、取代）赋予一定的权值来模拟基因序列的变异过程，找出序列间的最小距离或最大相似度。序列比对的数学模型大概可以分为两类一类是从序列的整体性出发，考虑序列的整体相似性，即全局比对，最有代表性的便是由SAULBNEEDLEMAN和CHRISTIANDWUNSCH于1970年提出的NEEDLEMANWUNSCH算法3，算法首次将动态规划（DYNAMICPROGRAMMING）引入到序列分析中，用来解决两条序列之间的全局比对问题；另一类则是考虑序列的局部相似性，即重点在于找出局部相似的部分，称为局部比对。最经典的方法便是在全局比对算法NEEDLEMANWUNSCH上改进得来的SMITHWATERMAN算法4。然而，这些方法都是利用动态规划的思想，其复杂度很高，随着生物序列数据的爆炸式增长，其运算速度已经很难达到实际应用的需求，而传统的计算机速度一方面受到加工工艺限制，计算机器件的尺寸不可能无限小，另一方面受到光速极限和量子效应的限制，所以计算速度的瓶颈一直不能突破，因此越来越多的人将目光转向并行算法。相对于串行计算，并行计算划分为时间并行和空间并行，随着近年来多核CPU和GPU（GRAPHICPROCESSINGUNIT）的快速发展，空间并行成为人们的研究的主要方向，即研究多处理器执行的高性能计算。例如，基于多核CPU的OPENMP（OPENMULTIPROCESSING）技术充分利用CPU的多核，对硬件进行抽象描述，使得开发人员能快速地掌握相关技术NVIDIA公司在GPU上推出了CUDA（COMPUTEUNIFIEDDEVICEARCHITECTURE）开发工具，该开发工具使程序员能在C语言环境下就可完成对GPU并行处理架构的调用，大大方便了对GPU的使用。近年来，高性能计算技术已经广泛应用于诸如天文学，气象学，生物科学，空间科学，图像处理，人工智能等各个研究领域。另一方面，由于基因数据库过于庞大并且还在以爆炸式的速度逐年增长，因此将未知序列与基因数据库中的序列逐一进行比较不可能用于实际应用。数据库相似性查询技术在近年来被广泛用于解决类似问题，其目的不在于搜索完全匹配的序列，而是在庞大的数据库中检索出与之有关联的序列，近而缩小比对范围，以提高整体查询效率。近年来，随着计算机硬件和软件的不断改进，针对基因数据库的相似性查询产生了很多可用软件，如BLAST5、FASTA6等，并且有相对满意的精度和效率，其核心思想是采用启发式的方法利用近似值加快序列比较。然而，这些软件的核心思想都是基于局部比对算法，将其应用于全局比对产生的效果并不是很理想，而QGRAM算法的提出则使得全局比对得到了快速的发展，QGRAM算法的实现比较简单，并且高效，因此受到越来越多人的关注，其主要依据是越相似的序列其局部相同的部分就越多。针对基因序列数据库的海量性，本文利用了QGRAM与编辑距离上下界的关系，采用QGRAM过滤算法对基因序列数据库进行预处理，并实现了多种过滤器，为了进一步的加快过滤速度，本文利用CPU的多核心硬件，使用OPENMP来加快过滤阶段的时间。在比对阶段，本文选用了生物序列比对的经典算法之一，针对其动态规划实现方法需要的NEEDLEMANWUNSCH算法2时间复杂度和空间复杂度，本文利用GPU上的CUDA开发工具将该算法移植到GPU上，提升其效率，并为后续的序列功能的分析提供比对结果。12国内外研究现状121序列的相似性查询研究现状序列是生物信息学中最普遍的存在，随着近年来科技的发展以及随之而来的海量基因序列，对序列进行相似性查询已经成为一个研究热点。序列的相似性查询是指在海量的基因序列数据库中根据用户给出的查询序列以及所要求的误差范围，找出满足条件的基因序列。对未知序列进行相似性查询不仅可以知道其是否已经被基因序列数据库收入，还可以查找出是否有与它相似的序列，进而通过与查找到的基因序列进行比对来预测所查询序列的功能。如何界定序列之间的相似程度是序列相似性查询最核心的问题。对于给定的两个序列A和B，度量它们相似程度最常用的方法就是编辑距离，记为ED（A，B）。编辑距离度量了将一个序列转化为另一个序列所需的最少操作（包括插入、删除、替换），求两个序列的编辑距离一般都会用到动态规划的思想7，其时间复杂度为。由于编辑距离能很好的反映序列中元（|）素的位置信息，因此很多算法都采用编辑距离作为基础的度量方法。可以看出，求编辑距离最大的难点在于其高昂的时间复杂度，因此相似性查询需要解决的核心难点即是怎么降低计算编辑距离的时间复杂度。1211序列相似性查询种类假定给出一个查询序列T，需要在给定的基因序列数据库D中查找与它在一定误差范围内的基因序列，序列之间的相似性用编辑距离ED来度量，根据查询精度的不同，主要分为以下两类1精确查询顾名思义，即在基因序列数据库D中找出所有与T等价的序列，其主要应用于早期的小型基因序列数据库中。经典的算法有BOYERMOORE算法8和KMP算法9以及KARPRABIN指纹算法10等。2近似查询在生物数据库中，生物序列由于变异的选择性，使人们必须对其进行近似查询才能了解各基因序列间的关系。主要分为两大类范围查询（RANGEQUERY）在数据库D中找到所有满足查询半径R的序列结果集R。由于范围查询需要指|,定查询半径，所以需对序列相似性程序有先验的领域知识。本文主要围绕范围查询展开研究。KNN查询（KTOPNEARESTNEIGHBORSQUERY）当不具备应用领域的先验知识时，可以使用KNN查询，即找出数据库D中前K条与Q最相似的序列，当其取值为1时，就变为最近邻查询，即在D中找出与所查询序列最相似的基因序列。实现KNN查询最常用的方法便是基于空间分割的KD树法11，它能高效插入和删除节点，可以解决动态环境下的最近邻搜索，其查找的时间复杂度为。1212近似查询精确查询的适用范围比较局限，而近似查询可以用误差参数来满足用户的需求，所以近似查询相比精确查询适用范围更广，而且近似查询中序列相似程序可以用以编辑距离为基础的各种距离函数来度量，因此更加灵活，但同时其实现比起精确查询也更加复杂。近似查询根据其实现方法的不同可分为以下四类1基于动态规划的方法在近似查询中，使用动态规划方法计算出得分矩阵，近而根据得分矩阵得出序列的最佳比对是计算编辑距离最基本的方法，该基本方法所需时间复杂度为。2EUKKONEN在12中提出，回溯路径只需两序列得分矩阵中主对角线附近区域的值便可，假设该区域的宽度为R，则该算法的时间和空间复杂度都为，该算法仅适合于较小查询误差的情况，其适用范围并不是很广。ON由上可以看出，基于动态规划的方法核心在于怎么高效的计算序列间的得分矩阵，其效率还是没有在根本上得到改善。FASTA算法6的提出使得序列相似性搜索进入新的局面。该算法于1985年由LIPMAN以及PEARSON提出，在随后的1988年，二者又对其做了进一步的改进，该算法以动态规划的算法作为其基本思想，并在其中引入了启发式算法，虽然其在精确性上不如纯粹的动态规划算法，但在速度上却是动态规划算法的几十倍，其基本思想是一个能揭示出真实序列关系的比对至少包含一个两序列都拥有的字片段，把查询序列中的所有字编成HASH表，然后在数据库搜索时查询这个HASH表，以检索出可能的匹配，这样那些命中的字就能很快地被鉴定出来。在随后的1990年，ALTSCHUL等人在FASTA算法的思想上提出了BLAST算法5，该算法采用了一种短片段匹配算法和一种有效的统计模型来找出目的序列和数据库之间的最佳局部比对效果。其基本思想是通过产生数量更少但质量更好的增强点来提高速度，利用片断匹配算法和有效的统计模型来找出目的序列和数据库之间的最佳局部比对效果。总体而言，这两种算法是目前使用最广泛的序列数据库搜索算法，它们都是基于局部相似的基础，用来在序列数据库中查找与查询序列局部完全匹配的片段，并在此基础上将其延伸。2基于编辑距离的过滤方法近年来，一种比较常见的方法是基于编辑距离上下界来实现序列数据库的过滤。当给定一个查询序列以及查询半径后，可以快速得出编辑距离的上下界，近而对数据库进行过滤。编辑距离满足三角不等式的性质，因此可以利用此性质形成编辑距离的上下界来过滤数据库，13中实现的算法即是依据了编辑距离的该性质，其思想是为每个数据库中的元素分配一个参考子集。当给定查询序列和查询半径后，由数据库D中的每个元素可以计算出EDQ,的上界和下界。3基于QGRAM的过滤方法基于QGRAM的过滤方法是近年来使用最广泛的一种算法，因其算法的灵活性和有效性而受到越来越多的关注，QGRAM算法的核心依据是越相似的基因序列，其局部相同的部分就越多，所以当用户给出一定的误差范围后，便能通过该算法计算出局部相同的部分是否满足给定误差范围，进而可以用此条件过滤数据库。最早提出QGRAM这个概念是在14中，文中首次提出可将序列分割成小片段，并且提出了数量过滤器的定理若，则其至少共享,个QGRAM。在文章15中，作者基于QGRAM提出了MAX|,|1一种查询算法QUASAR，其索引结构用后缀数组实现。16中，文章作者充分观察了序列长度之间的关系，在原有数量过滤器的基础上添加了限制序列之间长度的长度过滤器若，则，以及位置过滤器若,|，则序列A中的任意QGRAM与序列B中该QGRAM的位置之差不能,大于T，并对这些基本过滤器进行比较和组合，分析其过滤能力。然而上述过滤器要对序列的每个QGRAM都进行比较分析，一些研究人员发现了这一瓶颈，近而提出了前缀过滤1718。各种过滤器被不断提出，但是却没有文章去分析其过滤顺序对过滤能力有何影响，19中作者对各种过滤器的先后顺序进行分析，得出过滤器的先后顺序对过滤能力的影响。上述所有过滤器都只注意到去研究序列匹配应该满足什么样的条件，而20中的作者从相反的方向去考虑，研究了基于误配的QGRAM，并得出两个基于误配的过滤器基于内容误配（CONTENTMISMATCH）的过滤器以及基于位置误配LOCATIONMISMATCH的过滤器，在前人基础上近一步改进了算法，提高了效率。前面有关过滤器的研究都是基于连续的QGRAM，而21中的作者在选取QGRAM时采用了抽样的方法，每隔一定距离抽选一个QGRAM，使得过滤效率得到极大的提升，然而以此方法得出的QGRAM之间没有相关性，因此，虽然在效率上获得了提升，而其结果却并不理想，文献22的作者也就该方法与连续采样的方法相比较，证明了连续采样效果更加。前面介绍的算法都是基于定长的QGRAM，23中提出一种新颖的思想在数据库中选取变长的高品质GRAM来构建索引，文章中具体给出选取高品质GRAM的方法。在随后的文献24中，作者进一步对自己所提出的变长QGRAM索引进行研究，主要是研究其与查询性能之间的关系。4随机方法随机方法即利用各种哈希函数，使相似的序列在哈希函数的映射下分配在同一个桶内。随机方法可以会产生假阳性数据和假阴性数据，所以不能保证能找到全部的查询结果。如文章25中提出了一种处理高维向量空间中最近邻查询的随机工具位置敏感哈希函数（LOCALITYSENSITIVEHASHINGLSH），并有效利用该函数实现了高维图像中的相似性查询。LSH虽然在查询性能上有优势，但是其主要适用于序列在海明距离下的相似性查询，对于编辑距离而言，随机投影无法对应插入和删除操作引起的序列变化。122基因序列比对研究现状基因序列比对是生物信息学的基础，根据同时进行序列比对的序列个数，序列比对可分为双序列比对PAIRWISESEQUENCEALIGNMENT和多序列比对MULTIPLESEQUENCEALIGNMENT。下面分别介绍两类算法的实现细节。1221双序列比对算法介绍最早实现双序列比对所用的思想为动态规划思想，最典型的算法包括全局比对算法NEEDLEMANWUNSCH和局部比对算法SMITHWATERMAN。以下对两种经典算法做简单介绍1）NEEDLEMANWUNSCH算法NEEDLEMANWUNSCH算法3是由SAULBNEEDLEMAN以及CHRISTIANDWUNSCH于1970年提出的，该算法是关于序列比对最早的算法，是生物信息处理中最基本的算法。其算法使用动态规划的思想迭代计算两个序列，将其得分存储于一个得分矩阵，然后根据得分矩阵回溯找到序列的最佳比对。NEEDLEMANWUNSCH算法考察两个序列的全局相似性，适用于全局水平上相似性程序较高的两个序列的比对。2）SMITHWATERMAN算法SMITHWATERMAN算法4是在NEEDLEMANWUNSCH算法的基础上发展而来的，该算法于1981年被SMITH和WATERMAN提出的，它不考虑两序列的整体相似性，而将重点放在查找局部相似性最高的片段上，适用于亲源关系较远、整体上不具有相似性而在一些较小的区域上存在局部相似性的两个序列。1222多序列比对算法介绍虽然双序列比对是序列分析的基础，然而，对于由一个序列组构成的基因家族来说，需要建立多个序列之间的比对关系才能显示整个基因家族的特征。所以多序列比对也是生物信息学的重要课题之一。现有的多序列比对大都采用启发式思想，其中最常见的就是渐近式比对算法，在双序列比对的基础上逐步优化多序列比对的结果，然后对比对结果做进一步的加工处理。最广泛使用的CLUSTALW算法26便是基于多序列比对的，该算法于1987年由FENGDF，DOOLITTLERF提出，基本思想是在双序列比对的基础上通过建“树”的思路来进行多序列比对，具体步骤如下1）将序列组或者整个家族基因序列各自两两比对2）将两两比对的数据进行聚类分析，产生比对等级，该等级可用树的形式表示3）先从所有的比对结果中找出相似性最好的两条序列，然后从剩余的比对中找出相似性最好的，依次类推，直到所有的序列都被加入本文主要讨论了双序列比对算法的一些基本实现及其优化。传统的双序列比对算法在空间和时间上的复杂度都为，然而在实际的应用中，通常会2有海量的基因序列，这显然不满足应用需求，针对这一问题，近年来出现了很多对其加速的方法，按原理可组织为两类一是在算法本身方面的改进，主要是通过降低算法复杂度来提高速度；二则是通过硬件加速平台来提高速度。1223基于算法的改进HIRSCHBERGDS27提出了将得分矩阵压缩的方法，降低了空间复杂度，进而改进了整个算法，然而该算法以较高的时间复杂度来换取较低的空间复杂度，其时间复杂度增至原有算法时间复杂度的两倍。UKKONEN算法28中提出，基于最佳比对的获得并不需要得分矩阵中的所有元素，并且该算法针对对角线方向的值单调非递减这一原理进行改进，算法的空间复杂度为，时间复杂度2为N为序列长度，D为序列得分值。DIVIDEANDCONQUER算法29是上2述两种算法的结合，其算法将CHECKPOINT技术应用于UKKONEN算法来寻找CHECKPOINT点（即最佳路径经过的点），递归找出所有的CHECKPOINT点进而构成最佳的比对路径。DIVIDEANDCONQUER算法的空间复杂度为线性的，但是其时间复杂度相对于UKKONEN算法并没有很大的改进。1224基于硬件平台的改进在算法方面的改进虽然能在一定程度上提高序列比对效率，但是在基因序列爆炸式增长的今天，仅仅在算法层面上的改进已经满足不了实际应用的需求，随着硬件技术的不断改进，人们开始将之利用于生物信息学中。起初，人们利用特殊的硬件平台对序列比对进行加速，如BISP30，PNAC31以及SAMBA32等，它们利用现场可编辑逻辑门阵列硬件（FPGAFIELDPROGRAMMABLEGATEARRAY）来完成对核心算法的实现，虽然加快了速度，但是其不易进行扩展，灵活性很差，而且价格昂贵。随后人们尝试在多处理器上实现并行3336，但是商用的多核处理器价位较高，并且单纯依靠CPU进行计算已经无法满足人们的需求。GPU的出现引起了开发人员的关注，其不仅价格低廉，而且相对于CPU有更多的计算单元（图形渲染的需要），因此低价的GPU逐渐成为新的加速平台。LIU利用GPGPU平台实现了基因序列比对算法SMITHWATERMAN，使得性能得到提升37。JUNGS使用两个GPU完成SMITHWATERMAN算法的并行加速化38。陈波综合应用在CPGPU上的各种优化方法实现了并行的SMITHWATERMAN算法，获得了几十倍的加速比39。张林等综合分析了各种比对算法，并将准确性最高的动态规划算法在GPU上实现，结果证明其算法在效率上有很大的提升40。韦刚、马海晨等则通过在GPGPU上对读写延迟、数据传输以及重组函数进行优化，实现SMITHWATERMAN算法的加速，与传统的CPU算法相比较，能达到约100倍的性能提升41。MANAVSKI42第一次将算法移植到CUDA（NVIDIA的统一架构计算平台）上，在此之后，CUDA因其便捷的编程环境而成为研究人员关注最广泛的平台4345。林江、唐敏等人使用菱形数据布局来利用GPU的并行处理能力，使用查询串分批技术支持大规模的序列比对，使用树形算法优化最大匹配值的计算，并将该算法在NVIDIAGEFORCEGTX285显卡上测试，结果表明，其改进与CPU上的串行算法相比，能获得120倍以上的性能提升46。13本文主要工作首先，针对庞大的数据库，本文采用了QGRAM（即用长度为Q的窗口从序列的开始向后移动所得到的子序列）过滤技术，QGRAM算法的核心依据是两条序列越相似，其局部相同的地方就越多，近而拥有越多共同的QGRAM。所以，当用户给出一定的误差范围后，便能通过该算法计算出局部相同的部分是否满足给定的误差范围，进而可以用此条件过滤数据库。本文选择合适的Q值建立索引结构，构造了各种不同的过滤器，并且将各种过滤器组合，过滤掉不符合条件的序列，得到在一定误差范围内的序列，本文还对各种过滤器及其组合进行了分析研究，得出最佳的过滤流程，可以看到过滤器的合理选用明显提高了对数据库中数据的筛选能力和效率。其次，针对得到的满足条件的序列，本文选用了生物序列比对经典算法之一NEEDLEMANWUNSCH算法3，并且针对海量数据的问题，使用NVIDIA公司推出的基于GPU的CUDA开发工具，将算法移植到GPU上，并且在针对GPU各种内存大小限制的情况提出了分块的解决办法，极大的提升了算法效率，并且对其实现进行了分析研究。14本文组织结构本文的内容安排如下第一章介绍了生物信息学及它的两个重要研究方向基因序列比对和序列的相似性查询，并引出了高性能计算的内容，还介绍了国内外的研究现状以及本文的主要工作和组织结构。第二章介绍了QGRAM过滤算法的基本内容，介绍了基于动态规划的经典NEEDLEMANWUNSCH算法，还对本文中所用并行技术的基本概念做了介绍。第三章作为本文的主要工作，首先，介绍了基于QGRAM的各种不同过滤器的选用和实现的细节，其次介绍了NEEDLEMANWUNSCH算法在GPU上的并行实现。第四章对在过滤阶段的过滤器性能进行了分析，在比对阶段对其NEEDLEMANWUNSCH算法的并行实现的细节进行了分析，并且将两者结合，得出最终的实验结果。第五章总结了本文工作，并对将来算法的改进作了分析。2序列比对相关算法概述序列比对是本论文的核心内容，在上一章中已经介绍过，本文在实现序列比对之前，首先得到与查询序列相接近的序列数据库，然后再进行序列比对，因此，在本章中，首先介绍序列相似性查询的一些相关概念，然后介绍双序列比对算法的基础知识，并且还介绍了在具体的算法实现中用到的两种并行技术OPENMP技术和CUDA技术。21序列相似性查询在生物信息学中，序列的相似性查询即是指当用户给出一定的误差范围后，在生物序列数据库中找出满足该误差范围的生物序列，度量相似程度最常用的方法的就是进行编辑距离计算，然而上一章中提到计算编辑距离的时间复杂度较高，所以可以通过其他距离函数与编辑距离之间上下界的关系来过滤掉大部分不满足条件的序列，QGRAM算法的实现简单而高效，因此受到越来越多人的关注。其核心思想是越相似的基因序列，其局部相同的部分就越多。本文选择了QGRAM距离来实现序列的相似性查询。211QGRAM定义定义1给定序列，序列T1,2（其中为序列T的第I个元素）长度记为。对于任意的称为序列T的Q|1,GRAM。QGRAM体现了序列的局部特征，两条序列越相似，其局部相同的地方就越多，近而拥有越多共同的QGRAM。在序列T中拆分QGRAM的方法是用长为Q的窗口从序列T的开始向后移动，每次移动的字符个数称为滑动距离，记为D，这种拆分方式记为拆分方式。根据D的不同可以分为不同的拆分方式，当D1时，重叠Q1个字符，记为连续重叠拆分，如图21A所示，文献47即是采用这种拆分方式。当1Q时，相邻QGRAM之间没有重叠且不连接，称为不重叠不连接拆分，如图21B所示。序列提取QGRAM非常简单，并且其对错误的容忍度高。当在一个序列中有错误字符时，仅会影响到包含该错误字符的QGRAM，其余的不会受到影响，尤其当将滑动距离设为1时，能够通过相邻的QGRAM正确找到该段字符，因此本文选用D1的拆分方式对序列进行拆分。212QGRAM拆分及其与编辑距离之间的关系例如序列ATCGATCGAT以Q3来划分为,当将序列的其中一个字母替换ATCGTTCGAT,当在任意位置随机插入一个字母ATCGACTCGAT,当随机删除任意位置的一个字母ATCGATCGAT,例子中红色的部分表示编辑操作导致发生变化的QGRAM，由上面的例子可得出，对于一个序列来说，一个编辑距离（替换、插入、删除）最多能影响到Q个QGRAM的变化。213倒排索引结构为了要提高相似性查询的匹配速度，需要首先为序列数据库建立适当的索引，进而充分利用构建的索引信息快速进行相似性查找。相似性搜索中最常用到的索引结构便是倒排索引48（INVERTEDINDEX），倒排索引可以用来存储某一单词在文档中的存储位置，其结构如下图22所示，在倒排索引中最基本的结ATCGATCGATATCGATCGATA连续重叠拆分（Q4,D1）B不重叠不连接拆分（Q4,D5）图21QGRAM拆分QGRAMQGRAMQGRAMQGRAM构是倒排索引项，一个倒排索引项包含了某个单词的文档编号DOCID，在该文档中出现的次数TF以及该单词在哪些位置出现POS等，而由包含该单词的一系列倒排索引项组成了倒排列表，而所有文档集合中出现过的单词及其倒排列表则组成了倒排索引。在第三章中就基于QGRAM的倒排索引给出详细的构建过程。22双序列比对算法研究双序列比对算法是基因序列比对最基本的比对算法，其根本任务是通过比较两个基因序列，发现它们的相似性，找出序列之间的共同区域，同时辨别序列之间的差异。而比对算法的基础是编辑距离49的计算，编辑距离可以直观评价两个序列的相似程度，距离越小，两个序列的相似度就越高。我们知道，在基因序列的长期进化过程中，由于变异，原本相同的基因序列可能因为一些片段的缺失或者增加，或者发生变化，而导致它们的不同，为了模拟这种进化过程，引入字符串编辑操作的概念。其中，用字符“”代表空缺，定义了如下编辑操作MATCH（A，A）表示字符匹配DELETE（A，）表示从第一个序列删除一个字符，或在第二个序列相应的位置插入空白字符INSERT（，B）表示在第一个序列插入空白字符，或删除第二个序列中的对应字符单词DOCID,TFDOCID,TFDOCID,TF倒排索引项倒排索引项倒排索引项倒排列表词典项图22倒排列表REPLACE（A，B）用第二个序列中的字符B取代前一个序列中相应位置的字符A有了如上编辑操作的定义，便可进行编辑距离的计算。计算编辑距离最经典的方法便是使用动态规划思想，其具体算法如下公式2149（公式21）初始化ED,00,递归方程,1,1,111,110如下图23，第二个序列通过两次编辑操作可以变为第一个，那么两个序列的编辑距离为2。而序列比对就是将两序列相应的位置一一进行比较，对于匹配与否设置相应的值。定义2如果和是两个任意的字符，记为其比较所得的分值，也XY,称为记分函数，记分函数包含了或为空的情况，如果序列在某位置缺失一个字符，则可用“”来表示。如，可以这样定义记分函数，,1。比对分值SCORE可以用公式22来表示,0（公式22）1,如果我们设两序列的对应位置一致，就得1分，如果是空格对碱基，就不得分，如果两个碱基不匹配，就得1分，如图24所示。ATGGCGTATGAGT图23两个序列的比较其中最高分值SCORE所对应的排列称为全局最优比对，为了得到全局上的最优比对，NEEDLE和WUNSCH发明了一种快捷的算法NEEDLEMANWUNSCH算法。221NEEDLEMANWUNSCH算法NEEDLEMANWUNSCH算法3是由SAULBNEEDLEMAN以及CHRISTIANDWUNSCH于1970年提出的，该算法使用动态规划的思想来获得两条序列之间的全局最优比对。由于其考察的是两个序列在全局上的相似性，因此适用于在全局水平上相似程序比较高的序列。其基本思想是迭代计算由两个序列组成的得分矩阵（SCOREMATRIX），然后根据得分矩阵的值选用一定的算法回溯找到两个序列的最佳比对。其得分矩阵中的元素的值可通过公式233迭代得出,01,00,1,0ATGGCGTATGAGT1110111SCORE11101114图24序列比对ATGGCGTATGAGT1011111SCORE10111112,1,1,1,1,公式23得分矩阵中每个单元格的计算如图25，从图中可以看出，矩阵元素（I，J）的值主要依据是对角线方向元素、同一行或同一列的元素。与此相应，到达位置的某一元素也就有三种可能路径通过位置（，）的对角线方向，此种路径无空位罚分；通过列J的垂直方向以（1，1）及通过行I的水平方向，空位罚分值为。（，）或者，其算法描述如下算法21第一步计算得分矩阵ALGORITHM21计算得分矩阵1BEGIN2FORI0TO|X|DO,01,3FORJ0TO|Y|DO0,1,4FORI1TO|X|DO5FORJ1TO|Y|DO6,1,1,1,1,7RETURNM第二步回溯寻找最优比对ALGORITHM22回溯寻找最优比对10100213BEGIN40,1,05161，0，1171，1，1，0，8END而针对第二方面的改进，由于本文在实现过程中并没有涉及，因此不做详细介绍。23并行技术介绍本文为了加速实现基因序列比对，应用了并行技术。相对于串行计算，并行计算划分为时间并行和空间并行，时间并行指时间重叠，就是让多个处理过程在时间上相互错开，轮流重叠地使用同一套硬件设备的各个部分，以加快硬件周转而赢得速度。时间并行的实现方式采用流水处理部件，是一种非常经济实用的并行技术，能保证计算机系统具有较高的性能价格比。然而，随着近年来多核CPU和GPU的快速发展，空间并行成为人们的研究的主要方向，即研究多处理器执行的高性能计算。本文对论文中用到的相关并行技术做了简单的介绍。231OPENMP技术介绍我们知道，近年来，诸如提升时钟速度和指令吞吐量这些改善CPU性能的传统方法由于摩尔定律的限制，已经很难有突破性的进展，所以人们将提升单个处理器性能转而到研究多个处理器的方向上来，即在一个处理器芯片上集成多个处理器核心，并且每个核心都是一个独立的微处理器。随着多核心处理器的广泛使用，人们将更多的注意力集中到如何利用多核心处理器提高计算能力上来。传统程序都是基于单