CN107403075B - 比对方法、装置及系统 - Google Patents

比对方法、装置及系统 Download PDF

Info

Publication number
CN107403075B
CN107403075B CN201710652379.8A CN201710652379A CN107403075B CN 107403075 B CN107403075 B CN 107403075B CN 201710652379 A CN201710652379 A CN 201710652379A CN 107403075 B CN107403075 B CN 107403075B
Authority
CN
China
Prior art keywords
sequence
extended
seed
positioning result
window
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201710652379.8A
Other languages
English (en)
Other versions
CN107403075A (zh
Inventor
徐伟彬
金欢
颜钦
姜泽飞
周志良
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Genemind Biosciences Co Ltd
Original Assignee
Genemind Biosciences Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Genemind Biosciences Co Ltd filed Critical Genemind Biosciences Co Ltd
Priority to CN201710652379.8A priority Critical patent/CN107403075B/zh
Publication of CN107403075A publication Critical patent/CN107403075A/zh
Application granted granted Critical
Publication of CN107403075B publication Critical patent/CN107403075B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B25/00—ICT specially adapted for hybridisation; ICT specially adapted for gene or protein expression

Landscapes

  • Physics & Mathematics (AREA)
  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Medical Informatics (AREA)
  • Biophysics (AREA)
  • Theoretical Computer Science (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • General Health & Medical Sciences (AREA)
  • Evolutionary Biology (AREA)
  • Biotechnology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Data Mining & Analysis (AREA)
  • Software Systems (AREA)
  • Artificial Intelligence (AREA)
  • Evolutionary Computation (AREA)
  • Epidemiology (AREA)
  • Databases & Information Systems (AREA)
  • Public Health (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Bioethics (AREA)
  • Genetics & Genomics (AREA)
  • Molecular Biology (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Information Retrieval, Db Structures And Fs Structures Therefor (AREA)

Abstract

本发明公开了一种比对方法、装置及系统，比对方法包括：将每条读段转化成与该读段对应的一组短片段，获得多组短片段；确定短片段在参考库的对应位置，以获得第一定位结果，所称的参考库为基于参考序列构建的哈希表，参考库包含多个条目，参考库的一个条目对应一条种子序列，种子序列能够与参考序列上的至少一段序列匹配，参考库的相邻两个条目对应的两条种子序列在参考序列上的距离小于短片段的长度；去除第一定位结果中定位到参考库相邻条目中的任一条目上的短片段，获得第二定位结果；基于第二定位结果中来自相同读段的短片段的公共部分进行延伸，以获得读段的比对结果。该比对方法能够对测序数据进行高效准确的处理及定位。

Description

比对方法、装置及系统
技术领域
本发明涉及数据处理领域，具体地，涉及一种序列比对方法、一种序列比对装置及一种序列比对系统。
背景技术
在生物数据处理和分析中，比对作为生物信息分析的基础环节，对时间和效率提出了一定的要求，不同的比对模型和算法具有不同的敏感性和准确性。
目前主流的比对软件如bwa、bowtie、soap主要是针对二代测序平台产生的序列而开发的，而如blast和mapq等软件，主要是针对长序列而开发的。
对特定类型数据的敏感性、比对时间和/或比对效率，现有的比对方法有待提高。
发明内容
本发明实施方式旨在至少解决现有技术中存在的技术问题之一或者至少提供一种商业手段。
依据本发明的一方面提供的一种比对方法，包括如下步骤：将每条读段转化成与该读段对应的一组短片段，获得多组短片段；确定短片段在参考库的对应位置，以获得第一定位结果，所称的参考库为基于参考序列构建的哈希表，参考库包含多个条目，参考库的一个条目对应一条种子序列，所称的种子序列能够与参考序列上的至少一段序列匹配，参考库的相邻两个条目对应的两条种子序列在参考序列上的距离小于短片段的长度；去除第一定位结果中定位到参考库相邻条目中的任一条目上的短片段，获得第二定位结果；基于所述第二定位结果中来自相同读段的短片段进行延伸，以获得读段的比对结果。
依据本发明的另一方面提供的一种比对装置，用以实施上述比对方法的全部或部分步骤，该装置包括：转化模块，用于将每条读段转化成与该读段对应的一组短片段，获得多组短片段；查找模块，用于确定短片段在参考库的对应位置，以获得第一定位结果，所称的参考库为基于参考序列构建的哈希表，参考库包含多个条目，参考库的一个条目对应一条种子序列，所称的种子序列能够与参考序列上的至少一段序列匹配，参考库的相邻两个条目对应的两条种子序列在参考序列上的距离小于短片段的长度；剔除模块，用于去除第一定位结果中定位到参考库相邻条目中的任一条目上的短片段，获得第二定位结果；生长模块，用于基于所述第二定位结果中来自相同读段的短片段进行延伸，以获得读段的比对结果。
依据本发明的又一方面提供的一种计算机可读介质，用于存储/承载计算机可执行程序，本领域普通技术人员可以理解，在执行该程序时，通过指令相关硬件可完成上述比对方法的全部或部分步骤。所称介质可以包括：只读存储器、随机存储器、磁盘或光盘等。
依据本发明的再一方面提供的一种比对系统，包括：输入装置，用于输入数据；输出装置，用于输出数据；处理器，用于执行计算机可执行程序，执行所述计算机可执行程序包括完成上述比对方法；存储装置，用于存储数据，其中包括所述计算机可执行程序。
本发明的比对方法、装置和/或系统通过将读段转化成短片段以及将序列信息转化成位置信息，具有较高的灵敏性和准确性，能够高效准确地处理各种测序平台的下机数据，包括长读长和短读长的下机数据。特别是对于包含未能识别的碱基的读段，即包含gap的读段的高效精确的定位，尤其适用。
本发明的附加方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。
附图说明
本发明的上述和/或附加的方面和优点从结合下面附图对实施例的描述中将变得明显和容易理解，其中：
图1是本发明实施方式中的比对方法的流程示意图。
图2是本发明实施方式中的参考库的相邻两个条目的距离的示意图。
图3是本发明实施方式中的连通长度示意图。
图4是本发明实施方式中的比对装置的结构示意图。
图5是本发明实施方式中的比对装置的结构示意图。
图6是本发明实施方式中的比对装置的结构示意图。
图7是本发明实施方式中的比对装置的结构示意图。
图8是本发明实施方式中的比对装置的结构示意图。
具体实施方式
下面详细描述本发明的实施例，所述实施例的示例在附图中示出，其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附图描述的实施例是示例性的，仅用于解释本发明，而不能理解为对本发明的限制。
在本发明的描述中，所称的“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重要性和/或具有先后顺序。
所称的“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以具体情况理解上述术语在本发明中的具体含义。
所称的“读段”指测定DNA/RNA/蛋白序列所获得的序列片段，包括利用测序平台对DNA/RNA/蛋白序列的至少一部分进行测定识别所获得的序列。测序平台可选择但不限于Illumina公司的Hisq/Miseq/Nextseq测序平台、Thermo Fisher(Life Technologies)公司的Ion Torrent平台、BGI的BGISEQ平台和单分子测序平台，测序方式可以选择单端测序，也可以选择双末端测序，获得的下机数据是测读出来的片段，称为读段(reads)。
所称的“比对”指序列比对，包括将读段定位到参考序列上的过程，也包括获得读段定位/匹配结果的过程。
依据本发明的实施方式提供的一种比对方法，请参考图1，该方法包括如下步骤：S110将每条读段转化成与该读段对应的一组短片段，获得多组短片段；S120确定短片段在参考库的对应位置，以获得第一定位结果，参考库为基于参考序列构建的哈希表，参考库包含多个条目，参考库的一个条目对应一条种子序列，种子序列能够与参考序列上的至少一段序列匹配，参考库的相邻两个条目对应的两条种子序列在参考序列上的距离小于短片段的长度；S130去除第一定位结果中定位到参考库相邻条目中的任一条目上的短片段，获得第二定位结果；S140基于所述第二定位结果中来自相同读段的短片段进行延伸，以获得读段的比对结果。该比对方法通过将读段转化成短片段以及将读段序列信息转化成位置信息，即将序列形态转成数字形态，利于快速准确地实现各种测序平台的下机数据的比对定位。特别是对于包含未能识别的碱基的读段，即包含gap或者N的读段的快速准确比对，比如由于测序质量不佳、碱基识别不佳等所获得的读段的比对分析，尤其适用。
所称的参考序列(reference，ref)为预先确定的序列，可以是自己预先测定组装的DNA和/或RNA序列，也可以是他人测定公开的DNA和/或RNA序列，可以是预先获得的目标个体所属生物类别中的任意的参考模板，例如，同一生物类别的已公开的基因组组装序列的全部或者至少一部分，若目标个体是人类，其基因组参考序列(也称为参考基因组或者参考染色体组)可选择NCBI数据库提供的HG19；进一步地，也可以预先配置包含更多参考序列的资源库，在进行序列比对前，先依据目标个体的性别、人种、地域等因素选择或测定组装出更接近的序列来作为参考序列，有助于后续获得更准确的序列分析结果。参考序列包含染色体编号以及各个位点在染色体上的位置信息。所称的参考库实质为哈希表(hashtable)，可以直接以所称的种子序列为键(键名)、以所称的种子序列在参考序列上的位置(position)为值(键值)构建该参考库；也可以先将所称的种子序列转成数字或者整数字符串，以该数字或者整数字符串为键、以种子序列在参考序列上的位置为值建立该参考库。所称的以种子序列在参考序列上的位置为值，可以是该种子序列在参考序列/染色体上对应的一个或多个位置，位置可直接以真实数值或者数值范围表示，也可以重新编码以自定义的字符和/或数字表示。根据本发明的一个实施例，利用C++的向量vector实现哈希表的构建，可表示为：Hash(seed)＝Vector(position)，所称的向量vector是一种对象实体，能够容纳许多其他类型相同的元素，因此也被称为容器。可用二进制保存，以此建成该参考库。另外，也可以将哈希表分成块(block)存储，在block头设置块头键和块尾键，例如，对于顺序序列块{5,6,7,8...,19,20}，设置块头和块尾(或者说表头和表尾)5和20，若有个数为3，因3<5，可知3不属于该顺序序列块，若有个数为10，因5<10<20，可知10属于这个序列块。如此在查询的时候可以选择全局索引，也可以通过比较块头键和块尾键快速定位到所在的block，可不需要全局索引。
所称的参考库可以在要进行序列比对时构建，也可以预先构建保存。根据本发明的具体实施方式，预先构建参考库保存备用，参考库的构建包括：依据参考序列的碱基总数totalBase，确定种子序列的长度L，L＝μ*log(totalBase),
Figure BDA0001368381680000031
且L小于待比对分析的读段的长度(读长)；基于所述种子序列的长度，生成所有可能的种子序列，获得种子序列集；确定种子序列集中能够匹配到参考序列的种子序列以及该种子序列的匹配位置，以获得该参考库。该方法基于发明人多次假设试验验证而建立的种子序列长度与参考序列的关系，能够使构建得的参考库包含全面的种子序列以各种子序列在参考序列上对应的位置的关联信息，该参考库结构紧凑，内存占用小且能够用于序列定位分析中的高速访问查询。根据该实施方式获得的参考库的一个条目只包含一个键，一个键对应至少一个值。
本发明的该实施方式，对生成所有可能的种子序列、获得种子序列集的方法不作限制，对于输入一个集合，可遍历该集合中的元素来获得特定长度的、所有可能的元素组合，例如可以利用递归算法和/或循环算法来实现。
在一个例子中，参考序列为人基因组，包含大约30亿个碱基，待处理的读段的长度为不小于25bp，L取[11，15]中的整数，利于高效比对。
在一个例子中，参考序列为人cDNA参考基因组，统计该参考基因组的碱基总数totalBase，基于碱基总数设定种子序列(seed)的长度L，L(seed)＝log(totalBase)*μ,
Figure BDA0001368381680000032
基于L以及DNA序列的碱基种类包括A、T、C和G四种，利用递归算法，生成所有可能的种子序列的集合，获得种子序列集，该过程可表示为seed＝B1B2...BL,B∈{ATCG}；确定种子序列集中能够匹配到该参考基因组的种子序列以及该种子序列的匹配位置，以能够匹配到该参考基因组的种子序列为键、以该种子序列在参考基因组上的位置位为值来构建获得该参考库。
在一个例子中，参考序列为某物种的DNA基因组和转录组，统计该参考序列的碱基总数totalBase，基于碱基总数设定种子序列(seed)的长度L，L(seed)＝log(totalBase)*μ,
Figure BDA0001368381680000033
基于L、组成DNA序列的碱基种类包括A、T、C和G四种以及组成RNA序列的碱基种类包括A、U、C和G四种，利用递归算法，生成所有可能的种子序列的集合，获得种子序列集，该过程可表示为seed＝B1B2...BL,B∈{ATCG}∪{AUCG}；确定种子序列集中能够匹配到该参考序列的种子序列以及该种子序列的匹配位置，以能够匹配到该参考序列的种子序列为键、以该种子序列在参考序列上的位置位为值来构建获得该参考库。
在一个例子中，可将种子序列转化成由数字字符组成的字符串，以该字符串为键来建库，能够提高访问查询所建参考库的速度。例如，在获得能够匹配到参考序列的种子序列后，对种子序列进行如下编码：
Figure BDA0001368381680000041
又例如，在获得种子序列集后，对种子序列集中的种子序列进行编码，碱基编码规则可与上面相同，且对参考序列也可以进行同样规则的编码转换，利于快速获得种子序列在参考序列上对应的位置信息，也利于提高所建参考库的访问查询速度。
根据本发明的具体实施方式，确定种子序列集中能够匹配到参考序列的种子序列以及该种子序列的匹配位置，包括：利用大小为L的窗口对参考序列进行滑窗，将种子序列集中的种子序列与滑窗得的窗口序列进行匹配，以确定种子序列集中能够匹配到参考序列的种子序列以及该种子的匹配位置，进行匹配的容错率为ε1。如此，能够快速获得种子序列在参考序列上的对应位置信息，利于快速构建获得参考库。所称的容错率为允许的错配碱基所占的比例，错配选自置换、插入和缺失中的至少一种。
在一个例子中，所称的匹配为严格匹配，即容错率ε1为零，当种子序列与一条或多条滑窗序列完全一致时，滑窗序列的位置为该种子序列在参考序列上对应的位置。在另一个例子中，所称的匹配为容错匹配，容错率ε1大于零，当种子序列与一条或多条滑窗序列的相同位置的碱基不一致的比例小于容错率时，滑窗序列的位置为该种子序列在参考序列上对应的位置。在一个例子中，对种子序列在参考序列上对应的位置进行编码，以编码后的字符例如数字字符为值进行参考库的构建。
换个角度，容错率ε1为不为零，相当于将一条种子序列转变成ε1允许下的一组种子模板序列(seed template)，例如seed＝ATCG，ε1为0.25即四个碱基中允许一个错误，则seed template可以为ATCG、TTCG、CTCG、GTCG、AACG、ACCG、AGCG等等。在ε1为0.25下确定seed＝ATCG在参考序列上的位置时，相当于确定该seed对应的所有seed template在参考序列的位置，例如ref＝ATCG，前面所示的所有seed template都可匹配到该位置，ref＝TTCG，seed template为ATCG、TTCG、CTCG或者GTCG可以匹配到该位置。进而，构建得的参考库可以以一条seed为键，也可以以该条seed对应的所有seed template中的每一条为键，键与键不同，一个键至少对应一个值。
根据本发明的具体实施方式，在确定种子序列在参考序列上的对应位置时，对参考序列进行滑窗的步长依据L和ε1来确定。在一个例子中，进行滑窗的步长不小于L*ε1。在一个具体例子中，参考序列为人基因组，包含大约30亿个碱基，待处理的读段的长度为不小于25bp，L为14bp，ε1取0.2-0.3，进行滑窗的步长取3bp-5bp，使滑窗定位过程中相邻两个窗口能够跨过ε1条件下的连续错误组合，利于快速定位。在一个例子中，构建得的参考库的相邻两个条目之间的距离为滑窗的步长。
根据本发明的具体实施方式，S110包括：利用大小为L的窗口对读段进行滑窗，以获得与该读段对应的一组短片段，该滑窗的步长为1bp。如此，对于一条长度为K的reads，获得(K-L+1)条长度为L的短片段，将reads转成短片段，利用高速访问查询参考库，确定各短片段在参考库的对应位置，进而获得短片段对应的reads在参考库的信息。
根据本发明的具体实施方式，S120包括：将短片段与参考库的条目对应的种子序列进行匹配，以确定短片段在参考库的位置，进行匹配的容错率为ε2。
在一个例子中，所称的匹配为严格匹配，即容错率ε2为零，当一条短片段与参考库的一个条目对应的seed或者seed template完全一致时，获得该短序列在参考库上的位置信息。在另一个例子中，所称的匹配为容错匹配，容错率ε2大于零，当短序列与参考库的一个或多个条目对应的seed或者seed template不匹配的碱基的比例小于容错率ε2时，获得该短序列在参考库上的位置信息。在一个具体例子中，ε2＝ε1且不为零，使能够获得尽可能多的有效数据。
根据本发明的具体实施方式，参考图2，S120中，所称的参考库的相邻两个条目对应的两条种子序列X1和X2在参考序列ref上的距离，可分为以下两种情形：当参考库的两个条目的键和值均为唯一，即一个条目对应一[键，值]，参考图2a，相当于该X1和X2与参考序列均为唯一匹配时(X1和X2都只匹配到参考序列一个位置)，所称的距离为X1和X2在参考序列上对应的这两个位置之间的距离，加粗黑线显示这两个位置；当参考库的两个条目中至少一个条目的键对应多个值，参考图2b，相当于该两条种子序列X1和X2中的至少一条与参考序列为非唯一匹配即X1和X2中至少有一条匹配到参考序列的多个位置，所称的距离为该X1和X2在参考序列上对应的相距最近的两个位置之间的距离，加粗黑线显示这两个位置。该实施方式对两条序列之间的距离的表示方法不作限制，例如，可以表示为一条序列的两个末端中的任一末端到另一条序列的任一末端的距离，也可以表示为一条序列的中心到另一条序列的中心的距离。
根据本发明的具体实施方式，在获得第二定位结果之后，S130还包括：去除连通长度小于预定阈值的短片段，以去除后的结果替代第二定位结果，连通长度为第二比对结果中的来自相同读段且定位到参考库不同条目的短片段映射到参考序列的总长度。该处理有利于去除一些过渡冗余的和/或相对低质量的数据，利于提高比对速度。
连通长度可表示为来自相同读段且定位到参考库不同条目的短片段的长度总和减去映射到参考序列上短片段之间的重叠部分的长度。在一个例子中，来自一条读段且定位到参考库不同条目的短片段有4条，表示为Y1、Y2、Y3和Y4，各自的长度分别为L1、L2、L3和L4，定位位置示意如图3，其中的X1和X2映射到参考序列的位置有重叠，重叠部分的长度为J，连通长度为(L1+L2+L3+L4-J)。在一个例子中，不同短片段的长度均为L，所称的预定阈值为L，如此，可在允许丢失部分有效但质量相对低的数据的情况下，提高比对速度。
根据本发明的具体实施方式，在获得第二定位结果之后，S130还包括：依据第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行评判，去除评判结果不符合预定要求的读段。去除读段同时也是去除了该读段对应的短片段。如此，在满足一定的敏感性和准确性的前提下，基于第二定位结果，直接进行精确匹配/局部快速比对，能够加速比对。
该实施方式对评判的方法不作限定，例如可以利用量化打分的方式。在一个例子中，对来自相同读段的短片短的定位结果进行打分，打分规则是：与参考序列匹配的位点作减分，与参考序列不匹配的位点作加分；在获得第二定位结果之后，依据第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行计分，去除得分不大于第一预设值的读段。根据一个具体示例，读长为25bp，对来自相同读段的短片段进行序列构建，以获得重构序列，例如，可以根据具有更多短序列支持来确定某位点的碱基类型，若某位点没有支持的短片段即没有短片段比对到该位点，则该位点碱基类型不确定可以N来表示，以此来获得重构序列，可看出重构序列与读段是对应的，重构序列的长度为读长；重构序列与参考序列(ref)匹配的位点作减一分，与参考序列不匹配的位点作加一分，比对容错率即一条读段/重构序列允许的错配比例为0.12，比对容许错误的长度为3bp(25*0.12)，初始分数Scoreinit为读长，第一预设值为22(25-3)，如此，去除得分小于22即匹配不上参考序列的位点占比超过比对容错率的重构序列，利于在允许丢失部分有效但质量相对低的数据的情况下，加速比对。根据一个具体例子，使用位运算和动态规划算法[G.Myers.A fast bit-vector algorithm for approximate string matching based on dynamicprogamming.Journal of the ACM,46(3):395-415,1999]，对于每条重构序列，读入每个位点i的位置，利用64位的二进制掩模进行快速匹配计分，每个位点一分，初始分数Scoreinit为读长，可表示为Scoreinit＝length(read)，匹配计分获得分数Score，可表示为：
Figure BDA0001368381680000051
在一个例子中，对来自相同读段的短片短的定位结果进行打分，打分规则是：与参考序列匹配的位点作加分，与参考序列不匹配的位点作减分；在获得第二定位结果之后，依据第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行计分，去除得分不小于第二预设值的读段对应的短片段。根据一个具体示例，读长为25bp，对来自相同读段的短片段进行序列构建，以获得重构序列，例如，可以根据具有更多短序列支持来确定某位点的碱基类型，若某位点没有支持的短片段即没有短片段比对到该位点，则该位点碱基类型不确定可以N来表示，以此来获得重构序列，可看出重构序列与读段是对应的，重构序列的长度为读长；重构序列与参考序列(ref)匹配的位点作加一分，与参考序列不匹配的位点作减一分，比对容错率即一条读段/重构序列允许的错配比例为0.12，比对容许错误的长度为3bp(25*0.12)，初始分数Scoreinit为-25，第二预设值为-22(-25-3)，如此，去除得分大于-22的重构序列，在允许丢失部分有效但相对低质量的数据的情况下，加速比对。
根据本发明的具体实施方式，S140中的基于第二定位结果中来自相同读段的短片段进行延伸，包括：基于来自相同读段的短片段进行序列构建，获得重构序列；基于重构序列与该重构序列对应的参考序列的公共部分进行延伸，以获得延伸序列。如此，将短片段及短片段定位信息转化成短片段对应的读段(在此称为重构序列)的定位信息，利于后续比对处理快速准确的进行。
所称的公共部分，为多条序列共有的部分。根据本发明的具体实施方式，所称的公共部分为公共子串和/或公共子序列。公共子串指多条序列中共有的连续部分，公共子序列则不须连续。例如，对于ABCBDAB和BDCABA，公共子序列是BCBA，公共子串是AB。
所称的基于来自相同读段的短片段进行序列构建，获得重构序列，在一个例子中，可根据具有更多短片段支持来确定重构序列上某位点的碱基类型，若某位点没有支持的短片段即没有短片段比对到参考序列该位点，则该位点碱基类型不确定可以以N来表示，以此来获得所称的重构序列。可看出，重构序列与读段是对应的，重构序列的长度为读长。
所称的重构序列对应的参考序列，为与重构序列匹配的一段参考序列，该段参考序列的长度不小于读长。在一个例子中，重构序列对应的参考序列的长度与重构序列相同，均为读长。在另外一个例子中，允许重构序列与对应的参考序列容错匹配，重构序列对应的参考序列的长度为重构序列的长度加上两倍的容错匹配长度，例如，重构序列长度即读长为25bp，重构序列与参考序列的匹配允许错配12％，可以以重构序列对比上的那段参考序列以及该段参考序列两端各3bp(25*12％)序列来作为重构序列对应的参考序列。
根据本发明的一个具体例子，所称的公共部分为公共子串。S140中的基于第二定位结果中来自相同读段的短片段进行延伸，包括：查找所述重构序列与所述重构序列对应的参考序列的公共子串，确定所述重构序列和所述重构序列对应的参考序列的最长公共子串；基于编辑距离，延伸所述最长公共子串以获得延伸序列。如此，能够更加准确获得包含更长匹配序列的比对结果。
根据本发明的一个具体例子，所称的公共部分为公共子序列。S140中的基于第二定位结果中来自相同读段的短片段进行延伸，包括：查找所述重构序列与所述重构序列对应的参考序列的公共子序列，确定所述重构序列和所述重构序列对应的参考序列的最长公共子序列；基于编辑距离，延伸所述最长公共子序列以获得延伸序列。
所称的编辑距离，也叫Levenshtein距离，是指两个字符串之间，由一个转成另一个所需的最少编辑操作次数。编辑操作包括将一个字符替换成另一个字符、插入一个字符以及删除一个字符。一般来说，编辑距离越小，两个串的相似度越大。
在一个例子中，对于一条重构序列/读段，查找该重构序列与该重构序列对应的参考序列的最长公共子串，可表示为求两个字符串x1x2...xi和y1y2...yj的公共子串，字符串的长度分别为m和n，计算这两字符串的公共子串的长度c[i,j]，可以得到转移方程：
Figure BDA0001368381680000061
解方程可得这两条序列的最长公共子串的长度为max(c[i,j]),i∈{1,...,m},j∈{1,...,n}；接着利用编辑距离，将最长公共子串转化成对应的参考序列，可使最长公共子串两端不断生长，找出两个字符串之间需要的最小字符操作(替换，删除，插入)。可以使用动态规划算法确定编辑距离，该问题具备最优子结构，编辑距离d[i,j]的计算可表示为下列公式：
Figure BDA0001368381680000071
其中，洞/空缺(gap)表示插入或者删除一个字符，公式中的gap表示插入或者删除一个字符(对应序列中的位点)所需的罚分，匹配(match)表示两个字符一样，公式中的match表示两个字符一样时的得分，错配(mismatch)表示两个字符不相等/不同，公式中的mismatch表示表示两个字符不相等/不同时的阀分。d[i,j]取三者中最小的一项。在一个具体例子中，一gap罚3分，连续gap增加阀1分，一个位点错配罚2分，位点匹配得0分。如此，利于含gap序列的高效比对。
根据本发明的具体实施方式，所称的公共部分为公共子序列。根据本发明的具体实施方式，S140包括：查找第二定位结果中定位到参考库的相同条目的短片段的公共子序列，确定每条读段对应的最长公共子序列；基于编辑距离，延伸最长公共子序列以获得延伸序列。
在一个例子中，对于一条重构序列/读段，查找重构序列与该重构序列对应的参考序列的最长公共子序列，基于最长公共子序列，将最长公共子序列对应的那段重构序列转化为最长公共子序列对应的那段参考序列，可利用Smith Waterman算法找出这两段序列的编辑距离，对两个字符串x1x2...xi和y1y2...yj，可以通过以下公式求得：
Figure BDA0001368381680000072
其中，
σ表示记分函数，σ(i,j)表示字符(位点)xi和yj错配或者匹配的得分，σ(-,j)表示xi空缺(删除)或者yj插入的得分，σ(i,-)表示yj删除或者xi插入的得分；接着，可利用前面示例中的计算编辑距离的方法，将最长公共子序列对应的那段重构序列转化成重构序列对应的参考序列，可在最长公共子序列对应的那段重构序列的两端不断生长，找出最小字符操作(替换，删除，插入)。在一个具体例子中，一gap罚3分，连续gap增加罚1分，一个位点错配罚2分，位点匹配得4分。如此，能够实现含gap的序列的高效比对且能够保留既含gap而其它位点准确度高的序列。
根据本发明的具体实施方式，S140还包括：从延伸序列的至少一端对延伸序列进行截断，计算截断后的延伸序列的错误定位位点的比例，满足以下条件停止截断：截断后的延伸序列的错误定位位点的比例小于第三预设值。如此，采用截断和剔除的方式，可以较好的保留匹配良好的局部序列，利于提高数据的有效率。
具体地，根据本发明的具体实施方式，基于以下对延伸序列进行截断：i、计算第一错误率和第二错误率，若第一错误率小于第二错误率，则从延伸序列的第一端对延伸序列进行截断，若第一错误率大于第二错误率，则从延伸序列的第二端对延伸序列进行截断，以获得截断后的延伸序列，所称的第一错误率为从延伸序列的第一端对延伸序列进行截断获得的截断后的延伸序列的错误定位位点的比例，所称的第二错误率为从延伸序列的第二端对延伸序列进行截断、获得的截断后的延伸序列的错误定位位点的比例；ii、以截断后的延伸序列替代延伸序列进行i，直至截断后的延伸序列的错误定位位点的比例小于第四预设值。如此，采用双端截断和剔除的方式，可以较好的保留匹配良好的局部序列，利于提高数据的有效率。根据一个具体例子，延伸序列的长度为25bp，第四预设值为第三预设置为0.12。
根据本发明的具体实施方式，S140还包括：从延伸序列的至少一端对延伸序列进行滑窗，计算滑窗得的窗口序列的错误定位位点的比例，根据窗口序列的错误定位位点的比例对延伸序列进行截断，满足以下条件停止截断：滑窗得的窗口序列的错误定位位点的比例大于第五预设值。如此，采用截断和剔除的方式，可以较好的保留匹配良好的局部序列，利于提高数据的有效率。
具体地，根据本发明的具体实施方式，基于以下对延伸序列进行截断：i、计算第三错误率和第四错误率，若第三错误率小于第四错误率，则从延伸序列的第二端对延伸序列进行截断，若第三错误率大于第四错误率，则从延伸序列的第一端对延伸序列进行截断，以获得截断后的延伸序列，所称的第三错误率为从延伸序列的第一端对延伸序列进行滑窗、获得的窗口序列的错误定位位点的比例，所称的第四错误率为从延伸序列的第二端对延伸序列进行滑窗、获得的窗口序列的错误定位位点的比例；ii、以截断后的延伸序列替代延伸序列进行i，直至窗口序列的错误定位位点的比例大于第六预设值。如此，采用双端截断和剔除的方式，可以较好的保留匹配良好的局部序列，利于提高数据的有效率。
根据本发明的具体实施方式，滑窗的窗口不大于延伸序列的长度。根据一个具体例子，延伸序列的长度为25bp，滑窗的窗口大小为10bp，第六预设值为第五预设值为0.12。
根据本发明的具体实施方式，截断的大小为1bp，即一次截断为去掉1个碱基。如此，能够高效的获得包含更多长序列的比对结果。
在一个具体例子中，利用Bowtie(http://bowtie-bio.sourceforge.net/index.shtml)、BWA(http://bio-bwa.sourceforge.net/)以及上述比对方法对同一批模拟数据进行序列比对，模拟数据基于人类参考基因组设置、包含100K条长度为100bp的序列。各软件/方法运行所需空间、时间、比对上参考序列的比例(Map rate)以及准确性相当。在该示例中，利用该实施方式中的比对方法所需的时间和内存较Bowtie或BWA的稍长和大，但利用该实施方式的比对方法的序列比对上的比例以及比对准确性达到98.9％和99.9％，均较利用Bowtie和BWA的稍高。
依据本发明的实施方式提供的一种比对装置，参考图4，该装置用以实现上述任一实施例/实施方式中的方法，该装置100包括：转化模块10，用于将每条读段转化成与该读段对应的一组短片段，获得多组短片段；查找模块20，用于确定所述短片段在参考库的对应位置，以获得第一定位结果，所述参考库为基于参考序列构建的哈希表，所述参考库包含多个条目，所述参考库的一个条目对应一条种子序列，所述种子序列能够与所述参考序列上的至少一段序列匹配，所述参考库的相邻两个条目对应的两条种子序列在所述参考序列上的距离小于所述短片段的长度；剔除模块30，用于去除所述第一定位结果中定位到所述参考库相邻条目中的任一条目上的短片段，获得第二定位结果；生长模块40，用于基于所述第二定位结果中来自相同读段的短片段进行延伸，以获得所述读段的比对结果。
上述对本发明任一实施方式中的比对方法的技术特征和效果的描述，同样适用本发明这一实施方式中的比对装置，在此不再赘述。
例如，根据本发明的具体实施方式，参考图5，该装置100还包括建库模块12，用于构建所述参考库，利用所述建库模块12进行以下：依据所述参考序列的碱基总数totalBase，确定种子序列的长度L，L＝μ*log(totalBase),
Figure BDA0001368381680000081
基于所述种子序列的长度，生成所有可能的种子序列，获得种子序列集；确定所述种子序列集中能够匹配到所述参考序列的种子序列以及该种子序列的匹配位置，以获得所述参考库。
根据本发明的具体实施方式，所称的确定种子序列集中能够匹配到参考序列的种子序列以及该种子序列的匹配位置，包括：利用大小为L的窗口对所述参考序列进行滑窗，将所述种子序列集中的种子序列与滑窗得的窗口序列进行匹配，以确定所述种子序列集中能够匹配到所述参考序列的种子序列以及该种子的匹配位置，进行所述匹配的容错率为ε1。
根据本发明的具体实施方式，进行所述滑窗的步长依据L和ε1来确定。
根据本发明的具体实施方式，进行所述滑窗的步长不小于L*ε1。
根据本发明的具体实施方式，所述参考库的相邻两个条目之间的距离大于所述滑窗的步长。
根据本发明的具体实施方式，利用所述转化模块10进行以下：利用大小为L的窗口对所述读段进行滑窗，以获得与该读段对应的一组短片段，进行所述滑窗的步长为1bp。
根据本发明的具体实施方式，参考图6，还包括第一筛选模块32，所述第一筛选模块32与所述剔除模块30相连，用于对来自剔除模块30的第二定位结果进行以下：去除连通长度小于预定阈值的短片段，以去除后的结果替代所述第二定位结果，所述连通长度为所述第二比对结果中的来自相同读段且定位到所述哈希表不同条目的短片段映射到参考序列的总长度。
根据本发明的具体实施方式，参考图7，还包括第二筛选模块34，所述第二筛选模块34与所述剔除模块30相连，用于：依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行评判，去除评判结果不符合预定要求的读段对应的短片段。
根据本发明的具体实施方式，所述第二筛选模块34用于：依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行计分，去除得分不大于第一预设值的读段。
根据本发明的具体实施方式，所述第二筛选模块34用于：依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行计分，去除得分不小于第二预设值的读段。
根据本发明的具体实施方式，所述生长模块40用于：基于所述来自相同读段的短片段进行序列构建，获得重构序列；基于所述重构序列与所述重构序列对应的参考序列的公共部分进行延伸，以获得延伸序列。
根据本发明的具体实施方式，所述公共部分为公共子串。所述生长模块40用于：查找所述重构序列与所述重构序列对应的参考序列的公共子串，确定所述重构序列和所述重构序列对应的参考序列的最长公共子串；基于编辑距离，延伸所述最长公共子串以获得延伸序列。
根据本发明的具体实施方式，所述公共部分为公共子序列。所述生长模块40用于：查找所述重构序列与所述重构序列对应的参考序列的公共子序列，确定所述重构序列和所述重构序列对应的参考序列的最长公共子序列；基于编辑距离，延伸所述最长公共子序列以获得延伸序列。
根据本发明的具体实施方式，参考图8，还包括截断模块50，用于：从来自所述生长模块40的延伸序列的至少一端对所述延伸序列进行截断，计算截断后的延伸序列的错误定位位点的比例，满足以下条件停止所述截断：截断后的延伸序列的错误定位位点的比例小于第三预设值。
根据本发明的具体实施方式，还包括截断模块50，用于：i、计算第一错误率和第二错误率，若所述第一错误率小于所述第二错误率，则从所述延伸序列的第一端对所述延伸序列进行截断，获得截断后的延伸序列，若所述第一错误率大于所述第二错误率，则从所述延伸序列的第二端对所述延伸序列进行截断，获得截断后的延伸序列，所述第一错误率为从所述延伸序列的第一端对所述延伸序列进行截断获得的截断后的延伸序列的错误定位位点的比例，所述第二错误率为从所述延伸序列的第二端对所述延伸序列进行截断、获得的截断后的延伸序列的错误定位位点的比例；ii、以截断后的延伸序列替代所述延伸序列进行i，直至所述截断后的延伸序列的错误定位位点的比例小于第四预设值。
根据本发明的具体实施方式，还包括截断模块50，用于：从所述延伸序列的至少一端对所述延伸序列进行滑窗，计算滑窗得的窗口序列的错误定位位点的比例，根据所述窗口序列的错误定位位点的比例对所述延伸序列进行截断，满足以下条件停止所述截断：滑窗得的窗口序列的错误定位位点的比例小于第五预设值。
根据本发明的具体实施方式，还包括截断模块50，用于：i、计算第三错误率和第四错误率，若所述第三错误率小于所述第四错误率，则从所述延伸序列的第二端对所述延伸序列进行截断，获得截断后的延伸序列，若所述第三错误率大于所述第四错误率，则从所述延伸序列的第一端对所述延伸序列进行截断，获得截断后的延伸序列，所述第三错误率为从所述延伸序列的第一端对所述延伸序列进行滑窗、获得的窗口序列的错误定位位点的比例，所述第四错误率为从所述延伸序列的第二端对所述延伸序列进行滑窗、获得的窗口序列的错误定位位点的比例；ii、以截断后的延伸序列替代所述延伸序列进行i，直至所述窗口序列的错误定位位点的比例大于第六预设值。
根据本发明的具体实施方式，所述截断的大小为1bp。
根据本发明的具体实施方式，所述滑窗的窗口不大于所述延伸序列的长度。
依据本发明的一种实施方式提供的一种计算机可读介质，该介质用以承载上述任一实施方式中的比对方法的一部分或者全部步骤。所称介质包括但不限于只读存储器、随机存储器、磁盘和光盘等。
依据本发明的一种实施方式提供的一种比对系统，该系统1000包括：输入装置100，用于输入数据；输出装置200，用于输出数据；处理器300，用于执行计算机可执行程序，执行所述计算机可执行程序包括完成上述任一实施方式的比对方法；存储装置400，用于存储数据，其中包括所述计算机可执行程序。
在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不一定指的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。
尽管已经示出和描述了本发明的实施例，本领域的普通技术人员可以理解：在不脱离本发明的原理和宗旨的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由权利要求及其等同限定。

Claims (51)

1.一种比对方法，其特征在于，包括：
将每条读段转化成与该读段对应的一组短片段，获得多组短片段；
确定所述短片段在参考库的对应位置，以获得第一定位结果，
所述参考库为基于参考序列构建的哈希表，所述参考库包含多个条目，所述参考库的一个条目对应一条种子序列，所述种子序列能够与所述参考序列上的至少一段序列匹配，
所述参考库的相邻两个条目对应的两条种子序列在所述参考序列上的距离小于所述短片段的长度；
去除所述第一定位结果中定位到所述参考库相邻条目中的任一条目上的短片段，获得第二定位结果；
基于所述第二定位结果中来自相同读段的短片段进行延伸，以获得所述读段的比对结果。
2.权利要求1的方法，其特征在于，所述参考库的构建包括：
依据所述参考序列的碱基总数totalBase，确定种子序列的长度L，L＝μ*log(totalBase),
Figure FDA0002968858020000011
基于所述种子序列的长度，生成所有可能的种子序列，获得种子序列集；
确定所述种子序列集中能够匹配到所述参考序列的种子序列以及该种子序列的匹配位置，以获得所述参考库。
3.权利要求2的方法，其特征在于，所述确定种子序列集中能够匹配到参考序列的种子序列以及该种子序列的匹配位置，包括：
利用大小为L的窗口对所述参考序列进行滑窗，将所述种子序列集中的种子序列与滑窗得的窗口序列进行匹配，以确定所述种子序列集中能够匹配到所述参考序列的种子序列以及该种子序列的匹配位置，进行所述匹配的容错率为ε1。
4.权利要求3的方法，其特征在于，进行所述滑窗的步长依据L和ε1来确定。
5.权利要求3的方法，其特征在于，进行所述滑窗的步长不小于L*ε1。
6.权利要求3的方法，其特征在于，所述参考库的相邻两个条目之间的距离大于或者等于所述滑窗的步长。
7.权利要求1-6任一方法，其特征在于，所述将每条读段转化成与该读段对应的一组短片段，获得多组短片段，包括：
利用大小为L的窗口对所述读段进行滑窗，以获得与该读段对应的一组短片段，进行所述滑窗的步长为1bp。
8.权利要求1-6任一方法，其特征在于，所述确定短片段在参考库的对应位置，以获得第一定位结果，包括：
将所述短片段与所述参考库的条目对应的种子序列进行匹配，以确定所述短片段在所述参考库的位置，进行所述匹配的容错率为ε2。
9.权利要求1-6任一方法，其特征在于，在获得所述第二定位结果之后，
去除连通长度小于预定阈值的短片段，以去除后的结果替代所述第二定位结果，所述连通长度为所述第二定位结果中的来自相同读段且定位到所述参考库不同条目的短片段映射到参考序列的总长度。
10.权利要求1-6任一方法，其特征在于，在获得所述第二定位结果之后，
依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行评判，去除评判结果不符合预定要求的读段。
11.权利要求10的方法，其特征在于，在获得所述第二定位结果之后，
依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行计分，去除得分不大于第一预设值的读段。
12.权利要求10的方法，其特征在于，在获得所述第二定位结果之后，
依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行计分，去除得分不小于第二预设值的读段。
13.权利要求1-6任一方法，其特征在于，所述基于第二定位结果中来自相同读段的短片段进行延伸，包括：
基于所述来自相同读段的短片段进行序列构建，获得重构序列；
基于所述重构序列与所述重构序列对应的参考序列的公共部分进行延伸，以获得延伸序列。
14.权利要求13的方法，其特征在于，所述公共部分为公共子串。
15.权利要求14的方法，其特征在于，所述基于第二定位结果中来自相同读段的短片段进行延伸，包括：
查找所述重构序列与所述重构序列对应的参考序列的公共子串，确定所述重构序列和所述重构序列对应的参考序列的最长公共子串；
基于编辑距离，延伸所述最长公共子串以获得延伸序列。
16.权利要求13的方法，其特征在于，所述公共部分为公共子序列。
17.权利要求16的方法，其特征在于，所述基于第二定位结果中来自相同读段的短片段进行延伸，包括：
查找所述重构序列与所述重构序列对应的参考序列的公共子序列，确定所述重构序列和所述重构序列对应的参考序列的最长公共子序列；
基于编辑距离，延伸所述最长公共子序列以获得延伸序列。
18.权利要求15或17的方法，其特征在于，还包括：
从所述延伸序列的至少一端对所述延伸序列进行截断，计算截断后的延伸序列的错误定位位点的比例，满足以下条件停止所述截断：截断后的延伸序列的错误定位位点的比例小于第三预设值。
19.权利要求15或17的方法，还包括：
i、计算第一错误率和第二错误率，若所述第一错误率小于所述第二错误率，则从所述延伸序列的第一端对所述延伸序列进行截断，若所述第一错误率大于所述第二错误率，则从所述延伸序列的第二端对所述延伸序列进行截断，以获得截断后的延伸序列，
所述第一错误率为从所述延伸序列的第一端对所述延伸序列进行截断获得的截断后的延伸序列的错误定位位点的比例，
所述第二错误率为从所述延伸序列的第二端对所述延伸序列进行截断、获得的截断后的延伸序列的错误定位位点的比例；
ii、以截断后的延伸序列替代所述延伸序列进行i，直至所述截断后的延伸序列的错误定位位点的比例小于第四预设值。
20.权利要求15或17的方法，其特征在于，还包括：
从所述延伸序列的至少一端对所述延伸序列进行滑窗，计算滑窗得的窗口序列的错误定位位点的比例，
根据所述窗口序列的错误定位位点的比例对所述延伸序列进行截断，满足以下条件停止所述截断：滑窗得的窗口序列的错误定位位点的比例大于第五预设值。
21.权利要求15或17的方法，其特征在于，还包括：
i、计算第三错误率和第四错误率，若所述第三错误率小于所述第四错误率，则从所述延伸序列的第二端对所述延伸序列进行截断，
若所述第三错误率大于所述第四错误率，则从所述延伸序列的第一端对所述延伸序列进行截断，以获得截断后的延伸序列，
所述第三错误率为从所述延伸序列的第一端对所述延伸序列进行滑窗、获得的窗口序列的错误定位位点的比例，
所述第四错误率为从所述延伸序列的第二端对所述延伸序列进行滑窗、获得的窗口序列的错误定位位点的比例；
ii、以截断后的延伸序列替代所述延伸序列进行i，直至所述窗口序列的错误定位位点的比例大于第六预设值。
22.权利要求18的方法，其特征在于，所述截断的大小为1bp。
23.权利要求20的方法，其特征在于，所述截断的大小为1bp。
24.权利要求20的方法，其特征在于，所述滑窗的窗口不大于所述延伸序列的长度。
25.权利要求21的方法，其特征在于，所述滑窗的窗口不大于所述延伸序列的长度。
26.一种比对装置，其特征在于，包括：
转化模块，用于将每条读段转化成与该读段对应的一组短片段，获得多组短片段；
查找模块，用于确定所述短片段在参考库的对应位置，以获得第一定位结果，
所述参考库为基于参考序列构建的哈希表，所述参考库包含多个条目，所述参考库的一个条目对应一条种子序列，所述种子序列能够与所述参考序列上的至少一段序列匹配，
所述参考库的相邻两个条目对应的两条种子序列在所述参考序列上的距离小于所述短片段的长度；
剔除模块，用于去除所述第一定位结果中定位到所述参考库相邻条目中的任一条目上的短片段，获得第二定位结果；
生长模块，用于基于所述第二定位结果中来自相同读段的短片段进行延伸，以获得所述读段的比对结果。
27.权利要求26的装置，其特征在于，还包括建库模块，用于构建所述参考库，利用所述建库模块进行以下处理：
依据所述参考序列的碱基总数totalBase，确定种子序列的长度L，L＝μ*log(totalBase),
Figure FDA0002968858020000041
基于所述种子序列的长度，生成所有可能的种子序列，获得种子序列集；
确定所述种子序列集中能够匹配到所述参考序列的种子序列以及该种子序列的匹配位置，以获得所述参考库。
28.权利要求27的装置，其特征在于，所述确定种子序列集中能够匹配到参考序列的种子序列以及该种子序列的匹配位置，包括：
利用大小为L的窗口对所述参考序列进行滑窗，将所述种子序列集中的种子序列与滑窗得的窗口序列进行匹配，以确定所述种子序列集中能够匹配到所述参考序列的种子序列以及该种子序列的匹配位置，进行所述匹配的容错率为ε1。
29.权利要求28的装置，其特征在于，进行所述滑窗的步长依据L和ε1来确定。
30.权利要求28的装置，其特征在于，进行所述滑窗的步长不小于L*ε1。
31.权利要求28的装置，其特征在于，所述参考库的相邻两个条目之间的距离大于所述滑窗的步长。
32.权利要求26-31任一装置，其特征在于，利用所述转化模块进行以下处理：
利用大小为L的窗口对所述读段进行滑窗，以获得与该读段对应的一组短片段，进行所述滑窗的步长为1bp。
33.权利要求26-31任一装置，其特征在于，还包括第一筛选模块，所述第一筛选模块与所述剔除模块相连，用于对来自剔除模块的第二定位结果进行以下处理：
去除连通长度小于预定阈值的短片段，以去除后的结果替代所述第二定位结果，所述连通长度为所述第二定位结果中的来自相同读段且定位到所述哈希表不同条目的短片段映射到参考序列的总长度。
34.权利要求26-31任一装置，其特征在于，还包括第二筛选模块，所述第二筛选模块与所述剔除模块相连，用于：
依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行评判，去除评判结果不符合预定要求的读段对应的短片段。
35.权利要求34的装置，其特征在于，所述第二筛选模块用于，
依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行计分，去除得分不大于第一预设值的读段。
36.权利要求34的装置，其特征在于，所述第二筛选模块用于，
依据所述第二定位结果中来自相同读段的短片段的定位结果，对该读段的定位结果进行计分，去除得分不小于第二预设值的读段。
37.权利要求26-31任一装置，其特征在于，所述生长模块用于：
基于所述来自相同读段的短片段进行序列构建，获得重构序列；
基于所述重构序列与所述重构序列对应的参考序列的公共部分进行延伸，以获得延伸序列。
38.权利要求37的装置，其特征在于，所述公共部分为公共子串。
39.权利要求38的装置，其特征在于，所述生长模块用于：
查找所述重构序列与所述重构序列对应的参考序列的公共子串，确定所述重构序列和所述重构序列对应的参考序列的最长公共子串；
基于编辑距离，延伸所述最长公共子串以获得延伸序列。
40.权利要求37装置，其特征在于，所述公共部分为公共子序列。
41.权利要求40的装置，其特征在于，所述生长模块用于：
查找所述重构序列与所述重构序列对应的参考序列的公共子序列，确定所述重构序列和所述重构序列对应的参考序列的最长公共子序列；
基于编辑距离，延伸所述最长公共子序列以获得延伸序列。
42.权利要求39或41的装置，其特征在于，还包括截断模块，用于：
从来自所述生长模块的延伸序列的至少一端对所述延伸序列进行截断，计算截断后的延伸序列的错误定位位点的比例，满足以下条件停止所述截断：截断后的延伸序列的错误定位位点的比例小于第三预设值。
43.权利要求39或41的装置，其特征在于，还包括截断模块，用于：
i、计算第一错误率和第二错误率，若所述第一错误率小于所述第二错误率，则从所述延伸序列的第一端对所述延伸序列进行截断，获得截断后的延伸序列，
若所述第一错误率大于所述第二错误率，则从所述延伸序列的第二端对所述延伸序列进行截断，获得截断后的延伸序列，
所述第一错误率为从所述延伸序列的第一端对所述延伸序列进行截断获得的截断后的延伸序列的错误定位位点的比例，
所述第二错误率为从所述延伸序列的第二端对所述延伸序列进行截断、获得的截断后的延伸序列的错误定位位点的比例；
ii、以截断后的延伸序列替代所述延伸序列进行i，直至所述截断后的延伸序列的错误定位位点的比例小于第四预设值。
44.权利要求39或41的装置，其特征在于，还包括截断模块，用于：
从所述延伸序列的至少一端对所述延伸序列进行滑窗，计算滑窗得的窗口序列的错误定位位点的比例，
根据所述窗口序列的错误定位位点的比例对所述延伸序列进行截断，满足以下条件停止所述截断：滑窗得的窗口序列的错误定位位点的比例小于第五预设值。
45.权利要求39或41的装置，其特征在于，还包括截断模块，用于：
i、计算第三错误率和第四错误率，若所述第三错误率小于所述第四错误率，则从所述延伸序列的第二端对所述延伸序列进行截断，获得截断后的延伸序列，
若所述第三错误率大于所述第四错误率，则从所述延伸序列的第一端对所述延伸序列进行截断，获得截断后的延伸序列，
所述第三错误率为从所述延伸序列的第一端对所述延伸序列进行滑窗、获得的窗口序列的错误定位位点的比例，
所述第四错误率为从所述延伸序列的第二端对所述延伸序列进行滑窗、获得的窗口序列的错误定位位点的比例；
ii、以截断后的延伸序列替代所述延伸序列进行i，直至所述窗口序列的错误定位位点的比例大于第六预设值。
46.权利要求42的装置，其特征在于，所述截断的大小为1bp。
47.权利要求44的装置，其特征在于，所述截断的大小为1bp。
48.权利要求44的装置，其特征在于，所述滑窗的窗口不大于所述延伸序列的长度。
49.权利要求45的装置，其特征在于，所述滑窗的窗口不大于所述延伸序列的长度。
50.一种比对系统，其特征在于，包括：
输入装置，用于输入数据；
输出装置，用于输出数据；
处理器，用于执行计算机可执行程序，执行所述计算机可执行程序包括完成权利要求1-25任一方法；
存储装置，用于存储数据，其中包括所述计算机可执行程序。
51.一种计算机可读存储介质，用于存储计算机可执行程序，执行该计算机可执行程序包括完成权利要求1-25任一方法。
CN201710652379.8A 2017-08-02 2017-08-02 比对方法、装置及系统 Active CN107403075B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710652379.8A CN107403075B (zh) 2017-08-02 2017-08-02 比对方法、装置及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710652379.8A CN107403075B (zh) 2017-08-02 2017-08-02 比对方法、装置及系统

Publications (2)

Publication Number Publication Date
CN107403075A CN107403075A (zh) 2017-11-28
CN107403075B true CN107403075B (zh) 2021-04-27

Family

ID=60402372

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710652379.8A Active CN107403075B (zh) 2017-08-02 2017-08-02 比对方法、装置及系统

Country Status (1)

Country Link
CN (1) CN107403075B (zh)

Families Citing this family (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108182348B (zh) * 2018-01-12 2020-04-24 广州医科大学附属第三医院（广州重症孕产妇救治中心、广州柔济医院） 基于种子序列信息的dna甲基化数据检测方法及其装置
CN111128303B (zh) * 2018-10-31 2023-09-15 深圳华大生命科学研究院 基于已知序列确定目标物种中对应序列的方法和系统
CN109887547B (zh) * 2019-03-06 2020-10-02 苏州浪潮智能科技有限公司 一种基因序列比对滤波加速处理方法、系统及装置
CN110016498B (zh) * 2019-04-24 2020-05-08 北京诺赛基因组研究中心有限公司 Sanger法测序中确定单核苷酸多态性的方法
CN110362823B (zh) * 2019-06-21 2023-07-28 北京百度网讯科技有限公司 描述文本生成模型的训练方法和装置
CN110517727B (zh) * 2019-08-23 2022-03-08 苏州浪潮智能科技有限公司 序列比对方法及系统
CN110517728B (zh) * 2019-08-29 2022-04-29 苏州浪潮智能科技有限公司 一种基因序列比对方法及装置
CN112825268B (zh) * 2019-11-21 2024-05-14 深圳华大基因科技服务有限公司 测序结果比对方法及其应用
CN112712850A (zh) * 2020-12-29 2021-04-27 中南大学 一种可应用于传染病病原体测序读段映射的种子序列定位方法
CN115862735B (zh) * 2022-12-28 2024-02-27 郑州思昆生物工程有限公司 一种核酸序列检测方法、装置、计算机设备及存储介质

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103793627A (zh) * 2012-10-29 2014-05-14 三星Sds株式会社 碱基序列比对系统及方法
CN106096332A (zh) * 2016-06-28 2016-11-09 深圳大学 面向存储的dna序列的并行快速匹配方法及其系统
WO2017123864A1 (en) * 2016-01-13 2017-07-20 Seven Bridges Genomics Inc. Systems and methods for analyzing circulating tumor dna

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103793627A (zh) * 2012-10-29 2014-05-14 三星Sds株式会社 碱基序列比对系统及方法
WO2017123864A1 (en) * 2016-01-13 2017-07-20 Seven Bridges Genomics Inc. Systems and methods for analyzing circulating tumor dna
CN106096332A (zh) * 2016-06-28 2016-11-09 深圳大学 面向存储的dna序列的并行快速匹配方法及其系统

Also Published As

Publication number Publication date
CN107403075A (zh) 2017-11-28

Similar Documents

Publication Publication Date Title
CN107403075B (zh) 比对方法、装置及系统
CN107609350B (zh) 一种二代测序数据分析平台的数据处理方法
CN110692101B (zh) 用于比对靶向的核酸测序数据的方法
WO2018218788A1 (zh) 一种基于全局种子打分优选的三代测序序列比对方法
US9489414B2 (en) Prefix burrows-wheeler transformations for creating and searching a merged lexeme set
US10810239B2 (en) Sequence data analyzer, DNA analysis system and sequence data analysis method
CN105760706A (zh) 一种二代测序数据的压缩方法
KR20070083641A (ko) 전사 맵핑을 위한 유전자 식별 기호 분석방법
CN115631789A (zh) 一种基于泛基因组的群体联合变异检测方法
CN117766020A (zh) 检测染色体非整倍性的方法、装置及系统
CN111919256A (zh) 检测染色体非整倍性的方法、装置及系统
CN109658981B (zh) 一种单细胞测序的数据分类方法
US11482304B2 (en) Alignment methods, devices and systems
Vaddadi et al. Read mapping on genome variation graphs
EP3938932B1 (en) Method and system for mapping read sequences using a pangenome reference
US20220171815A1 (en) System and method for generating filters for k-mismatch search
US9600599B2 (en) Prefix burrows-wheeler transformation with fast operations on compressed data
CN112397148A (zh) 序列比对方法、序列校正方法及其装置
US10867134B2 (en) Method for generating text string dictionary, method for searching text string dictionary, and system for processing text string dictionary
CN110866088B (zh) 一种语料库之间的快速全文检索方法及系统
Chen et al. CGAP-align: a high performance DNA short read alignment tool
US11929150B2 (en) Methods and apparatuses for performing character matching for short read alignment
CN115497567A (zh) 核酸序列聚类方法、装置、计算机可读存储介质、终端
Zhang MST Based Ab Initio Assembler of Expressed Sequence Tags
Khiste Efficient Alignment Algorithms for DNA Sequencing Data

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
REG Reference to a national code

Ref country code: HK

Ref legal event code: DE

Ref document number: 1239873

Country of ref document: HK

CB02 Change of applicant information
CB02 Change of applicant information

Address after: 518000 Shenye Jinyuan Building, No. 116 Qingshuihe Road, Qingshuihe Street, Luohu District, Shenzhen City, Guangdong Province, 2 5th and 6th floors

Applicant after: Shenzhen Zhenmai Biotechnology Co., Ltd.

Address before: 518000 First Floor of 111 High-tech Industrial Park, No. 72 Guowei Road, Liantang Street, Luohu District, Shenzhen City, Guangdong Province

Applicant before: SHENZHEN HANHAI GENE BIOTECHNOLOGY CO., LTD.

GR01 Patent grant
GR01 Patent grant