一种基因组测序数据序列组装方法.pdf

本发明实施例提供了一种基因组测序数据序列组装方法，可以简便而准确地实现基因组测序数据序列的组装。该方法包括：根据序列片段之间的重叠关系，构建重叠关系图G以及反向互补图G’；任取重叠关系图G中的未被检查过的节点nx，以任意方向D遍历G和G’；判断任意方向D中是否存在一节点ny与节点nx有连接关系；如果存在，则继续判断在方向D中的节点ny与节点nx是否为双向唯一关系；当所述节点ny与节点nx为双向唯一关系时，将所述节点ny与nx的关系判定为可信连接关系；否则，将所述节点ny与nx的所有关系判定为不确认连接关系。

1.  一种基因组测序数据序列组装方法，其特征在于，包括：
根据序列片段之间的重叠关系，构建重叠关系图G以及反向互补图G’；
任取重叠关系图G中的未被检查过的节点n_x，以任意方向D遍历G和G’；
判断任意方向D中是否存在一节点n_y与节点n_x有连接关系；如果存在，则继续判断在方向D中的节点n_y与节点n_x是否为双向唯一关系；
当所述节点n_y与节点n_x为双向唯一关系时，将所述节点n_y与n_x的关系判定为可信连接关系；否则，将所述节点n_y与n_x的所有关系判定为不确认连接关系。

2.  如权利要求1所述的方法，其特征在于，进一步包括：
将已经检查过的关系从G和G’中删除；并转向任取重叠关系图G中的未被检查过的节点的步骤，直至重叠关系图G中每一个节点均被检查过。

3.  如权利要求1所述的方法，其特征在于，所述判断节点n_y与节点n_x是否为双向唯一关系包括：
当且仅当在G中，找到一节点n_y为该节点n_x下游的唯一路径；并且在G’中，节点n_x是节点n_y的下游唯一路径时，认定所述节点n_y与节点n_x为双向唯一关系。

4.  如权利要求1所述的方法，其特征在于，任意方向D为从所述节点n_x出发的出度方向或指向所述节点n_x的入度方向。

5.  如权利要求1所述的方法，其特征在于，进一步包括：
将所述可信连接关系、节点n_x和节点n_y放入可靠序列连接关系集合；和、或将所述不确认连接关系放入到可疑序列片段关系集合中。

6.  如权利要求1至5任一所述的方法，其特征在于，进一步包括：
将所有不确认连接关系与参考序列进行比对，判断所述参考序列中是否存在与所述不确认连接关系相同的连接关系；如果存在，则将所述出现在参 考序列中的不确认连接关系重新判定为可信连接关系。

7.  如权利要求6所述的方法，其特征在于，判断所述参考序列中是否存在与所述不确认连接关系相同的连接关系包括：
根据已有的参考序列，建构参考序列的重叠关系图R和反向互补图R’；
判断任意一出发节点M的多个不确认连接关系中，是否有一个在R或者R’中存在。

8.  如权利要求6所述的方法，其特征在于，若在所述参考序列中存在涉及同一节点的多种摆放可能时，将比对得分最高的摆放方式用于构建R和R’。

9.  如权利要求6所述的方法，其特征在于，所述根据已有的参考序列，建构参考序列的重叠关系图R和反向互补图R’包括：
使用参考序列与测序得到的读序进行比对，将读序按照比对的位置进行排序，得到一个参考序列重叠关系图R和反向互补图R’。

10.  如权利要求9所述的方法，其特征在于，采取容忍大片段重排和缺失的比对程序将参考序列与测序得到的读序进行比对。

11.  如权利要求9所述的方法，其特征在于，在将参考序列与测序得到的读序进行比对之前，进一步包括：
对测序的所有读序进行聚类分析，并将聚类化简后的结果与参考序列进行比对。

12.  如权利要求1至6任一所述的方法，其特征在于，进一步包括：
根据可信连接关系对序列进行组装，得到组装结果序列，供用户参考。

13.  如权利要求12所述的方法，其特征在于，所述得到的组装结果序列为唯一的结果序列或一系列的多个结果序列。

一种基因组测序数据序列组装方法
技术领域
本发明涉及基因技术，尤其涉及一种基因组测序数据序列组装方法。
技术背景
随着测序技术的不断进步，已经有大量的微生物基因组被完成并提交到数据库。具有工业用途的微生物，其工业用菌株大多是对已有菌株进行不断筛选和改造得到的。
为了获取一株工业用菌株的基因组全图，目前常用的分析方案有De novo(从头测序)。De novo是指在没有任何背景信息的基础上使用测序和常规分子生物学实验手段对目的物种基因组进行测序，组装，框架搭建和空白(Gap)填补的技术流程。这一方案在面对较为复杂或者较大的基因组时，会耗费大量时间和成本，但其获得的结果是最为可靠的，可以获得一个物种的质粒、质体、特异性序列和突变结果，进而能够获取一个物种的全部功能序列并分析得到其生理生化能力，重构其生活史。
目前De novo最常用的是“overlap-layout-consensus”方案，通过序列比对的方法对所有测序所得的读序(又称为read，是指测序出的一条序列)的边缘序列进行侦测，寻找可能存在的重叠区域。而后，根据这些读序的重叠关系对这些读序进行合并，构成重叠群从而完成组装。
图1所示为现有技术De novo中“overlap-layout-consensus”算法的示意图。如图1所示，假设在基因组上存在两个高度相似序列REP1和REP2。Read1和Read2分别位于REP1两侧，并且其重叠区域坐落于REP1中，两者重叠区域长度为L1；Read3和Read4位于REP2两侧，并且其重叠区域位于REP2内，这两者重叠长度为L2，并且L2>L1。
如果组装程序使用贪婪算法，假设在遍历序列过程中首先遍历到Read1，由于其具有最好重叠关系的序列为Read4，因此，Read1->Read4这个连接结果被带入到最终结果中，导致组装结果发生错误。因此，只有Read3或者Read4的遍历顺序大于Read1、Read2时，才能得到正确的结果。
如果组装程序使用图论理论，虽然能够识别出贪婪算法可能导致的错误，但也是有局限性的。假设在读序重叠侦测(overlap)步骤中设定一个参数L，如果两个序列的重叠长度大于L才认为两个序列具有重叠关系，假设L1<L<L2；那么在读序重叠侦测过程中，由于参数选择的原因，Read1->Read2这个关系就不会被识别到，导致在读序摆放(layout)过程中，如果Read1或者Read2首先被遍历到，由于其只有一种连接方式(Read1->Read4，Read2->Read3)被识别到，该部分会被认为是可信的，进而被带入到最终结果中导致组装错误。
发明内容
有鉴于此，本发明实施例提供一种新的基因组测序数据序列组装方法，可以简便而准确地实现基因组测序数据序列的准确还原。
为了达到上述目的，本发明实施例提供的一种基因组测序数据序列组装方法，包括：
根据序列片段之间的重叠关系，构建重叠关系图G以及反向互补图G’；
任取重叠关系图G中的未被检查过的节点n_x，以任意方向D遍历G和G’；
判断任意方向D中是否存在一节点n_y与节点n_x有连接关系；如果存在，则继续判断在方向D中的节点n_y与节点n_x是否为双向唯一关系；
当所述节点n_y与节点n_x为双向唯一关系时，将所述节点n_y与n_x的关系判定为可信连接关系；否则，将所述节点n_y与n_x的所有关系判定为不确认连接关系。
其中，进一步包括：
将已经检查过的关系从G和G’中删除；并转向任取重叠关系图G中的未被检查过的节点的步骤，直至重叠关系图G中每一个节点均被检查过。
其中，所述判断节点n_y与节点n_x是否为双向唯一关系包括：
当且仅当在G中，找到一节点n_y为该节点n_x下游的唯一路径；并且在G’中，节点n_x是节点n_y的下游唯一路径时，认定所述节点n_y与节点n_x为双向唯一关系。
其中，任意方向D为从所述节点n_x出发的出度方向或指向所述节点n_x的入度方向。
其中，进一步包括：
将所述可信连接关系、节点n_x和节点n_y放入可靠序列片段重叠群中；和、或
将所述不确认连接关系放入到可疑序列片段关系集合中。
其中，其特征在于，进一步包括：
将所有不确认连接关系与参考序列进行比对，判断所述参考序列中是否存在与所述不确认连接关系相同的连接关系；如果存在，则将所述出现在参考序列中的不确认连接关系重新判定为可信连接关系。
其中，判断所述参考序列中是否存在与所述不确认连接关系相同的连接关系包括：
根据已有的参考序列，建构参考序列的重叠关系图R和反向互补图R’；
判断任意一出发节点M的多个不确认连接关系中，是否有一个在R或者R’中存在。
其中，若在所述参考序列中存在涉及同一节点的多种摆放可能时，将比对得分最高的摆放方式用于构建R和R’。
其中，所述根据已有的参考序列，建构参考序列的重叠关系图R和反向互补图R’包括：
使用参考序列与测序得到的读序进行比对，将读序按照比对的位置进行排序，得到一个参考序列重叠关系图R和反向互补图R’。
其中，采取容忍大片段重排和缺失的比对程序将参考序列与测序得到的读序进行比对。
其中，在将参考序列与测序得到的读序进行比对之前，进一步包括：
对测序的所有读序进行聚类分析，并将聚类化简后的结果与参考序列进行比对。
其中，进一步包括：
根据可信连接关系对序列进行组装，得到组装结果序列，供用户参考。
其中，所述得到的组装结果序列为唯一的结果序列或一系列的多个结果序列。
其中，如果结果是一系列多个结果序列时，所有无法选定的不确定连接关系也会作为结果的一部分，反馈给用户。
在用户提供重叠关系图的情况下，利用本发明实施例提供的技术方案，可对重叠关系图进行检查，并自动生成初步组装结果，该结果可以修正现有技术中贪婪算法或图论算法的错误。
附图说明
图1是现有技术中De novo中“overlap-layout-consensus”算法的示意图。
图2是本发明实施例中基因组测序数据序列组装方法的流程示意图。
图3是本发明另一实施例中基因组测序数据序列组装方法的流程示意图。
图4是本发明一实施例提供的参考序列重叠关系图构建示意图。
图5是本发明实施例中参考序列定量化利用方法的示例示意图。
图6是本发明实施例中基因组测序数据序列组装方法的示例示意图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚，下面结合附图对本发明作进一步的详细描述。
图2是本发明实施例中基因组测序数据序列组装方法的流程示意图。如图2所示，该方法包括：
步骤201：根据测序得到的读序之间的重叠关系，构建重叠关系图以及反向互补图。重叠关系图及其反向互补图中的所有的对应节点彼此之间都是反向互补等效的关系。由于我们通常只知道两个序列是否具有重叠关系，但不确定组装结果中序列组的最终摆放顺序，因此，我们需要同时构建两个图，重叠关系图G以及其反互补序列图G’。只要两个序列片段之间有重叠关系，即可在重叠关系图中给予标记。
步骤202：判断重叠关系图G中是否每一个节点均被检查过；如果是，则结束整个流程；否则转向步骤203。
步骤203：任取重叠关系图G中的未被检查过的节点n_x，以方向D遍历G和G’。其中方向D可以是出度方向(从该节点出发)或入度方向(指向该节点)。
步骤204：判断任意方向D中是否存在某一节点n_y与节点n_x有连接关系；如果存在，则进入步骤205；否则转入步骤206。
步骤205：如果节点n_y与节点n_x为双向唯一关系，则进入步骤208；否则进入步骤209。
这里，当且仅当在G中，找到一节点n_y为该节点n_x在D方向下游的唯一路径；并且在G’中，节点n_x是节点n_y在D方向下游唯一路径时，认定节 点n_y与节点n_x为双向唯一关系。
步骤206：判断在方向D’中是否存在某一节点n_z与节点n_x有连接关系；如果存在，则进入步骤207。这里，由于n_x必然会与一个节点相连，所以如果执行到步骤206，则此处的判断结果一定为存在节点n_z有连接关系。
步骤207：如果节点n_z与节点n_x为双向唯一关系，则进入步骤208；否则进入步骤209。
步骤208：将该n_y与n_x的关系，和/或该n_z与n_x的关系确认为可信连接关系，并将该可信连接关系n_x->n_y，和/或n_z->n_y放入可靠序列片段重叠群中；并转入步骤210。
步骤209：将该n_y与n_x的关系，和/或该n_z与n_x的关系确认为不确认连接关系，将该不确认连接关系n_x->n_y，和/或n_z->n_y放入到可疑序列片段关系集合中；并转入步骤210。
步骤210：将已经检查过的关系，包括关系n_x->n_y、关系n_y->n_x，和/或关系n_x->n_z、关系n_z->n_x分别从G和G’中删除，同时将该节点n_x设置为“已检查”。
这样，通过双向检察方法对重叠关系图进行修剪，获得可靠序列片段重叠群以及可疑序列片段关系集合。
在本发明一实施例中，以上步骤204和步骤205与步骤206和步骤207可以平行执行，这样进一步提高运算效率。
在本发明一实施例中，还提供了一种方法，进一步将可疑序列片段关系集合中的不确认连接关系确认为可信连接关系。即将所有可疑序列片段关系集合中的不确认连接与参考序列进行比对，如果某一个不确认连接关系可以在参考序列中找到，则认为该不确认连接关系为可信连接关系。假设涉及同一节点的多个不确认关系在参考序列中都可以找到(即一个读序可以被比对到参考序列的多个位置上)，则使用对比得分最高的结果作为其唯一的比对位置。这里由于常常会采用一定的衡量算法来衡量两个读序的相近或相远程序，这样一个衡量算法常常可以数字化地评估对比得分。
此时，如图3所示，本发明实施例提供的参考序列定量化利用方法包括：
步骤301：根据已有的参考序列以及初始测序得到的读序，建构参考序列的重叠关系图R和反向互补图R’。这里参考序列可以是已经为现有技术公开的序列，这些序列的组装关系已经确定，在具有参考序列情况下，使用参考序列与测序得到的读序进行比对，将读序按照比对的位置进行排序，进 而得到一个参考序列重叠关系图R和反向互补图R’。
在本发明一实施例中，在将初始测序得到的读序与参考序列进行比对时，可以使用blat等可以容忍大片段重排和缺失的比对程序进行比对。另外由于测序数据太大，并且其长度并不均一，因此，在进行比对之前要对测序的所有读序进行聚类分析，使用95％的相似度进行序列聚类，只选取每个类中最长的一个或者几个序列用于后续分析，并将聚类化简后的结果与参考序列进行比对，而后对比对后的结果的方向和其比对位置进行排序。
图4为本发明一实施例提供的参考序列重叠关系图构建示意图。如图4a所示，假设读序1、2、3、4为初始测序得到的，由于读序2和3完全被读序1所包含，在建构参考序列重叠关系图R前先通过聚类分析，将读序1、2和3聚成一类，把读序1选取出来，读序2和3舍弃。将聚类化简后的结果(读序1和读序4)与参考序列比对，如图4b所示，而后将读序1和读序4在参考序列上的比对位置进行排序。发现读序1和读序4是相邻的，因此将1->4的关系作为参考序列重叠关系图的一部分。
步骤302：在可疑序列片段关系集合中，如果任意一出发节点M的多个不确定关系(M->X或M->Y)中的某一个连接关系(M->X)在R或者R’中存在，则将该连接关系存入可靠序列片段重叠群中，其余涉及M的连接关系在可疑序列片段关系集合中被清除。
步骤303：根据可靠序列片段重叠群中的可信连接关系进行组装，得到组装结果序列。在本发明一实施例中，得到的组装结果序列有可能是唯一的结果序列，也有可能是一系列的多个结果序列。
如果结果是一系列多个结果序列时，所有无法选定的不确定连接关系也会作为结果的一部分，反馈给用户。
图5是本发明实施例中利用参考序列组装的示例示意图。如图5所示，根据一个已有的参考序列(图5第一列)构成的一个参考重叠图(图5的第二列)，将该参考重叠图与可疑序列片段关系组合中的关系(图5的第三列)进行对比，将可疑序列片段关系组合中在参考重叠图中出现的关系判定为可信连接关系，将其他关系删除(图5中的第四列中的关系B、C、E为可信连接关系，关系D和F删除)，最后得到无误的结果(图5中的第五列)。
在用户提供重叠关系图或者提供原始数据的情况下，利用本发明实施例提供的技术方案，可对重叠关系图进行构建和检查，并自动生成初步组装结果，该结果可以修正现有技术中贪婪算法或图形算法的错误。
下面以一个具体事例为例来说明以上的基因组测序数据序列组装方法。如图6a为一段序列组的重叠关系图G以及其反互补序列图G’的示例。图中每一个结点代表一个读序，记为n，{n₁,n₂,n₃,…n_n}∈G。图中每一条边记为e，{e₁,e₂,e₃,…e_n}∈G。以节点5为例，每一条入度方向(指向节点5)代表节点5代表的读序能够和节点4所代表的读序重叠，每一条出度方向(从节点5出发)代表节点5代表的读序能够和节点7代表的读序进行连接。
以任何一个节点如节点4为起始节点，以任意方向如节点4的出度方向根据图G和G’进行检查，在G中节点4的出度方向上有两个关系：与节点5的关系和与节点6的关系。因此可以确定该节点延长关系并不唯一，将其关系记录到可疑序列片段关系集合中。由于通过图G已经可以确定节点4的出度方向上并不具有唯一的可信连接关系，因此，不需要再检查G’。该过程如图6b所示。
再对节点4的入度方向进行检查，在图G中发现节点3是节点4入度方向的唯一节点，而后，在G’中检查节点3，发现节点4是节点3上游的唯一节点。因此确定其在G和G’中节点4和节点3均为彼此的唯一延长节点。因此，认为该3->4关系是可信的，放入到可靠序列片段重叠群中。该过程如图6c所示。
将已经检查过的n_x节点相关的关系从G和G’中删除。如图6d所示，从图G和G’中删除以节点4为中心的边。
在未被检查过的节点中再任选一个节点，如节点1，在G中检查节点1的出度方向，发现节点3为其最优解，在图G’中检查节点3，发现节点3在出度方向上并不唯一，因此将1->3，2->3放入到可疑序列片段关系集合中。该过程如图6e所示。
将已经检查过的关系从G和G’中删除。如图6f所示，从图G和G’中删除以节点1、节点2和节点3为中心的边。
在剩余的节点中再任取一个节点，如节点5，发现在G和G’中，其与节点7为唯一关系，因此，将5->7放入到可靠序列片段重叠群中。该过程如图6g所示。
将已经检查过的关系从G和G’中删除。如图6h所示，从图G和G’中删除以节点5、节点7为中心的边。
重复之前的过程，发现节点6和节点8关系为双向唯一，因此，将6->8放入到可靠序列片段重叠群中。该过程如图6i所示。
将已经检查过的关系从G和G’中删除。如图6j所示，从图G和G’中删除以节点6、节点8为中心的边。
对所有的节点进行检查之后，发现图中已经没有任何的边，整个双边检查过程结束。
假设存在一个参考序列r，并且我们发现在该参考序列r中，有一段序列为2->3以及3->4。根据对这一段序列的分析，我们可以选择将原可疑关系2->3判定为可信连接关系，同时删除关系1->3，如图6k所示。
本发明实施例提供的一种基因组测序数据序列组装方法，利用双向最优算法构建序列重叠图，很好地解决了De novo的错拼问题。结合参考序列、使用重测序相关方法对De novo重叠图进行简化，使得组装结果和实际序列最相近。该方法对参考序列进行了充分的挖掘，并对De novo方法的局限性和复杂性进行了弥补和化简。根据优化后的重叠图，本发明实施例的方法还可以自动地对已有的序列片段重叠群进行错误检查和错误拆分，防止错拼发生。
另外，此方法还有助于细胞内质粒序列的组装。在微生物中，大部分质粒都是环状序列存在，在重叠关系图中的表现形式是n₁->n₂->n₃……n₁，即在重叠关系图中存在一个欧拉回路。因此，在构建结果序列过程中，该方法能够自动筛选出较长的在参考序列上比对不到的序列作为质粒寻找的种子，使用弗勒里算法，寻找能够通过这些种子节点的欧拉回路作为可能存在的质粒候选；而后通过组装结果的长度(不超过1MB)，对质粒序列进行筛选。
因此，对无法进行优化的重叠图，本发明实施例的方法还能够进行质粒拆分和质粒猜测，最大可能地将可能存在的质粒序列进行筛选、分离、并进行自动化成环，方便进一步的分析。
对大肠杆菌K₁₂MG1655菌株的基因组、大肠杆菌某苏氨酸产生菌株、Meiothermus Ruber DSM1299、Pedobacter heparinus DSM2366、黄色短杆Zl5进行测序，并利用本专利发明的基因组测序数据组装方法进行测试，测试结果将从以下几个方面进行评价：
组装得到的基因组的完整性，得到的组装结果是否能还原所有的基因组区域，是否能够完成一个微生物基因组；
组装得到的基因组的整体可靠性，组装结果与参考序列相比是否具有结构性差异，是否有区域被遗漏或者错误带入；
组装需要使用测序数据的获取成本，以及完成同一个样本基因组测序数 据组装所需要的数据量下限。
根据测试结果显示，利用本发明实施例提供的方法，得到的初步组装结果与真实情况相比，几乎不存在结构性错误，其结果明显优于已经公布的其他工具。在具有较为近缘的参考序列的情况下，本发明实施例提供的方法只需要其他方法所需约2/3的数据量，就可以得到相同或者更优的分析结果。
这里，本领域技术人员可以理解，上文中出现的可疑、不可信、不确认可以为同一意思；同样，可信、可靠、确认也可以为同一意思；同时，序列关系组合、序列集合等根据上下文理解也可以等同，本发明对这些用词并不给予严格区分。
以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换等，均应包含在本发明的保护范围之内。