用于组装基因组序列的方法、系统及装置的制造方法
【专利摘要】本申请提供一种利用锚点序列对待组装序列进行压缩并进行基因序列组装的方法,包括获取一待组装序列集的数据和多个锚点序列的数据;利用所述多个锚点序列,对所述待组装序列集中的每一个待组装序列进行压缩以获得该待组装序列的索引表,其中所述索引表以锚点序列的标号的序列来表示;以及利用所述待组装序列集中各待组装序列压缩得到的所述索引表进行比对并构建组装骨架。
【专利说明】
用于组装基因组序列的方法、系统及装置
技术领域
[0001] 本申请设及基因测序领域,具体而言,设及组装基因序列的方法、系统及装置,其 可W提高基因序列组装的速度及精准度,尤其是基因组序列组装的速度及精准度。
【背景技术】
[0002] 人类基因组计划是迄今为止最大的生物医药项目之一,其大大加快了 DNA测序技 术的发展。在过去的Ξ十年中,已经发展了Ξ代DNA测序技术,目前正处于在第二代测序技 术和第Ξ代测序技术更替的十字路口。
[0003] 其中,第一代Sanger测序技术所产生的序列长度一般小于80化P。根据第 一代测序技术,利用重叠图(Overlapping Graph)在测序得到的序列集之中,捜索共 有区域W组装基因序列,但是运种基于重叠图的捜索和比对算法在时间和存储上耗费 较高。而第二代测序技术(也称为下一代测序技术,NG巧所产生的序列较短,一般在 30-300bp的范围内。在第二代测序技术的情况下,利用基于核巧k联体化-mer)的德 布鲁因图(de Bruijn Graph,DB(i)算法(Conway, T.C.&Bromage,A.J. Succinct data structures for assembling large genomes. Bioinformatics 27,479-486(2011); Melsted,P. &Pritchard,J. Efficient counting of k_mers in DNA sequences using a bloom filter,BMC Bio informatics 12, 333 (2011) ; Ye,C.,Ma,Z. S.,Cannon, C. H.,Pop,M.&Yu, D. W. Exploiting sparseness in de novo genome assembly, BMC Bioinformatics 13 Suppl 6,S1(2012))对测序得到的序列集进行组装D与重叠图相 关联的时间和存储问题在处理第一代测序技术生成的相对少量的序列的时候还可^承 受,但在处理由第二代测序仪生成的巨大数量的短序列时,就成为很大的问题。虽然采 用DBG模式避免了成对对比的问题,但还是不得不面对巨大的内存占用的问题。在DBG 成为标准方法之前,原本用于组装第一代的序列的基于重叠图模型或重叠-排布-整 合(Overlap-Layout-Consensus,0LC)的软件,如 Celera(Venter, J. C. et al. The sequence of the human genome. Science 291,1304-1351 (2001)),AMOS(Treangen, T. J.,Sommer, D. D.,Angly, F. E.,Koren, S. &Pop, M. Next generation sequence assembly with AMOS, Current protocols in bioinformatics Chapter 11,Unit 11.18(2011)) 和 ARACHNE(Batzoglou, S. et al. ARACHNE: a whole-genome shotgun assembler. Genome research 12, 177-189(2002)),也被应用于第二代测序技术所得序列的组装。字串图认及 最优重叠图是0LC方法的特定形式,其通过简化整体重叠图谱在组装长序列方面非常有 效。
[0004] 目前,已经开发了第兰代测序技术(其也被称为单分子测序技术HLevene,M. J.et al.Zero-mode waveguides for single-molecule analysis at high concentrat ions,Science299,682-686(2003);加,L·Q·,C?eley,S·&Bayley,比Cap化reofasinglΘ molecule in a nanocavity,Science291, 636-640 (2001)),从而可生成的最长序列范围为 500bp-120化P之间。第兰代测序技术最大的优点是能够实现高通量长序列测序。一方面, 增加测序长度具有协助辨别基因序列重复区域的显著优势,从而有助于人类能更加直接并 有效地探索生物学问题,例如基因分型和突变等问题。然而,另一方面,增加测序长度就会 提高错误率,使得在错误的检测和修正方面存在极大问题。 阳0化]第Ξ代测序的高错误率使得最初开发用于相对准确的第一代测序技术的基于重 叠组装的OLC方法不适用。同样,第Ξ代测序得到的带错长序列集使得最初设计用于第二 代测序技术的短序列组装的DBG方法也不再适用。
[0006] 因此,在现阶段,基因测序技术所面临的很大障碍是缺少一种有效的基因序列组 装方法。
【发明内容】
[0007] 因此,本申请提供了用于组装基因序列的方法、系统及装置,其适用于任意长度测 序序列,尤其适用于组装通过第Ξ代测序技术测得的长基因序列,从而能够大大降低基因 序列组装在计算时间和内存上的高度复杂性,并且提高基因序列组装的速度和精度。
[0008] 本申请的一个方面提供了一种用于组装基因序列的方法,包括:获取一待组装序 列集的数据和多个错点序列的数据;利用所述多个错点序列,对所述待组装序列集中的每 一个待组装序列进行压缩W获得该待组装序列的索引表,其中所述索引表W错点序列的标 号来表示;W及利用所述待组装序列集中各待组装序列压缩得到的所述索引表构建组装骨 架。
[0009] 在一些实施方式中,上述压缩W获得待组装序列的索引表包括:确定每个错点序 列在待组装序列中的排列位置,并按照错点序列在待组装序列上的排列位置,获得对应于 该待组装序列的索引表。在一些实施方式中,上述确定每个错点序列在待组装序列中的排 列位置包括:将所述错点序列和每一个待组装序列分别打断成一系列核巧k联体,将所述 错点序列的核巧k联体分别与所述待组装序列的核巧k联体进行比对,当所述待组装序列 的核巧k联体中与错点序列的核巧k联体相匹配的个数超过预定阔值,则认为所述错点序 列可W匹配至该待组装序列;通过每个所述匹配的核巧k联体分别在所述错点序列及所述 待组装序列中的位置的对应关系,计算得到所述匹配的错点序列在待组装序列中的排列位 置。对于有多个匹配的核巧k联体的错点序列,所述错点序列在待组装序列中的排列位置 可W是根据每个匹配的核巧k联体计算得到的所述错点序列在待组装序列中的排列位置 的平均值。在一些实施方式中,所述k-mer的长度为:巧bp-5化P。在一些实施方式中,所述 阔值为所述错点序列的长度的0. 1 %至2%。
[0010] 在一些实施方式中,上述构建组装骨架包括:将所述待组装序列集中各待组装序 列压缩得到的所述索引表进行一一比对,筛选出具有最优重叠的待组装序列子集,W构建 待组装序列间的最优重叠图谱,然后根据所述子集中的待组装序列索引表中错点序列标号 的位置对应关系对所述子集中的各待组装序列进行排布,从而形成组装骨架。在一些实施 方式中,上述构建组装骨架包括:选取所述待组装序列集中的一条序列作为当前序列,将 所述当前序列的索引表与待组装序列集中的所有其它序列的索引表进行一一比对,找出与 所述当前序列具有重叠关系的备选序列,将备选序列与当前序列进行索引表比对或碱基比 对W找出与当前序列具有最优重叠和延伸的序列,并将其相对于当前序列延伸的部分加到 当前序列上即获得当前序列的扩展序列,通过将所述待组装序列集中的每条序列作为当前 序列重复上述步骤w得到一组扩展序列,然后根据所述扩展序列中错点序列标号的位置对 应关系对各扩展序列进行排布,从而形成组装骨架。在一些实施方式中,本申请提供的用 于组装基因序列的方法进一步包括:在获得每条待组装序列的索引表后,利用错点序列对 待组装序列的索引表进行反向检索,得到与所述错点序列相对应的待组装序列的索引表子 集,通过将所述索引表子集中所有其他序列的索引表比对到当前待组装序列的索引表上来 修正所述待组装序列的索引表子集中的待组装序列的索引表。在一些实施方式中,上述通 过比对来修正被压缩的待组装序列包括:对待组装序列的索引表子集中每条待组装序列的 索引表,将其作为当前序列索引表与所述索引表子集中所有其他待组装序列的索引表进行 一一比对,根据每个位置上不同的错点序列标号出现的频次进行投票,利用每个位置的投 票结果对每条当前序列的索引表进行修正,从而获得修正后的待组装序列的索引表。在一 些实施方式中,上述对待组装序列的索引表进行修正包括:将比对结果中在待组装序列的 索引表子集中出现不多于5次、或不多于4次、或不多于3次、或不多于2次、或不多于1次 的错点从包含其的待组装序列的索引表中删除;将存在断点的杂合序列删除;W及将比对 结果中的被包含的序列删除。在一些实施方式中,构建组装骨架是在对待组装序列的索引 表进行修正后进行的,其包括:将所有上述修正后的待组装序列的索引表进行一一比对,筛 选出具有最优重叠的修正后的待组装序列的子集,然后根据所述子集中的修正后的待组装 序列索引表中错点序列标号的位置对应关系对待组装序列进行排布,从而形成组装骨架。 在一些实施方式中,本申请提供的用于组装基因序列的方法进一步包括对构建的组装骨架 进行如下所述的再修正:在构建组装骨架后,将所有待组装序列与组装骨架进行比对,对于 组装骨架上的任一位点,通过计算所有可比对至该位点的各条序列中在该位点的每种碱基 出现的频次进行投票,利用每个位点的投票结果对组装骨架进行再修正,从而得到再修正 后的基因序列组装结果。在一些实施方式中,上述投票是通过利用跳跃式核巧k联体实现 的。在一些实施方式中,上述对组装骨架进行再修正包括:获得覆盖待组装基因序列不同区 域的组装骨架,利用多个计算机处理器同时分别修正所述覆盖待组装基因序列不同区域的 组装骨架,W及将再修正后的组装骨架进行整合得到完整的基因序列组装结果。
[0011] 在一些实施方式中,本申请所述的待组装序列集包括通过第二代基因组测序 技术测序得到的序列。在一些实施方式中,本申请所述的待组装序列集中的序列长度 在l(K)bp-120化P之间、或者在5(K)bp-120化P之间、或者在5(K)bp-100化P之间、或者在 化bp-100化P之间、或者在l(K)bp-化bp之间、或者在l(K)bp-300bp之间。在一些实施方式 中,本申请所述的待组装序列集包括通过第Ξ代基因组测序技术所得到的序列,并且所述 待组装序列集中的序列长度在5(K)bp-100化P之间。在一些实施方式中,本申请所述的错点 序列的长度在化p-30bp之间、或者在30bp-100bp之间、或者在l(K)bp-300bp之间、或者在 30bp-100化P之间、或者在3(K)bp-100化P之间。在一些实施方式中,本申请所述的错点序 列选自下组:精度大于95%的序列长度在l(K)bp-300bp之间的短序列、基因组高保守区域 内的精度大于95%的序列长度在10化口-3006口之间的短序列和待组装序列集中序列长度 在化bp-120化P之间的长序列的局部截取。在一些实施方式中,本申请所述的错点序列的 长度在2-30bp之间。
[0012] 本申请的另一个方面提供了一种用于组装基因序列的装置,包括:获取模块,所述 获取模块用于获取一待组装序列集的数据和多个错点序列的数据;压缩模块,所述压缩模 块用于利用所述多个错点序列,对所述待组装序列集中的每一个待组装序列进行压缩w获 得该待组装序列的索引表,其中所述索引表W错点序列的标号来表示;W及组装模块,所述 组装模块用于利用所述待组装序列集中各待组装序列压缩得到的所述索引表构建组装骨 架。
[0013] 在一些实施方式中,本申请提供的用于组装基因序列的装置,进一步包括:待组装 序列的索引表修正模块,所述待组装序列的索引表修正模块用于利用所述错点序列对待组 装序列的索引表进行反向检索,得到与所述错点序列相对应的待组装序列的索引表子集, 通过将所述索引表子集中所有其他序列的索引表比对到当前待组装序列的索引表上来修 正所述待组装序列的索引表。
[0014] 在一些实施方式中,本申请提供的用于组装基因序列的装置,进一步包括:基因序 列组装结果修正模块,所述基因序列组装结果修正模块用于对构建的组装骨架进行如下所 述的再修正:在构建组装骨架后,将所有待组装序列与组装骨架进行比对,对于组装骨架上 的任一位点,通过计算所有可比对至该位点的各条序列中在该位点的每种碱基出现的频次 进行投票,利用每个位点的投票结果对组装骨架进行再修正,从而得到再修正后的基因序 列组装结果。
[0015] 在一些实施方式中,本申请提供的用于组装基因序列的装置,进一步包括:基因序 列组装结果整合模块,所述基因序列组装结果整合模块用于获得覆盖待组装基因序列不同 区域的组装骨架,利用多个计算机处理器同时分别修正所述覆盖待组装基因序列不同区域 的组装骨架,W及将再修正后的组装骨架进行整合得到完整的基因序列组装结果。
[0016] 本申请的再一个方面提供了一种用于组装基因序列的系统,包括:输入组件,其用 于获取一待组装序列集的数据和多个错点序列的数据;存储器,其用于存储所述待组装序 列集的数据和多个错点序列的数据;W及处理器,其被配置于;利用所述多个错点序列,对 所述待组装序列集中的每一个待组装序列进行压缩W获得该待组装序列的索引表,其中所 述索引表W错点序列的标号来表示;W及利用所述待组装序列集中各待组装序列压缩得到 的所述索引表构建组装骨架。
[0017] 在一些实施方式中,本申请提供的用于组装基因序列的系统中,处理器被进一步 配置于:在获得所述索引表后,利用所述错点序列对待组装序列的索引表进行反向检索,得 到与所述错点序列相对应的待组装序列的索引表子集,通过将所述索引表子集中所有其他 序列的索引表比对到当前待组装序列的索引表上来修正所述待组装序列的索引表子集中 每条待组装序列的索引表。
[0018] 在一些实施方式中,本申请提供的用于组装基因序列的系统中,处理器被进一步 配置于对构建的组装骨架进行如下所述的再修正:在构建组装骨架后,将所有待组装序列 与组装骨架进行比对,对于组装骨架上的任一位点,通过计算所有可比对至该位点的各条 序列中在该位点的每种碱基出现的频次进行投票,利用每个位点的投票结果对组装骨架进 行再修正,从而得到再修正后的基因序列组装结果。
[0019] 在一些实施方式中,本申请提供的用于组装基因序列的系统中,处理器对组装骨 架进行的再修正包括:获得覆盖待组装基因序列不同区域的组装骨架,利用多个计算机处 理器同时分别修正所述覆盖待组装基因序列不同区域的组装骨架,W及将再修正后的组装 骨架进行整合得到完整的基因序列组装结果。
[0020] 相应地,通过应用本申请所述的方法、系统和装置,由高精度测序短序列或者高精 度测序短序列的组装片段组成错点序列,并利用错点序列对待组装序列进行标记,并将待 组装序列压缩成对应的索引表,然后再进行比对和组装,从而大大加快了组装速度并减少 了计算量。
[0021] 并且,考虑到通过第Ξ代测序技术获得的待组装序列的错误率比较高,因此进一 步对被压缩的待组装序列进行了修正或移除,从而使得待组装序列的错误被修正,提高了 组装骨架的正确率。最后,对于组装骨架上的任一位点,通过计算所有可比对至该位点的各 条序列中在该位点的每种碱基出现的频次进行投票,利用每个位点的投票结果对组装骨架 进行修正从而得到了准确率更高的基因序列组装结果。
[0022] 因此,本申请所述的方法、系统和装置只需要较低的基因序列覆盖率就能提供更 高的连续性,并且只需要最小的计算机内存,而且在处理大的基因组测序序列集时,速度有 多个数量级的提升,因此适用于任意长度测序序列,尤其适用于组装通过第Ξ代测序技术 测得的长基因序列。
【附图说明】
[0023] 本申请的上述及其他特征将通过下面结合附图及其详细描述作进一步说明。应当 理解的是,运些附图仅示出了根据本申请的若干示例性的实施方式,因此不应被视为是对 本申请保护范围的限制。除非特别说明,附图不必是成比例的,并且其中类似的标号表示类 似的部件。
[0024] 图1是根据一示例性实施例示出的一种基因序列组装方法的流程图。
[00巧]图2是根据另一示例性实施例示出的一种基因序列组装方法的流程图。
[0026] 图3是根据又一示例性实施例示出的一种基因序列组装方法的流程图。
[0027] 图4A至图4G是根据图1至图3所示的方法进行序列组装的示意图。
[0028] 图5A至图5D是用于说明图3所示的整合步骤的示例性图示。
[0029] 图6是用于实现上述图1中所述基因序列组装方法的装置框图。
[0030] 图7是用于实现上述图2中所述基因序列组装方法的装置框图。
[0031] 图8是用于实现上述图3中所述基因序列组装方法的装置框图。
[0032] 图9是根据一示例性实施例示出的一种基因序列组装系统的框图。
【具体实施方式】
[0033] W下的详细描述中引用了构成本说明书一部分的附图。说明书和附图所提及的示 意性实施方式仅仅出于是说明性的目的,并非意图限制本申请的保护范围。本领域技术人 员可W理解,也可采用许多其他的实施方式,并且可对所描述实施方式做出各种改变,而不 背离本申请的主旨和保护范围。应当理解的是,在此说明并图示的本申请的各个方面可W 按照很多不同的配置来布置、替换、组合、分离和设计,运些不同配置都包含在本申请中。
[0034] 定义
[0035] 本申请中对基因序列的"组装"是指对多个待组装的序列进行比对,并利用序列间 的重叠关系构建基于所述待组装序列的基因图谱。
[0036] 本申请中的术语"比对"当用于对核酸序列的比对时可W通过本领域技术人员 公知的任何方法完成,例如,可^通过现有的软件BLASR、DALIGN、BLAST、BLAST-2、ALIGN、 ALIGN-2BWA、B0WTIE或者MegalignONASTAR)等。本申请中的术语"比对"当用于对索引 表之间的比对时可W通过本领域技术人员公知的任何其他方法完成。(如Smith, Temple F. and Waterman,Michael S. (1981),"Identification of Common Molecular Subsequences", Journal of Molecular Biology 147:195 - 197; 或者 Lipman, DJ ; Pearson, W民(1985),"Rapid and sensitive protein similarity searches",Science 227(4693):1435 - 41:或者 Pearson,WR;Lipman,DJ (1988)'"Improved tools for biological sequence comparison". Proceedings of the National Academy of Sciences of 化e United States of America 85 巧):2444 - 8)。
[0037] 本申请中的"OLC算法"是指由Staden开发的算法(Staden R. A new computer method for the storage and manipulation of DNA gel reading data. Nucleic Acids Res.8(16):3673 - 3694. (1980)) W及在其基础上开发出的相关算法,主要用于 数据量较少的长序列(特别是第一代测序序列)的拼接。基于0LC算法的常见软件 有 Arachne (参见 Batzoglou S, Jaffe DB, Stanley Κ, et al. ARACHNE:a whole-genome shotgun assembler. Genome Resl2:177 - 89. (2002)),Celera Assembler(参见 Myers EW, Sutton GG, Delcher AL,etal.A whole-genome assembly of Drosophila. Science 287:2196 - 204. (2000)),CAP3(参见 Huang X,Madan A. CAP3:A DNA sequence assembly program. GenomeRes 9:868 - 77. (1999)),PCAP(参见加 ang X,Yang SP.Generating a genome assembly with PCAP. Curr Protoc Bioinformatics. Chapter 11 :Unitll. 3. (2005)),Phrap (参见 de la Bastide M,McCombie WR. Assembling genomic DNA sequences with PHRAP. Curr Protoc Bioinformatics. Chapter11:Unitil.4. (2007)),Phusion(参见 Mullikin JC, Ning Z.The phusion assembler. Genome Res 13:81 - 90. (2003))及 Newbler (参见 Marcel MargulieslME, William EA, Said A, et al. Genome sequencing in open microfabricated high density picoliter reactors. NaUire437:376 - 80. (2005))。
[003引 本申请中的"DBG算法,,是指由Idu巧和Waterman开发的算法(Idu巧RM,Waterman MS. A new algorithm for DNA sequence assembly. JComput Biol 2:291 - 306. (19%))臥 及在其基础上开发出的相关算法,后来被广泛应用于生物信息学中基因片段的组装中,主 要适合高覆盖率的短序列的拼接。其相应的软件有Euler-USR(参见化aisson MJ,Brinza D,Pevzner PA.De novo fragment assembly with short mate-paired reads: does 化e read length matter ? Genome Res 19:336 - 46. (2009)),Velvet (参见 Zerbino D民,Birney E.Velvet: algorithms for de novo short read assembly using de Bruijn graphs. Genome Resl8:821 - 9. (2008)),ABySS(参见 Simpson JT, Wong K, Jackman SD, et al. ABySS:a parallel assembler for short read sequence data. GenomeRes 19:1117 - 23. (2009)),AllPa化-LG(参见 Gnerre S,Maccallum I,Przybylski D,et al. High-quality draft assemblies of mammalian genomes from massively parallel sequence data.Rroc Natl Acad Sci USA108:1513 - 18.(2011)),S0APdenovo(参见 Li 民,Zhu H, Ruan J, et al. De novo assembly of human genomes with massively parallel short read sequencing. Genome Res20:265 - 72. (2009))0
[0039] 本申请中的"k-mer"即核巧k联体,是指将基因序列W-定长度k连续切割,按照 每个碱基顺次移动一位切割相同序列长度为k的核巧酸序列,例如对于序列为:ATCGTTGCT TAATGACGTCAGTCGAATGCGATGACGTGACTGACTG 的核巧酸序列,如果 k-mer 为 13-mer 的话,则可 W得到如下的一系列长度为13的核巧酸序列:
[0040] ATCGTTGCTTAAT
[0041] TCGTTGCTTAATG
[0042] CGTTGCTTAATGA
[0043] GTTGCTTAATGAC
[0044] …
[0045] ACGTGACTGACTG。
[0046] 本申请中的"Sparse k-mer"即跳跃式核巧k联体,是指在分析比对结果时,对于 基因序列W k长度进行跳跃性的分割,而不是顺移一位依次切割。
[0047] 本申请中的术语"覆盖率"或者"深度"当用于测序时是指测序得到的碱基总量 化P)与被测基因序列大小的比值,它是评价测序量的指标之一。测序带来的错误率或假阳 性结果会随着测序深度的提升而下降,因此基因组测序中通常会对原始样品进行多次重复 的扩增W增加测序深度。 引实施方式
[0049] 如【背景技术】中所述,第Ξ代测序的高错误率使得最初开发用于相对准确的第一代 测序技术的基于重叠组装的0LC算法不适用。同样,第Ξ代测序得到的带错长序列集使得 最初设计用于第二代测序技术的短序列组装的DBG算法也不再适用。因此,本申请提出了 一种采用混合策略进行基因序列组装的方法。并且,由于第Ξ代测序方法获得的序列集所 需要的数据容量很大,因此即使是在小的微生物的基因组对于运些长序列的纠错都需要很 高的计算资源。为此,本申请进一步通过引入"错点序列"来降低基因序列组装在计算时间 和内存上的高度复杂性,并且确保基因序列组装的精度。
[0050] W下,将参考图1对于本申请进行示例性地说明。图1所示为本申请的一示例性 实施例的一种基因序列组装方法的流程图。首先,本申请所述的基因序列组装方法可W在 任何形式的数字平台系统上实现,该数字平台包括台式计算机、笔记本、便携式掌上电脑、 服务器、计算机集群、网络工作站等,后面将结合图9对实现本申请的方法的一种示例性实 施例系统900加 W详细说明,运里先暂不详述。
[0051] 如图1所示,本申请的基因序列组装方法主要包括S101、S102和S103 Ξ个步骤。 W下将结合附图分别介绍该Ξ个步骤。
[0052] 首先,在步骤S101中,获取一待组装序列集的数据和多个错点序列的数据。
[0053] 具体而言,在步骤S101中,由系统900(例如计算机)通过输入组件(例如数据传 输接口)获取基因测序的数据。该基因测序数据既包括由基因组测序仪器得到的原始序列 的数据,也可W包括基于原始序列的数据进行组装而得到的经过预处理之后的数据。所述 经过预处理之后的数据包括:例如,通过第二代基因组测序技术可W得到高精度的短序列, 然后利用现有的组装软件将其组装成片段的序列(contig)。在一些实施方式中,待组装序 列集可W是任何序列集。在一些实施方式中,待组装序列集是第二代基因组测序技术所得 到的序列集。在另一些实施方式中,待组装序列集是第Ξ代基因组测序技术所得到的序列 集。W下,将W待组装序列集是第Ξ代基因组测序技术所得到的序列集为例来进行说明。
[0054] 错点序列为对应于待测基因组中的某一片段区域的序列,错点序列与其对应的待 测基因组的片段具有相同的或者相似的核酸序列。错点序列可W对应于待测基因组有代表 性的区域,例如待测基因组的高保守区域。在一些实施方式中,所述错点序列为由现有数 据库中获取的序列。在另一些实施方式中,所述错点序列是第二代基因组测序技术所得到 的短序列,或者是通过利用现有的组装软件(例如华大基因公司提供的SOAPdenovo组装软 件)对第二代基因组测序技术所得到的短序列进行组装而得到的序列,或者上述两者的组 合。在另一些实施方式中,所述错点序列为待组装序列本身的局部截取,例如通过将待组装 序列集中最长的一段序列截短为一系列的短序列,之后将上述短序列与待组装序列集中所 有其他序列进行比对,对于那些与上述截短得到的短序列没有重叠的序列的部分进行进一 步的截取,再加入到上述短序列的集合中,直到得到的短序列集可W与所有序列均有重叠。
[0055] 本领域一般技术人员可W根据待组装序列的长度选择错点序列的长度,对于第Ξ 代测序,错点序列的长度一般设为待组装序列集中序列的平均长度的10-30%。对于第二代 测序技术,错点长度可W长于平均测序长度。对于一组待组装序列集可W选择固定的错点 序列长度,也可W选择不同长度的错点序列。应当理解,错点序列的平均长度越长,则压缩 后的待组装序列的索引表的平均长度越短,数据处理量越小;同时错点序列数量越多则对 其整个待组装基因序列的覆盖范围越高,一般选择的错点序列数量乘W错点序列的平均长 度应该不低于待组装基因序列的总长度的30%。在一些实施方式中,所述错点序列长度在 2bp-30bp之间、或者在30bp-100bp之间、或者在100bp-300bp之间、或者在30bp-1001Ap之 间、或者在3(K)bp-100化P之间。不同的错点序列可W由不同的标号表示,所述标号可W是 数字或者字母,或者数字与字母的组合,或者其他适合的表示形式。例如,可W用若干位数 字为错点序列标号,例如多个错点序列可分别标记为数字001、002、003等等。又例如可W 用代数标号,如XI、X2、X3等分别表示不同的错点序列。
[0056] 然后,在步骤S102中,将利用所述多个错点序列,对所述待组装序列集中的每一 个待组装序列进行压缩,W获得该待组装序列的W错点序列的标号来表示的索引表。
[0057] 具体而言,将错点序列与待组装序列进行比对,找出错点序列在待组装序列上 的排放位置,并按照错点序列在待组装序列上的排列位置,获得对应于该待组装序列的 索引表,由此将每条待组装序列W错点序列的标号表示,即,将基因序列压缩为W标号形 式表示的一串数字或字符,完成序列压缩,并获得对应于多个错点序列的多个索引表, 所述多个索引表形成对应于该待组装序列集的压缩的数据集。例如:如果错点序列为 contig_331, contig_101, contig_17,当运些contig被依次比对到一条待组装序列上,且 错点序列在待组装序列上的排列位置依序为contig_331, contig_101, contig_17,则该待 组装序列对应的索引表为巧31,101,17]。
[0058] 在某些实施方式中,通过W下的k-mer方式压缩每一个待组装序列,并获得每一 个待组装序列的索引表:将错点序列和待组装序列打断成相同长度的k-mer顺位集合(例 如,设k = 13,将错点序列和待组装序列分别从其起始位点起打断成13-mer的顺位集合)。 将获得的错点序列和待组装序列的k-mer分别进行一一比对,并且根据错点序列的长度选 择合适的阔值,当所述待组装序列的核巧k联体中与错点序列的核巧k联体相匹配的个数 超过预定阔值,则认为所述错点序列可W匹配至该待组装序列;通过每个所述匹配的核巧 k联体分别在所述错点序列及所述待组装序列中的相对位置,计算得到该匹配的核巧k联 体所对应的错点序列在待组装序列中的排列位置。对于有多个匹配的核巧k联体的错点 序列,所述错点序列在待组装序列中的排列位置可W是根据每个匹配的核巧k联体计算得 到的所述错点序列在待组装序列中的排列位置的平均值。其中,在本发明的一些实施方 式中,所述预定阔值为所述错点序列的长度的0. 1%至2%。本领域技术人员可W根据计 算得到的每个错点序列在待组装序列中的排列位置,获得待组装序列的索引表。其中所 述1<-11161'的长度为:269-5化9、1569-1%9、1569-2化9、269-1569、269-1%9、3化9-5化口或 15bp-51bp。
[0059] 另外,可W将索引表被全部包含在另一待组装序列的索引表之内的待组装序列删 除,W简化待组装序列数据集。并且,W-待组装序列为当前序列,将所述当前序列的索引 表与待组装序列集中的所有其它序列的索引表进行一一比对,找出与当前序列索引表具有 重叠关系的序列(也可称为"备选序列"),之后将备选序列与当前序列进行索引表比对或 碱基比对,找到其中相对于当前序列包含在某一方向上索引表的延伸的被压缩的待组装序 列,并将完全被当前序列所包含的序列移除;在所剩的备选序列中挑选出与当前序列具有 向左或向右最优重叠的序列,移除不具有在某一方向上的最大重叠度的序列,仅保留与当 前序列具有向左或向右最优重叠的序列(其也可W称为当前序列的最优延伸序列),并将 所述最优延伸序列相对于当前序列延伸的部分加到当前序列上即获得当前序列的扩展序 列。在上述比对过程中,若当前序列被其他与其具有重叠关系的序列完全包含,则删除当前 序列,不生成当前序列的扩展序列;若没有找到相对于当前序列具有向左或向右延伸的序 列,则W所述当前序列作为该轮比对所得到的扩展序列。通过将所述待组装序列集中的每 条序列作为当前序列重复上述步骤W得到一组扩展序列,然后根据所述扩展序列中错点序 列标号的位置对应关系对各扩展序列进行排布,从而完成构建组装骨架。可W使用本领域 技术人员所知的任何合适的方法和软件来计算序列之间的重叠,例如通过BLASR、DALIGN、 BLAST、BLAST-2、ALIGN、ALIGN-2BWA、B0WTIE、SOAP 或者 Megali即值NASTAR)等比对软件对 当前序列与待组装序列集中的所有其他序列进行一一比对,基于比对结果找出与当前序列 具有重叠关系的序列。
[0060] 在某些实施方式中,在最优重叠图的构建步骤中或者在对最优延伸序列的筛选步 骤中,每个节点代表一条待组装序列,其为一个或多个连续的错点序列标号的组合。为获得 最优重叠图可W进行两轮运算,第1轮,删除所有被全部包含的节点,例如,如果存在节点 [又1,而,&],则删除被全部包含在上述节点中的节点[Xi,X2]及节点[而,&],运样做的目的 是为了避免不必要的重复和包含的节点的比对。第2轮,计算不相互包含的节点之间的所 有后缀-前缀重叠。然后对最优重叠图进行简化,当从节点A到节点B的延伸的重叠度高 于所有从节点A向同一方向延伸的其他重叠时,节点A至节点B的延伸存在优势。在某些 实施方式中,重叠度可由对应的重叠区域的长度或重叠区域所设及多个错点的总长度或重 叠区域中待组装序列与错点匹配的k-mer数量来衡量。一般来说,将每条序列作为当前序 列与其他序列进行比对,都可能找到一个向左具有与当前序列最大重叠的最优延伸序列, W及一个向右具有与当前序列最大重叠的最优延伸序列,将所述向左的最优延伸序列和向 右的最优延伸序列相对于当前序列延伸的部分加到所述当前序列上即获得该序列的扩展 序列。但由于待组装序列长度有限,因此可能出现无法处理的重复(比如多条待组装序列 均为某一当前序列的最优重叠,称为并列最优)此时会导致在图谱上出现分支(Miller, J. R. et al. Aggressive assembly of pyrosequencing reads with mates. Bioinformatics 24, 2818-2824 (2008))。当出现分支时,组装算法在出现分支处断开,或选择其中一条序列 间相互最优(即A至B的延伸是A在右向的最优重叠延伸,B至A的延伸是B在左向的最 优重叠延伸)的重叠进行延伸。将最优重叠图中不出现分叉(即序列间不存在并列最优) 的每个线性区域单独输出为组装结果。将W索引表形式表示的压缩序列串联起来,之后重 新转化为W碱基形式表示的未压缩的待组装序列作为基因组轮廓的组装骨架。
[0061] 因此,通过利用如图1所述的组装基因序列方法,由错点序列(例如,高精度测序 短序列或者高精度测序短序列的组装片段)对待组装序列(例如,由第二代基因组测序技 术或者第Ξ代基因组测序技术所得到的长序列)进行标记,并将待组装序列压缩成对应的 索引表,然后再进行比对和组装,从而大大加快了组装速度并减少了计算量。
[0062] 在某些实施方式中,对待组装序列的组装可W依靠上述待组装序列的索引表,具 体地包括:选取所述待组装序列集中的一条序列作为当前序列,将所述当前序列的索引表 与待组装序列集中的所有其它序列的索引表进行一一比对,找出与所述当前序列具有重叠 关系的序列,将所述当前序列与所述与当前序列具有重叠关系的序列进行索引表比对或碱 基比对W找出当前序列的最优延伸序列,并将其相对于当前序列延伸的部分加到当前序列 上即获得当前序列的扩展序列(即向左或向右最优重叠延伸),并通过将所述待组装序列 集中的每条序列作为当前序列重复上述步骤W得到一组扩展序列,然后根据所述扩展序列 中错点序列标号的位置对应关系对各扩展序列进行排布,从而形成组装骨架。
[0063] 图2所示为根据本申请另一示例性实施例的一种基因序列组装方法的流程图。图 2中所示基因序列组装方法主要包括S201、S202、S203 W及S204四个步骤。其中S201步骤 用于获取一待组装序列集和多个错点序列,其与上文所述的图1中的S101基本相似,在此 就不再重复。同样,S202步骤(压缩待组装序列W获得由错点序列的标号构成的索引表) 与上文的S102类似,在此也不再重复。
[0064] W下将详细描述S203步骤。通过S202步骤获得W错点序列的标号来表示的索引 表后,在S203步骤中进一步通过比对来修正被压缩的待组装序列。 阳0化]具体地,利用所述错点序列对待组装序列的索引表进行反向检索,通过比对来修 正被压缩的待组装序列包括:根据S201中得到的错点序列对待组装序列的索引表进行反 向检索,得到S202中所述索引表内与所述的一条错点序列或者所述错点序列中的多条连 续排列的错点序列的组合相对应的待组装序列的集合(即,待组装序列的索引表子集);对 利用所述反向检索得到的待组装序列集合中的每条待组装序列的索引表,将其作为当前序 列索引表与所述索引表子集中所有其他待组装序列进行一一比对(即将每两条待组装序 列中对应位置的错点序列的标号进行一一比对),根据每个位置上不同的错点序列标号出 现的频次进行投票,修正每条当前序列中出现的错误的错点序列标号,并且若当前序列的 索引表中具有重复的或被包含的错点序列排列关系,则将该当前序列删除,从而获得更为 精准且精简的待组装序列的索引表集。
[0066] W错点序列为contig-Xi、-X2、-X3。。。-Xi为例,上述通过比对来修正被压缩的待 组装序列的步骤包括:利用contig的标号Xi、X2、&。。。Xi构建索引表,使得对每条错点序 列,都可借助索引表迅速寻找所有含有该错点序列的待组装序列。对每个contig错点序列 标号Xl,记录包含该Xl的待组装序列的序列索引表,得到相对应的待组装序列的索引表的 集合,运个过程可W被看作为一个序列分组的过程。之后,用动态规划算法或者其他更高效 的算法对得到的待组装序列集合内的各条待组装序列进行组内一一比对。上述一一比对可 在被压缩的待组装序列上进行也可在未压缩的待组装序列上进行。
[0067] 对于被压缩的待组装序列,组内的一一比对是指将每条压缩序列(当前序列) 的索引表与组内其他序列的索引表进行比对:一方面,将当前序列中在多序列比对仅仅 出现很少次数(例如,仅仅出现1-5次、或不多于5次、或不多于4次、不多于3次、不多 于2次、不多于1次)的错点序列标号从该序列的索引表中删除,例如当前序列索引表为 巧31,101,21,17],而在多序列比对中标号为21的错点序列仅在当前序列的索引表中出现 了一次,而未在任何其他序列的索引表中出现过,则删除当前序列索引表巧31,101,21,17] 中的错点序列标号21,即原序列的索引表由[331,101,21,17]变为巧31,101,17];另一方 面,如果在多序列比对中,发现当前序列中存在断点(即所有与断点左端重叠的序列均与 断点右边不重叠,反之也是如此,则判定该点为断点),则将此整条序列称为杂合序列,并将 它删除。同时在一一比对过程中,若当前序列被其他序列包含,则删除该被包含的当前序列 (例如,如果当前序列索引表为[XI,而]或者[而,X3],当存在索引表为[XI,枯X3]的序列,贝U 删除当前序列)。被删除的序列将不再参与后续的组装骨架的构建。
[0068] 考虑到通过第Ξ代测序技术获得的待组装序列的错误率比较高,因此通过S203 步骤对被压缩的待组装序列进行了修正或移除,可使得后续的组装骨架的正确率得到了很 大提局。
[0069] 接着,在步骤S204中,构建组装骨架。具体而言,利用S203中修正后的被压缩的 待组装序列构建最优重叠图谱,从而完成构建组装骨架。因为该步骤S204和上文所述的步 骤S104相类似,运里就不再重复。
[0070] 并且,本申请还提供了又一种基因序列组装方法的实施方式。图3是本申请又一 示例性实施例的一种基因序列组装方法的流程图。图3中所示的基因序列组装方法主要包 括 S301、S302、S303、S304 W及 S305 五个步骤。
[0071] 其中在步骤S301中,获取一待组装序列集的数据和多个错点序列的数据;然后在 步骤S302中,压缩待组装序列为W错点序列的标号表示的索引表;在S303步骤中,通过比 对删除被压缩的待组装序列中的杂合序列及被包含序列,W及将待组装序列的多序列比对 中仅仅出现很少次数(例如,仅仅出现1-5次、或不多于5次、或不多于4次、不多于3次、 不多于2次、不多于1次)的错点序列的标号从包含其的序列的索引表中删除;在步骤S304 中,构建组装骨架;并且在步骤S305中,通过整合得到准确率更高的基因序列组装结果。在 本实施方式中,步骤S301、S302、S303和S304分别与上文中的步骤S101、S102、S103和S104 相类似,在此就不再重复。W下将详细阐述步骤S305。
[0072] 具体而言,在步骤S305中,对构建的组装骨架进行如下所述的再修正:将所有待 组装序列与步骤S304所得到的组装骨架进行比对,对于骨架上的任一位点,通过计算所有 可比对至该位点的各条序列中每种碱基出现的频次进行投票,利用每个位点的投票结果对 组装骨架进行再修正,从而最终整合得到准确率较高的基因序列组装结果。
[0073] 在某些实施方式中,上述投票整合是通过将所有待组装序列比对到组装骨架 上,并根据比对结果构建跳跃式 k-mer 图(Ye, C. , Ma, Z. S. , Cannon, C. H. , Pop, M. &Yu, D. W. Exploiting sparseness in de novo genome assembly. BMC BioinformaticslSSuppl 6, SI (2012)),利用同一位点的多条序列进行投票整合,再修正骨架中的错误,并利用跳跃 式k-mer图的结构获取最终准确率更高的基因序列组装结果。跳跃式k-mer的用意是在 大数据量情况下通过跳跃方式来节省计算时占用的内存。k值一般用1~5,用稍微大的 k值可W绕开比对算法中的错误比对。跳跃步长g-般也是1~5,稍微大的步长可W显 著减少计算时占用的内存。在整合精度相当的情况下,与现有的PBagcon算法(化in, C. S. et al. Nonhybrid,finished microbial genome assemblies from long-read SM民T sequencing data. Na 化 re me 化 odslO, 563-569(2013))相比,使用 1^下所述的 Sparc 整合算 法速度为其5倍左右,计算时占用的内存为其1/2-1/5。
[0074] Sparc整合算法是基于跳跃式k-mer的一种方法,具体可W包括Ξ个步骤,第一 步:构建骨架k-mer图;第二步进行序列比对并在比对过程中根据比对结果记录新的节点 和边W及更新各边的权重;W及第Ξ步:捜索最大路径。W下将结合附图5A-5D对上述各步 骤作出详细说明。 阳0巧]如图5A所示,组装骨架序列为ACTGGACTAAA,W k-mer长度为2,跳跃步长g为3 构建k-mer图,则:所述k-mer的起始位置分别为1, 4, 7, 10,其相对应的k-mer分别为:AC、 GG、CT、AA,被选为节点。记录下k-mer之间的连接关系:例如第一个k-mer AC往右延伸经 过TGGS个碱基后,到达第二个k-mer GG的结束位置。TGG则被记录为连接AC节点与GG 节点的边(edge),此时因为此TGG仅出现过一次,所W边的权重是1,依此类推。
[0076] W下将示例性地分析几条序列与骨架的比对结果。
[0077] 如图5B所示,第一条原始序列(即待组装序列或其片段)为ACTGGACCAAA,将其与 组装骨架比对。在比对结果中,起始位点1、4的AC、GG节点能够完全比对上,则对该节点不 需进行额外操作,仅仅将边的权重增加,即路径TGG的权重变为2。到达第7位点的时候, 由于CC与骨架中的CT不完全匹配,则将CC作为新的节点,并记录下一从第4位节点GG到 第7位节点CC的新路径ACC,作为新的边,权重为1。由于第10位节点的AA与骨架匹配, 所W不需要生成新的节点,仅仅记录下从新的第屯位节点CC到之前已有的第十位点的新 的边AAA,权重为1。 阳07引如图5C所示,第二条原始序列为ACGGACCAAA,将其与骨架比对。在比对结果中,起 始位点1的节点AC匹配,经过两个碱基GG到达第四位节点GG,由于GG节点与前面生成的 边TGG不同,所W生成新的边缘GG,权重为1。接下来,从第4位节点GG到第7位节点的时 候,由于节点CC在之前序列中已经生成,路径也已经被记录,则增加对应的ACC边的权重, 权重为2。
[0079] 该过程可逐序列进行,然后可在最终结果中捜索权重最大的路径并将所述路径相 组合(如图5D所示),从而得到整合序列:ACTGGACCAAA。在某些实施方式中,图中的所有 权重需要经过调节,比如减去某个小的数值(例如该区域测序深度的10% -30% ),W避免 由测序区域中的长插入性错误导致的较大的权重总和的错误路径。
[0080] 因此,通过将所有待组装序列与组装骨架比对,对多条序列在同一位点的核巧酸 进行投票整合(即找出出现频率最高的核巧酸并选取该核巧酸作为该位点的核巧酸),优 选的,通过Sparc整合算法再修正所述组装骨架,从而得到准确率更高的基因序列组装结 果。
[0081] 图4A至图4G是根据图1至图3所示方法进行序列组装的示意图。首先,如图4A 所示,在步骤S10US201或S301中,获取待组装序列集和错点序列。其中图4A左方的各种 类型的细线条表示带有数字标号的错点序列,而右方的各种粗线条则表示各个待组装序列 (例如图4A所示的原始长序列)。如上文所述,所述待组装序列集可W是由基因组测序仪 器得到的原始序列的数据或者是基于原始序列的数据进行组装而得到的预处理过的数据; 所述错点序列可W是由现有数据库中筛选获取的,或者是通过利用现有的组装软件(例如 华大基因公司提供的组装软件)对第二代基因组测序技术所得到的短序列进行组装而得 到的,或者是测序结果中待组装序列本身的局部截取。
[0082] 图4B示出了与S102、S202、S302步骤相对应的序列标记和压缩过程:即将错点序 列与待组装序列进行基因序列比对,计算错点序列在待组装序列上的排放位置,并如图4B 所示,用错点序列的标号标记待组装序列;被标记的各条待组装序列即被压缩为由各个错 点序列的标号W其在长序列上出现次序排列的列表,完成序列压缩并获得W错点序列为标 记的待组装序列的索引表。
[0083] 如图4C和4D所示为与S203、S303相对应的步骤:利用从S102、S202、S302中获 得的索引表,反向检索得到与某一错点序列相关的所有待组装序列的分组(如图4C所示, 反向检索所有包含错点序列1的序列得到的序列分组);之后将反向检索得到的组内序列 进行一一比对,对于每条序列(当前序列),将它的索引表与所有组内其他序列的索引表进 行比对:一方面,将当前序列中在多序列比对中仅仅出现极少次数(例如,仅仅出现1-5次、 或不多于5次、或不多于4次、不多于3次、不多于2次、不多于1次)的错点序列标号从该 序列的索引表中删除;另一方面,如果在多序列比对中,发现当前序列中存在断点(即所有 与断点左端重叠的序列均与断点右边不重叠,反之也是如此,则判定该点为断点),则将该 当前序列作为杂合序列删除;再一方面,在一一比对过程中,若当前序列被其他序列包含, 则删除该被包含的当前序列(图4D)。被删除的错点及杂合序列或被包含的序列将不再参 与后续的组装骨架的构建。上述反向检索后对待组装序列进行分组,并对每条待组装序列 的索引表进行比对修正的过程,最终得到的结果是更为精准且精简的待组装序列的索引表 集,而该精准且精简的索引表集中所有的序列将被用于后续的组装骨架的构建过程中。
[0084] 而图4E-图4F则示意了与S104、S204、S304相对应的构建组装骨架的过程:即基 于所述索引表,将每一条被压缩的待组装序列(当前序列)的索引表与用于组装骨架的其 他被压缩序列的索引表进行一一比对,找出与所述当前序列具有重叠关系的备选序列,将 备选序列与当前序列进行索引表比对或碱基比对,首先移除比对结果中索引表整个被其他 序列的索引表完全包含的序列,之后寻找到与它具有重叠的错点序列并且能够在某一方向 上延展的多条被压缩的待组装序列,将比对结果中索引表不具有在某一方向上的最大重叠 的序列移除,仅保留当前序列的最优延伸序列,并将其相对于当前序列延伸的部分加到当 前序列上即获得当前序列的扩展序列,并通过将所述待组装序列集中的每条序列作为当前 序列重复上述步骤W得到一组扩展序列,然后根据扩展序列中错点序列标号的位置对应关 系对各扩展序列进行排布(如图4E所示);从而完成构建组装骨架(如图4F所示)。 阳0化]如图4G所示为与S305相对应的步骤:在构建组装骨架后,将所有待组装序列比对 到组装骨架上,对于组装骨架上的任一位点,通过计算所有可比对至该位点的各条序列中 在该位点的每种碱基出现的频次进行投票,利用每个位点的投票结果对组装骨架进行再修 正,w得到准确的基因序列组装结果。
[0086] 图6-8是用于实现上述图1-3中所述基因序列组装方法的装置框图。参照图6, 该装置包括获取模块601、压缩模块602 W及组装模块604。其中获取模块601被配置为获 取一待组装序列集的数据和多个错点序列的数据;压缩模块602被配置为利用所述多个错 点序列,对所述待组装序列集中的每一个待组装序列进行压缩W获得该待组装序列的索引 表,其中所述索引表W错点序列的标号来表示;组装模块603被配置为利用所述待组装序 列集中各待组装序列压缩得到的所述索引表构建组装骨架。
[0087] 如图7所示,上述实现基因序列组装方法的装置除了包括与图6中601、602、604 具有相同功能的获取模块701、压缩模块702 W及组装模块704 W外,还可W进一步包括待 组装序列的索引表修正模块703。其中该待组装序列的索引表修正模块703被配置为利用 所述错点序列对待组装序列的索引表进行反向检索,得到与所述错点序列相对应的待组装 序列的索引表子集,通过将所述索引表子集中所有其他序列的索引表比对到当前待组装序 列的索引表上来修正所述待组装序列的索引表子集中每条待组装序列的索引表。
[0088] 如图8所示,上述实现基因序列组装方法的装置除了包括与图7中701、702、703、 704具有相同功能的获取模块801、压缩模块802、待组装序列的索引表修正模块803 W及组 装模块804 W外,还可W进一步包括基因序列组装结果修正模块805。其中该基因序列组 装结果修正模块805被配置为对构建的组装骨架进行如下所述的再修正:在构建组装骨架 后,将所有待组装序列与组装骨架进行比对,对于骨架上的任一位点,通过计算所有可比对 至该位点的各条序列中在该位点的每种碱基出现的频次进行投票,利用每个位点的投票结 果对组装骨架进行再修正,从而得到再修正后的基因序列组装结果。在一些实施方式中,如 图8所示的装置还可W进一步包括:基因序列组装结果整合模块806,所述基因序列组装结 果整合模块用于获得覆盖待组装基因序列不同区域的组装骨架,利用多个计算机处理器同 时分别修正所述覆盖待组装基因序列不同区域的组装骨架,W及将再修正后的组装骨架进 行整合得到完整的基因序列组装结果。
[0089] 本申请的实现基因序列组装方法的装置可W包括和上述方法步骤相对应的模块 和/或模块组合,并且随着上述方法步骤的改变,上述装置也可W进行相应的变化W实现 本申请的基因序列组装方法。
[0090] 并且,本申请还提供了一种用于实现上述图1-3中所述方法的系统。图9是根据 一示例性实施例示出的一种本申请的系统900。系统900可W是任何形式的数字平台,包括 台式计算机、笔记本、便携式掌上电脑、数字平台、计算机集群、网络工作站等。
[0091] 如图9所示的系统900包括:通过总线950内部连接的一个或多个输入组件910 ; 一个或多个输出组件920 ; -个或多个中央处理器(CPU)930 及一个或多个存储组件 940,所述存储组件940在其中存储一个或多个计算机程序944、一个或多个操作系统946、 W及可选的一个或多个数据库948。
[0092] 其中,输入组件910使得可W向系统900提供数据输入。常用的输入组件910包 括数据输入接口、网络传输接口,或者其它类型的输入部件。在本申请中,输入组件910用 于获取一待组装序列集的数据和多个错点序列的数据。
[0093] 并且,输出组件920用于向用户提供计算结果,常见的输出组件920包括用户图形 界面(GUI)、Ξ维显示界面、输出数据接口、网络传输接口,或者其它类型的输出部件等。
[0094] 中央处理器(CPU) 930通常控制系统980的整体操作,诸如与显示、数据处理、数据 存储、数据通信W及记录操作相关联的操作等。中央处理器(CPU)930可W包括一个或多个 处理器来执行指令,W完成上述的方法的全部或部分步骤。此外,中央处理器(CPU)930可 W包括一个或多个模块,便于处理中央处理器(CPU)930和其他组件之间的交互。例如,中 央处理器(CPU) 930可W包括输入/输出模块,W方便输入组件910、输出组件920和中央处 理器(CPU)930之间的交互。
[0095] 存储器940可W由任何类型的随机存取存储设备RAM 941或只读存储设备ROM 942或者它们的组合实现,如静态随机存取存储器(SRAM)、电可擦除可编程只读存储器 巧EPROM)、可擦除可编程只读存储器巧PROM)、可编程只读存储器(PROM)、只读存储器 (ROM)、CD-ROM、磁带、软盘和光数据存储设备等。存储器940可W被配置为存储各种类型的 数据W支持在系统900的操作。运些数据的示例包括用于在系统900上操作的任何操作系 统946、计算机程序944、数据库948等。其中可W在所述系统900中执行所述计算机程序 944,计算机程序944通过指令指导中央处理器930的相应模块执行本申请中所述的基因序 列组装的方法。并且,在本发明中,存储器940还被配置于存储所述待组装序列集的数据和 多个错点序列的数据。
[0096] 其中,如图1-3所示的方法可W用例如计算机软件、硬件或者其组合来在计算机 可读介质中实现。对于硬件实现而言,运里所述的实施例可W通过专用集成电路(ASIC)、数 字信号处理器值SP)、数字信号处理设备值SPD)、可编程逻辑器件(PLD)、现场可编程口阵 列(FPGA)、处理器、控制器、微控制器、微处理器、设计成执行运里所述功能的其它电子单元 或者其选择性组合内的一个或多个来实现。
[0097] 对于软件实现而言,运里所述的实施例可W用单独的软件模块,诸如过程模块和 功能模块来实现,其每一个都执行一个或多个运里所述的功能和操作。软件代码可W用任 何适当的编程语言编写的软件应用来实现,并且可W存储在专用的计算机系统的存储器或 者其他计算机可读介质中,并且由计算机系统的处理器来执行,也可W安装在具备数据存 储和处理功能的其他电子设备中,如平板电脑、服务器等。
[0098] 在一些实施方式中,可W通过系统900实现:利用所述多个错点序列,对所述待组 装序列集中的每一个待组装序列进行压缩W获得该待组装序列的索引表,其中所述索引表 W错点序列的标号来表示;W及利用所述待组装序列集中各待组装序列压缩得到的所述索 引表构建组装骨架。在一些实施方式中,可W通过系统900进一步实现:在获得每条待组装 序列的索引表后,利用错点序列对待组装序列的索引表进行反向检索,得到与所述错点序 列相对应的待组装序列的索引表子集,通过将所述索引表子集中所有其他序列的索引表比 对到当前待组装序列的索引表上来修正待组装序列的索引表子集中每条待组装序列的索 引表。在另一些实施方式中,可W通过系统900进一步实现对构建的组装骨架进行如下所 述的再修正:在构建组装骨架后,将所有待组装序列与组装骨架进行比对,对于组装骨架上 的任一位点,通过计算所有可比对至该位点的各条序列中在该位点的每种碱基出现的频次 进行投票,利用每个位点的投票结果对组装骨架进行再修正,从而得到再修正后的基因序 列组装结果。在一些实施方式中,可W通过系统900进一步实现:获得覆盖待组装基因序列 不同区域的组装骨架,利用多个计算机处理器同时分别修正所述覆盖待组装基因序列不同 区域的组装骨架,W及将再修正后的组装骨架进行整合得到完整的基因序列组装结果。
[0099] 系统900可W被配置于实现本申请所记载和包含的方法的任一步骤W及任一步 骤的组合,在此不再重复。
[0100] 综上所述,本申请提供了一种高效的并支持各种各样的测序技术的从头开始的 (de novo)组装算法。与针对第Ξ代测序的一些现有的算法工具相比,本申请的方法需要 较低的序列覆盖率、提供更高的连续性、需要最小的内存并且在处理大的基因组测序序列 集时速度有多个数量级的提升。例如,根据本申请的方法所组装得到的基因组大小可W在 任意范围,例如在30化P-300抓P之间、或者在lMbp-4. 6Mbp之间、或者在lMbp-12Mbp之间、 或者在lMbp-120Mbp之间、或者在lMbp-3抓P之间。并且测序深度可W在5X-200X之间、 5x-20x之间、或者在5x-30x之间、或者在5x-40x之间、或者在5x-50x之间、或者在5X-128X 之间、或者在5X-237X之间。而且通过利用本申请的方法、系统W及装置所进行的基因组 组装耗时为现有组装工具的10 1-10 5倍、或者10 1-10 3倍、或者10 2-1〇 3倍、或者10 2-1〇 5 倍、或者10 3-10 5倍。 阳101] 本申请所公开的方法、系统W及装置与进行第Ξ代测序技术的公司(例如太平洋 生物科学公司(Pacific Biosciences)和牛津纳米孔公司(Oxford Nanopore))的测序仪 均可兼容,也与第二代测序技术,如Illumina公司的各种测序技术兼容。并且,本申请所公 开的方法、系统W及装置可W用于大动物、植物、微生物或者人的基因组的研究。进一步,本 申请所公开的方法、系统W及装置允许在工作站级别的电脑上甚至可能在普通台式电脑上 高效组装大型基因组,而现有的组装工具通常需要超级计算机或计算机集群。因此,本申请 所公开的方法、系统W及装置将能够大大加快第Ξ代测序技术扩展至设及植物、动物、微生 物和人类基因组等的研究的更加广泛的应用领域。 阳102] W下将通过多个实施例将本申请所述的方法与其他现有的组装方法进行比较。 阳103] 实施例1 :
[0104] 对第Ξ代测序序列集进行组装
[0105] 在本实施例中,比较了基于本申请的方法、系统W及装置所设计的程序(W下简 称为"DBG20LC"程序)与其他组装程序化cBio2CA(参见Koren,S.et al. Hybrid error correction and de novo assembly of single-molecule sequencing reads. Nature biotechnology 30, 693-700 (2012)), HGAP(参见 Qiin,C.S.et al. Nonhybrid, finished microbial genome assemblies from long-read SMRT sequencing data. Nature methods 10, 563-569(2013) ;Berlin, K. et al. Assembling Large Genomes with Single-Molecule Sequencing and Locality Sensitive Hashing, (2014))和 ECtools (参见 Lee, H. et al. Error correction and assembly complexity of single molecule sequencing reads, (2014))得到的结果,其中上述其他组装程序均是为第Ξ代测序在不同的测序序列 集中设计的基因组组装程序。
[0106] 在本实施例中选择的待组装序列是基于化cBio测序平台测得的S. cerevisiae w303基因组(基因组大小为12Mbp)的测序序列集,该测序序列集中待组装序列平均长度 为51Λρ。错点序列来自:利用Illumina Miseq测序序列集,WSparseAssembler的跳跃 k-mer方法进行组装,得到约2500条错点序列,错点序列平均长度约为4. 3化P,错点序列总 长度约为11. 4Mbp,对待测基因组序列的覆盖率约为95%。基于本申请的DBG20LC对第二 代基因测序巧Ox深度/覆盖率)和第Ξ代基因测序(10X/20X深度/覆盖率)得到的序 列集进行组装,W k-mer长度为17bp、阔值为(0. 1 % *每个错点序列的长度)的参数将错 点序列比对到待组装序列从而得到压缩的待组装序列索引表的集合,之后修正压缩的待组 装序列中的错误,并且完成组装骨架的构建的步骤一共用了大概10CPU分钟,然后用Sparc 运行了两次对基因序列组装结果的修正步骤W获得合适的序列组装结果。在本实施例中, Sparc中将所有待组装序列与组装骨架W修正基因序列组装结果的步骤通过现有的Blasr 矛呈序(Chaisson, Μ. feTesler, G. Mapping sin邑le molecule sequencin邑 reads usin邑 basic local alignment with successive refinement(BLASR):application and theory. BMC 61〇111化^11日*山313,238(2012))实现的,其中利用61日3'程序进行序列比对^修正基因序 列组装结果的过程占据了整个组装过程中大部分的运算时间(约2个CPU小时)。
[0107] 将本申请的方法与其他Ξ个现有的软件进行比对,并且比对了修正后的 W303化cBio组装序列结果。 阳10引表1. S. cerevisiae基因组组装性能比较(基因组大小:12M bp) 阳 109]
[0110] 表1中的术语"同一性"是指两条或多条核酸序列进行比对时序列相同的部分的 碱基数占整条序列的全部碱基数的百分比。表1中的术语"NG50"是指在基因序列组装时 当50%的待测基因原长被包含在内的组装序列的最小大小。表1中的术语"修正NG50"是 指在基因序列组装中当50%的待测基因原长被包含在内的修正后的组装序列的最小大小。 除非另有说明,上述术语在W下各实施例中或本申请说明书其他部分均表示相同的含义。 阳111] 如表1所示的结果表明DBG20LC需要最少的覆盖率(深度)和计算时间,并且能 够得到与其他组装程序类似质量的组装结果。在一般情况下,DBG20LC可W用比现有的大 多数组装程序低一到两个数量级的时间和内存就产生比较好的结果,并且可W在较低的 (10X-20X)覆盖率(深度)的化cBio测序序列集中实现较好的组装质量。在本实施例中, 注意到大部分时间是用在使用Blasr将所有测序得到的原始序列比对至组装骨架的组装 结果修正部分(2CPU小时)。因此,表1中分开标注了运部分需要的CPU时间。由于运项对 基因序列组装结果的比对修正步骤可W通过将组装骨架分割成多个不同区域,利用多个计 算机处理器同时分别修正组装骨架的不同区域,之后将修正后的不同区域进行整合得到完 整的基因序列组装结果,所W其经历时间可W更短。 阳11引实施例2:
[0113] 在本实施例中选择的待组装序列是基于化cBio测序平台测得的A. thaliana ler-0基因组(基因组大小为llSMbp)的测序序列集,该测序序列集中待组装序列平均长 度为51Λρ。错点序列来自:利用Illumina Miseq测序序列集,W Sparse Assembler的 跳跃k-mer方法进行组装,得到约83000条错点序列,错点序列平均长度约为1. 3化p,错 点序列总长度约为115Mbp,对待测基因组序列的覆盖率约为97 %。在本实施例中使用 Sparse Assembler 组装软件(Ye, C. , Ma, Z. S. , Cannon, C. H. , Pop, M. &Yu, D. W. Exploiting sparseness in de novo genome assembly. BMC BioinformaticslSSuppl 6,SI (2012)) 化-mer长度为5化p ;sparse步长g = 15)在2小时内预组装了 50x深度的Illumina (第二 代基因组测序平台)测序序列集。基于本申请的DBG2化C,利用上述组装结果作为错点,对 20x/40x化cBio测序序列集进行组装,W k-mer长度为17bp、阔值为(0. 5% *每条错点序 列长度)的参数将错点序列比对到待组装序列从而得到压缩的待组装序列索引表的集合, 之后修正压缩的待组装序列中的错误,并且完成组装骨架的构建的步骤一共用了大概1个 CPU小时。而之后用Sparc (k-mer长度1,步长g = 1)对基因序列组装骨架进行修正。又 花费了 18个CPU小时来完成最终的对组装结果的修正过程。在整个过程中,最大内存使用 量为6GB。
[0114] 相比较而言,由于现有的软件在处理相同规模的问题时需要耗时超过一千个CPU /j、日寸(Berlin, K. et al. Assembling Large Genomes with Sin邑le_Molecule Sequencing and Locality Sensitive Hashing. (2014)),?下为收集的数据对比(Xee,H.et al.Rrror correction and assembly complexity of single molecule sequencing reads. (2014))。 为了评估序列的组装,将本申请的方法得到的组装结果与使用高覆盖率的化cBio测序序 列集并且经如iver修正后的HGAP组装结果进行比对,计算纠错修正序列。 阳11引表2. A. thaliana基因组组装性能比较(基因组大小:118M bp) 阳116]
[0117] 实施例3:
[0118] 在本实施例中,待组装序列集是基于化cBio测序平台测得的54x覆盖率的人类 化.sapiens)(基因组大小:3. OGbp)测序序列集中的较长的30x覆盖率序列的序列集,该测 序序列集中待组装序列平均长度为14. 5化P。错点序列来自:利用80x覆盖率的Illumina Miseq测序序列集,W Sparse Assembler的跳跃k-mer方法化-mer长度为5化P ;sparse步 长g = 20)进行组装,得到约二百二十万条错点序列,错点序列平均长度约为1.化bp,错点 序列总长度约为2. 6Gbp,对待测基因组序列的覆盖率约为85%。上述预组装80x覆盖率的 IIlumina (第二代测序)的测序序列集得到错点序列的步骤花费34个CPU小时。
[0119] 之后基于本申请的DBG2化C,利用第二代Illumina预组装的错点序列,W k-mer 长度为17bp、阔值为(1%*每条错点序列长度)的参数将错点序列比对到上述较长的 30x覆盖率序列的序列集中(平均长度14.5化P),从而得到压缩的待组装序列索引表的集 合,运个步骤花费了 3个CPU小时,而17-mer的存储共占用了 70GB的内存。实验结果证 明,即使是在54x覆盖率的完全序列集中将错点序列比对到所有待组装序列上的过程也 只需花费少于6小时(数据未在表3中显示)。最后化-mer长度1,步长g= 1)的对基 因序列组装骨架进行修正的步骤花费了大约2000个CPU小时。在最初的报导化ttp:// blog. P曰cificbiosciences. com/2014/02/d3t3-rele3se-54x-long-re3d-cover3ge-for. html)中,构建重叠图谱的过程花了 405000个CPU小时,而运个过程可W使用局部敏感性 拼接(locality sensitive h曰sh) (Berlin, Κ. et 曰1. Assembling Large Genomes with Single-Molecule Sequencing and Locality Sensitive Hashing. (2014))的方法减少 大约一半的时间。而在另一个直接而高度优化的实施方式中(Myers, G. in Algorithms in Bioinformatics, Vol. 8701. (eds. D. Brown&B. Morgenstern), p52-67(Springer Berlin 化ide化erg, 2014)),能够在利用更多的内存的基础上进一步减少所需时间至15600CPU小 时。
[0120] DBG20LC最终的组装结果的质量可W与目前先进的组装工具使用更大量资源得到 的结果相媳美。为了评估结果的准确性,W下表格中将组装结果与人类参考基因组染色体 14进行比对,从而得到同一性和修正NG50的比对结果。 阳12U 表3. Η. sapiens基因组组装性能比较(基因组大小:3. 0G bp) 阳1。]
[012引由于测序序列集规模的问题,表中同一性及修正NG50仅使用染色体14的序列集 进行计算。 阳124] 实施例4: 阳1巧]对基于化化rd Nanopore平台的第Ξ代测序序列集进行组装
[0126] DBG20LC可W正常用于计算错误率30% W上的序列。然而,为了保证准确的整合 结果则需要更高的测序深度。为了避开运个问题,在实验中给高精度序列(例如第二代测 序中的contig)分配了更高的权重。
[0127] 在本实施例中,待组装序列集是基于Nanopore测序平台的30x覆盖率的大肠杆 菌K12基因组(基因组大小:4.6Mbp)测序序列集,该测序序列集中待组装序列平均长度 为6. 5化P。错点序列来自:利用30x覆盖率的Illumina Miseq测序序列集,^Sparse Assembler的跳跃k-mer方法进行组装,得到约2500条错点序列,错点序列平均长度约为 1. 7化P,错点序列总长度约为115Mbp,对待测基因组序列的覆盖率约为97%。基于本申请 的DBG2化C,利用上述组装结果作为错点,对30x覆盖率的Nanopore测序序列集进行组装, Wk-mer长度为15bp、阔值为(2%*每条错点序列长度)的参数将错点序列比对到待组装 序列从而得到压缩的待组装序列索引表的集合,之后修正压缩的待组装序列中的错误,并 且完成组装骨架的构建。而之后用Sparc (k-mer长度1,步长g = 1)运行了 3轮对基因序 列组装骨架进行修正的步骤,结果列于表4中。 阳12引表4. E. coli K12基因组组装性能比较(基因组大小:4, 6M bp.) 阳 129]
阳130] 实施例5: 阳13U 对基于Illumina平台的第二代测序数据进行组装
[0132] 为了评估短序列集组装工具在处理通过精确的第二代测序技术获得的长序列时 的性能,下载了针对大肠杆菌E. coli K-12菌株MG1655克隆值NAnexus数据库(具体可 参考 ht1:p://sra. anexus.com)登录号:SRA073308)的 50x 覆盖率的 Illumina Miseq 的 测序序列集。在本实验中,将本申请所述的组装工具对较长的Illumina测序序列集进行 组装的结果与其他几个常用组装工具(包括SGA(Simpson, J. T. &Durbin, R. Efficient de novo assembly of large genomes using compressed data structures. Genome researc h22, 549-556(2012)), S0APdenovo2(Luo, R. et al. S0APdenovo2:an empirically improved memory-efficient short-read de novo assembler. GigaSciencel, 18(2012)), SPAdesS ( Bankevich, A. et al. SPAdes:a new genome assembly algorithm and its applications to single-cell sequencing. Journal of computational biology:a journal of computational molecular cell biologyl9, 455-477(2012))的组装结果进行了 比较。
[0133] 有意思的是,对于在第Ξ代测序技术中生成的相对较长的测序序列集,传统的具 有固定k-mer长度的DBG的组装结果并非最优的,其已经暴露了难W分辨重复区域并且当 使用较大的k值时需要更高的精准度及更高覆盖率的测序序列集的缺点。运对矛盾的要 求,使得它很难获得有限的计算资源的最优组合。迭代DBG W更多的计算时间、费用和更复 杂的算法设计与实施为代价,部分解决了上述问题。例如,必须使用一个误差校正程序W生 成长且正确的k-mer字串。SGA利用FM-index索引方法寻找准确匹配,运也带来了对序列 质量的限制。
[0134] 相反的,本申请的算法是稳定的并且对质量的限制更松散,因此其能够更有效及 准确地寻找重叠区域。本申请的DBG20LC使用Sparse Assembler组装程序设置k = 31来 组装最初的错点序列。使用所述错点序列及待组装序列集计算压缩序列。
[0135] 其中错点序列来自:利用50x覆盖率的长度为150bp的Illumina Miseq测序序列 集,W SparseAssembler的跳跃k-mer方法进行组装,得到约400条错点序列,错点序列平 均长度约为13化P,错点序列总长度约为4. 6Mbp,对待测基因组序列的覆盖率约为100%。 基于本申请的DBG2化C,利用上述组装结果作为错点序列,对50x覆盖率的11 lumina测序序 列集进行组装,Wk-mer长度为3化P、阔值为0的参数将错点序列比对到待组装序列从而得 到压缩的待组装序列索引表的集合,并且完成组装骨架的构建。对于Illumina数据,由于 精度高,组装骨架即可作为最终结果。在运个例子里错点序列和待组装序列来自同一数据 集。主要原因是传统的k-mer法仅仅能解决比k短的重复区域,并不能解决基因组中比k 长,但比读长短的所有重复区域。所W虽然错点序列平均长度很长,但仍有大量重复区域不 能解开。DBG20LC则用错点压缩待组装序列,并最大程度化解了运些比读长短的重复区域。 阳136] 在此序列集中,本申请所述方法中的重叠图谱的构建只花费了大概1秒的时间, 而其中序列的压缩花费了大部分的计算时间-大约两分钟。 阳137] 表5.对使用Miseq读取的E. coli基因组序列的组装性能比较(基因组大小:4. 6M bp) 阳13引
阳139] 实施例6:
[0140] 在本实施例中,本申请中描述的Sparc整合算法在各种从细菌基因组到哺乳动物 基因组大小的数据集中被验证。Sparc整合算法在化cBio数据集(http://schatzl油. cshl.edu/data/ectools/)中及在 Oxford Nanopore 数据集化ttp://gigadb.org/ dataset/100102)中的验证结果如下。
[0141] 与Sparc整合算法最为相似的现有算法是PBdagcon,W下将显示在上述两个 数据集中PBdagcon算法与Sparc整合算法的比较结果。基于完全相同的数据输入,使 用如本申请所述的方法首先利用DBG20LC生成组装骨架,最后使用MUMmer 3(Kudz,S. et al. Versatile and open software for comparing large genomes. Genome BiologyS,R12 (2004))的化adiff功能来计算最终的整合同一性。所有的实验均在具有AMD 化teron 2425肥CPUs处理器(在800MHz的频率下)的工作站中进行。
[0142] Sparc被设计为允许使用一种(如第二代或第Ξ代测序)或多种混合测序数据 (如同时使用第二代测序W及第Ξ代测序序列)并且赋予相应的路径更高的权重。将整合 的结果再作为输入进行迭代再运行整合运算有助于提高同一性。在混合数据集的整合实验 中,使用Sparse Assembler组装工具组装的Illumina组装错点序列并且将边的权重增加 至5。在化cBio数据集中,Wk= l,g= 1运行两轮整合算法,第一轮及第二轮后的逐碱 基准确率在表6中W同一性1及同一性2来表示。在第一个实验中,使用化cBio中得到的 E. coli测序数据集并且在不同的测序深度情况下对同一性进行检测。在10X/30X深度的测 序数据中使用DBG20LC得到的最长组装骨架分别为1. 3Mbp和4. 6Mbp。在现有数据库中可 W找到登记号为NC_000913,长度为4. 6Mbp的E. coli参考基因组。在混合数据集中,仅使 用lOx深度的数据Sparc整合算法就能够达到99. 91 %的同一性,而当使用30x深度的数据 时其能达到更好的整合质量巧9. 98%的同一性)。而在化cBio数据集中得到的同一性只 略低于在混合数据集中得到的结果。
[0143] 表6.在化cBio中得到的E. coli测序数据集中的计算结果
[0144]
[0145] 使用在化cBio中得到的测序深度为20x的A. thaliana(基因组大小为120Mbp) 测序数据集,检测Sparc整合算法在计算更大数据容量时的性能。使用DBG20LC得到的 最长组装骨架为 7. IMbp。WHGAP(Qiin,C.S. et al.Nonhybrid, finished microbial genome assemblies from long-read SMRT sequencing data. Nature methods 10, 563-569(2013))方法得到的纯化cBio数据集基因组组装被作为参照W计算同一性。
[0146] 表7.在化cBio中得到的A. thaliana测序数据集中的计算结果 阳 147]
[0149] 在化化rd Nanopore数据集中,基于其更高的错误率,W k = 2, g = 2运行四轮整 合算法,第一轮及第四轮后的逐碱基准确率在表8中W同一性1及同一性2来表示。在混 合数据集中与在化化rd Nanopore数据集中使用Sparc与PBdagcon整合算法的结果列于 表8中。运用Sparc整合算法在具有约40%原始错误率的化化rd Nanopore数据集中能够 得到高于99. 5%的同一性,其比在混合数据集中得到的结果同一性更高。在此实验中,得到 的最长组装骨架为4.6Mbp。
[0150] 表8.在化化rd Nanopore中得到的E. coli测序数据集中的计算结果
[0151]
阳152] 在30x PacBio E.coli的混合数据集中使用不同的参数对内存及质量进行比较, 如表9中所示,除了在设定更高的g值能够降低内存的消耗W外,其他参数的变化并未能引 起结果的较大改变。 阳153] 表9.使用不同的参数对内存及质量进行比较 阳 154]
[01巧]上述实验的结果显示在使用Sparc整合算法时,一般可设定k = 1-3, g = 1-5。使 用混合数据集可W减少对第Ξ代测序数据深度的需求,运将大大降低测序成本。由于第二 代测序具有更高的准确率,将第二代测序的结果用于整合将有助于提高基因组装的质量。
[0156] 尽管在此公开了本申请的各个方面和实施例,但其他方面和实施例对于本领域技 术人员而言也是显而易见的。在此公开的各个方面和实施例仅用于说明目的,而非限制目 的。本申请的保护范围和主旨仅通过后附的权利要求书来确定。
[0157] 同样,各个图表可W示出所公开的方法和系统的示例性架构或其他配置,其有助 于理解可包含在所公开的方法和系统中的特征和功能。要求保护的发明并不限于所示的示 例性架构或配置,而所希望的特征可W用各种替代架构和配置来实现。除此之外,对于流程 图、功能性描述和方法权利要求,运里所给出的方框顺序不应限于W同样的顺序实施W执 行所述功能的各种实施例,除非在上下文中明确指出。
[0158] 除非另外明确指出,本文中所使用的术语和短语及其变体均应解释为开放式的, 而不是限制性的。在一些实例中,诸如"一个或多个"、"至少"、"但不限于"运样的扩展性词 汇和短语或者其他类似用语的出现不应理解为在可能没有运种扩展性用语的示例中意图 或者需要表示缩窄的情况。
【主权项】
1. 一种用于组装基因序列的方法,包括: 获取一待组装序列集的数据和多个锚点序列的数据; 利用所述多个锚点序列,对所述待组装序列集中的每一个待组装序列进行压缩以获得 该待组装序列的索引表,其中所述索引表以锚点序列的标号来表示;以及 利用所述待组装序列集中各待组装序列压缩得到的所述索引表构建组装骨架。2. 如权利要求1所述的方法,其中,所述压缩以获得待组装序列的索引表包括: 确定每个锚点序列在待组装序列中的排列位置,并按照锚点序列在待组装序列上的排 列位置,获得对应于该待组装序列的索引表。3. 如权利要求2所述的方法,其中,所述确定每个锚点序列在待组装序列中的排列位 置包括: 将所述锚点序列和每一个待组装序列分别打断成一系列核苷k联体,将所述锚点序列 的核苷k联体分别与所述待组装序列的核苷k联体进行比对,当所述待组装序列的核苷k 联体中与锚点序列的核苷k联体的匹配率超过预定阈值,则认为所述锚点序列可以匹配至 该待组装序列; 通过所述匹配的核苷k联体在所述锚点序列上的位置及所述匹配的核苷k联体在所述 待组装序列上的位置之间的相互对应关系,确定匹配的锚点序列在待组装序列中的排列位 置。4. 如权利要求3所述的方法,其中,所述核苷k联体的长度为:15bp-51bp、15bp-25bp、 或 31bp_51bp。5. 如权利要求3所述的方法,其中,所述阈值为所述锚点序列的长度的0. 1%至2%。6. 如权利要求1所述的方法,其中,所述构建组装骨架包括: 将所述待组装序列集中各待组装序列压缩得到的所述索引表进行一一比对,筛选出具 有最优重叠的待组装序列子集,以构建待组装序列间的最优重叠图谱,然后根据所述子集 中的待组装序列中锚点序列标号的位置对应关系对所述子集中的各待组装序列进行排布, 从而形成组装骨架。7. 如权利要求1所述的方法,其中,所述构建组装骨架包括: 选取所述待组装序列集中的一条序列作为当前序列,将所述当前序列的索引表与待 组装序列集中的所有其它序列的索引表进行一一比对,找出与所述当前序列具有重叠关系 的备选序列,将备选序列与当前序列进行索引表比对或碱基比对以找出与当前序列具有最 优重叠和延伸的序列,并将其相对于当前序列延伸的部分加到当前序列上即获得当前序列 的扩展序列,通过将所述待组装序列集中的每条序列作为当前序列重复上述步骤以得到一 组扩展序列,然后根据所述扩展序列中锚点序列标号的位置对应关系对各扩展序列进行排 布,从而形成组装骨架。8. 如权利要求1所述的方法,进一步包括:在获得每条待组装序列的索引表后,利用所 述锚点序列对待组装序列的索引表进行反向检索,得到与所述锚点序列相对应的待组装序 列的索引表子集,通过将所述索引表子集中所有其他序列的索引表比对到当前待组装序列 的索引表上来修正所述待组装序列的索引表子集中的每条待组装序列的索引表。9. 如权利要求8所述的方法,其中,所述通过比对来修正所述待组装序列的索引表子 集中的每条待组装序列的索引表包括: 对待组装序列的索引表子集中每条待组装序列的索引表,将其作为当前序列索引表与 所述索引表子集中所有其他待组装序列的索引表进行一一比对,根据每个位置上不同的锚 点序列标号出现的频次进行投票,利用每个位置的投票结果对每条所述当前序列索引表进 行修正,从而获得修正后的待组装序列的索引表。10. 如权利要求9所述的方法,其中,所述对待组装序列的索引表进行修正包括: 将比对结果中在待组装序列的索引表子集中出现不多于5次、或不多于4次、或不多于 3次、或不多于2次、或不多于1次的锚点序列标号从包含其的序列的索引表中删除; 将存在断点的杂合序列删除;以及 将比对结果中的被包含的序列删除。11. 如权利要求9所述的方法,其中,所述构建组装骨架包括:将所有所述修正后的待 组装序列的索引表子集中的索引表进行一一比对,筛选出具有最优重叠的修正后的待组装 序列的子集,然后根据所述子集中的修正后的待组装序列索引表中锚点序列标号的位置对 应关系对待组装序列进行排布,从而形成组装骨架。12. 如权利要求1-11中任一项所述的方法,进一步包括对构建的组装骨架进行如下所 述的再修正:在构建组装骨架后,将所有待组装序列与组装骨架进行比对,对于组装骨架上 的任一位点,通过计算所有可比对至该位点的各条序列中在该位点的每种碱基出现的频次 进行投票,利用每个位点的投票结果对组装骨架进行再修正,从而得到再修正后的基因序 列组装结果。13. 如权利要求12所述的方法,其中,所述投票是通过利用跳跃式核苷k联体实现的。14. 如权利要求12-13中任一项所述的方法,其中,所述对组装骨架进行再修正包括: 获得覆盖待组装基因序列不同区域的组装骨架,利用多个计算机处理器同时分别修正 所述覆盖待组装基因序列不同区域的组装骨架,以及将再修正后的组装骨架进行整合得到 完整的基因序列组装结果。15. 如权利要求1-14中任一项所述的方法,其中,所述待组装序列集包括通过第二代 基因组测序技术测序得到的序列。16. 如权利要求1-15中任一项所述的方法,其中,所述待组装序列集中的序列长度 在100bp-120kbp之间、或者在500bp-120kbp之间、或者在500bp-100kbp之间、或者在 lkbp-100kbp之间、或者在100bp-lkbp之间、或者在100bp_300bp之间。17. 如权利要求1-16中任一项所述的方法,其中,所述待组装序列集包括通过第三代 基因组测序技术所得到的序列,并且所述待组装序列集中的序列长度在500bp-100kbp之 间。18. 如权利要求1-17中任一项所述的方法,其中,所述锚点序列的长度在2bp-30bp之 间、或者在30bp-100bp之间、或者在100bp-300bp之间、或者在30bp-100kbp之间、或者在 300bp-100kbp 之间。19. 如权利要求18所述的方法,其中,所述锚点序列选自下组:精度大于95%的序列 长度在100bp-300bp之间的短序列、基因组高保守区域内的精度大于95%的序列长度在 100bp-300bp之间的短序列和待组装序列集中序列长度在lkbp-120kbp之间的长序列的局 部截取。20. -种用于组装基因序列的系统,包括: 输入组件,其用于获取一待组装序列集的数据和多个锚点序列的数据; 存储器,其用于存储所述待组装序列集的数据和多个锚点序列的数据;以及 处理器,其被配置于: 利用所述多个锚点序列,对所述待组装序列集中的每一个待组装序列进行压缩以获得 该待组装序列的索引表,其中所述索引表以锚点序列的标号来表示;以及利用所述待组装 序列集中各待组装序列压缩得到的所述索引表构建组装骨架。21. 如权利要求20所述的系统,其中,所述处理器被进一步配置于: 在获得每条待组装序列的索引表后,利用所述错点序列对待组装序列的索引表进彳丁反 向检索,得到与所述锚点序列相对应的待组装序列的索引表子集,通过将所述索引表子集 中所有其他序列的索引表比对到当前待组装序列的索引表上来修正所述待组装序列的索 引表子集中每条待组装序列的索引表。22. 如权利要求20所述的系统,其中,所述处理器被进一步配置于对构建的组装骨架 进行如下所述的再修正: 在构建组装骨架后,将所有待组装序列与组装骨架进行比对,对于骨架上的任一位点, 通过计算所有可比对至该位点的各条序列中在该位点的每种碱基出现的频次进行投票,利 用每个位点的投票结果对组装骨架进行再修正,从而得到再修正后的基因序列组装结果。23. 如权利要求22所述的系统,其中,所述处理器被进一步配置于: 获得覆盖待组装基因序列不同区域的组装骨架,利用多个计算机处理器同时分别修正 所述覆盖待组装基因序列不同区域的组装骨架,以及将再修正后的组装骨架进行整合得到 完整的基因序列组装结果。
【文档编号】G06F19/22GK105989249SQ201510084489
【公开日】2016年10月5日
【申请日】2015年2月16日
【发明人】叶承羲
【申请人】叶承羲