CN112786110A - 一种序列组装方法及系统 - Google Patents

Abstract

本发明公开了一种序列组装方法及系统，对待组序列进行比对，获得首尾重叠比对序列；根据首尾重叠比对序列构建第一字符串，对第一字符串的目标特征进行处理，获得第二字符串图，将第二字符串图转换为序列组装结果。目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种。本发明待组装序列可为经矫正的序列，将其作为字符串构图的数据基础，并且对字符串图进行处理，使得处理后的字符串图更加满足序列组装的需求，提升组装结果的准确性、连续性以及实现快速组装的目的。

Description

一种序列组装方法及系统

技术领域

本发明涉及生物信息技术领域，特别是涉及一种序列组装方法及装置。

背景技术

生物的遗传物质主要包括DNA(脱氧核糖核酸)、RNA(核糖核酸)，这些遗传物质是生物生存、繁殖、演化中的关键。生物物种遗传物质的测定和分析，是了解生命存在和演化的必要手段。目前，主要是通过各种测序技术得到大量的遗传物质片段信息，而根据所使用的测序技术的不同，所获得的片段大小从几百bp到几Mb不等。相对于几kb到几Gb不等的遗传物质大小，这些片段是零碎的，需要通过生物信息的方法进行拼接、组装，重构出物种的完整遗传信息。

DNA测序技术经过了基于终止法的第一代测序技术，边合成边测序的第二代测序技术，以及单分子测序的第三代测序技术。第三代测序技术以长度长为显著技术特点，现主流商业化平台出自PacBio和Oxford Nanopore Technologies。但是，在对长度长序列进行组装时，存在数据处理的缺陷，无法实现高准确性、高连续性的组装，并且需要消耗大量的计算资源。

发明内容

针对于上述问题，本发明提供一种序列组装方法及系统，实现了提升组装结果准确性、连续性以及快速组装的目的。

为了实现上述目的，本发明提供了如下技术方案：

一种序列组装方法，包括：

S101、对待组装序列进行比对，获得首尾重叠比对序列；

S102、根据所述首尾重叠比对序列构建第一字符串图；

S103、对所述第一字符串图的目标特征进行处理，获得第二字符串图，其中，所述目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种；

S104、将所述第二字符串图转换为序列组装结果。

可选地，步骤S101中所述待组装序列为经矫正的序列。

可选地，步骤S102中的所述首尾重叠比对序列为经过过滤处理后的首尾重叠比对序列，所述过滤处理同时包括对嵌合体序列过滤、对被包含序列过滤、对冗余首尾重叠比对序列过滤。

可选地，步骤S103中的所述目标特征包括满足预设条件的边和/或节点，其中，对所述满足预设条件的边和/或节点进行处理，包括：

根据覆盖度、出度，删除高深度边和/或节点；

删除传递简化对应边及反向边；

根据比对长度，删除低数值边；

根据覆盖度、出度、比对长度，保留高数值边和/或节点；

删除冗余边；

删除嵌合体边。

可选地，在S104步骤前，所述方法还包括：

对所述第一字符串图或第二字符串图所含的子图进行处理；

所述对所述第一字符串图或第二字符串图所含的子图进行处理，包括：

在所述第一字符串图或所述第二字符串图所含的子图内进行搜索，根据匹配碱基数对搜索路径评分，保留评分最高路径。

可选地，所述方法还包括：

S105、根据比对深度，对所述序列组装结果连接点进行验证，以使打断组装错误的连接点。

可选地，所述方法还包括：

对所述序列组装结果进行纠错，得到纠错后的序列组装结果。

一种序列组装系统，包括：

比对单元，用于对待组装序列进行比对，获得首尾重叠比对序列；

构建单元，用于根据所述首尾重叠比对序列构建第一字符串图；

处理单元，用于对所述第一字符串图的目标特征进行处理，获得第二字符串图，其中，所述目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种；

转换单元，用于将所述第二字符串图转换为序列组装结果。

相较于现有技术，本发明提供了一种序列组装方法及系统，对待组序列进行比对，获得首尾重叠比对序列；根据首尾重叠比对序列构建第一字符串，对第一字符串的目标特征进行处理，获得第二字符串图，将第二字符串图转换为序列组装结果。目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种。本发明待组装序列为经矫正的序列，将其作为字符串构图的数据基础，并且对字符串图进行冗余处理，使得处理后的字符串图更加满足序列组装的需求，提升组装结果的准确性、连续性以及快速组装的目的。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据提供的附图获得其他的附图。

图1为本发明实施例一提供的一种序列组装方法的流程示意图；

图2为本发明实施例二提供的一种序列矫正的方法的流程示意图；

图3为本发明实施例三提供的一种序列比对方法的流程示意图；

图4为本发明实施例六提供的一种序列组装结果纠错方法的流程示意图；

图5为本发明实施例九提供的一种序列组装系统的结构示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明的说明书和权利要求书及上述附图中的术语“第一”和“第二”等是用于区别不同的对象，而不是用于描述特定的顺序。此外术语“包括”和“具有”以及他们任何变形，意图在于覆盖不排他的包含。例如包含了一系列步骤或单元的过程、方法、系统、产品或设备没有设定于已列出的步骤或单元，而是可包括没有列出的步骤或单元。

为了便于对本发明实施例进行说明，对本发明实施例中应用到的术语进行说明。

测序技术：测定遗传物质序列的技术，例如测定组成DNA分子的碱基序列。

第三代测序：是一种测序技术。在测序时，不经过PCR扩增，直接对每条分子进行单分子从头测序。代表性的技术有PacBio的SMRT技术和Oxford Nanopore Technologies的纳米孔测序技术。

基因组：指生物体所有遗传物质的总和。

测序序列(read，一般也可称为读段)：通过测序技术，得到的由碱基组成的一段序列信息。

序列比对：将两段序列根据具体的碱基进行排列，标明匹配、错配、间隔等。两段序列匹配的碱基数达到一定阈值则认为两序列能比对上。

重叠关系(overlap)：如果两序列能够比对上，存在重叠区域，则称两序列存在重叠关系。

被包含关系：如果两条序列存在重叠关系，且一序列(即包含序列)覆盖另一序列(即被包含序列)，当被包含序列两端未比对上的区域长度都小于阈值，则两序列的重叠比对关系称为被包含关系。

嵌合体(chimera)：测序过程中，两段来源于基因组不同位置的序列由于外界因素(例如实验建库)连接在一起，形成的非自然存在的序列，即为嵌合体。

首尾重叠关系(dovetail)：如果两序列的重叠关系位于这两条序列不同末端，且末端未比对上区域长度均小于阈值，则两序列的重叠关系称为首尾重叠关系。含有首尾重叠关系的比对称为首尾重叠比对。本领域通常是以5'-3'方向书写碱基序列，存在首尾重叠关系的序列具体而言是，上游序列的3’端区域与下游序列的5’端区域存在重叠关系。其中，5'端指的是碱基序列中带有游离5′-羟基或其磷酸酯的末端。3'端指的是碱基序列中带有游离3′-羟基或其磷酸酯的末端。

有向图：图由一个非空有限顶点集和一个非空有限边集组成，以集合G表示。图G的顶点集和边集分别用V和E表示，则图G可表示成G＝(V,E)。如果图的每条边都用箭头指明了方向，则为有向图。

出边：有向图中，从某个顶点起始的所有边，为该顶点的出边。

入边：有向图中，到某个顶点终止的所有边，为该顶点的入边。

出度：有向图的顶点的出边条数称为该顶点的出度。

入度：有向图的顶点的入边条数称为该项点的入度。

出度和入度统称为度。

字符串图(StringGraph)：是根据序列的首尾重叠关系构建的有向图。顶点为序列的某一端，也称为节点。边表示测序序列间的首尾重叠关系。

反向节点：字符串图中，存在首尾重叠关系的每条序列的5'端和3'端各有一个节点，同一条序列的两个节点互为反向节点。

反向边：字符串图中，每一个首尾重叠关系会产生两条边，且边的方向相反，这两条边互为反向边。

直接前任节点：字符串图中，某个节点的所有入边的另一端的节点，为直接前任节点。

直接后继节点：字符串图中，某个节点的所有出边的另一端的节点，为直接后继节点。

前任节点：字符串图中，某个节点的直接前任节点，及继续往前逐步搜索的所有直接前任节点的并集，为该节点的前任节点。

后继节点：字符串图中，某个节点的直接后继节点，及继续往后逐步搜索的所有直接后继节点的并集，为该节点的后继节点。

宽度优先搜索(BFS，BreadthFirstSearch)：是从某节点开始，沿字符串图的宽度遍历图的节点，如果发现目标，则搜索终止。

重叠群(contig)：根据序列间的重叠关系拼接得到的连续的、长序列。

传递简化(Transitive reduction，也称为传递归约)：字符串图的简化方法。例如，如果有向图中节点A根据首尾重叠关系往后延伸最长可延伸至节点C，且存在A→B→C、A→C的不同路径，若A→B→C路径的总延伸长度小于等于A→C路径的总延伸长度与长度波动阈值之和，则可将A→C边删除，实现简化。

末端节点：如果某个节点只有出边，且所有直接后继节点都有其他直接前任节点，则认为该节点为末端节点。反过来，如果某个节点只有入边，且所有直接前任节点都有其他直接后继节点，则也认为该节点为末端节点。

简单路径：从某个末端节点开始往后搜索，直至搜索到的节点入度或者出度不为1为止(搜索停止的节点即为终止节点)，认为从末端节点到终止节点间形成的路径为简单路径。

支线路径：从只有一条出边的某个节点开始，往后搜索形成一条简单路径，搜索深度小于阈值，且终止节点入度大于1，则认为该路径为支线路径。

Z型路径：从出度大于1的末端节点开始，往某一个直接后继节点搜索形成一条简单路径，搜索深度小于阈值，且终止节点的入度大于1,则认为从末端节点到终止节点间的简单路径为Z型路径。

泡状结构：从出度大于1的末端节点开始，往每一个直接后继节点搜索形成简单路径，搜索深度小于阈值。如果有两条或以上简单路径的终止节点为同一个节点，则认为这些简单路径形成的图形为泡状结构。

复合路径：从出度大于1的末端节点开始往后搜索，直到搜索到的节点出度为0或者深度大于阈值，则终止。如果终止节点的个数为1，则所有遍历到的边和节点形成的路径为复合路径。

实施例一

参见图1，其示出了本发明实施例一提供的一种序列组装方法的流程示意图，该方法可以包括以下步骤：

S101、对待组装序列进行比对，获得首尾重叠比对序列。

待组装序列可以是测序序列(read)，也可以是经矫正的序列，对序列的矫正处理会在本发明的后续实施例二、三中进行详细说明。待组装序列可以是全基因组的序列，也可以是经筛选后的部分基因组序列，例如某条或某几条染色体、染色体的某段或某几段区域。

S102、根据所述首尾重叠比对序列构建第一字符串图。

比对结果中包括了一定的错误信息、冗余信息，因此对其进行过滤处理，挑选高质量、有限数量的首尾重叠比对用于构图，有助于节约计算空间、缩短计算时间并提高构图的准确性。具体过滤处理在本发明的后续实施例四中进行详细说明。

在获得了首尾重叠比对后，根据首尾重叠比对构建序列末端间的第一字符串图，所述第一字符串图的节点为序列的5'端或者3'端，边表示序列末端间的首尾重叠关系。所述字符串图边的属性信息包括边的长度(length属性)、比对长度(score属性)、比对相似性(identity属性)和匹配碱基数(match属性)。所述边的长度为基于首尾重叠比对序列往前延伸的长度，比对长度为首尾重叠比对中重叠区域的长度，比对相似性为首尾重叠比对中重叠区域的匹配碱基数占比，匹配碱基数为首尾重叠比对中重叠区域的匹配碱基对数目。

S103、对所述第一字符串图的目标特征进行处理，获得第二字符串图。其中，所述目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种。

构建的第一字符串图只是一个初始的字符串图，其中可能会包括过多的错误、冗余信息，因此，在本发明实施例中需要对第一字符串图进行处理，得到第二字符串图，该处理包括对特定边、节点的处理以及特殊路径的处理中的一种或多种，可以基于实际情况进行选定。具体的，可以包括对边/节点进行筛选，对支线路径、Z型路径、泡状结构、复合路径进行处理。具体处理方法在本发明的后续实施例五中进行详细说明。

若所构建的第一字符串图或第二字符串图含有子图，可进一步对子图进行处理，包括：子图内进行搜索，根据匹配碱基数对搜索路径评分，保留评分最高路径。

S104、将所述第二字符串图转换为序列组装结果。

在本发明实施例中还包括对序列组装结果的验证，即所述方法还包括：

S105、根据比对深度，对所述序列组装结果连接点进行验证，以使打断组装错误的连接点。

在一种更优的实施方式中，还包括：

对所述序列组装结果进行纠错，得到纠错后的序列组装结果。

纠错的具体方法在本发明的后续实施例六中进行详细说明。

本发明实施例一提供了一种序列组装方法，对待组序列进行比对，获得首尾重叠比对序列；根据首尾重叠比对序列构建第一字符串，对第一字符串的目标特征进行处理，获得第二字符串图，将第二字符串图转换为序列组装结果。目标特征包括满足预设条件的边和/或点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种。本发明待组装序列可为经矫正的序列，将其作为比对的数据基础，进一步对比对结果进行过滤处理，使得过滤后的比对结果更加满足构图的需求，配合对字符串图的处理，实现准确、连续以及快速的组装。

实施例二

在实施例一种的序列组装方法步骤S101中的待组装序列可以为经矫正的序列，即通过序列矫正(即Correct)的方法实现对待组装序列的矫正。矫正后的序列的准确性可明显提高，有助于进一步提升组装结果的准确性。参见图2，其示出了本发明实施例二提供的一种序列矫正的方法的流程示意图，该序列矫正的方法可以包括以下步骤：

S001、对多个待组装序列片段之间进行比对，得到参考序列、比对序列；

S002、以k-mer为单位，对所述参考序列进行划分，对所述比对序列进行对应的划分，得到多个k-mer片段；

S003、从正序，以k-mer片段为单位，计算每个位置上所有碱基的得分；

S004、从反序，确定最末碱基为最末位置上最高得分碱基，根据所述最末碱基对应的所述k-mer片段，确定倒数第二位碱基；根据所述倒数第二位碱基对应的所述k-mer片段，确定倒数第三位碱基；依次类推，获得k-mer得分链；

S005、使用k-mer得分链替换所述参考序列，获得第一矫正结果；

S006、对所述第一矫正结果中的低质量序列，使用步骤S001比对结果中覆盖所述低质量序列的所述比对序列的比对区间序列作为种子序列，基于所述种子序列使用有向图算法进行矫正；使用矫正后的序列替换所述低质量序列，获得第二矫正结果；其中，所述低质量序列为满足以下所有条件的序列，以连续碱基数目不小于预设第一数值的低质量碱基为起始、区间平均质量值小于预设第二数值的最长序列、区间长度超过预设第三数值、区间内连续高质量碱基数目不大于预设第一数值。

可选的，步骤S002所述k值为3。

可选的，步骤S003所述计算每个位置上所有碱基得分的方法包括，根据公式score(p,b)＝max{score(p-1,b)+count_k-mer}-C×0.3计算每个位置上所有碱基的得分；其中，score(p,b)为所述参考序列p位置处碱基的得分，b为碱基A、T、G、C或缺失，count_k-mer为所述k-mer片段在所述比对序列中的出现次数；C为所述参考序列p位置处所述比对序列对所述参考序列的覆盖度。

可选地，步骤S006所述有向图算法进行矫正的方法，包括：

S061、任选一条所述种子序列作为基序绘制有向图；

S062、根据得分矩阵计算所述种子序列中碱基得分，更新有向图；

S063、对有向图进行拓扑排序；

S064、根据拓扑排序结果，并确定每个节点的碱基为该节点上出现概率最高的碱基，获得矫正后的序列。

可选的，序列矫正方法还包括：

S007、以种子序列为比对序列，以第二矫正结果为参考序列，使用前述步骤S001至S005的序列矫正方法，再次进行矫正，得到第三矫正结果。

本实施例所述矫正方法使用k-mer片段的方式。k-mer，monomeric unit(mer)，相当于nt或者bp，100mer DNA相当于每一条链有100nt，那么整条链就是100bp。一般长度为m的reads可以分成m-k+1个k-mers。以k-mer为3举例说明，以k-mer为3对上述待矫正序列进行划分，得到多个3-mer片段，比如位点1、2、3为第一个3-mer，位点2、3、4为第二个3-mer，依次类推。

本实施例二通过上述序列矫正方式先对待组装序列片段进行比对，然后再基于k-mer片段对序列获得矫正结果，可以避免因只统计单碱基比对所导致的矫正结果存在偏好性，从而提高序列矫正的准确性，还能够提高序列矫正的处理效率。此外，本实施例提供的上述方式是基于序列比对结果进行的矫正，不对获得比对结果的比对方法进行限制，可以与各种序列比对方法进行组合。

进一步的，考虑到长读长测序技术存在准确度分布不均匀的情况。以三代测序结果为例，某些区间准确度极低(<70％)，该区间的比对结果可靠性差，因此基于比对结果进行矫正获得的矫正序列准确度低。通过查找低准确度的区间、并针对该区间进行再次矫正(即步骤S006)，可进一步提高准确性。

在一种具体的实现方式中，步骤S006中低质量序列为，以连续碱基数目不小于预设第一数值(如4)的低质量碱基为起始、区间平均质量值小于预设第二数值(如40)的最长序列、区间长度超过预设第三数值(如16)、区间内连续高质量碱基数目不大于预设第一数值(如4)。

在一种具体的实现方式中，定义碱基质量值Q为Q＝count_k-mer÷C×100；将碱基质量值大于预设的第一数值4的碱基确定为低质量碱基，将碱基质量值不大于第一数值4的碱基确定为高质量碱基；区间平均质量值为该区间内所有碱基质量值的平均值。

可参照如下步骤(a)至(c)，基于低质量碱基和高质量碱基查找第一矫正结果中的低质量序列：

(a)对第一矫正结果按反序进行低质量碱基查找，记录第一个低质量碱基位置lq_s。

(b)以lq_s位置为起始位置，按反序继续低质量碱基查找，记录查找的连续低质量碱基数目lq_c，并判断lq_c与第一数值4的关系。

若lq_c不小于第一数值4，继续按反序查找，若查找到高质量碱基则记录连续的高质量碱基数目hq_c。若hq_c大于第一数值4，初始化lq_c、hq_c为0，按反序重新查找低质量碱基，重复步骤a、b。若hq_c不大于第一数值4，初始化hq_c为0，继续按反序查找，连续记录lq_c，执行步骤(c)；

若lq_c小于第一数值4，初始化lq_c为0，按反序重新查找低质量碱基，重复步骤a、b。

(c)计算区间平均质量值，当区间平均质量值为小于等于第二数值40的最大值时，判断区间长度lq_l是否不小于第三数值16。

若大于等于第三数值，记录反序查找的最后一位碱基位置lq_e为终止位置，标记lq_s位置与lq_e位置之间的序列为低质量序列。同时，初始化低质量区间的长度lq_l、低质量碱基数目lq_c及高质量碱基数目hq_c为0，按反序重新查找低质量碱基，重复步骤a、b、c，直至完成所有参考序列的查找。

若小于第三数值，初始化低质量区间的长度lq_l、低质量碱基数目lq_c及高质量碱基数目hq_c为0，按反序重新查找低质量碱基，重复步骤a、b、c，直至完成所有参考序列的查找。

在步骤S006中，有向图算法进行矫正的具体步骤为：

S061、任选一条种子序列作为基序绘制有向图。

在本实施例中，为了减少无效数据的比对，可以在进行有向图算法矫正前，按照以下任一条件对第一矫正结果进行过滤，即对满足以下任一条件的第一矫正结果不进行有向图算法矫正：a.低质量序列累计长度超过预设第四数值(如10000)的第一矫正结果；b.种子序列数目小于预设第一数值(如4)的第一矫正结果；c.长度超过预设第四数值(如10000)的比对序列占比超过第四预设倍数(如1/3)的第一矫正结果。过滤掉的数据以低可信度、超长为特征，可以解决有向图算法中因超长序列造成计算速度极慢的瓶颈问题。

具体的，遍历每一个低质量序列，利用该低质量序列的起点位置和终点位置，查找出所有覆盖该低质量序列的比对序列。如果查找到的比对序列大于60条，则停止查找低质量区间序列的比对序列。

进一步的，将查找到的比对序列按序列长度从长到短排序，选取种子序列。在一种较佳实施例中，种子序列数目不超过预设第五数值(如30)，优选为长度排序位于中间预设第五数值(如30)的比对序列包含的种子序列。若查找到的比对序列不足30条，则将所有的比对序列均作为种子序列S。

在种子序列S中任选一条作为基序绘制有向图G；其中，有向图的节点Node为序列中的碱基，有向图的有向边Edge为任意两个连续碱基之间的连线，有向图的起始节点为序列的起始碱基。

S062、根据得分矩阵计算种子序列中碱基得分，更新有向图，即将全部种子序列均加入有向图。

在一种具体的实现方式中，按以下步骤(1)至(4)执行：

(1)根据公式score(x)＝max{score(x-1)}-GP对节点得分进行初始化，其中score(x)为节点得分，x为节点索引，GP为空位罚分(如2)；设定无前任节点的节点得分为0。

(2)本实施例可以遍历有向图中的节点，基于每一个节点，遍历该节点的所有前任节点。对于任意一个前任节点，遍历S中的每一个碱基，并按照如下计算公式更新矩阵得分：

其中，x为节点索引，y为种子序列中碱基索引，GP为空位罚分(如2)，M为错配得分(如-2)或匹配得分(如1)。

(3)以出度为0的最高得分节点为起点，选择最高得分的前任节点，进行回溯，以节点索引为0或碱基索引为0为终点。

具体的，遍历有向图中的节点，如果节点的出度为0，则定义该节点为终节点，记录该终节点的得分和索引。找到所有的终结点后，查找出得分最大的终节点，将其确定为回溯起点，并记录该起始的节点索引和对应的碱基索引。从起点开始循环往回追溯，查找任意节点x的得分最高的前任节点和前任碱基，直到节点索引为0或者序列索引为0终止，将其确定为回溯终点，记录终点的节点索引，以及对应的碱基索引。

(4)若终点的碱基索引大于0，将种子序列中位于终点的碱基索引之前的碱基加入有向图；若起点的碱基索引小于种子序列长度，将种子序列中位于起点的碱基索引之后的碱基加入有向图。

具体加入有向图的操作为：如果当前碱基在当前节点下为匹配，则跳过该碱基；如果为错配，则将当前碱基加入到当前节点包含的碱基列表里；如果为插入，则以新节点和新边加入到有向图中。

S063、对有向图进行拓扑排序。

对于更新后的有向图还可以进行拓扑排序，使得更新后的有向图中任意一对顶点u和v，若满足边(u,v)∈E(G)，则u在排序后的有向图中出现在v之前。

S064、根据拓扑排序结果，并确定每个节点的碱基为该节点上出现概率最高的碱基，获得矫正后的序列。

在一种具体的实现方式中，按照拓扑排序结果，在每个节点上均选取出现次数最多的碱基，由此确定矫正后的低质量序列，替换第一矫正结果的对应序列，获得第二矫正结果。

为了进一步提高准确性，序列矫正方法还可以包括：

S007、以种子序列为比对序列，以第二矫正结果为参考序列，使用前述步骤S001至S005的序列矫正方法，再次进行矫正，得到第三矫正结果。具体细节不再赘述。

采用上述方法进行矫正时，考虑前向序列与当前碱基的组合情况，即考虑最可能存在的连续小片段作为矫正结果，有效地避免了只统计单碱基比对次数对于矫正结果的偏好性影响。

综上，上述实施例中所提供的序列矫正方法，可以实现序列高效、高准确度的矫正，有助于确定准确的首尾重叠比对用于后续构图。

实施例三

在本发明实施例三提供了一种序列比对的方法，该序列比对方法可以对实施例二中提供的序列矫正的方法中步骤S001中的对多个待组装序列之间进行比对。即序列比对可以对参考序列和比对序列进行多次筛选，可以减少无效比对数量、提高比对准确性，减少数据存储和运算空间。此外，本实施例三通过该比对方式提供待组装序列之间的对比结果，可以与各种基于比对结果进行矫正的矫正方法进行组合，对序列矫正方法不进行限制。

参见图3，其示出了实施例三提供的一种序列比对方法的流程示意图，所述序列比对方法可以包括以下步骤：

S011、获取待组装序列，根据第一预设数据长度在所述待组装序列中筛选第一序列，根据第二预设长度在所述待组装序列中筛选第二序列；

S012、将所述第一序列与第二序列进行比对，得到第一参考序列、第一比对序列；

S013、计算所述第一比对序列对所述第一参考序列的覆盖度，基于所述覆盖度筛选所述第一参考序列、第一比对序列，获得第二参考序列、第二比对序列；

S014、计算所述第二参考序列与第二比对序列之间的比对路径，基于编辑距离筛选所述第二参考序列、第二比对序列，获得第三参考序列、第三比对序列。

可选的，步骤S013所述覆盖度包括窗口覆盖度、整体覆盖度；其中，所述窗口覆盖度为，按照预设第一碱基数将所述第一参考序列划分为多个窗口，计算所述第一比对序列对每个窗口的覆盖度；所述整体覆盖度为，所述第一参考序列的平均覆盖度。

可选的，步骤S013所述基于所述覆盖度筛选的方法为，将满足预设第一覆盖条件的所述第一比对序列确定为所述第二比对序列，且所述第二比对序列对所述第一参考序列的所述整体覆盖度达到预设第一覆盖度；其中，所述第一覆盖条件为，所述窗口覆盖度小于预设第二覆盖度，且小于所述整体覆盖度的第一预设倍数。

可选的，步骤S013按照所述第一比对序列与所述第一参考序列重叠长度由长至短的顺序，进行所述覆盖度计算及所述筛选。

可选的，步骤S014所述基于编辑距离的筛选方法，包括：

I、使用贪婪算法计算所述第二参考序列与第二比对序列之间的比对路径，选择最优比对路径并确定其编辑距离；

优选的，贪婪算法计算时约束线范围为，在所述比对路径的编辑距离的约束下，最大延伸距离减去指定距离值时最小约束线与最大约束线之间；

II、从所述第二比对序列中，选择不大于预设编辑距离条件的第二比对序列为第三比对序列，其对应的参考序列为第三参考序列；其中，所述预设编辑距离条件为，所述比对区间内第一参考序列与第二比对序列长度之和的第三预设倍数；

III、循环步骤I至II，直至第三比对序列对第三参考序列的覆盖度达到预设第六覆盖度；

优选的，计算所述覆盖度的区间为，根据所述最优比对路径，回溯确定第三比对序列与第三参考序列的比对区间，查找所述比对区间内的第一个预设第三碱基数完全匹配的起始位置、最后一个预设第三碱基数完全匹配的终止位置，将所述起始位置、所述终止位置确定为所述区间的起始位置、终止位置。

参考序列与比对序列因序列相似性存在比对关系，具有比对关系的区间即为比对区间，即体现比对结果；其中，所述参考序列包括第一参考序列、第二参考序列和第三参考序列中的一种或多种，所述比对序列包括第一比对序列、第二比对序列和第三比对序列中的一种或多种。

在本实施例步骤S011中，具体实现可参考如下内容。在待组装序列中选去数据长度大于第一预设长度(如最长的40×序列长度)的序列作为第一序列，并将第一序列输出至第一序列文件；在待组装序列中选取第二预设长度(如高于1k)的序列作为第二序列，并将第二序列输出至第二序列文件。

为便于后续的多任务并行运算，可以将第一序列文件和第二序列文件进行切割，本实施例对所切割的份数不做限定。

在本实施例步骤S012中，基于上述第一序列和第二序列，提供了一种区间比对方法，将第一序列中所有序列与第二序列中所有序列进行两两比对，得到第一参考序列、第一比对序列。此步骤可以使用现有比对软件实现，例如minimap2。现有技术中通常仅过滤掉短序列(如小于1K)，对过滤后的序列之间全部进行比对，即第二序列与第二序列进行比对，这样就包括了第二序列中非第一序列之间的比对，此类序列较短(如最长的40×序列长度、1K之间)且数量多，其比对结果相对冗余。本实施例仅将第一序列与第二序列进行比对，可以避免冗余比对、提高比对效率，并避免冗余比对对后续矫正准确性的负影响。

在本实施例步骤S013中，为了有效减少用于矫正的序列的数据集，同时降低因重复序列、比对错误导致的多重比对以及比对噪音，基于覆盖度筛选第一参考序列、第一比对序列，获得比对结果更准确的第二参考序列、第二比对序列。同时，经覆盖度筛选后，获得的各第二参考序列的覆盖度水平相对一致，有利于后续的矫正。

在一种可能的实施方式中，步骤S013具体实现可参考如下内容：

①：计算第一比对序列与第一参考序列比对区间长度，将第一比对序列按比对区间由长至短排列；对于比对区间长度相等的第一比对序列，按照其编号由大到小进行排序。

②：对第一参考序列，按预设第一碱基数(如64)划分为多个窗口，并将其窗口覆盖度、整体覆盖度均初始化为0。

③：对当前第一参考序列的每个窗口，按照①排列顺序依次查找第一比对序列。如果该窗口的窗口覆盖度小于预设第二覆盖度(如50)，且低于当前第一参考序列整体覆盖度的第一预设倍数(如1.3)，则认为该第一比对序列有效，输出该第一比对序列，将与该第一比对序列对应的窗口的窗口覆盖度加1；反之，表示该第一比对序列冗余，舍弃。直至，当前第一参考序列整体覆盖度达到预设第一覆盖度(如40)，即表示该第一参考序列有足够的第一比对序列，舍弃后续所有的第一比对序列。

在本实施例步骤S014中，基于步骤S012的区间比对、S013的覆盖度筛选，进一步的在第二参考序列、第二比对序列的比对区间内进行更为精确的碱基比对。在区间比对结果上进行碱基比对，而非直接进行碱基比对，其原因在于，可基于更少但更准的数据进行更高效的比对。通过碱基比对，能够进一步保留可信的、剔除不可信的、去掉冗余的比对信息，进一步缩减了比对结果、提高了比对精度。

在一种可能的实施方式中，步骤S014具体实现可参考如下步骤(1)至(4)：

(1)利用Diff贪婪算法，在k∈[k_min，k_max]的约束范围内，计算各路径的编辑距离，选择最优比对路径并确定其编辑距离。其中k_min约束线为在当前编辑距离d的约束下，找到最大延伸距离的最大值减去指定距离值(如150)时最小的k线。k_max为约束线为在当前编辑距离d的约束下，找到最大延伸距离的最大值减去指定距离值(如150)时最大的k线。最优编辑距离为编辑操作次数最少时所对应的编辑距离。

(2)判断最优编辑距离是否大于预设编辑距离条件。预设编辑距离条件为(L1+L2)×0.4，其中，L1是第二参考序列的数据长度，L2是第二比对序列的数据长度，0.4为第三预设倍数。

若是，表示第二参考序列与第二比对序列为错误比对，舍弃当前最优编辑距离对应的第二比对序列；若否，表示第二参考序列与第二比对序列为正确比对，保留当前最优编辑距离对应的第二比对序列，即为第三比对序列，其对应的为第三参考序列。

(3)根据最优编辑距离回溯至比对区间，在比对区间内，从前往后搜索，将第一个第三参考序列与第三比对序列长度为预设第三碱基数(如8)完全匹配的起点重新设置为比对区间起点；从后往前搜索，将最后一个第三参考序列与第三比对序列长度为预设第三碱基数(如8)完全匹配的终点设置为比对区间终点。本步骤重新确定计算覆盖度的比对区间，可有效解决比对区间边界(也即端点)错误率较高的问题，有助于进一步准确的计算覆盖度。

(4)计算比对区间的覆盖度，直至第三比对序列对第三参考序列的覆盖度达到预设第六覆盖度(如90)。

综上，本实施例提供的序列比对方法，通过两次比对、多次筛选，有效的减少无效比对、剔除错误比对，实现快速、准确、高效的比对。

实施例四

在本发明实施例四中提供了一种构建第一字符串图的方法以及该方法应用到的子流程。

在本发明实施例中可以直接对待组装序列进行比对，来获得首尾重叠比对序列，也可以是对待组装序列进行切分，获得切分结果，然后基于切分结果进行比对，来获得首尾重叠比对序列。切分处理后可采取并行的方式进行计算，有助于计算速度的进一步提升。

在一种更佳实施方式中，所述首尾重叠比对序列为经过过滤处理后的首尾重叠比对序列，所述过滤处理同时包括对嵌合体序列进行过滤、对被包含序列过滤、对冗余首尾重叠比对序列过滤。经过滤处理后，实现了使用数目更少、比对结果更准确、比对质量更高的数据用于构图，减少构图所需的计算资源、并可提高图的准确性，避免非准确信息的干扰。

其中，被包含序列包括末端未比对上区域长度小于等于第一阈值且被包含次数大于等于第二阈值的序列。第一阈值为序列两端未比对上区域的长度阈值，可以设置为300bp；第二阈值为被包含次数阈值，可以设置为2。

在优选的实施方式中，还可以对被包含序列进行重判定，所述对所述被包含序列进行重判定，包括：

去除序列末端未比对上区域，进行重比对，若重比对结果中末端未比对上区域长度小于等于第一阈值且被包含次数大于等于第二阈值的序列确定为被包含序列。

具体的，重比对可以是基于在先比对结果进行计算，也可以利用比对软件自动完成重比对过程。

对应的，所述对冗余首尾重叠比对序列过滤，包括：

对不满足第一比对条件和第二比对条件的首尾重叠比对序列进行过滤。

其中，所述第一比对条件为比对相似性或者比对长度占最大比对相似性或最大比对长度的比例大于等于第三阈值，所述第二比对条件为首尾重叠区域内末端未比对上区域长度小于等于第四阈值。例如，第三阈值可以是0.7，即第一条件可以是比对相似性或者比对长度占最大比对相似性或最大比对长度的70％；第四阈值可以是500bp。

对应的，嵌合体序列表征被其他待组装序列覆盖的区域分成至少两端的序列。

需要说明的是，以上述具体阈值的数值为例，对本发明实施例中获取首尾重叠比对序列进行说明。设定序列两端未比对上区域长度阈值为300bp，进行判断，并统计被包含的次数。如果被包含次数大于等于2，判定为可信被包含关系，为冗余信息会被舍弃，即不会用于后续字符串构图应用。如果被包含次数小于2，说明这个被包含不可信，该比对信息需要保留，用于后续分析，记录比对结果。比对结果包括末端覆盖深度，该末端覆盖深度是指序列的5'端和3'端未比对上区域小于300bp的比对个数；比对结果还包括比对的最大相似性、最大比对长度，统计除被包含的比对关系以外的最大长度的比对的起始坐标和终止坐标，统计被其他序列覆盖的区域坐标。被包含次数小于2的序列进行重比对，得到重比对结果。进行重比对的目的是只针对有限数目的、高质量首尾比对进行重比对来精确计算，避免了无效重比对，加速比对时间。

在获得比对结果后可以进行进一步过滤，去除低质量比对、嵌合体序列及被包含序列。

依次读取比对结果，储存每条序列的比对的最大相似性、最大比对长度，如果当前比对结果的比对相似性或者比对长度低于最大值的一定比例(如0.7)，则忽略该条记录，即不使用该比对用于后续构图。

依次读取比对结果，若某条序列被其他序列覆盖的区域分成两段或以上区域，则认为是嵌合体序列，忽略该序列，不用于后续分析。

依次读取比对结果，计算每条序列两端未比对上区域的长度，去除未比对上区域后，根据已有比对结果重新计算该序列是否被包含(阈值为300bp)，并统计被包含次数。如果被包含次数大于等于2,判定为可信被包含关系，为冗余信息，舍弃，不用于后续作图。需要说明的是，因为去除了末端未比对上的序列，相当于序列有变化，所以比对结果会与之前不一样。在一种更优的实施方式中，此处并不使用比对软件进行重比对，而是在之前比对结果的基础上进行。之前的比对结果中记录有位置的相关信息，可以直接进行判断。

根据比对结果，若形成比对的两条序列都不是嵌合体或者被包含关系，选择具有首尾重叠比对(根据定义阈值可以设为500bp)。在获得了首尾重叠比对后，需要进行字符串图的构建。将首尾重叠比对转换为第一字符串图。

实施例五

在本发明实施例五中主要提供了对第一字符串图的目标特征进行处理的实施方式。由于初始生成的第一字符串图中可能会存在过多冗余及错误信息，需要对第一字符串图进行进一步处理，其中，目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种。

其中，在该实施例中还包括对满足预设条件的边和/或节点进行处理，包括：

根据覆盖度、出度，删除高深度边和/或节点；

删除传递简化对应边及反向边；

根据比对长度，删除低数值边；

根据覆盖度、出度、比对长度，保留高数值边和/或节点；

删除冗余边；

删除嵌合体边。

对应的，所述根据覆盖度、出度的高深度边和/或节点，至少满足如下一种：序列末端覆盖深度大于等于第五阈值(如平均覆盖深度的2000倍)对应的节点；出度大于等于第六阈值(如平均出度的2000倍)的节点；出度大于等于第七阈值(如平均出度的十倍)的节点且不满足第三条件的出边，所述第三条件表征按比对长度递减排序前第一数量(如10)的出边，即第三条件可以是score属性值最大的10条边。

所述传递简化对应边及反向边，包括：所述第一字符串图中的大于等于长度波动阈值(如1000)的边及其相反边。

所述根据比对长度的低数值边，包括：比对长度小于等于第八阈值(如该节点所有出边的score属性值最大值的50％)的出边及其反向边。

所述根据覆盖度、出度、比对长度的高数值边和/或节点，包括：满足第四条件的节点，所述第四条件节点对应的序列末端覆盖深度大于等于第九阈值(如平均深度的1.5倍)或者出度大于等于第十阈值(如平均出度的1.5倍)；以及不满足第四条件的节点的比对长度最大的出边，第四条件可以是节点对应的序列末端覆盖深度大于平均深度的1.5倍，或者出度大于平均出度的1.5倍。

对应的，所述目标特征包括支线路径，其中，对所述支线路径进行处理，包括：删除所含边的数目小于等于第十一阈值(如15)的支线路径。

对应的，所述目标特征包括Z型路径，其中，对所述Z型路径进行处理，包括：

当Z型路径包含嵌合体边时，根据比对相似性和比对长度对所含边评分，从低分值依次删除对应边，直至不为Z型路径；

删除不满足第五条件边所在的Z型路径，所述第五条件表征边的匹配碱基数小于第十二阈值或节点支持数小于第十三阈值；

根据匹配碱基数对Z型路径评分，从低分值依次删除路径并重新搜索路径，直至遍历全部Z型路径。

其中，第十二阈值和第十三阈值可以相同，也可以不同，如都可以设置为60，即第五条件为匹配碱基数小于60或节点支持数小于60。

对应的，所述目标特征包括泡状结构，其中，对所述泡状结构进行处理，包括：根据所含边的匹配碱基数对泡状结构评分，保留评分最高路径。

所述目标特征包括复合路径，其中，对所述复合路径进行处理，包括：当复合路径以末端节点为起始时，根据边的匹配碱基数对复合路径评分，保留评分最高路径。

在本发明实施例中，为了提升结果的准确性，在对第二字符串图转换为序列组装结果之前，还可以针对第一字符串图或第二字符串图所含的子图进行处理；所述对所述第一字符串图或第二字符串图所含的子图进行处理，包括：在所述第一字符串图或所述第二字符串图所含的子图内进行搜索，根据匹配碱基数对搜索路径评分，保留评分最高路径。

在本发明实施例中，S105、根据比对深度，对所述序列组装结果连接点进行验证，以使打断组装错误的连接点。

下面以具体的阈值数值对本发明实施例中对目标特征进行处理的过程进行说明，其主要包括以下内容：

构建字符串图：根据首尾重叠比对构建序列末端间的字符串图，节点表示序列的5'端或者3'端，边表示序列末端间的首尾重叠关系。length属性为边的长度，score属性为比对长度，identity属性为比对相似性，match属性为匹配碱基数，edge为节点支持数。

删除高深度节点和边：如果节点对应的序列末端覆盖深度大于平均覆盖深度的2000倍，或者节点出度大于平均出度的2000倍，则删除该高深度节点。如果节点出度大于平均出度的10倍，则仅保留该节点score属性值最大的10条出边，删除其他出边。

标记嵌合体和候选嵌合体边：对于某个节点，如果节点的直接前任节点与直接后继节点之间都没有边，则认为该节点为候选嵌合体节点，该节点的所有入边和出边标记为候选嵌合体边。如果从某个节点的直接后继节点往前BFS(即宽度优先搜索)搜索2层(即2层搜索宽度)，从该节点的直接前任节点往后BFS搜索2层，找不到公共节点，则认为该节点为嵌合体节点，该节点的所有入边和出边标记为嵌合体边。

计算边的支持度：设边的末端节点、终止节点分别为s、e，s的直接前任节点与e的直接前任节点的公共节点个数为n1，s的直接后继节点与e的直接前任节点的公共节点个数为n2，s的直接后继节点与e的直接后继节点的公共节点个数为n3，则n1、n2、n3之和为该边的支持度值，记录为边的edge属性。该计算是针对所有边进行。

对字符串图进行传递简化，长度波动阈值为1000，删除对应的边及其反向边。

删除低分值的边：对每个节点的出边，如果某条出边的score属性小于该节点所有出边score属性最大值的50％，则将边及其反向边删除。

保留高分值的边：按节点遍历字符串图。如果read末端覆盖深度大于平均深度的1.5倍，或者节点出度大于平均出度的1.5倍，则不作删除处理。否则仅保留score属性最大的出边，删除其他边及其反向边。

搜索并删除边的数目小于等于阈值(如该阈值可以设置为15)的支线路径。搜索并删除末端节点入度为0、且终止节点出度为0、边的数目小于等于阈值(如该阈值可以设置为16)的简单路径。

按嵌合体、候选嵌合体标记做Z型剪切。从某个出度大于1的节点的直接后继节点开始，找出所有Z型路径(例如，路径中包含的边的数目阈值为8；即小于等于8的为Z型路径)。对包含嵌合体标记的边的Z型路径，先按路径第一条边的identity属性从小到大排序，再按路径第一条边的score属性从小到大排序。按排序后的顺序依次删除，如果路径仍为Z型路径且边都存在，则删除路径中的边及其反向边。需要说明的是，每删除一条边，就判断一次是非仍为Z型路径，若是就继续按排列顺序删除下一条；直到不为Z型路径。按节点遍历字符串图，并重复以上操作。之后再按候选嵌合体标记对整个字符串图进行一次以上操作。该过程针对之前标记为嵌合体边和候选嵌合体边组合的Z型路线进行处理；若Z型路径不包括嵌合体边和候选嵌合体边，就需要用通过其他方式进行处理。

在本发明实施例中处理复合路径的一种方式是，如果存在从某个末端节点开始的复合路径(例如，搜索深度阈值为500)，遍历复合路径时找出并记录从末端节点开始到每个遍历到的节点的匹配数评分最高的路径(即最优路径)。路径的匹配数评分为路径中所有边的匹配数评分之和。边的匹配数评分为“边的match属性-边的末端节点所有出边match属性最大值×0.9”。从末端节点到某个节点的最优路径，可根据从末端节点到该节点所有直接前任节点的最优路径的评分，以及该节点的所有入边的匹配数评分计算得到，满足最优子结构条件(如，最优是指最大值)。保留复合路径中从末端节点到终止节点的最优路径中的边和节点，删除复合路径中的其他边和节点，及其反向边和反向节点。按节点遍历字符串图，并重复以上操作。

找共有前任节点，处理冗余边。从某个出度大于1的末端节点往后做BFS搜索(如，搜索深度阈值为500)。搜索时记录从末端节点到当前节点的最优路径，路径评分同前述步骤。直至搜索达到最大深度限制或者所有节点都无法继续往下搜索为止。找出所有终止节点的搜索深度最大的已遍历了的共有的前任节点。保留从末端节点到该前任节点、从该前任节点到每个终止节点的最优路径中的边和节点，删除遍历到的其他边和节点，及其反向边和反向节点。按节点遍历字符串图，并重复以上操作。

按阈值进行Z型剪切。先计算每条边的in_out_match值(in_match特指入边的match属性，out_match特指出边的match属性)。设边的末端节点、终止节点分别为s、e，边的match属性值为m1。s的所有出边最大match属性值为m2。e的所有入边最大match属性值为m3，则该边的in_out_match属性值为“m1/m2×50+m1/m3×50”。如果m2或者m3为0，则设置m1/m2为1或者m1/m3为1再计算。按同样方式计算每条边的in_out_edge值(其中，s的直接前任节点与e的直接前任节点的公共节点个数为n1，s的直接后继节点与e的直接前任节点的公共节点个数为n2，s的直接后继节点与e的直接后继节点的公共节点个数为n3，则n1、n2、n3之和为该边的支持度值)。从某个出度大于1的节点的直接后继节点开始，找出所有Z型路径(如，路径中包含的边的数目阈值为8；小于等于8的为Z型路径)。对Z型路径，先按路径第一条边的in_out_match属性从小到大排序，in_out_match属性大于等于60则按路径第一条边的in_out_edge属性从小到大排序，如果in_out_edge属性也大于等于60，则保留这个路径。对排序后且未略过的路径，如果路径仍为Z型路径且边都存在，则删除路径中的边及其反向边。按节点遍历字符串图，并重复以上操作。

在本发明实施例中Z型剪切的一种处理方式是：按前述方式重新计算每条边的in_out_match值。找出字符串图中的所有Z型路径(如，路径中包含的边的数目阈值为16)。计算路径评分，路径评分为第一条边的in_out_match值与最后一条边的in_out_match值之和，按路径评分从小到大排序、删除。找出路径列表中评分最小的路径，如果路径仍为Z型路径且边都存在，则删除路径中的边及其反向边。从删除的路径两端节点开始搜索，如果有新增Z型路径，则计算评分并加入新的Z型路径列表。重新找出两个列表中评分最小的Z型路径，重复上述操作，直到所有Z型路径都处理完为止。

删除短的环状路径：找出字符串图中的首尾相连的简单路径(环状路径)，如果路径中边的数目小于10，则删除路径中的边和节点，及其反向边和反向节点。

处理泡状结构：如果存在从末端节点开始的泡状结构图(如，路径中包含的边的数目阈值为500)，按上述描述的处理复合路径的一种实施方式中所述方法计算泡状结构的多条路径的评分，保留评分最高的路径，删除其他路径中的边和节点，及其反向边和反向节点。按节点遍历字符串图，并重复以上操作。

在本发明实施例中处理复合路径的另一实施方式可以为：按照上述处理复合路径的实施方式从末端节点开始往后做BFS搜索，找出最优路径，搜索深度为30，忽略出度为0的终止节点。如果只有一个终止节点，则只保留从末端节点到终止节点的最优路径中的边和节点，删除其他遍历到的边和节点，及其反向边和节点。按节点遍历字符串图，并重复以上操作。

还可以按照嵌合体标记、低分值标记处理复合路径：按照上述处理复合路径的实施方式从末端节点开始往后做BFS搜索，找出最优路径。搜索深度为20，忽略出度为0的终止节点，忽略带有嵌合体标记的边(即这条路径的搜索停止)。如果只有一个终止节点，则只保留从末端节点到终止节点的最优路径中的边和节点，删除其他遍历到的边和节点，及其反向边和节点。按节点遍历字符串图，并重复以上操作。按同样的方式，忽略带有低分值标记的边，处理字符串图。

删除contig末端符合条件的候选嵌合体边：从带有候选嵌合体标记的边的末端节点开始往后搜索50层。如果碰到入度不为1或者出度不为1的节点，则提前终止搜索，否则继续处理下一条带有候选嵌合体标记的边。如果候选嵌合体标记的边的edge属性为0，则删除该边及其反向边。如果边的identity属性小于0.85或者边带有低分值标记，并且边的score属性小于0.5，则删除该边及其反向边。继续处理下一条带有候选嵌合体标记的边，直至处理完所有边。

对以上处理得到的字符串图，从分叉节点处切割字符串图，生成路径，并根据路径和read序列生成最终的组装contig序列，并保留构成contig的每条reads的连接点位置信息。

在上述实施例的基础上，本发明实施例四还提供了对第一字符串图对应的子图进行处理的方法，该方法包括：

获取所述第一字符串图对应的子图；确定所述子图的目标节点，依据所述目标节点对所述子图进行搜索，获得包含起始分叉节点的子图；基于所述起始分叉节点的子图，确定目标分叉节点；利用所述目标分叉节点获取目标子图；对所述目标子图进行搜索，保留所有末端节点和终止节点对间的最优路径的边和节点。其中子图是指交集少的个图。

从某个出度大于1或者入度大于1的节点(即分叉节点)起始，往前和往后同时进行BFS搜索，搜索深度不超过8层。如果搜索中碰到新的分叉节点，则继续往前和往后同时进行BFS搜索，搜索深度不超过8层。搜索时每个节点只遍历一次。则可得到包含起始分叉节点的子图。如果子图中分叉节点数目大于等于1600或者节点数目大于等于4000，则略过该子图(即直接完全保留)，从另一个不在任一子图中的分叉节点开始重新获取子图。将所有反向节点和反向边也加入子图中，找出子图中的所有末端节点(入度为0或者至少有一个直接前任节点都不在子图中)。从末端节点开始往后在子图中做BFS搜索，除末端节点外，只有当节点的所有入边都遍历到了，才继续往后搜索。搜索时记录从末端节点到当前节点的最优路径，路径评分同前述步骤。搜索深度不做限制，但只搜索子图中的节点。搜索完成后，保留所有末端节点-终止节点对间的最优路径中的边和节点，删除遍历到的其他边和节点。遍历字符串图，从另一个不在任一子图中的分叉节点开始，获取子图并重复以上操作。

在上述实施例的基础上，本发明实施例四还提供了对组装结果进行验证的方法，即步骤S105。具体实施方式可为：将待组装比对至序列组装结果上，并计算每条组成conitg的reads的连接点处的比对深度，如果比对深度显著低于平均深度，则将该contig从该位置区域切断。

实施例六

本发明实施例六主要提供了在上述序列组装的方法上，还包括：对所述序列组装结果进行纠错(即Polish)，得到纠错后的序列组装结果，参见图4，所述纠错方法包括：

S201、将多个测序片段与所述序列组装结果进行比对，将多个所述测序片段排列至所述序列组装结果上相对应的位置，获得第一排列结果；

S202、以k-mer方式对所述序列组装结果进行划分，得到多个k-mer片段；基于所述第一排列结果，对所述测序片段进行对应的k-mer划分；

S203、从正序，以k-mer片段为单位，计算每个位置上所有碱基的得分；

S204、从反序，确定最末碱基为最末位置上最高得分碱基，根据所述最末碱基对应的所述k-mer片段，确定倒数第二位碱基；根据所述倒数第二位碱基对应的所述k-mer片段，确定倒数第三位碱基；依次类推，获得k-mer得分链；

S205、使用k-mer得分链替换所述参考序列，获得第一纠错结果。

可选的，步骤S202所述k值为3。

可选的，步骤S203所述计算每个位置上所有碱基得分的方法包括，根据公式score(p,b)＝max{score(p-1,b)+count_k-mer}-C计算每个位置上所有碱基的得分；其中，score(p,b)为所述参考序列p位置处碱基的得分，b为碱基A、T、G、C或缺失，count_k-mer为所述k-mer片段在所述比对序列中的出现次数；C为所述参考序列p位置处所述比对序列对所述参考序列的覆盖度。

可选的，纠错方法还包括：将所述第一纠错序列作为所述序列组装结果，使用如步骤S201至S205所述方法进行迭代纠错。

可选的，纠错方法还包括：

S206、将低质量区间的序列替换为该区间内出现概率最高的测序片段的序列，得到第二纠错序列；

其中，所述低质量区间为，所述序列组装结果上，出现两个以上低质量位点，且两个以上所述低质量位点之间最大间隔长度小于或等于低质量间隔预设值的区间；所述低质量位点为，在所述序列组装结果上，根据第一排列结果，将每个所述最后位点上，获得最高所述概率的k-mer片段对应的碱基占该位点上k-mer片段总数的比值小于低质量占比预设值的位点。

可选的，纠错方法还包括：

S207、当前后相邻所述低质量位点的间隔距离大于低质量间隔预设值时，将所述前后相邻所述低质量位点的间隔距离大于低质量间隔预设值的低质量位点的碱基串联成低质量长序列，使用如步骤S201至S205的方法进行纠错，其中以第二纠错序列为所述序列组装结果，获得第三纠错结果。

可选的，纠错方法还包括：

S208、对低覆盖区域使用如步骤S201至S205的方法进行纠错，其中以第一纠错序列或第二纠错序列或第三纠错序列为所述序列组装结果，获得第四纠错结果；

其中，所述低覆盖区域为，所述序列组装结果上，对应的所述测序片段的数量低于覆盖阈值的区域。

可选的，纠错方法还包括，使用如步骤S206或S207或S208所述方法进行迭代纠错。

在一种具体的实现方式中，步骤S201所述多个测序序列可以通过各种测序手段获得，例如二代测序技术、三代测序技术，其来源可以是序列组装结果的原始测序数据，所述序列组装结果即为待纠错序列。

在一种具体的实现方式中，步骤S206的具体实施方式如下：

将第一纠错序列作为序列组装结果，在序列组装结果上，根据第一排列结果，将每个最后位点上，获得最高概率的k-mer片段对应的碱基占该位点上k-mer片段总数的比值小于低质量占比预设值(如80％)的位点，标记为低质量位点；

将序列组装结果上出现两个以上低质量位点，且两个以上低质量位点之间最大间隔长度小于或等于低质量间隔预设值(如50碱基)的区间划分为低质量区间；

将低质量区间的序列替换为该区间内出现概率最高的测序片段的序列，具体包括：将多个测序片段排列至第一纠错序列上相对应的位置，获得第二排列结果；基于第二排列结果，确定低质量区间内重复次数最多的测序片段，将序列组装结果上低质量区间内的序列纠错为重复次数最多的测序片段对应的序列，得到第二纠错序列。

在另一种具体实施方式中，步骤S207具体实施方法包括：

将与前后相邻低质量位点的间隔距离大于低质量间隔预设值(如50碱基)的低质量位点的碱基串联成低质量长序列；

将多个测序片段排列至第二纠错序列上相对应的位置，挑选出与低质量长序列对应的排列，获得第三排列结果；

以k-mer方式(k-mer优选为3)对低质量长序列进行划分，得到多个k-mer片段，基于第三排列结果，对测序片段进行对应的k-mer划分；

确认测序片段每个k-mer片段的最后位点，基于第三排列结果，从正序，以k-mer片段为单位，计算每个位置上所有碱基的得分；计算每个位置上所有碱基得分的方法为，根据公式score(p,b)＝max{score(p-1,b)+count_k-mer}-C计算每个位置上所有碱基的得分；其中，score(p,b)为所述参考序列p位置处碱基的得分，b为碱基A、T、G、C或缺失，count_k-mer为所述k-mer片段在所述比对序列中的出现次数；C为所述参考序列p位置处所述比对序列对所述参考序列的覆盖度；

从反序，确定最末碱基为最末位置上最高得分碱基，根据所述最末碱基对应的所述k-mer片段，确定倒数第二位碱基；根据所述倒数第二位碱基对应的所述k-mer片段，确定倒数第三位碱基；依次类推，获得k-mer得分链；

使用k-mer得分链替换所述参考序列，获得第三纠错序列。

为了获得更准确的纠错序列，在本实施例步骤S207中，采用优选方案，使用长读长的测序片段进行纠错，长读长的测序片段由三代测序获得。

在另一种具体实施方式中，步骤S208包括：

将多个测序片段排列至纠错序列上相对应的位置，挑选出低覆盖区域对应的排列，获得第四排列结果；其中，地覆盖区域为，在序列组装结果上，对应的测序片段的数量低于覆盖阈值(如3)的区域，覆盖阈值为位点上对应排列的测序片段的数量；

以k-mer方式(k-mer优选为3)对低覆盖区域进行划分，得到多个k-mer片段，基于第四排列结果，对测序片段进行对应的k-mer划分；

基于第四排列结果，从正序，以k-mer片段为单位，计算每个位置上所有碱基的得分；计算每个位置上所有碱基得分的方法为，根据公式score(p,b)＝max{score(p-1,b)+count_k-mer}-C计算每个位置上所有碱基的得分；其中，score(p,b)为所述参考序列p位置处碱基的得分，b为碱基A、T、G、C或缺失，count_k-mer为所述k-mer片段在所述比对序列中的出现次数；C为所述参考序列p位置处所述比对序列对所述参考序列的覆盖度；

从反序，确定最末碱基为最末位置上最高得分碱基，根据所述最末碱基对应的所述k-mer片段，确定倒数第二位碱基；根据所述倒数第二位碱基对应的所述k-mer片段，确定倒数第三位碱基；依次类推，获得k-mer得分链；

使用k-mer得分链替换所述参考序列，获得第四纠错序列。

为了获得更准确的纠错序列，在本实施例S208步骤中，采用优选方案，使用高准确度的测序片段进行纠错，高准确度的测序片段由二代测序获得。

综上，上述实施例中所提供的序列纠错方法，可以对组装结果中SNP、InDel等错误的纠错，快速获得高准确度的组装结果。

实施例七

本发明所述序列组装方法还包括：基于二进制对待组装序列、对待组装序列进行比对获得的比对结果以及首尾重叠比对序列进行存储。

为了便于后续处理，本发明所涉及的序列信息、比对结果均可以采用二进制形式进行存储。在实际应用中，序列的数据量非常庞大，特别是对于大型及超大型的基因组，需占用较大的硬盘储存空间。将具体的序列转换为二进制形式，可以有效的进行压缩，节约存储空间。

本实施例所提供的基于二进制的存储方式，适用于待组装序列、首尾重叠比对序列等序列信息，同时也适用于对待组装序列进行比对获得的比对结果。

对于序列信息的二进制存储具体方法包括：(1)按照预设分组将序列划分为多个碱基组合；其中，预设分组策略包括将相邻的四个碱基确定为一个碱基组合；(2)根据碱基与二进制数码之间预设的转换关系，将序列的碱基转换为二进制数码；(3)对序列中少于四个碱基的组合，采用指定的二进制数码对少于四个碱基的组合进行补位扩充，得到满足分组的二进制序列；(4)按照划分的碱基组合对二进制序列进行存储。

具体的，测序序列包含腺嘌呤(A)、胞嘧啶(C)、鸟嘌呤(G)和胸腺嘧啶(T)四种形式，每一个碱基都需要1个字节(8个比特位)来存储数据。在本实施例中可以首先将上述碱基采用二进制形式进行转换，诸如：A＝00，C＝01，G＝10，T＝11，即将每种碱基转换为占用两个比特位的二进制表现形式。然后将相邻的四个碱基作为一组构成一个字节存储；因此采用二进制形式压缩后的测序序列约为原始数据大小的四分之一，有效的减少了所需占用的内存空间。

对于序列中序列长度是四的倍数，每四个碱基组合均可表示为完整的一个字节形式。对于序列长度不是四的倍数，最后一组的碱基组合不能占用八个比特位，可以采用00对最后一组的碱基组合进行补位扩充，使其满足八个比特位即一个字节的形式。

对于比对结果的二进制存储与序列不同。通常比对结果是一连串数字，其固定了每一列所记录的指标以及具体数值的含义，在基于二进制数码序列进行转换时可以根据实际情况中每列的具体含义、顺序、具体数值含义进行定义。还可以对比对结果进行压缩，例如，通过行或列之间做减法，实现节省存储空间。

比对结果可以通过记录比对区间的特征进行存储，为了尽可能的缩小存储空间。本实施例所提供的序列比对方法还包括：纪录比对序列编号、比对方向、比对序列比对区间起始、比对序列比对区间终止、参考序列编号、参考序列比对区间起始、参考序列比对区间终止存储比对结果；优选的，纪录比对序列编号差值、参考序列编号差值、参考序列比对区间长度与比对序列比对区间长度差值存储比对结果。此外，按照4字节进行存储；并且，每字节使用7比特位纪录数值，不足7比特位的使用特定二进制数码进行填充；剩余1比特位使用特定二进制数码标识比对结果是否终止。

由于当前储存数值是按照4个字节(32个比特位)来储存的，而输出结果中的值一般都很小，只需要2个字节或者更小，高位的比特位都为零，为无效位。因此为了进一步减小存储空间，可以对第一比对序列中各列的数值做如下压缩处理：按每七个比特位作为一组进行分割，丢弃前端分组中值为0的组合，对保留的分组除最后一组最高位置为0其他组最高位均置为1，最后一组不足七比特位的使用0填充存储数值压缩后的字节组合。通过采用上述方式压缩处理后的第一比对结果的大小通常为原始大小的1/3，明显占用存储空间更小。

实施例八

本发明实施例八为提供实际应用场景的示例。需要说明的是，因基因组组装所需数据量大且复杂，能更好的体现本发明的技术效果，本实施例仅作为示例并不用作限制本发明的应用范围。

(1)对6种不同的参考基因组已报道的植物样品，使用三代测序技术获取长度长序列数据，对相同的数据分别使用不同的组装方式进行组装，具体如下：

对比例1：现有技术中广泛使用的、针对长度长序列的组装软件WTDBG(https://github.com/ruanjue/wtdbg)；

效果例1：测序获得的长度长序列使用如实施例2及实施例3所述的方法进行矫正，对矫正数据使用WTDBG进行组装；与对比例1相比，组装效果明显提升；

效果例2：测序获得的长度长序列使用如实施例2及实施例3所述的方法进行矫正，使用本发明所述方法进行组装；组装效果最佳。

通过测序获得的若干条读段(reads)，根据重叠关系对读段进行拼接，所获得的无间隔(gap)的序列为重叠群(contig)。将组装出的重叠群(contig)从大到小排列，当累计长度超过组装序列总长度的50％时，最后一个重叠群或骨架序列的长度即为对应的N50的大小。各对比例及效果例的详细组装结果如下表，其中，N50数值越大、重叠群数目越小说明组装的连续性、完整性越好。

表1

(2)水稻基因组的组装：水稻是世界上最终要的农作物之一，基因组大小约为400Mb-430 Mb左右，使用三代测序技术获取长度长序列数据，共获得6,100,295条reads，共92,930Mb的碱基，约230倍的覆盖，平均长度为17,905bp。对相同的数据分别使用不同的组装方式进行组装，具体如下：

对比例2：组装所使用软件CANU(https://github.com/marbl/canu)；

对比例3：组装所使用软件WTDBG(https://github.com/ruanjue/wtdbg)；

效果例3：使用如实施例2及实施例3所述的方法进行矫正，使用本发明所述方法进行组装；经过过滤处理的、用于构建第一字符串图的首尾重叠比对序列为71,940条、共7,459Mb；特别的，使用本发明所述方法对字符串图的进行处理(如步骤S103、实施例五所示)仅需1s即可完成；

对比例4：组装所使用软件SHASTA(https://github.com/chanzuckerberg/shasta)；

对比例5：组装所使用软件NECAT(https://github.com/xiaochuanle/NECAT)。

各对比例及效果例的详细组装结果如表2：

表2

从上述结果得出本发明实施例提供的组装结果长度最长，此外该结果与参考序列R498做共性比较，发现组装的序列与参考序列匹配性更好，几乎没有组装错误。

本发明所述方案极大的提升了组装序列连续性：(1)对构图所用首尾重叠比对处理时，通过重比对使得序列间比对的末端边界更准确，同时也使匹配数计算更准确，后续依赖于匹配数的构图算法能得到更好的结果；(2)对构图所用首尾重叠比对进行了过滤，使得错误的首尾重叠关系更少，构建的字符串图更简化；(3)对首尾重叠比对进行嵌合体序列过滤，以及构图时标记嵌合体边并用于Z型剪切，极大减小了嵌合体序列对构图的影响，使得字符串图更简化；(4)标记的候选嵌合体边很多是重复导致的边，将含有这些边的Z型路径删除，使得字符串图更简化；(5)对边和/或节点、支线路径、Z型路径、泡状结构和复合路径的分步处理，也在尽量保证准确性的情况下使得字符串图更简化。这些简化字符串图的方法，使得最终的字符串图中分叉较少，按分叉节点切割字符串图生成的路径和序列才更长，组装序列连续性得到了很大提升。

本发明所述方案对处理高杂合基因组也有较好的技术效果。主要在复合路径处理时，删除了冗余结构，使得组装基因组大小尽量接近预估基因组大小。此外，通过找共有前任节点并处理冗余边、获取并处理子图、处理泡状结构的步骤，使杂合导致的泡状结构和两条路径交联在一起的结构仅保留了一条路径，删除了杂合导致的冗余路径，减小了杂合对基因组组装大小和组装连续性的影响。

实施例九

参见图5，在本发明实施例九中提供了一种序列组装系统，该系统包括：

比对单元10，用于对待组装序列进行比对，获得首尾重叠比对序列；

构建单元20，用于根据所述首尾重叠比对序列构建第一字符串图；

处理单元30，用于对所述第一字符串图的目标特征进行处理，获得第二字符串图，其中，所述目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种；

转换单元40，用于将所述第二字符串图转换为序列组装结果。

可选地，比对单元中所述待组装序列为经矫正的序列。

可选地，所述比对单元还可具体用于：

对所述待组装序列进行切分，获得切分结果；

基于所述切分结果进行比对，获得首尾重叠比对序列。

可选地，所述系统还包括：

存储单元，用于基于二进制对待组装序列、对待组装序列进行比对获得的比对结果以及首尾重叠比对序列进行存储。

可选地，构建单元中所述首尾重叠比对序列为经过滤处理后的首尾重叠比对序列，所述过滤处理同时包括对嵌合体序列过滤、对被包含序列过滤、对冗余首尾重叠比对序列过滤。

可选地，所述被包含序列包括末端未比对上区域长度小于等于第一阈值且被包含次数大于等于第二阈值的序列；

优选的，所述方法还包括：对所述被包含序列进行重判定，所述对所述被包含序列进行重判定，包括：

去除序列末端未比对上区域，进行重比对，若重比对结果中末端未比对上区域长度小于等于第一阈值且被包含次数大于等于第二阈值的序列确定为被包含序列。

其中，所述对冗余首尾重叠比对序列过滤，包括：

对不满足第一比对条件和第二比对条件的首尾重叠比对序列进行过滤；其中，所述第一比对条件为比对相似性或者比对长度占最大比对相似性或最大比对长度的比例大于等于第三阈值，所述第二比对条件为首尾重叠区域内末端未比对上区域长度小于等于第四阈值。

对应的，所述嵌合体序列表征被其他待组装序列覆盖的区域分成至少两段的序列。

可选地，所述处理单元中所述目标特征包括满足预设条件的边和/或节点，其中，所述处理单元具体用于：

根据覆盖度、出度，删除高深度边和/或节点；

删除传递简化对应边及反向边；

根据比对长度，删除低数值边；

根据覆盖度、出度、比对长度，保留高数值边和/或节点；

删除冗余边；

删除嵌合体边。

可选地，所述根据覆盖度、出度的高深度边和/或节点，至少满足如下一种：序列末端覆盖深度大于等于第五阈值对应的节点；出度大于等于第六阈值的节点；出度大于等于第七阈值的节点且不满足第三条件的出边，所述第三条件表征按比对长度递减排序前第一数量的出边；

所述传递简化对应边及反向边，包括：所述第一字符串图中的大于等于长度波动阈值的边及其相反边；

所述根据比对长度的低数值边，包括：比对长度小于等于第八阈值的出边及其反向边；

所述根据覆盖度、出度、比对长度的高数值边和/或节点，包括：满足第四条件的节点，所述第四条件节点对应的序列末端覆盖深度大于等于第九阈值或者出度大于等于第十阈值；以及不满足第四条件的节点的比对长度最大的出边。

可选地，所述处理单元中的所述目标特征包括支线路径，其中，所述处理单元具体用于：

删除所含边的数目小于等于第十一阈值的支线路径。

可选地，所述处理单元中所述目标特征包括Z型路径，其中，所述处理单元具体用于：

当Z型路径包含嵌合体边时，根据比对相似性和比对长度对所含边评分，从低分值依次删除对应边，直至不为Z型路径；

删除不满足第五条件边所在的Z型路径，所述第五条件表征边的匹配碱基数小于第十二阈值或节点支持数小于第十三阈值；

根据匹配碱基数对Z型路径评分，从低分值依次删除路径并重新搜索路径，直至遍历全部Z型路径。

可选地，所述处理单元中的所述目标特征包括泡状结构，其中，所述处理单元具体用于：

根据所含边的匹配碱基数对泡状结构评分，保留评分最高路径。

可选地，所述处理单元中的所述目标特征包括复合路径，其中，所述处理单元具体用于：

当复合路径以末端节点为起始时，根据边的匹配碱基数对复合路径评分，保留评分最高路径。

可选地，所述系统还包括：

字符串图处理单元，用于对所述第一字符串图或第二字符串图所含的子图进行处理；

所述字符串图处理单元具体用于：

在所述第一字符串图或所述第二字符串图所含的子图内进行搜索，根据匹配碱基数对搜索路径评分，保留评分最高路径。

可选地，所述系统还包括：

验证单元，用于根据比对深度，对所述序列组装结果连接点进行验证，以使打断组装错误的连接点。

可选地，所述系统还包括：

纠错单元，用于对所述序列组装结果进行纠错，得到纠错后的序列组装结果。

本发明实施例九提供了一种序列组装系统，比对单元对待组序列进行比对，获得首尾重叠比对序列；构建单元根据首尾重叠比对序列构建第一字符串，处理单元对第一字符串的目标特征进行处理，获得第二字符串图，转换单元将第二字符串图转换为序列组装结果。目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种。本发明待组装序列可为经矫正的序列，将其作为字符串构图的数据基础，并且对字符串图进行处理，使得处理后的字符串图更加满足序列组装的需求，提升组装结果的准确性、连续性以及快速组装的目的。

本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的装置而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。

对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其它实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。

Claims (15)

1.一种序列组装方法，其特征在于，包括：

S101、对待组装序列进行比对，获得首尾重叠比对序列；

S102、根据所述首尾重叠比对序列构建第一字符串图；

S103、对所述第一字符串图的目标特征进行处理，获得第二字符串图，其中，所述目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种；

S104、将所述第二字符串图转换为序列组装结果。

2.根据权利要求1所述的方法，其特征在于，步骤S101中所述待组装序列为经矫正的序列。

3.根据权利要求1所述的方法，其特征在于，步骤S102中的所述首尾重叠比对序列为经过滤处理后的首尾重叠比对序列，所述过滤处理同时包括对嵌合体序列过滤、对被包含序列过滤、对冗余首尾重叠比对序列过滤。

4.根据权利要求3所述的方法，其特征在于，所述被包含序列包括末端未比对上区域长度小于等于第一阈值且被包含次数大于等于第二阈值的序列；

优选的，所述方法还包括：对所述被包含序列进行重判定，所述对所述被包含序列进行重判定，包括：

去除序列末端未比对上区域，进行重比对，若重比对结果中末端未比对上区域长度小于等于第一阈值且被包含次数大于等于第二阈值的序列确定为被包含序列。

5.根据权利要求3所述的方法，其特征在于，所述对冗余首尾重叠比对序列过滤，包括：

对不满足第一比对条件和第二比对条件的首尾重叠比对序列进行过滤；其中，所述第一比对条件为比对相似性或者比对长度占最大比对相似性或最大比对长度的比例大于等于第三阈值，所述第二比对条件为首尾重叠区域内末端未比对上区域长度小于等于第四阈值。

6.根据权利要求1所述的方法，步骤S103中的所述目标特征包括满足预设条件的边和/或节点，其中，对所述满足预设条件的边和/或节点进行处理，包括：

根据覆盖度、出度，删除高深度边和/或节点；

删除传递简化对应边及反向边；

根据比对长度，删除低数值边；

根据覆盖度、出度、比对长度，保留高数值边和/或节点；

删除冗余边；

删除嵌合体边。

7.根据权利要求6所述的方法，其特征在于，

所述根据覆盖度、出度的高深度边和/或节点，至少满足如下一种：序列末端覆盖深度大于等于第五阈值对应的节点；出度大于等于第六阈值的节点；出度大于等于第七阈值的节点且不满足第三条件的出边，所述第三条件表征按比对长度递减排序前第一数量的出边；

所述传递简化对应边及反向边，包括：所述第一字符串图中的大于等于长度波动阈值的边及其相反边；

所述根据比对长度的低数值边，包括：比对长度小于等于第八阈值的出边及其反向边；

所述根据覆盖度、出度、比对长度的高数值边和/或节点，包括：满足第四条件的节点，所述第四条件节点对应的序列末端覆盖深度大于等于第九阈值或者出度大于等于第十阈值；以及不满足第四条件的节点的比对长度最大的出边。

8.根据权利要求1所述的方法，其特征在于，步骤S103中的所述目标特征包括支线路径，其中，对所述支线路径进行处理，包括：

删除所含边的数目小于等于第十一阈值的支线路径。

9.根据权利要求1所述的方法，其特征在于，步骤S103中的所述目标特征包括Z型路径，其中，对所述Z型路径进行处理，包括：

当Z型路径包含嵌合体边时，根据比对相似性和比对长度对所含边评分，从低分值依次删除对应边，直至不为Z型路径；

删除不满足第五条件边所在的Z型路径，所述第五条件表征边的匹配碱基数小于第十二阈值或节点支持数小于第十三阈值；

根据匹配碱基数对Z型路径评分，从低分值依次删除路径并重新搜索路径，直至遍历全部Z型路径。

10.根据权利要求1所述的方法，其特征在于，步骤S103中的所述目标特征包括泡状结构，其中，对所述泡状结构进行处理，包括：

根据所含边的匹配碱基数对泡状结构评分，保留评分最高路径。

11.根据权利要求1所述的方法，其特征在于，步骤S103中的所述目标特征包括复合路径，其中，对所述复合路径进行处理，包括：

当复合路径以末端节点为起始时，根据边的匹配碱基数对复合路径评分，保留评分最高路径。

12.根据权利要求1所述的方法，其特征在于，在S104步骤前，所述方法还包括：

对所述第一字符串图或第二字符串图所含的子图进行处理；

所述对所述第一字符串图或第二字符串图所含的子图进行处理，包括：

在所述第一字符串图或所述第二字符串图所含的子图内进行搜索，根据匹配碱基数对搜索路径评分，保留评分最高路径。

13.根据权利要求1或12所述的方法，其特征在于，所述方法还包括：

S105、根据比对深度，对所述序列组装结果连接点进行验证，以使打断组装错误的连接点。

14.根据权利要求1或12所述的方法，其特征在于，所述方法还包括：

对所述序列组装结果进行纠错，得到纠错后的序列组装结果。

15.一种序列组装系统，其特征在于，包括：

比对单元，用于对待组装序列进行比对，获得首尾重叠比对序列；

构建单元，用于根据所述首尾重叠比对序列构建第一字符串图；

处理单元，用于对所述第一字符串图的目标特征进行处理，获得第二字符串图，其中，所述目标特征包括满足预设条件的边和/或节点、支线路径、Z型路径、泡状结构和复合路径中的一种或多种；

转换单元，用于将所述第二字符串图转换为序列组装结果。

Images