《生物信息学序列拼接精品文稿.ppt》由会员分享,可在线阅读,更多相关《生物信息学序列拼接精品文稿.ppt(56页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、生物信息学 序列拼接第1页,本讲稿共56页序列拼接序列拼接 序列拼接任务即将测序生成的reads短片段拼接起来,恢复出原始的序列。该问题是序列分析的最基本任务,是基因组研究成功与失败的关键,拼接结果直接影响到序列标注,基因预测、基因组比较等后续任务。基因组序列的拼接也是基因组研究必须解决的首要难题。其困难不仅来自它的海量数据(以人类基因组序列为例,从数量为10兆级的片断恢复出长度为亿级的原始序列),而且源于它含有高度重复的序列。第2页,本讲稿共56页第3页,本讲稿共56页拼接问题的难点拼接问题的难点 DNA测序数据有其固有的四个的特点,他们也正是解决实际的序列拼接问题的难点所在:1测序有误差
2、2不完全覆盖性 3序列所在链不确定 4重复序列的干扰第4页,本讲稿共56页1测序有误差测序有误差 由于测序技术的局限,难免会出现测序错误,尤其是在序列的末端,一般错误率可控制在1以下。所以对每个碱基一般有一个正确概率,以质量打分的形式给出。因此每个ri都有个可信度。而read与read之间有不同程度的重叠,由此导致有的重叠可信度高,有的重叠可信度低。第5页,本讲稿共56页2不完全覆盖性不完全覆盖性 不是所有的碱基被测序的次数都等于平均测序覆盖度。极端的情况,可能会出现源基因组序列上部分区域未被测序的情况(这段区域称为gap)。即,测序的reads集合不是原始基因组序列一个完整覆盖。此时需要借助
3、于各种图谱如:基因组指纹图谱(genome fingerprint map),基因组级物理图谱(genome-wide physical map),细胞发生图谱(cytogenetic maps)等协助对reads进行定位第6页,本讲稿共56页3序列所在链不确定序列所在链不确定 由于测序过程中无法确定特定片断属于DNA双链中的哪一条链上,所以我们在拼接过程中并不清楚使用的是read的正义链,还是其互补链。4重复序列的干扰重复序列的干扰 DNA序列自身含有高度重复的子序列,它们一种表现为短序列的串级重复,比如:(GGAA)n。或AmTn等。另一种表现为大量相似序列(其拷贝数可达几十万)散布在基因
4、组的各个地方。Repeat的存在,将导致fragments间overlap的不真实性,进而产生错拼的结果。因此在拼接过程中耍确定这些序列的形式及大小,才能保证以高概率恢复出其在原始真实序列中的位置第7页,本讲稿共56页拼接算法评价拼接算法评价 以上拼接问题的四个难点不仅极大的增加了解决实际拼接问题的难度,而且从某种程度上说无法完整地恢复出原始DNA序列来。即实际上仅能构建出若干个contig(重建的fragments的一种排列形式,它覆盖基因组上一段连续区域)这些contig将指导测序项目finishing阶段的实验方法最终构建DNA完整序列。第8页,本讲稿共56页 目前,国际上对拼接软件的公
5、认评价标准包括两方面,即重建出的contig的数目和准确度。我们发展的基因组序列拼接新算法的目标是在确保准确性的前提下,构建尽量少的contig,以减少测序后期大量的人力和财力的投入。第9页,本讲稿共56页基因组序列拼接算法研究现状基因组序列拼接算法研究现状 现在最常用的拼接程序使用的拼接算法可分成两类,一类是将拼接问题转化为在图中寻找的Hamilton路径的问题;另一类是将拼接问题在某种特殊情况下转化成寻求图中的Euler路径的问题。他们均有其成功的典型算法。第10页,本讲稿共56页1.转化为转化为Hamilton Path问题问题 每个DNA片段(read)相当于图中一个结点,如果两个片段
6、之间存在着重叠(overlap)关系,则在两个结点之间定义一条边,而沿着DNA原始序列从头到尾,则必然经过每个结点一次且仅一次,即是一条Hamilton路径。一条contig表示图中一条简单路,此类算法以Phrap,TIGR Assembler,CAP3,GigAssemble等为代表。第11页,本讲稿共56页 他们都是遵循“overlap-layout-consensus”的框架。首先,为了构建图。计算任意两个read间可能的比对情况。其次,通过去除歧义的或者不确信的边得到较为准确的图,并在其上寻找非交叉的简单路的集合,该集合对应于contig的集合。最终,通过对包含在一个简单路上的所有re
7、ad进行多序列比对,为每一个contig构建一个一致性序列(consensus sequence)。第12页,本讲稿共56页2.转化为转化为Euler Path问题问题 EULER是这类算法的代表。与传统方法沿着“OverlapLayoutConsensus”路线不同,它不计算各个read之间的Overlap,即没有Overlap步骤。它的大致想法如下:为了排除read中的错误,获得Error-Free的read,将所有的read切割成小片n-mers。第13页,本讲稿共56页 将每个read和Gk的近似进行比对,寻求read的最小改变能够使得read的所有n-mers包含在Gk的近似集合中。
8、从而构建了高质量序列,而对于Poor read,直接抛弃,对Chimeric read(两端在n-mers中但整体不在的reads)进行特殊处理。第14页,本讲稿共56页 初始的想法是要实现去除reads中的测序错误的目的,如果知道原始序列G,那么直接使用测序获得的read和G进行比较即可。但是实际上G并不可知,那么退而求其次,G的序列片断Gk亦可,事实上Gk亦不可知。所以将所有的read切割成小片n-mers,所有Solid的n-mers形成的集合称为Gk的近似。最后,构造De Bruijn图。第15页,本讲稿共56页现有算法的主要问题现有算法的主要问题 虽然已经开发了以上的算法,基因组序列
9、拼接问题尚未彻底解决,以上两类算法都存在着各自的缺陷。第16页,本讲稿共56页 对于第一类算法来说,实际上是在图中寻找一条使得评价函数值最优的Hamilton路径,这是一个NP完全问题。一般都采用greedy-merging的算法近似求解。由于这种step-by-step的局部贪心算法,其明显的局部特性忽略了reads间“长距离”或者整体性的联系,从而导致了拼接错误,即拼接结果和真实的DNA原始序列不同。最近研究指出,在对已知序列的流行性感冒嗜血杆菌基因组的拼接过程中,无论是Phrap,TIGR Assembler,还是CAP3,都发生了拼接错误的现象。第17页,本讲稿共56页 对于第二类算法
10、来说,它只能在特殊的情况下,才能将问题简化成寻求一条Euler路径,最终的结果是从多条候选的Euler超路中选择出来的。EULER算法依然存在拼接错误,且结果选择的过程没有理论依据。EULER软件在实际数据集的运行速度上和第一类算法相当。更重要的是,EULER采用的算法过于独立,很难利用其他辅助生物信息,导致其实用性和流行性大打折扣。第18页,本讲稿共56页局部搜索局部搜索(Local Search)方法方法胡杰胡杰第19页,本讲稿共56页 将局部搜索方法运用于一个具体的问将局部搜索方法运用于一个具体的问题,需要对以下四项进行明确的定义。题,需要对以下四项进行明确的定义。1将原同题表示成一个优
11、化问题,即定义一个可行域以及在该可行域上的一个目标函数。2定义可行域中的邻域结构,即说明满足什么条件的两个点相邻。第20页,本讲稿共56页3确定在邻域中的搜索方式。4局部极值点的处理。如果当前解点邻域中的所有点的目标函数值都比当前点大,那么这点就称为局部极值点。在一些问题中,局部极值点就是全局最优点。而对另一些问题而言,局部极值点已经满足求解实际问题的需求。第21页,本讲稿共56页基于局部搜索的序列拼接算法框架基于局部搜索的序列拼接算法框架 主要目标主要目标是在layout阶段尽可能准确的前提下,获得更长的contig。具体的,使用局部搜索算法求得数据集上近似全局的最短超串;然后根据求得的最短
12、公共超串对应的fragment的排列关系为基础获得“consensus segment”第22页,本讲稿共56页1.Overlap定义定义 如果一个串的前缀是另一个串的后缀则认为这两个串之间存在overlap,并根据over-lap构建超串。对给定的串f和g,存在多个可能的overlap关系 比如说,若f=ACTGGGAGCAGC,g=AGCAGCTTTTACT,那么他们之间至少存在两个overlap形式。第23页,本讲稿共56页第24页,本讲稿共56页 在我们的算法中,仅考虑两个串之间最大的overlap情况,并定义overlap(f,g)表示f和g之间存在的多个overlap关系中最长的一
13、个overlap所包含的字符个数。在上面的例子中overlap(f,g)=6。如果f和g之间overlap区域长度小于M(M是一个足够小的正整数),则overlap(f,g)=0。第25页,本讲稿共56页2.优化目标定义优化目标定义 我们对reads集合S中的每个元素按照其左端在超串t中首次出现的位置进行排序,并沿超串t从左向右依次读入,并将其对应为序S=sl,s2,Sn)。我们用P(S)表示这个由字符串为元素构成的序。在序在序P(S)中,对于每一个连续的字符串元素对中,对于每一个连续的字符串元素对si,si+1,都存在都存在overlap(si,Si+1)。第26页,本讲稿共56页 因此,字
14、符串的一个排列等价于一个超串t及作用在其上的函数overlap,超串t的长度length(t)=t =ssIs|-overlap(si,Si+1),为了求超串t具有最小的长度。因为在给定集合s中ssIs是确定的值,我们可以进一步把问题转化成寻找一个排列P使得 最大化。第27页,本讲稿共56页3.邻域定义邻域定义 在使用局部搜索方法前,我们首先定义排列的邻域这一概念。我们称一个排列为问题的一个解。设reads集合S是一个具有n个元素的字符串集合,P(S)是reads集合S所有解的集合。定义操作定义操作rshift(i,j),即,即对一个解对一个解PP(S),将,将P第第i位置的元素向右位置的元素
15、向右移至第移至第j个位置个位置(1ijn)。第28页,本讲稿共56页也就是说,如果我们应用这一操作于给定的解P=sl,s2,sn),则结果P=rshift(i,j)(P)定义为;第29页,本讲稿共56页 类似的,定义操作lshift(i,j)将P第i位置的元素向左移至第j个位置(1jin),定义P=lshift(i,j)(P)为;第30页,本讲稿共56页 这样,一个解PP(s)的邻域可以定义为N(P)=rshift(i,j)(P)|1ijn)lshift(i,j)(P)|1|1j0时,p比p更好,则由p转移到p第33页,本讲稿共56页 此外,对给定的解P=s1,s2Sn,内部的元素均包含前趋和
16、后继。但是,头元 素仅有后继而尾元素仅有前驱。在算法中,为了消除头尾元素差异,我们将排列的头尾 元素连结起来并使用局部搜索方法寻找最优的“Loop”超串。接着在overlap最小的地方切分该环状超串,最终还原成一条线性超串。第34页,本讲稿共56页5.局部极值点的处理局部极值点的处理 当经过以read为单元的搜索后,可以获得一个当前邻域内的局部最优解1,k1,k1+1,k2,k2+1,k3km。它对应集合s上的一个superstring。该解对应的reads的排列中,任意相邻两个reads间的overlop关系薄弱,即overlap(i,i+1)M(其中,M是一个足够小的正整数,1i n)。第
17、35页,本讲稿共56页例如:如果overlap(rkl,rkI+1)M,overlap(rk2,rk2+1)M,,overlap(rkm,rkm+1)M1,overlap(r1,r3)M1,但overlap(r2,r3)M2和overlap(r3,r2)50k)Contigs平均长度平均长度错误错误contigs数数目目错误错误contigs总长总长(bp)AE0078725428697.5CAP35414652132876300PHRAP5405881932845200B-LSA5410821442864800AE00065715513357.5CAP31539089777199886144
18、539PHRAP1528509717215003104164B-LSA15459535011309192101102INRUENSS18301388CAP317957526812264085227941PHRAP17986456314287484188705B-LSA181116944124116215219AE00078221784007.5CAP3216969310311215096279380PHRAP21506509411223814212060B-LSA2171203751328949339270AE00786928415817.5CAP32826991126142243613627
19、4PHRAP28155571161524272229589B-LSA283368499203148500AL00912642148147.5CAP3416548417110243586138549PHRAP413779915519256955168824B-LSA419066712321314849122422CAP3、PHRAP、B-LAS的性能比较的性能比较第44页,本讲稿共56页结果一contigs的个数contig平均长度B-LSA40833.661PHRAP53325.483CAP356623.035-结果:结果:产生更少的产生更少的contig,而且,而且contig更长更长有效的
20、减少了测序项目后续工作的工作量有效的减少了测序项目后续工作的工作量第45页,本讲稿共56页结果二总contigs数错误数目错拼率B-LSA408153.67CAP3745526.97PHRAP566244.24-结果:结果:包含更少的错误拼写现象包含更少的错误拼写现象第46页,本讲稿共56页 基于原型算法框架实现的基于原型算法框架实现的Basic LSA系系统在针对实际序列模拟的统在针对实际序列模拟的shotgun reads数数据集上的拼接结果表明:据集上的拼接结果表明:在不引入其他辅助信息的前提下,在不引入其他辅助信息的前提下,Basic LSA比比PHRAP(1999)及及CAP3都有都
21、有明显的改进这也证明了原型算法的有效明显的改进这也证明了原型算法的有效性和可行性性和可行性第47页,本讲稿共56页基于原型算法的优化基于原型算法的优化 原型算法的运行速度过慢,以至成为原型算法的运行速度过慢,以至成为整个系统运行的整个系统运行的瓶颈瓶颈,从而导致该系统无,从而导致该系统无法满足实验和可能的进一步应用的需要法满足实验和可能的进一步应用的需要这也成了将原型系统实用化必须解决的问这也成了将原型系统实用化必须解决的问题题.提速优化的目标是在不改变运行结果的提速优化的目标是在不改变运行结果的前提下,提高系统运行速度前提下,提高系统运行速度第48页,本讲稿共56页“邻域剪枝邻域剪枝”-(N
22、eighborhood Pruning)尽可能通过排除那些目标函数值低于尽可能通过排除那些目标函数值低于当前值的邻接点来缩小搜索空间这里的当前值的邻接点来缩小搜索空间这里的关键在于找到了一个关于能够增加目标函关键在于找到了一个关于能够增加目标函数值的变换的必要条件数值的变换的必要条件第49页,本讲稿共56页After transposition“领域剪枝领域剪枝”(Neighborhood Pruning)方法)方法Befor transposition第50页,本讲稿共56页 在搜索过程中,每一个在搜索过程中,每一个transposition将导将导致目标函数改变致目标函数改变overlap
23、,Aoverlap=Overlap4+Overlap5+Overlap6 Overlapl-Overlap2-Overlap3 若该若该transposition操作能改进当前解,必定操作能改进当前解,必定有有Aoverlap0 因此,也必然满足以下条件因此,也必然满足以下条件 (Overlap4-Overlapl-Overlap2)+(Overlap5+Overlap6)Overlap30第51页,本讲稿共56页两种可能性:两种可能性:1Overlap4Overlapl+Overlap2第52页,本讲稿共56页 基于上述基于上述Neighborhood Pruning方法,方法,我们优化了原
24、型算法,由于局部搜索我们优化了原型算法,由于局部搜索(Local Search)是一种启发式的算法,其算法复杂是一种启发式的算法,其算法复杂度难以精确估计,随着真实数据特性的不同度难以精确估计,随着真实数据特性的不同存在显著的差异存在显著的差异.我们算法优化后的改进程度也会有所波我们算法优化后的改进程度也会有所波动。我们选择了三个具有代表性的数据集,动。我们选择了三个具有代表性的数据集,描述在不改变算法运行结果的前提下,提速描述在不改变算法运行结果的前提下,提速优化前后的运行时间的比较优化前后的运行时间的比较第53页,本讲稿共56页原始序列原始序列序列长度序列长度覆盖度覆盖度系统系统Overl
25、ap阶段用时阶段用时Layout阶阶段用时段用时Conmmous阶段阶段用时用时系统总用系统总用时时LSA layout阶阶段优化后段优化后Artificial3M8PHRAP30047901390718B-LSA5571350643414497Im-LSA692195821293BA0000044.2M8PHRAP1379387182135191B-LSA138182346646898Im-LSA1899435171964B42.76M23PHRAP1073075728021428989B-LSA108041693852788182975Im-LSA108611903282115585LSA提速优化前后性能比较提速优化前后性能比较第54页,本讲稿共56页 经过提速优化之后,经过提速优化之后,Layout阶段不再阶段不再是整个系统的时间瓶颈因此整个系统的是整个系统的时间瓶颈因此整个系统的运行速度与运行速度与Phrap大致相当,从而朝着实大致相当,从而朝着实现全基因组序列拼接系统的终极目标走出现全基因组序列拼接系统的终极目标走出了关键一步了关键一步第55页,本讲稿共56页Thank You!第56页,本讲稿共56页