CN104424398A - 碱基序列对准系统及方法 - Google Patents
碱基序列对准系统及方法 Download PDFInfo
- Publication number
- CN104424398A CN104424398A CN201410196155.7A CN201410196155A CN104424398A CN 104424398 A CN104424398 A CN 104424398A CN 201410196155 A CN201410196155 A CN 201410196155A CN 104424398 A CN104424398 A CN 104424398A
- Authority
- CN
- China
- Prior art keywords
- seed
- short
- movie section
- mapping position
- mispairing
- 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.)
- Pending
Links
Classifications
-
- C—CHEMISTRY; METALLURGY
- C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
- C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
- C12Q1/00—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
- C12Q1/68—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving nucleic acids
- C12Q1/6869—Methods for sequencing
-
- 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
- G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
-
- 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
- G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
-
- 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
- G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
- G16B30/10—Sequence alignment; Homology search
Landscapes
- Life Sciences & Earth Sciences (AREA)
- Physics & Mathematics (AREA)
- Chemical & Material Sciences (AREA)
- Engineering & Computer Science (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Health & Medical Sciences (AREA)
- Bioinformatics & Cheminformatics (AREA)
- General Health & Medical Sciences (AREA)
- Biotechnology (AREA)
- Analytical Chemistry (AREA)
- Biophysics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Evolutionary Biology (AREA)
- Medical Informatics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Theoretical Computer Science (AREA)
- Organic Chemistry (AREA)
- Wood Science & Technology (AREA)
- Molecular Biology (AREA)
- Zoology (AREA)
- Genetics & Genomics (AREA)
- Microbiology (AREA)
- Immunology (AREA)
- Biochemistry (AREA)
- General Engineering & Computer Science (AREA)
- Image Analysis (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明公开一种碱基序列对准系统及方法。根据本发明的一个实施例的碱基序列对准系统包括:精确匹配模块,用于执行输入的短片段(read)与参考序列(reference sequence)的精确匹配(exact matching);二次匹配模块,当所述短片段无法精确匹配于所述参考序列时,通过考虑所述短片段与所述参考序列之间的错配(mismatch)而将所述短片段映射于所述参考序列;全局对准模块,当所述短片段无法通过所述二次匹配模块而得到映射时,执行所述短片段与所述参考序列的全局对准(global alignment)。
Description
技术领域
本发明的实施例涉及一种用于分析基因组的碱基序列的技术。
背景技术
用于生产大容量的短序列的下一代测序(NGS:Next GenerationSequencing)因其低廉的费用和快速数据生产方式而正在迅速地取代传统的桑格(Sanger)测序方式。并且,着重于准确率而开发出多种多样的NGS序列对准程序。
序列重组的第一个步骤为通过碱基序列对准(alignment)算法将短片段(read)映射(mapping)到参考序列的准确的位置。为此,普通的现有碱基序列对准算法首先将从短片段中选取的预定长度的种子(seed)映射到参考序列,然后在映射的位置上对其余短片段进行全局对准(Global Alignment)。
对于这种现有技术中的序列对准算法而言,需要在利用种子获得的参考序列内的所有候选位置上执行全局对准。然而全局对准的复杂度为O(N2),其执行时间很长,因此如果按照现有技术,则突出的问题在于碱基序列对准时间随着候选位置的增加而按照几何级数增长。
发明内容
本发明的实施例的目的在于提供一种用于在利用从测序仪输入的短片段进行序列对准时有效地减少需要很长的执行时间和很强的处理能力的全局对准的执行次数的技术方案。
根据本发明的一个实施例的一种碱基序列对准系统,包括:精确匹配模块,用于执行输入的短片段(read)对于参考序列(reference sequence)的精确匹配(exact matching);二次匹配模块,当所述短片段无法精确匹配到所述参考序列时,通过考虑所述短片段与所述参考序列之间的错配(mismatch)而将所述短片段映射到所述参考序列;全局对准模块,当通过所述二次匹配模块无法映射所述短片段时,执行所述短片段对于所述参考序列的全局对准(global alignment)。
所述碱基序列对准系统还可以包括:种子生成模块,当所述短片段无法精确匹配到所述参考序列时,从所述短片段生成多个种子(seed)。
所述种子生成模块可从所述短片段的整个区间生成所述多个种子。
所述种子生成模块从所述短片段的第一个碱基(base)开始移动设定的间距并读取与设定的大小一样大的所述短片段的值,从而可以生成所述多个种子。
所述种子生成模块可将所述多个种子生成为使生成的种子具有15bp~30bp的长度。
所述二次匹配模块可计算生成的各个种子在所述参考序列上的映射位置,并可以通过考虑在所述各个种子的映射位置上将所述短片段精确匹配到参考序列时发生的错配而确定所述短片段在所述参考序列上的映射位置。
所述二次匹配模块可将所述各个种子的映射位置当中的所述错配的数量最少的位置确定为所述短片段的映射位置。
所述二次匹配模块可将所述各个种子的映射位置当中的所述错配的质量分数之和最小的位置确定为所述短片段的映射位置。
所述二次匹配模块可将所述各个种子的映射位置当中的所述错配的数量为设定值以下的同时错配的质量分数之和最小的位置确定为所述短片段的映射位置。
所述全局对准模块可在所述各个种子的映射位置上执行所述短片段对于所述参考序列的全局对准。
所述全局对准模块可从所述各个种子的映射位置当中的所述错配的质量分数之和最小的映射位置开始依次执行所述全局对准。
所述全局对准模块可从所述各个种子的映射位置当中的所述错配的数量与间隙(gap)的数量之和为设定值以下的同时所述错配的数量与间隙的数量之和最小的映射位置开始依次执行所述全局对准。
所述全局对准模块可从所述各个种子的映射位置当中的所述错配的数量与间隙的数量之和为设定值以下的同时所述错配以及所述间隙的质量分数之和最小的映射位置开始依次执行所述全局对准。
另外,根据本发明的一个实施例的一种碱基序列对准方法,包括如下步骤:精确匹配步骤,由精确匹配模块执行输入的短片段对于参考序列的精确匹配;二次匹配步骤,当所述短片段无法精确匹配到所述参考序列时,由二次匹配模块通过考虑所述短片段与所述参考序列之间的错配而将所述短片段映射到所述参考序列;全局对准步骤,当通过所述二次匹配步骤无法映射所述短片段时,由全局对准模块执行所述短片段对于所述参考序列的全局对准。
所述的碱基序列对准方法在执行所述二次匹配步骤之前还可以包括如下步骤:种子生成步骤,当所述短片段无法精确匹配到所述参考序列时,从所述短片段生成多个种子。
在所述种子生成步骤中,可从所述短片段的整个区间生成所述多个种子。
在所述种子生成步骤中,可从所述短片段的第一个碱基开始移动设定的间距并读取与设定的大小一样大的所述短片段的值,从而生成所述多个种子。
在所述种子生成步骤中,可将所述多个种子生成为使生成的种子具有15bp~30bp的长度。
所述二次匹配步骤还可以包括如下步骤:计算生成的各个种子在所述参考序列上的映射位置;通过考虑在所述各个种子的映射位置上将所述短片段精确匹配到参考序列时发生的错配而确定所述短片段在所述参考序列上的映射位置。
在确定所述映射位置的步骤中,可将所述各个种子的映射位置当中的所述错配的数量最少的位置确定为所述短片段的映射位置。
在确定所述映射位置的步骤中,可将所述各个种子的映射位置当中的所述错配的质量分数之和最小的位置确定为所述短片段的映射位置。
在确定所述映射位置的步骤中,可将所述各个种子的映射位置当中所述错配的数量为设定值以下的同时错配的质量得分之和最小的位置确定为所述短片段的映射位置。
在所述全局对准步骤中,可在所述各个种子的映射位置上执行所述短片段对于所述参考序列的全局对准。
在所述全局对准步骤中,可从所述各个种子的映射位置当中的所述错配的质量分数之和最小的映射位置开始依次执行所述全局对准。
在所述全局对准步骤中,可从所述各个种子的映射位置当中的所述错配的数量与间隙的数量之和为设定值以下的同时所述错配的数量与间隙的数量之和最小的映射位置开始依次执行所述全局对准。
在所述全局对准步骤中,可从所述各个种子的映射位置当中的所述错配的数量与间隙的数量之和为设定值以下的同时所述错配以及所述间隙的质量分数之和最小的映射位置开始依次执行所述全局对准。
根据本发明的实施例,经过如下的步骤化的过程:首先尝试对由测序仪生成的整个短片段的精确匹配,并对没有精确匹配的短片段执行仅考虑错配(mismatch)的二次匹配,并只对所述二次匹配中也没有被映射的短片段进行选择性的全局对准,即,执行同时考虑错配和间隙(gap)的对准。此时,由于所述二次匹配基本上是精确匹配(exact matching)过程,因此处理速度明显快于具有O(N2)的复杂度的全局对准。即,根据本发明的实施例,可在执行全局对准之前通过精确匹配和二次匹配过程而事先筛选出与参考序列完全一致的短片段以及只存在一些错配的短片段,从而与简单地直接将短片段全局对准于参考序列的现有技术相比,可以有效地减少碱基序列对准时间。
而且,根据本发明的实施例,在进行所述二次匹配时通过考虑错配的质量得分而确定短片段的映射位置,从而具有可以减少碱基序列对准时间的同时能够维持碱基序列对准准确率的优点。
附图说明
图1为用于说明根据本发明的一个实施例的碱基序列对准方法100的图。
图2为用于举例表示在根据本发明的一个实施例的碱基序列对准方法100中计算mEB(误差数量最小值)的过程的图。
图3至图5为用于举例表示根据本发明的实施例的种子生成过程的图。
图6为用于举例说明根据本发明的一个实施例而将短片段精确匹配到参考序列时的错配(mismatch)的图。
图7为用于举例说明根据本发明的一个实施例的二次匹配过程的图。
图8为用于说明根据本发明的一个实施例的碱基序列对准系统800的模块图。
符号说明:
800:碱基序列对准系统 802:精确匹配模块
804:种子生成模块 806:二次匹配模块
808:全局对准模块
具体实施方式
以下,参照附图对本发明的具体实施方式进行说明。然而这仅仅是示例,本发明并不局限于此。
在对本发明进行说明时,如果认为对有关本发明的公知技术的具体说明有可能对本发明的主旨造成不必要的混乱,则省略其详细说明。另外,后述的术语为考虑到在本发明中的功能而定义的术语,其可能因使用者、运用者的意图或惯例等而不同。因此,要以整个说明书的内容为基础而对其进行定义。
本发明的技术思想由权利要求书确定,以下的实施例只是用于向本发明所属技术领域中具有普通知识的人员有效地说明本发明的技术思想的一种手段。
在对本发明的实施例进行详细说明之前,先对本发明中使用的术语进行如下说明。首先,“短片段(read)”是指由基因组测序仪(genome sequencer)输出的长度较短的碱基序列数据。短片段的长度通常根据测序仪的类型而大致构成为35~500bp(base pair,碱基对),通常对于DNA碱基而言是用字母A、C、G、T来表示。
“参考序列(reference sequence)”是指从所述短片段生成整个碱基序列时作为参照的碱基序列。在碱基序列分析中,通过参照参考序列对基因组测序仪中输出的大量短片段进行映射,来完成整个碱基序列。在本发明中,所述参考序列可以是在碱基序列分析时所预先设定的序列(例如,人类的整个碱基序列等),或者将基因组测序仪中制作出的碱基序列使用为参考序列。
“碱基(base)”为构成参考序列和短片段的最小单位。如前所述,对于DNA碱基而言可以由A、C、G、T这四种字母构成,将这些分别称为碱基。即,对于DNA碱基而言,通过四个碱基来表达,这对于短片段也一样。只是对于参考序列而言,由于多种多样的原因(测序错误、样本错误等),可能会出现无法确定应该用A、C、G或T中的哪种碱基来表示特定位置的碱基的情形,对于这种无法明确的碱基通常是用N等另外的文字进行标记。
“种子(seed)”是指为了短片段的映射而将短片段与参考序列进行比较时成为单位的序列。理论上,为了将短片段映射到参考序列,需要将整个短片段从参考序列的起始部分开始依次比较下去并计算短片段的映射位置。然而对于这种方法而言,映射一个短片段需要很长的时间和很强的计算能力,因此实际上,首先将作为由短片段的一部分构成的种子映射到参考序列,从而找出整个短片段的映射候选位置,并将整个短片段映射到对应的候选位置(Global Alignment,全局对准)。
图1为用于说明根据本发明的一个实施例的碱基序列对准方法100的图。在本发明的实施例中,碱基序列对准方法100是指通过将基因组测序仪输出的短片段与参考序列进行比较来确定短片段在所述参考序列上的映射(或对准)位置的一系列过程。如图所示,根据本发明的一个实施例的碱基序列对准方法100大体上分为包括如下三个步骤:对整个短片段的精确匹配步骤;对没有得到精确匹配的短片段进行二次匹配的步骤;对没有得到二次匹配的其余短片段进行全局对准的步骤。
首先,如果由基因组测序仪输入短片段(步骤102),则尝试整个短片段与所述参考序列的精确匹配(exact matching)(步骤104)。如果执行所述步骤104的结果针对整个短片段的精确匹配成功,则不再执行后续的对准步骤,并判断为对准成功(步骤106)。将人类的碱基序列作为对象而进行实验的结果显示,将基因组测序仪输出的100万个短片段精确匹配到人类的碱基序列时在一共为200万次的对准(正向序列100万次,反向互补(reverse complement)方向序列100万次)中发生了231,564次的精确匹配。因此,执行所述步骤104的结果大致可以减少11.6%的对准量。
然而,如果在所述步骤106中判断为对应短片段没有得到精确匹配,换言之,参考序列中不存在与短片段完全一致的区域,则接着从短片段中生成多个种子(步骤108),并考虑所述种子在参考序列的映射位置上的所述短片段与所述参考序列之间的错配(mismatch)而尝试将所述短片段映射到所述参考序列的二次匹配(步骤110)。如果所述步骤110的二次匹配结果存在至少一个满足二次匹配条件的映射位置,则将其中之一选为所述短片段的映射位置(步骤112)。即,在此情况下二次匹配成功。但如果在所述步骤112中不存在满足所述二次匹配条件的映射位置,则最后,在所述种子在所述参考序列中的映射位置上执行所述短片段对所述参考序列的全局对准(步骤114)。此时,在进行所述全局对准的结果,如果短片段的错误数量超过已设定的最大错误数量允许值(maxError)则判断为对准失败,否则判断为对准成功(步骤116)。
另外,虽然没有图示,然而如果在所述步骤106的判断结果是短片段没有精确匹配到参考序列,则在执行针对短片段的二次匹配之前,还可根据实施例而包括估计将相关短片段对准到所述参考序列时可能出现的错误数量的步骤。
在本发明的实施例中,所述错误数量的估计可以通过计算将所述短片段对准至所述参考序列时可能出现的错误数量最小值(mEB:minimum ErrorBound)而实现。图2为用于举例表示所述mEB计算过程的图。首先,如图2的(a)所示,将初始mEB设定为0,并从短片段的第一个碱基开始朝短片段的末端方向每次移动一个碱基来尝试精确匹配。此时,如图2的(b)所示,假定从短片段的特定碱基(图中以箭头表示的部分)处开始无法再进行精确匹配。这种情况表示从短片段的匹配起始位置到当前位置之间的区间某处发生了错误。因此,在此情况下将mEB增加1并从下一个位置开始新的精确匹配(图中表示为(c))。以后如果在特定位置上再次判断为无法精确匹配,则说明重新开始精确匹配的位置到当前位置之间的区间某处又发生了错误,因此将mEB再增加1并从下一个位置开始新的精确匹配(图中表示为(d))。经过这种过程而到达短片段的末尾(即,图中以(e)表示的情形)时的mEB成为对应短片段中可能存在的错误数量的最小值。
如果经过如上所述的过程而计算出短片段的mEB,接着便判断计算出的mEB是否超过已设定的最大错误数量允许值(MaxError),并在超过的情况下判断为针对相关短片段的对准失败并终止对准。在前述的以人类碱基序列为对象的实验中,将最大错误数量允许值取为3,并计算剩余短片段的mEB的结果显示,共有相当于844,891次的短片段超过所述最大错误数量允许值。即,执行所述mEB计算步骤的结果,可以减少大约42.2%的对准量。但如果所述判断的结果mEB为所述最大错误数量允许值以下,则依次执行前述步骤108以后的步骤。
以下详细说明所述步骤108至步骤116的具体过程。
从短片段生成多个种子
本步骤是为了正式执行短片段的对准而从短片段中生成作为多个小的片段的种子的步骤。在本步骤中,考虑所述短片段的一部分或全部而生成多个种子。
图3至图5为用于举例说明这种通过考虑短片段的整个区间而生成种子的方法的图。只是,在本发明中说明的种子生成方法仅仅是示例,本发明并不局限于特定的种子生成过程。例如,可通过将整个短片段或者其一部分特定区间分割为多个片段亦或是将分割的片段进行组合而生成种子。在此情况下,生成的种子可相互连续连接,然而并非一定要如此,也可以组合短片段内相互分离的片段来构成种子。而且,从一个短片段生成的种子不一定非要具有相同的长度,从一个短片段内也可生成具有多种长度的种子。简而言之,在本发明的实施例中,从短片段生成种子的方法并不特别受限,可不受限制地采用从短片段的一部分或整个短片段中提取种子的多种算法。
首先,图3为用于举例表示根据本发明的一个实施例的种子生成过程的图。如图所示,在本实施例中可通过将整个短片段分割为与设定的大小一样大的片段而生成种子。即,分割为预定长度的所述片段分别可以成为本发明的种子。在图中图示了将短片段分为6个片段的实施例,然而片段的数量以及每个片段的长度并不特别受限,其可以通过考虑参考序列的类型、短片段的长度、或者短片段的最大错误数量允许值等而适当地进行调整。并且,在图中只图示了将短片段分割为没有相互重叠部分(overlap)的示例,然而也可以将短片段分割为分割的各片段中存在一部分重叠的部分。
图4为用于举例表示根据本发明的另一实施例的种子生成过程的图。如图所示,在本实施例中,可以在将整个短片段分割为与设定大小一样大的片段之后,通过将分割的所述短片段的片段中的两个以上的片段进行组合而生成所述种子。例如,如图所示,在将短片段分割为4个片段(片段1~4)之后,将其两两组合的情况下可以生成总共6个种子。与前述的实施例相同,分割的片段的数量、各片段的长度、以及组合的片段数量等并不特别受限,其可以通过考虑参考序列的类型、短片段的长度、或者短片段的最大错误数量允许值等而适当地调节。
图5为用于举例表示根据本发明的又一实施例的种子生成过程的图。在本实施例中,从所述短片段的第一个碱基(base)开始移动设定的间距并读取与设定大小一样大的所述短片段的值,从而生成所述种子。在图示的实施例中,表示了短片段的长度为75bp(base pair,碱基对)、短片段的最大错误数量允许值为3bp、种子的大小(fragment size)为15bp、移动间距(shift size)为4bp的情况下的实施例。即,从短片段的第一个碱基开始以4bp为单位向右移动并生成种子。只是图示的实施例仅仅是示例性的,例如可通过考虑短片段的长度、短片段的最大错误数量允许值等而适当地确定所述移动间距、种子大小等。换言之,本发明的权利范围并不局限于特定的种子大小和移动间距。
另外,虽然在本发明的实施例中种子的长度并不特别受限,然而可以优选将所述种子的长度设定为所述短片段长度的20%~30%。通常,种子的长度越短,对应种子在参考序列上的映射数增加,而种子的长度越长,对应种子在参考序列上的映射数越少。通常情况下,在考虑由基因组测序仪生成的短片段的长度时,如果种子的长度为短片段长度的20%以下,则种子在参考序列上的映射数过多,因此出现后续的全局对准过程中全局对准次数不必要地增加的问题。相反,如果所述种子的长度为短片段长度的30%以上,则种子在参考序列上的映射数过少,从而使映射的准确率下降。因此在本发明中考虑到短片段的长度而将种子的长度构成为所述短片段长度的20%~30%,从而在确保映射的质量的同时使映射时可发生的复杂度最小化。
并且,如果所述参考序列为人类的碱基序列,则可以将所述种子生成为具有15bp~30bp的长度。如前所述,通常情况下种子的长度越短对应种子在参考序列上的映射数越多,而种子的长度越长对应种子在参考序列上的映射数越少。尤其,对于人类的碱基序列而言,在种子的长度为14以下的情况下参考序列内的映射位置的数量急剧增加。如下的表1表示基于种子长度的人类基因组内的种子平均出现频率。
[表1]
种子的长度 | 平均出现频率 |
10 | 2,726.1919 |
11 | 681.9731 |
12 | 170.9185 |
13 | 42.7099 |
14 | 10.6470 |
15 | 2.6617 |
16 | 0.6654 |
17 | 0.1664 |
从上面的表中可知,各个种子的出现频率在种子的长度为14以下时为10以上,然而在种子的长度为15时却减小为3以下。即,相比于将种子的长度构成为14以下的情形,将种子的长度构成为15以上时可以大幅度减少种子的重复。另外,如果所述种子的长度为30以上,则种子在参考序列上的映射数将过少,从而使映射的准确率下降。因此,在本发明中,在参考序列为人类的碱基序列的情况下将种子的长度构成为15~30,从而可以在确保映射的质量的同时使映射时可发生的复杂度最小化。
尝试短片段的二次匹配以及确定映射位置
如果按照以上所述从短片段生成种子,则接着,在将生成的所述种子分别映射到参考序列之后,在各个种子的映射位置上执行短片段的二次匹配。
在本发明的实施例中,短片段的二次匹配是指如下的一系列过程:从短片段生成种子,并在生成的各个种子在参考序列中的映射位置上将对应短片段的其余区间与参考序列互相比较而计算短片段的错配,并根据计算出的错配和已设定的二次匹配条件而确定短片段在参考序列上的映射位置。此时,由于对应短片段为通过执行所述步骤104的结果而被判断为无法完全匹配到参考序列的短片段,因此在种子的映射位置上将短片段与参考序列进行比较时必然会发生一个以上的不一致(错配,mismatch)。据此,在本发明的实施例中,利用在从短片段中生成的各个种子的映射位置上映射短片段时的错配而将所述种子的映射位置之一确定为短片段的映射位置。即,在本发明的实施例中,二次匹配相当于只考虑短片段的错配而不考虑短片段中可能发生的间隙(gap)的对准方法(ungapped alignment,非间隙对准)。
图6为用于举例说明根据本发明的一个实施例的将短片段精确匹配到参考序列时的错配的图。在图示的实施例中表示了将总长度为12bp的短片段精确匹配到参考序列的实施例,并将对应短片段中的前4bp部分假设为种子。对于图3所示的短片段而言,从第一个碱基开始的5个碱基与参考序列精确匹配,然而后面的第6个、第7个、第10个碱基却不是精确匹配于参考序列。即,对于图示的例而言,对应的映射位置上的短片段的错配数量为3。
在本发明的实施例中,可通过多种多样的方式考虑短片段的错配。在一个实施例中,所述步骤110可构成为在各个种子的映射位置上将短片段精确匹配于参考序列并计算这一过程中发生的错配数量。在此情况下,用于确定短片段的映射位置的二次匹配条件可以是错配的数量。例如,在所述种子的映射位置当中计算出的错配数量最少的位置可以成为所述短片段的映射位置。在这样通过考虑错配的数量而设定二次匹配条件时,可以获得将短片段映射到从算术角度上错误最少的位置的效果。
在另一实施例中,所述步骤110可以构成为考虑各个种子的每一映射位置的错配质量分数。即,在所述步骤110中,可以计算发生错配的位置上的质量分数(Quality Score)之和,并将计算出的质量分数之和最小的位置确定为所述短片段的映射位置。即,在此情况下二次匹配条件可以是错配的质量分数。短片段的质量分数是指将构成由基因组测序仪输出的短片段的各碱基的错误概率(error probability)换算为分数值而表示的值。计算短片段的质量得分的方法有多种,例如可以使用弗雷德分数(Phred Quality Score)等。然而本发明并不局限于特定的质量得分计算方法。与质量得分相关的详细事项为本领域普通技术人员所公知的,故在此省略其详细说明。
通常,在短片段中质量分数较低的位置表示那是发生错误的可能性相应地较高的位置。因此,越是错配的质量分数之和小的位置,短片段映射到对应位置的可能性越高,于是根据本实施例时,提高了将短片段映射到准确的位置的可能性。
另外,在又一实施例中,所述步骤110可被构成为一并考虑各个种子的每一映射位置的错配数量与错配的质量分数。在此情况下,在所述步骤110中可以首先从所述各个种子的映射位置中一次选取所述失配数量为设定值以下的位置,并将选取的位置当中失配的质量分数之和最小的位置确定为所述短片段的映射位置。
图7为用于举例表示这样的根据本发明的一个实施例的二次匹配过程的图。例如假设特定种子在参考序列的A、B、C位置上分别映射到参考序列,并假设对应位置上的短片段的错配数量(number of mismatches)以及错配的质量分数(QS)之和如图所示。如果二次匹配条件被设定为“错配数量为5以下的位置当中错配的质量分数之和最小的位置”,则由于满足这一条件的位置为A,因此可以确定对应短片段映射于参考序列的位置A。
短片段的全局对准
另外,即使进行所述二次匹配,也可能遇到不存在满足二次匹配条件的映射位置的情形。在此情况下,与一般的短片段映射方法相同,通过在所述各个种子的映射位置上执行所述短片段对所述参考序列的全局对准而将短片段映射于参考序列。在本发明的实施例中,全局对准为不仅考虑到短片段的错配而且还考虑了间隙的对准方法,例如可以使用史密斯-沃特曼(Smith-Waterman)算法或内德勒曼-文施(Needleman-Wunsch)算法等,然而本发明的实施例并不局限于特定算法。
在一个实施例中,在所述全局对准步骤中,可以从所述各个种子的映射位置当中所述错配质量分数之和最小的映射位置开始依次执行所述全局对准。这是由于质量分数之和越小,映射于对应位置的可能性越高。
而且,所述全局对准的顺序还可以通过考虑所述错配的同时一并考虑将短片段匹配到参考序列时出现的间隙而被确定。例如在所述全局对准步骤中,可以从所述各个种子的映射位置当中所述错配的数量与间隙的数量之和为设定值以下的同时所述错配的数量与间隙的数量之和最小的映射位置开始依次执行所述全局对准。或者,在所述全局对准步骤中,也可以从所述各个种子的映射位置当中所述错配的数量与间隙的数量之和为设定值以下的同时出现所述错配以及所述间隙的位置的所述短片段的质量分数之和最小的映射位置开始依次执行所述全局对准。
图8为用于说明根据本发明的一个实施例的碱基序列对准系统800的模块图。如图所示,根据本发明的一个实施例的碱基序列对准系统800包括精确匹配模块802、种子生成模块804、二次匹配模块806以及全局对准模块808。
精确匹配模块802执行输入的短片段与参考序列的精确匹配。
当所述短片段没有通过精确匹配模块802而精确匹配到所述参考序列时,种子生成模块804从所述短片段中生成多个种子。关于种子生成模块804中生成种子的具体方法已在前面叙述。
在所述短片段没有精确匹配于所述参考序列时,二次匹配模块806通过考虑所述短片段与所述参考序列之间的错配而将所述短片段映射于所述参考序列。二次匹配模块806计算生成的所述各个种子在所述参考序列上的映射位置,并通过考虑在所述各个种子的映射位置上将所述短片段精确匹配于参考序列时发生的错配而可以确定所述短片段在所述参考序列上的映射位置。
在一个实施例中,二次匹配模块806可将所述各个种子的映射位置当中所述错配数量最少的位置确定为所述短片段的映射位置。另外,二次匹配模块806可将所述各个种子的映射位置当中所述错配质量分数之和最小的位置确定为所述短片段的映射位置,或者可以将所述各个种子的映射位置当中所述错配数量为设定值以下的同时错配的质量分数之和最小的位置确定为所述短片段的映射位置。
在所述短片段甚至通过所述二次匹配模块也未得到映射时,全局对准模块808可以执行所述短片段对所述参考序列的全局对准。如前所述,全局对准模块808可在所述各个种子的映射位置上执行所述短片段对所述参考序列的全局对准,在此情况下,可以从所述各个种子的映射位置当中所述错配的质量分数之和最小的映射位置开始依次执行所述全局对准。当然,如前所述,在此情况下也可以只考虑间隙的数量与错配的数量之和为设定值以下的映射位置而执行全局对准。
以下的表2为通过将根据本发明实施例的碱基序列对准方法与现有技术(即,只进行全局对准的情形)进行比较而用于说明本发明的有益效果的表。为了进行比较,分别计算出了将分别具有75bp的长度的100万个短片段对准于参考序列时的映射时间、映射率以及误差率。
[表2]
指标 | 现有技术 | 本发明 |
映射时间 | 00:58:52 | 00:07:57 |
映射率 | 91.11% | 93.52% |
误差率 | 3.89% | 3.90% |
从上面的表中可知,在应用本发明时与现有技术相比映射时间从58分52秒显著减少到7分57秒。这是由于根据本发明的实施例时在执行全局对准之前能够在精确匹配步骤和二次匹配步骤中确定相当多的短片段的映射位置。即,根据本发明,可通过减少需要消耗大量时间的全局对准的执行次数而提高碱基序列对准速度。
而且,在映射率和误差率方面,本发明显示出比现有技术稍好或者具有相近的数值。即,根据上述对比实验结果可知本发明的实施例具有能够维持映射质量的同时可提高映射速度的效果。
另外,本发明的实施例可以包括记录有用于在计算机上执行本说明书中记载的方法的程序的计算机可读记录介质。所述计算机可读记录介质可将程序命令、本地数据文件、本地数据结构等单独或组合而包含在内。所述介质可以是为了本发明而特别设计并构成的,也可以是被计算机软件领域中具有普通知识的人员所公知而可以使用的。计算机可读记录介质的例中有硬盘、软盘、磁带之类的磁介质;CD-ROM、DVD之类的光记录介质;软盘之类的磁光介质以及ROM、RAM、闪存等为了存储并执行程序命令而特别构成的硬件装置。程序命令的例中不仅包括通过编译器制作的机器语言代码,而且还可以包括利用解释器并通过计算机而得以执行的高级语言代码。
以上已通过代表性实施例而对本发明进行了详细说明,然而本发明所属技术领域中具有普通知识的人员应该会理解可以在不脱离本发明范围的限度内对所述的实施例进行多种多样的变形。
因此,本发明的权利范围不应局限于所述的实施例而确定,而是要根据权利要求书及其等价内容来确定。
Claims (26)
1.一种碱基序列对准系统,包括:
精确匹配模块,用于执行输入的短片段对于参考序列的精确匹配;
二次匹配模块,当所述短片段无法精确匹配到所述参考序列时,通过考虑所述短片段与所述参考序列之间的错配而将所述短片段映射到所述参考序列;
全局对准模块,当通过所述二次匹配模块无法映射所述短片段时,执行所述短片段对于所述参考序列的全局对准。
2.如权利要求1所述的碱基序列对准系统,其中,还包括:
种子生成模块,当所述短片段无法精确匹配到所述参考序列时,从所述短片段中生成多个种子。
3.如权利要求2所述的碱基序列对准系统,其中,所述种子生成模块从所述短片段的整个区间生成所述多个种子。
4.如权利要求2所述的碱基序列对准系统,其中,所述种子生成模块从所述短片段的第一个碱基开始移动设定的间距并读取与设定的大小一样大的所述短片段的值,从而生成所述多个种子。
5.如权利要求2所述的碱基序列对准系统,其中,所述种子生成模块将所述多个种子生成为使生成的种子具有15bp~30bp的长度。
6.如权利要求2所述的碱基序列对准系统,其中,所述二次匹配模块计算生成的各个种子在所述参考序列上的映射位置,并通过考虑在所述各个种子的映射位置上将所述短片段精确匹配到参考序列时发生的错配而确定所述短片段在所述参考序列上的映射位置。
7.如权利要求6所述的碱基序列对准系统,其中,所述二次匹配模块将所述各个种子的映射位置当中的所述错配的数量最少的位置确定为所述短片段的映射位置。
8.如权利要求6所述的碱基序列对准系统,其中,所述二次匹配模块将所述各个种子的映射位置当中的所述错配的质量分数之和最小的位置确定为所述短片段的映射位置。
9.如权利要求6所述的碱基序列对准系统,其中,所述二次匹配模块将所述各个种子的映射位置当中的所述错配的数量为设定值以下的同时错配的质量分数之和最小的位置确定为所述短片段的映射位置。
10.如权利要求6所述的碱基序列对准系统,其中,所述全局对准模块在所述各个种子的映射位置上执行所述短片段对于所述参考序列的全局对准。
11.如权利要求10所述的碱基序列对准系统,其中,所述全局对准模块从所述各个种子的映射位置当中的所述错配的质量分数之和最小的映射位置开始依次执行所述全局对准。
12.如权利要求10所述的碱基序列对准系统,其中,所述全局对准模块从所述各个种子的映射位置当中的所述错配的数量与间隙的数量之和为设定值以下的同时所述错配的数量与间隙的数量之和最小的映射位置开始依次执行所述全局对准。
13.如权利要求10所述的碱基序列对准系统,其中,所述全局对准模块从所述各个种子的映射位置当中所述错配的数量与间隙的数量之和为设定值以下的同时所述错配以及所述间隙的质量得分之和最小的映射位置开始依次执行所述全局对准。
14.一种碱基序列对准方法,包括如下步骤:
精确匹配步骤,由精确匹配模块执行输入的短片段对于参考序列的精确匹配;
二次匹配步骤,当所述短片段无法精确匹配到所述参考序列时,由二次匹配模块通过考虑所述短片段与所述参考序列之间的错配而将所述短片段映射到所述参考序列;
全局对准步骤,当通过所述二次匹配步骤无法映射所述短片段时,由全局对准模块执行所述短片段对于所述参考序列的全局对准。
15.如权利要求14所述的碱基序列对准方法,其中,在执行所述二次匹配步骤之前,还包括如下步骤:
种子生成步骤,当所述短片段无法精确匹配到所述参考序列时,从所述短片段中生成多个种子。
16.如权利要求15所述的碱基序列对准方法,其中,在所述种子生成步骤中,从所述短片段的整个区间中生成所述多个种子。
17.如权利要求15所述的碱基序列对准方法,在所述种子生成步骤中,从所述短片段的第一个碱基开始移动设定的间距并读取与设定的大小一样大的所述短片段的值,从而生成所述多个种子。
18.如权利要求15所述的碱基序列对准方法,其中,在所述种子生成步骤中,将所述多个种子生成为使生成的种子具有15bp~30bp的长度。
19.如权利要求15所述的碱基序列对准方法,其中,所述二次匹配步骤还包括如下步骤:
计算生成的各个种子在所述参考序列上的映射位置;
通过考虑在所述各个种子的映射位置上将所述短片段精确匹配到参考序列时发生的错配而确定所述短片段在所述参考序列上的映射位置。
20.如权利要求19所述的碱基序列对准方法,其中,在确定所述映射位置的步骤中,将所述各个种子的映射位置当中的所述错配的数量最少的位置确定为所述短片段的映射位置。
21.如权利要求19所述的碱基序列对准方法,其中,在确定所述映射位置的步骤中,将所述各个种子的映射位置当中的所述错配的质量分数之和最小的位置确定为所述短片段的映射位置。
22.如权利要求19所述的碱基序列对准方法,其中,在确定所述映射位置的步骤中,将所述各个种子的映射位置当中的所述错配的数量为设定值以下的同时错配的质量分数之和最小的位置确定为所述短片段的映射位置。
23.如权利要求19所述的碱基序列对准方法,其中,在所述全局对准步骤中,在所述各个种子的映射位置上执行所述短片段对于所述参考序列的全局对准。
24.如权利要求23所述的碱基序列对准方法,其中,在所述全局对准步骤中,从所述各个种子的映射位置当中的所述错配的质量分数之和最小的映射位置开始依次执行所述全局对准。
25.如权利要求23所述的碱基序列对准方法,其中,在所述全局对准步骤中,从所述各个种子的映射位置当中的所述错配的数量与间隙的数量之和为设定值以下的同时所述错配的数量与间隙的数量之和最小的映射位置开始依次执行所述全局对准。
26.如权利要求23所述的碱基序列对准方法,其中,在所述全局对准步骤中,从所述各个种子的映射位置当中的所述错配的数量与间隙的数量之和为设定值以下的同时所述错配以及所述间隙的质量分数之和最小的映射位置开始依次执行所述全局对准。
Applications Claiming Priority (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
KR10-2013-0105529 | 2013-09-03 | ||
KR1020130105529A KR101584857B1 (ko) | 2013-09-03 | 2013-09-03 | 염기 서열 정렬 시스템 및 방법 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN104424398A true CN104424398A (zh) | 2015-03-18 |
Family
ID=52584376
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201410196155.7A Pending CN104424398A (zh) | 2013-09-03 | 2014-05-09 | 碱基序列对准系统及方法 |
Country Status (3)
Country | Link |
---|---|
US (1) | US20150066384A1 (zh) |
KR (1) | KR101584857B1 (zh) |
CN (1) | CN104424398A (zh) |
Families Citing this family (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
KR101508817B1 (ko) * | 2012-10-29 | 2015-04-08 | 삼성에스디에스 주식회사 | 염기 서열 정렬 시스템 및 방법 |
KR102377810B1 (ko) | 2020-03-11 | 2022-03-23 | 한국전자통신연구원 | 유전체 서열 정렬 방법 및 그 장치 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101539967A (zh) * | 2008-12-12 | 2009-09-23 | 深圳华大基因研究院 | 一种单核苷酸多态性检测方法 |
CN101751517A (zh) * | 2008-12-12 | 2010-06-23 | 深圳华大基因研究院 | 一种基因组短序列映射的快速处理方法及系统 |
CN101887493A (zh) * | 2010-08-05 | 2010-11-17 | 湖州瑞万思信息技术有限公司 | 一种基于单比特化的基因比对处理方法 |
US20110136686A1 (en) * | 2009-12-09 | 2011-06-09 | Oracle International Corporation | Sequence matching allowing for errors |
-
2013
- 2013-09-03 KR KR1020130105529A patent/KR101584857B1/ko active IP Right Grant
-
2014
- 2014-05-08 US US14/273,205 patent/US20150066384A1/en not_active Abandoned
- 2014-05-09 CN CN201410196155.7A patent/CN104424398A/zh active Pending
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101539967A (zh) * | 2008-12-12 | 2009-09-23 | 深圳华大基因研究院 | 一种单核苷酸多态性检测方法 |
CN101751517A (zh) * | 2008-12-12 | 2010-06-23 | 深圳华大基因研究院 | 一种基因组短序列映射的快速处理方法及系统 |
US20110136686A1 (en) * | 2009-12-09 | 2011-06-09 | Oracle International Corporation | Sequence matching allowing for errors |
CN101887493A (zh) * | 2010-08-05 | 2010-11-17 | 湖州瑞万思信息技术有限公司 | 一种基于单比特化的基因比对处理方法 |
Non-Patent Citations (3)
Title |
---|
GREGORY G.FAUST ET AL.: "YAHA:fast and flexible long-read alignment with optimal breakpoint detection", 《BIOINFORMATICS ORIGINAL PAPER》 * |
HENG LI ET AL.: "Mapping short DNA sequencing reads and calling variants using mapping quality scores", 《GENOME RESOURCE》 * |
M.KIELBASA ET AL.: "Adaptive seeds tame genomic sequence comparison", 《GENOME RESEARCH》 * |
Also Published As
Publication number | Publication date |
---|---|
KR101584857B1 (ko) | 2016-01-12 |
KR20150026542A (ko) | 2015-03-11 |
US20150066384A1 (en) | 2015-03-05 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
US10204207B2 (en) | Systems and methods for transcriptome analysis | |
WO2016141294A1 (en) | Systems and methods for genomic pattern analysis | |
AU2014340461A1 (en) | Systems and methods for using paired-end data in directed acyclic structure | |
CN103793627B (zh) | 碱基序列比对系统及方法 | |
JP2019512795A5 (zh) | ||
KR101481457B1 (ko) | 리드 전체를 고려한 염기 서열 정렬 시스템 및 방법 | |
CN107480470B (zh) | 基于贝叶斯与泊松分布检验的已知变异检出方法和装置 | |
CN108595915B (zh) | 一种基于dna变异检测的三代数据校正方法 | |
Voshall et al. | Next-generation transcriptome assembly: strategies and performance analysis | |
RU2018109529A (ru) | Многоуровневая архитектура распознавания паттернов в биологических данных | |
US20140121983A1 (en) | System and method for aligning genome sequence | |
CN104424398A (zh) | 碱基序列对准系统及方法 | |
US20150142328A1 (en) | Calculation method for interchromosomal translocation position | |
CN103793626B (zh) | 碱基序列比对系统及方法 | |
KR101394339B1 (ko) | 시드의 길이를 고려한 염기 서열 처리 시스템 및 방법 | |
CN104239748A (zh) | 考虑失配的碱基序列对准系统及方法 | |
CN104239749A (zh) | 碱基序列对准系统及方法 | |
KR101538852B1 (ko) | 정확도를 고려한 염기 서열 정렬 장치 및 방법 | |
Dutheil et al. | Optimization of sequence alignments according to the number of sequences vs. number of sites trade-off | |
US20140121992A1 (en) | System and method for aligning genome sequence | |
CN117238368B (zh) | 分子遗传标记分型方法和装置、生物个体识别方法和装置 | |
CN111583997B (zh) | 杂合变异下校正第三代测序数据中测序错误的混合方法 | |
CN103793624A (zh) | 考虑重复性的碱基序列重组系统及方法 | |
Agarwal et al. | CPDP: A connection based PDP algorithm | |
Jean et al. | The contig assembly problem and its algorithmic solutions |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
WD01 | Invention patent application deemed withdrawn after publication |
Application publication date: 20150318 |
|
WD01 | Invention patent application deemed withdrawn after publication |