一种基于长读数和contig分类的scaffolding方法技术

技术编号:19480655 阅读:434 留言:0更新日期:2018-11-17 10:36
本发明专利技术公开了一种基于长读数和contig分类的scaffolding方法。本方法首先把长读数比对到contig集合上,根据比对结果生成局部scaffold集合。一条局部scaffold是由比对到同一条长读数的contig构成。基于每条contig在局部scaffold中出现的位置信息,把所有的contig分成两类,一类是重复contig,另一类是非重复contig。构建只包含非重复contig的scaffold图,图中每一个节点代表一个非重复contig。接着利用线性规划方法消除scaffold图中的方向和顺序冲突,并使scaffold图中只包含简单路径,其中每条简单路径对应一条scaffold。然后把重复contig插入到scaffold中,形成最终的scaffolding结果。本发明专利技术简单易用,在不同的真实数据上表现出良好的scaffolding结果,较其它scaffolding方法具有更高的准确性和连续性。

【技术实现步骤摘要】
一种基于长读数和contig分类的scaffolding方法
本专利技术涉及生物信息学的序列组装领域,特别是一种基于长读数和contig分类的scaffolding方法。
技术介绍
基因组一般是指全部编码和非编码的脱氧核糖核酸(DNA)序列,它是由四种碱基:腺嘌呤(A)、胸腺嘧啶(T)、胞嘧啶(C)与鸟嘌呤(G)组成的序列,即基因组序列是一个字符串,这个字符串中只包含四个字符A,T,G,C。在实际基因组序列中也包含另一个字符N,代表该位置的碱基无法确定。基因组DNA序列包含了遗传和调控信息,引导生物发育与生命机能运作。在基础生物学研究和众多应用领域中,如诊断、生物技术、法医生物学、生物系统学中,完整和正确的基因组DNA序列已成为不可缺少的知识。通过基因组测序,可以获得大量基因组序列上碱基序列片段(读数或read)。序列组装是通过这些得到的序列片段还原整个基因组DNA序列的方法。而由于重复区、测序错误以及测序不均衡等问题,序列组装方法往往先生成一些比较独立和零散的序列片段,即contig,这些contig可能分布在基因组DNA序列的任意区域,并且由于DNA序列是双链结构,这些contigs可能处在双链上的任意一条链上。scaffolding方法就是确定这些contigs之间的方向和顺序关系,进而产生更长的scaffold。scaffolding会使序列组装结果更加连续和完整,这有助于后续基因识别,基因组比对,结构变异检测等研究,是序列组装研究中的热点之一。目前,以Illumina/Solexa以及AB/SOLid公司为代表的第二代测序技术在显著降低成本的同时,单次运行也能够产生海量和错误率较低的读数。因此,第二代测序技术在国内外得到了广泛的应用。由第二代测序技术得到的双端短读数(pairedreads)是来自一段较长原始基因组序列片段两端的两个序列片段。双端短读数的间距(insertsize)可以达到数千碱基,所以双端短读数能够跨过一段较长的区域并克服序列组装中的部分重复区问题,因此基于双端短读数的scaffolding方法获得了研究人员广泛的关注。其步骤一般是先利用已有的序列组装工具生成contig,然后把双端短读数比对到contig上,再通过比对信息构建scaffold图(scaffoldgraph或者bidierctedgraph),进而推断contigs之间的方向和顺序关系。随着测序技术的迅速发展,速度更快通量更高的第三代测序技术正在逐步完善成熟。第三代测序技术主要有太平洋生物科学公司(PacificBiosciences)的单分子实时测序技术和牛津纳米技术公司(OxfordNanoporeTechnology)的纳米孔单分子技术。第三代测序技术所产生的长读数长度可以达到数万碱基,这些长读数可以跨过基因组中大部分的重复区,进而帮助研究人员获得完整的基因组序列。由于长读数的长度较长,能够跨过大部分重复区,但是长读数的测序错误率较高,一般达到15%左右。由于第二代测序技术比较成熟,并具有测序数据正确率高、成本低和通量高等优势,所以在国内外得到了广泛的应用。虽然第二代测序技术产生的读数比较短,但是测序得到的双端短读数的间距(insertsize)可以达到数千碱基,能够克服部分重复区带来的问题。现有的基于双端短读数的scaffolding方法一般包含以下两个步骤:(1)scaffold图构建。在scaffold图中,一个节点往往代表一个contig,边代表相应的两个contigs在基因组序列上是紧邻的。(2)推断contigs之间的方向和顺序。基于scaffold图的拓扑结构,不同的scaffolding方法采取不同的策略抽取相应的路径。每条路径对应一条scaffold。SCARPA首先只保留双端短读数都是唯一比对的信息。在构建scaffold图时,把每个contig对应两个节点,分别代表该contig的3’端和5’端。SCARPA把确定contig之间的方向问题转化为最小奇数圈遍历问题,采用固定参数算法进行求解,通过删除一些边使contig的关联图中不存在奇数圈。再通过线性规划方法分配给每个contig一个起始坐标,删除不符合距离约束的边。最后基于scaffold图采用启发式方法决定contig之间的顺序关系。SSPACE在构建scaffold图时,以contig为节点,如果两个contigs之间能够匹配上的双端短读数大于一个阈值,则在它们之间添加一条边。SSPACE根据scaffold图的拓扑结构,从最长的contig开始扩展路径,当有多条边供后续扩展时,根据边的权重大小,设定了一种贪婪策略进行扩展。MIP根据连通性对scaffold图进行划分,然后在每个子图上进行scaffolding。针对scaffold图中的每条边,设计一种融合方向关系和顺序关系的综合约束,每条边的权重仍设置为两个contigs之间能够匹配上的双端短读数个数。SOPRA根据scaffold图中的边,设计contig之间的方向约束和距离约束。接着设计一种贪婪的启发式算法,删除最小权重的边集找到一个contig方向分配方案满足所有剩下的方向约束。最后,在确定contig之间的顺序关系时,仍然是删除最小权重的边集找到一个contig位置分配方案满足剩下的距离约束。BESST根据匹配到两个contig之间的双端短读数,计算这两个contig之间距离的理论标准差和实际标准差之间的差异,以及双端短读数在两个contig上的位置分布差异,推断两个contig之间是否应该添加一条边。上述方法中边的权重往往设为两个节点之间能够匹配上的双端短读数个数。BESST在构建的scaffold图中,选出能够匹配最多双端短读数的路径作为最终结果。ScaffMatch把整个方向和顺序推断问题转化为图的最大权重无环二分匹配问题,迭代的删除一些边使该图中不存在环。在最终构建的无环图中,线性化contig节点。GRASS提供了一种混合整数规划方法使contig之间的方向和顺序推断规约到一个单独的优化目标中。SLIQ设计了一组线性不等式用来约束contig之间方向和顺序关系。SILP2利用最大似然模型去除可能错误的匹配信息和发现读数覆盖度异常的区域,并采用整数线性规划方法求解。Briot等人则在scaffold图上采用迭代破圈的方法确定contig之间方向和顺序关系。WiseScaffolder采用了一种需要人工干预的方式优化scaffolding方法。Weller等人利用一种树分解方法对contig之间的方向和顺序关系进行求解。ScaffoldScaffolder将contig之间的方向确定问题转化为双向图中至少包含k条有向边的求解问题,并证明该问题是NP-hard的。然后根据最大生成树算法,设计了一种新的贪婪算法解决该问题,并和其他几种启发式算法在性能上进行比较。Bambus将contig方向推断问题转化为图着色问题。Bambus2通过删除图中最小权重的边,使图中不存在不一致的方向,接着利用最优线性组合模型求解contig之间的方向和顺序关系,并根据输出结果分析基因组变异可能发生的位点。第三代测序技术产生的长读数能够达到数万碱基,因此长读数可以跨过本文档来自技高网...

【技术保护点】
1.一种基于长读数和contig分类的scaffolding方法,其特征在于,包括以下步骤:1)首先将长读数比对到contig集合上,并生成局部scaffold集合;1.1)利用比对工具BWA,把长读数集合比对到contig集合上,生成比对结果。其中只考虑长度大于Lr的长读数和长度大于Lc的contig,Lr=500,Lc=3000。1.2)针对一条长读数,抽取出所有能够比对到它上的contig集合,并计算比对区间位置。如果没有或者只有一条contig比对上,则该条长读数不做后续处理。如果有两条或者更多条contig能够比对上,则根据该条长读数和这些contig之间的比对位置和方向信息,确定这些contig之间的先后顺序和方向,并生成一条局部scaffold。当处理完所有的长读数后,生成一个局部scaffold集合。2)contig分类;如果一条contig出现在两条或者更多条局部scaffold的中间位置(即在一个局部scaffold中,它既不是第一个,也不是最后一个contig),并且在不同的局部scaffold中紧邻它5’端(或者3’端)的contig并不全一样,则该条contig是重复contig。或者一条contig的长度小于MIN,MIN=2000,则也认为该条contig是重复contig。剩下的contig为非重复contig。当处理完所有的contig后,则所有的contig被分成两类:重复contig和非重复contig。3)构建和优化scaffold图;3.1)首先针对每个非重复contig构建一个节点;针对两个非重复contig,判断它们是否同时出现在同一条局部scaffold中,如果可以,则根据比对信息,确定这两个非重复contig之间的方向和顺序信息,并计算它们之间的距离。然后判断它们之间是否能够添加一条边,并确定边的权重。当处理完所有的两两节点后,则一个初始scaffold图构建完成。3.2)scaffold图中每条边约束了其相连接两个节点之间的方向、顺序和距离信息,因此根据scaffold图中所有的边,构建线性规划模型,检测和移除造成方向和顺序冲突的边,保证scaffold图中不存在方向和顺序冲突。3.3)消除冲突后,在scaffold图中,如果存在多个节点同时和某个节点的5’端(或者3’端)相连,则只保留权重最大的边,剩余的边进行移除。通过上述步骤的处理,scaffold图中只包含简单路径。4)生成scaffold集合;scaffold图中的一条简单路径包含了节点的顺序和方向信息,以及相邻节点之间的距离信息,因此每条简单路径对应一个scaffold,并生成一个scaffold集合。针对两个在scaffold中紧邻的非重复contig,如果一条局部scaffold包含它们,并且在该局部scaffold中它们之间包含的全部是重复contig,则这些方向和顺序确定的重复contig是一个插入候选项。如果两条非重复contig出现在多条局部scaffold中,则每条局部scaffold对应一个插入候选项,则选择具有最多频次的插入候选项插入到这两个非重复contig中间。如果最多频次的插入候选项有多个,则这两个非重复contig不进行插入操作。当处理完所有的相邻非重复contig后,生成最终的scaffold集合。其中contig为基因组序列片段;scaffold为基因组超长序列片段;scaffolding方法是指用于确定各条contig的方向,以及它们在基因组序列上的先后顺序,从而产生一些基因组超长序列片段,即scaffold的方法。一条序列的最左端是其5’端,最右端是其3’端。...

【技术特征摘要】
1.一种基于长读数和contig分类的scaffolding方法,其特征在于,包括以下步骤:1)首先将长读数比对到contig集合上,并生成局部scaffold集合;1.1)利用比对工具BWA,把长读数集合比对到contig集合上,生成比对结果。其中只考虑长度大于Lr的长读数和长度大于Lc的contig,Lr=500,Lc=3000。1.2)针对一条长读数,抽取出所有能够比对到它上的contig集合,并计算比对区间位置。如果没有或者只有一条contig比对上,则该条长读数不做后续处理。如果有两条或者更多条contig能够比对上,则根据该条长读数和这些contig之间的比对位置和方向信息,确定这些contig之间的先后顺序和方向,并生成一条局部scaffold。当处理完所有的长读数后,生成一个局部scaffold集合。2)contig分类;如果一条contig出现在两条或者更多条局部scaffold的中间位置(即在一个局部scaffold中,它既不是第一个,也不是最后一个contig),并且在不同的局部scaffold中紧邻它5’端(或者3’端)的contig并不全一样,则该条contig是重复contig。或者一条contig的长度小于MIN,MIN=2000,则也认为该条contig是重复contig。剩下的contig为非重复contig。当处理完所有的contig后,则所有的contig被分成两类:重复contig和非重复contig。3)构建和优化scaffold图;3.1)首先针对每个非重复contig构建一个节点;针对两个非重复contig,判断它们是否同时出现在同一条局部scaffold中,如果可以,则根据比对信息,确定这两个非重复contig之间的方向和顺序信息,并计算它们之间的距离。然后判断它们之间是否能够添加一条边,并确定边的权重。当处理完所有的两两节点后,则一个初始scaffold图构建完成。3.2)scaffold图中每条边约束了其相连接两个节点之间的方向、顺序和距离信息,因此根据scaffold图中所有的边,构建线性规划模型,检测和移除造成方向和顺序冲突的边,保证scaffold图中不存在方向和顺序冲突。3.3)消除冲突后,在scaffold图中,如果存在多个节点同时和某个节点的5’端(或者3’端)相连,则只保留权重最大的边,剩余的边进行移除。通过上述步骤的处理,scaffold图中只包含简单路径。4)生成scaffold集合;scaffold图中的一条简单路径包含了节点的顺序和方向信息,以及相邻节点之间的距离信息,因此每条简单路径对应一个scaffold,并生成一个scaffold集合。针对两个在scaffold中紧邻的非重复contig,如果一条局部scaffold包含它们,并且在该局部scaffold中它们之间包含的全部是重复contig,则这些方向和顺序确定的重复contig是一个插入候选项。如果两条非重复contig出现在多条局部scaffold中,则每条局部scaffold对应一个插入候选项,则选择具有最多频次的插入候选项插入到这两个非重复contig中间。如果最多频次的插入候选项有多个,则这两个非重复contig不进行插入操作。当处理完所有的相邻非重复contig后,生成最终的scaffold集合。其中contig为基因组序列片段;scaffold为基因组超长序列片段;scaffolding方法是指用于确定各条contig的方向,以及它们在基因组序列上的先后顺序,从而产生一些基因组超长序列片段,即scaffold的方法。一条序列的最左端是其5’端,最右端是其3’端。2.根据权利要求1所述的基于长读数和contig分类的scaffolding方法,其特征在于,所述步骤1.2)具体为:1.2.1)根据比对工具BWA产生的比对结果,如果一条长读数和一个contig能够比对上,可以得到其比对区间的位置和比对方向。假设第j条长读数(lrj)和第i个contig(ci)能够比对上,并且在lrj上的比对区间是[SPR(ci,lrj),EPR(ci,lrj)],在ci上的比对区间是[SPC(ci,lrj),EPC(ci,lrj)]。由于长读数的测序错误率比较高,所以比对工具给出的比对区间位置往往具有一些偏差。本方法采用下述步骤对比对区间进行修正。如果SPR(ci,lrj)<SPC(ci,lrj),则SPC’(ci,lrj)=SPC(ci,lrj)-SPR(ci,lrj),SPR’(ci,lrj)=0。如果SPR(ci,lrj)>=SPC(ci,lrj),则SPR’(ci,lrj)=SPR(ci,lrj)-SPC(ci,lrj),SPC’(ci,lrj)=0。如果LEN(lrj)-EPR(ci,lrj)>LEN(ci)-EPC(ci,lrj),则EPC’(ci,lrj)=LEN(ci)-1,EPR’(ci,lrj)=EPR(ci,lrj)+LEN(ci)-EPC(ci,lrj)。如果LEN(lrj)-EPR(ci,lrj)<=LEN(ci)-EPC(ci,lrj),则EPC’(ci,lrj)=EPC(ci,lrj)+LEN(lrj)-EPR(ci,lrj),EPR’(ci,lrj)=LEN(lrj)-1。如果|SPR(ci,lrj)-SPR’(ci,lrj)|<α,或者|SPC(ci,lrj)-SPC’(ci,lrj)|<α,或者|EPC(ci,lrj)-EPC’(ci,lrj)|<α,或者|EPR(ci,lrj)-EPR’(ci,lrj)|<α,则认lrj和ci无法比对上,本方法忽略该比对,其中α是一个参数(α=500)。经过上述修正后,本方法可以得到长读数和contig之间的比对区间,以及比对方向。1.2.2)如果一条长读数上只能够比对上一条contig,则该条长读数不做后续处理。如果有两个或者更多个contig能够比对到该条长读数上,则根据这些contig在该长读数上比对区间的起始位置,由小到大对它们进行排序。如果有多条contig能够同时比对到该条长读数的5’端(或者3’端),则只保留一条具有最长比对区间的conti...

【专利技术属性】
技术研发人员:罗军伟王俊峰张波张霄宏贾利琴
申请(专利权)人:河南理工大学
类型:发明
国别省市:河南,41

网友询问留言 已有0条评论
  • 还没有人留言评论。发表了对其他浏览者有用的留言会获得科技券。

1