CN105483244A - 一种基于超长基因组的变异检测算法及检测系统 - Google Patents
一种基于超长基因组的变异检测算法及检测系统 Download PDFInfo
- Publication number
- CN105483244A CN105483244A CN201511005316.0A CN201511005316A CN105483244A CN 105483244 A CN105483244 A CN 105483244A CN 201511005316 A CN201511005316 A CN 201511005316A CN 105483244 A CN105483244 A CN 105483244A
- Authority
- CN
- China
- Prior art keywords
- reference sequences
- local matching
- local
- sequenced fragments
- matching event
- 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.)
- Granted
Links
- 238000001514 detection method Methods 0.000 title claims abstract description 32
- 238000000034 method Methods 0.000 claims abstract description 6
- 239000012634 fragment Substances 0.000 claims description 98
- 238000012163 sequencing technique Methods 0.000 claims description 82
- 238000012217 deletion Methods 0.000 claims description 8
- 230000037430 deletion Effects 0.000 claims description 8
- 230000008878 coupling Effects 0.000 claims description 6
- 238000010168 coupling process Methods 0.000 claims description 6
- 238000005859 coupling reaction Methods 0.000 claims description 6
- 238000002865 local sequence alignment Methods 0.000 claims description 6
- 238000005516 engineering process Methods 0.000 description 7
- 241000894007 species Species 0.000 description 6
- 108020004414 DNA Proteins 0.000 description 4
- 241001282153 Scopelogadus mizolepis Species 0.000 description 4
- 230000009286 beneficial effect Effects 0.000 description 4
- 230000008034 disappearance Effects 0.000 description 4
- 238000003780 insertion Methods 0.000 description 4
- 230000037431 insertion Effects 0.000 description 4
- 230000035772 mutation Effects 0.000 description 4
- 239000002773 nucleotide Substances 0.000 description 4
- 125000003729 nucleotide group Chemical group 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 108090000623 proteins and genes Proteins 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- 238000001712 DNA sequencing Methods 0.000 description 2
- 241000196324 Embryophyta Species 0.000 description 2
- 108091028043 Nucleic acid sequence Proteins 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 2
- 150000001875 compounds Chemical class 0.000 description 2
- 230000006872 improvement Effects 0.000 description 2
- 230000035945 sensitivity Effects 0.000 description 2
- 241001269238 Data Species 0.000 description 1
- 241000206602 Eukaryota Species 0.000 description 1
- 230000004075 alteration Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000008901 benefit Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000006073 displacement reaction Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000013467 fragmentation Methods 0.000 description 1
- 238000006062 fragmentation reaction Methods 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 230000017105 transposition Effects 0.000 description 1
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
Landscapes
- Life Sciences & Earth Sciences (AREA)
- Chemical & Material Sciences (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Organic Chemistry (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Physics & Mathematics (AREA)
- Biotechnology (AREA)
- Analytical Chemistry (AREA)
- Genetics & Genomics (AREA)
- Zoology (AREA)
- Biophysics (AREA)
- Molecular Biology (AREA)
- General Health & Medical Sciences (AREA)
- Wood Science & Technology (AREA)
- Evolutionary Biology (AREA)
- Theoretical Computer Science (AREA)
- Immunology (AREA)
- Microbiology (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Medical Informatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Biochemistry (AREA)
- General Engineering & Computer Science (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明涉及一种基于超长基因组的变异检测算法，简称VariationBlast算法，在得到长序列的情况下，大尺度结构变异通常可以通过序列与参考基因组进行比对来检测到，因为跨越结构变异的序列会产生比对上参考序列的部分片段，然后通过序列部分片段与相对应的参考序列片段之间的比较，结构变异的精确点位可以被检测到，VariationBlast是通过逐次对比对方法来检测每一个序列与参考基因组序列之间的比对，然后把所有代表结构变异的序列进行分类与筛选最后从比对的位点与方向得到可能的结构变异以及它们各自的类型。
Description
技术领域
本发明涉及基因序列结构变异检测，具体涉及到一种基于超长基因组的变异检测算法及检测系统。
背景技术
上世纪九十年代启动的人类基因组计划(HGP)，加上随之而来的千人基因组计划，DNA元件百科全书计划(ENCODE)等的实施加快了基因组时代的发展。第二代、三代DNA促使许多物种的基因组测序计划得以成功完成，从而积累了大量的生物数据。这些生物大数据必须通过适当的分析手段才能挖掘出具有潜在理论价值和应用价值的信息。基因组序列多态性，是指物种群体内及群体间DNA序列和结构的差异。人类的这些基因组差异决定了不同种族、群体和个体间基因组的差异或多态性。基因组多样性包括单核苷酸多态性(Singlenucleotidepolymorphism,SNP)、小片断插入和缺失(smallindel)；基因组结构变异通常指DNA序列长度大于1kb的基因组差异,包括DNA片段插入(insertion)与缺失(deletion)、移位(transposition)、倒位(inversion)以及DNA拷贝数目变化(Copynumbervariations,CNVs)等。结构变异由于覆盖的核苷酸总数远超过SNP的总数,因此对个体表型的影响可能更大。科学家们开发了很多用于预测基因组结构变异的软件，而且这些软件已被广泛应用到人类基因组多态性研究中。但是已开发的绝大多数软件只能有效地预测基因组单核苷酸多态性，小片段结构变异，以及结构变异中某种或某几种类型。这些软件对某些大尺度的结构变异(例如，长片段插入、缺失序列等)不能很好的挖掘，以至不同软件预测结果之间的差异很大。常用的变异检测软件，例如Pindel、Breakdancer、CNVnator等都无法检测大型的基因组变异，少部分算法可以检测到基因组大型缺失，但是几乎所有的算法流程都无法鉴定出大型插入、移位、倒置等大型的基因组结构变异。本发明中开发的VariationBlast算法，有效地将具有高灵敏度的比对算法和图形理论相结合，实现了碱基水平检测基因组变异，变异检测范围可以从最小的单核苷酸改变到大型复杂的基因组结构变异。
发明内容
本发明提供一种基于超长基因组的变异检测算法，简称VariationBlast算法，解决了现有算法不能准确检测大型超长基因组的结构变异的问题。
本发明解决上述技术问题的技术方案如下：
一种基于超长基因组的变异检测算法，包括以下步骤：
S1、采用局部序列比对算法，检测出测序片段和参考序列间存在的所有匹配，得到局部匹配事件，每一个局部匹配事件均包括测序片段和参考序列上的参考片段；
S2、将所有局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，将测序片段比对至参考序列上的位置重叠或顺次相连的局部匹配事件分到一个组，有多少能够相连或重叠的测序片段，就有多少个组；
S3、对每一个组内的局部匹配事件按照不同的排序方式进行评分，建立一个用来表示每一个组内的所有局部匹配事件之间不同的排序方式的图形模型，根据评分情况确定每一组内的局部匹配事件的最佳排序方式；
S4、按照最佳排序方式将每一个组内的局部匹配事件中的测序片段进行串联，一个组得到一个测序序列read；
S5、将每一组的测序序列read和参考序列进行对比，在测序序列read上和参考序列不匹配的区域添加gap区域在参考序列上和测序序列read不匹配的区域添加gap区域，gap区域即空白区域，用于使测序序列read和参考序列的匹配区域位置完全对应，便于确认变异的类型；
S6、对添加了gap区域后的最终测序序列进行结构变异检测。
本发明的有益效果是：能够充分利用由不断提高(例如测序长度不断增加)的第二代以及第三代DNA测序技术所带来的可能性，利用测序结果检测各种尺度的结构变异，包括其它软件所不能够预测的大尺度插入以及所有现有软件都忽略的复合变异。本发明的VariationBlast将大大提高基因组变异和结构变异预测的灵敏度、效率和准确度，可以有效地应用于来自不同物种、不同复杂度的基因组序列，从而为各个领域的科学研究和应用提供强急需的强有力的支撑。
在上述技术方案的基础上，本发明还可以做如下改进。
进一步的，述步骤S1里采用的局部序列比对算法，为BLASTn、MegaBlast、BLASTz、cross_match、BLAT和wublast算法中的一种。
进一步的，所述步骤S2中需要同时考虑了测序片段两条链上的情况，来建立局部匹配事件之间的关系模型。
采用上述进一步方案的有益效果是：本发明可检测出倒置等结构变异。
进一步的，所述步骤S4中的具体实现为：
S4.1、根据打分情况确定每一组内的局部匹配事件的最佳排序方式后，
将每一组内的局部匹配事件按照最佳排序方式进行排序，并将所有局部匹配事件中的测序片段串联成一个基础测序序列；
S4.2、对基础测序序列中的相邻的测序片段的重叠部分进行剔除，得到测序序列read。
进一步的，所述步骤S4.2中，对基础测序序列中的相邻的测序片段的重叠部分进行剔除的方法是，若相邻两个测序片段之间存在重叠部分，分别将两个测序片段的重叠部分和对应的参考片段进行对比，选取匹配度更高的一个重叠部分进行保留，并把另一个重叠部分剔除。
进一步的，步骤S4.2中对基础测序序列中的相邻的测序片段的重叠部分进行剔除后，还包括将基础测序序列和参考序列进行对比，若存在缺失片段，检测剩余未归入组内的测序片段与缺失片段的匹配度，选取匹配度最高的测序片段放入相应的基础测序序列内，再得到每一组的测序序列read的步骤。
采用上述进一步方案的有益效果是：本发明可检测出转座结构变异。
进一步的，步骤S3中，将符合线性连接关系的局部匹配事件分到一个组时,每一个组内的局部匹配事件的重叠度不大于阈值t。
进一步的，所述步骤S3的具体实现为：
S3.1、将每个组内的局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，所有的局部匹配事件按照位置从前到后依次排列；
S3.2、按照所有的局部匹配事件均可以与位置排在其后的局部匹配事件连接的条件，按照所有的可能性建立排序方式，构建有向非循环图形模型；
S3.3、设图形模型中，在前局部匹配事件和在后局部匹配事件之间的关系为边e，则每个边e的减分Length(e)＝P(e)–R(e)，一个边e的减分越少，则评分越高；
其中，P(e)代表边e的减分、R(e)代表边e的得分；
R(e)＝PID*len，PID在后的局部匹配事件的匹配质量，len是在后局部匹配事件的长度；
P(e)＝α*P_skip(e)+β*(P_ext(e)+P_sep(e))；
α和β为权重系数，可根据所分析的物种设定不同的值；
跳过中间节点的边的减分P_skip(e)＝min(PIDi*leni),其中i为所有被跳过的局部对比序列的序号，
边e有线性连接和非线性连接两种情况，若两个局部匹配事件中测序片段比对至参考序列的位置相邻则属于线性连接，不相邻则属于非线性连接，
P_ext(e)表示边e是线性连接时，测序片段和参考序列无法比对的区域，P_ext(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
P_sep(e)表示边e是非线性连接时，测序片段和参考序列无法比对的区域，P_sep(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
不同的排序方式组合成的测序序列减分为Length(p)＝Σe∈pLength(e)，一种排序方式测序序列的减分越小，则得分越高。
进一步的，步骤S3的具体实现为：
对Length(p)进行判断，选取Length(p)值最小情况下的边e对应的连接关系，确定为最优排序方式。
采用上述进一步方案的有益效果是：本发明可以有效建立排序模型，得到最佳排序方式。
一种基于超长基因组的变异检测系统，包括以下模块：
局部序列比对模块，用于检测出测序片段和参考序列间存在的所有匹配，得到局部匹配事件，每一个局部匹配事件均包括参考序列上的参考片段和测序片段；
局部匹配事件分组模块，用于将所有局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，将测序片段比对至参考序列上的位置重叠或顺次相连的局部匹配事件分到一个组；
局部匹配事件排序模块，用于对每一个组内的局部匹配事件按照不同的排序方式进行评分，建立一个用来表示每一个组内的所有局部匹配事件之间不同的排序方式的图形模型，根据评分情况确定每一组内的局部匹配事件的最佳排序方式；
局部匹配事件串联模块，用于按照最佳排序方式将每一个组内的局部匹配事件中的测序片段进行串联，得到一个测序序列read；
gap添加模块，用于将每一组的测序序列read和参考序列进行对比，在测序序列read上和参考序列不匹配的区域添加gap区域，在参考序列上和测序序列read不匹配的区域添加gap区域，使测序序列read和参考序列的匹配区域的位置完全对应；
检测模块，用于对添加了gap区域后的最终测序序列进行结构变异检测。
进一步的，所述局部匹配事件分组排序模块包括以下单元：
位置排序单元，用于将每个组内的局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，每一个局部匹配事件按照位置从前到后依次排列；
图形建模单元，用于按照所有的局部匹配事件均可以与位置排在其后的局部匹配事件连接的条件，按照所有的可能性建立排序方式，构建有向非循环图形模型；
评分单元，用于对不同的排序关系进行评分，设图形模型中，任一个局部匹配事件和位置排在其后的局部匹配事件之间的关系为边e，则每个边e的减分Length(e)＝P(e)–R(e)，一个边e的减分越少，则评分越高；
其中，P(e)代表边e的减分、R(e)代表边e的得分；
R(e)＝PID*len，PID位置排在后的局部匹配事件的匹配质量，len是位置排在后的局部匹配事件的长度；
P(e)＝α*P_skip(e)+β*(P_ext(e)+P_sep(e))；
α和β为权重系数，可根据所分析的物种设定不同的值；
跳过中间节点的边的减分P_skip(e)＝min(PIDi*leni),其中i为所有被跳过的局部对比序列的序号，
边e有线性连接和非线性连接两种情况，若两个局部匹配事件中测序片段比对至参考序列的位置相邻则属于线性连接，不相邻则属于非线性连接，
P_ext(e)表示边e是线性连接时，测序片段和参考序列无法比对的区域，P_ext(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
P_sep(e)表示边e是非线性连接时，测序片段和参考序列无法比对的区域，P_sep(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
不同的排序方式组合成的测序序列减分为Length(p)＝Σe∈pLength(e)，一种排序方式测序序列的减分越小，则得分越高；
决定单元，用于选取Length(p)最小的一种排序方式，确定为最优的排序方式。
附图说明
图1为本发明方法流程图；
图2为局部对比图像化结果图；
图3为步骤S3.2的图形模型示意图；
图4为步骤S5测序序列read和参考序列进行对比时示意图；
图5为本发明系统示意图。
具体实施方式
以下结合附图对本发明的原理和特征进行描述，所举实例只用于解释本发明，并非用于限定本发明的范围。
如图1所示，一种基于超长基因组的变异检测算法，包括以下步骤：
S1、采用局部序列比对算法，检测出测序片段和参考序列间存在的所有匹配，得到局部匹配事件，每一个局部匹配事件均包括测序片段和参考序列上的参考片段；
S2、将所有局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，将测序片段比对至参考序列上的位置重叠或顺次相连的局部匹配事件分到一个组；
S3、对每一个组内的局部匹配事件按照不同的排序方式进行评分，建立一个用来表示每一个组内的所有局部匹配事件之间不同的排序方式的图形模型，根据评分情况确定每一组内的局部匹配事件的最佳排序方式；
S4、按照最佳排序方式将每一个组内的局部匹配事件中的测序片段进行串联，一个组得到一个测序序列read；
S5、如图4所示，将每一组的测序序列read和参考序列进行对比，在测序序列read上和参考序列不匹配的区域添加gap区域在参考序列上和测序序列read不匹配的区域添加gap区域，gap区域即空白区域，用于使测序序列read和参考序列的匹配区域位置完全对应，便于确认变异是属于插入还是缺失的类型；
S6、对添加了gap区域后的最终测序序列进行结构变异检测。
所述步骤S1里采用的局部序列比对算法，为BLASTn、MegaBlast、BLASTz、cross_match、BLAT和wublast算法中的一种。
所述步骤S2中需要同时考虑了测序片段两条链上的情况，来建立局部匹配事件之间的关系模型。
所述步骤S4中的具体实现为：
S4.1、根据打分情况确定每一组内的局部匹配事件的最佳排序方式后，
将每一组内的局部匹配事件按照最佳排序方式进行排序，并将所有局部匹配事件中的测序片段串联成一个基础测序序列；
S4.2、对基础测序序列中的相邻的测序片段的重叠部分进行剔除，得到测序序列read。
所述步骤S4.2中，对基础测序序列中的相邻的测序片段的重叠部分进行剔除的方法是，若相邻两个测序片段之间存在重叠部分，分别将两个测序片段的重叠部分和对应的参考片段进行对比，选取匹配度更高的一个重叠部分进行保留，并把另一个重叠部分剔除。
所述步骤S4.2中对基础测序序列中的相邻的测序片段的重叠部分进行剔除后，还包括将基础测序序列和参考序列进行对比，若存在缺失片段，检测剩余未归入组内的测序片段与缺失片段的匹配度，选取匹配度最高的测序片段放入相应的基础测序序列内，再得到每一组的测序序列read的步骤。
步骤S3中，将符合线性连接关系的局部匹配事件分到一个组时,每一个组内的局部匹配事件的重叠度不大于阈值t。
如图2所示，为局部对比图像化结果展示，每一项局部比对(或者高分片段对，HSP)包含了一对序列片段[Q,T]，Q代表测序片段序列，T代表参考序列上匹配上的序列，H代表局部比对事件，H1包括[T1,Q1],H2包括[T2,Q1],H3包括[T3,Q2],H4包括[T4,Q3],H5包括[T5,Q4],H6包括[T6,Q1],所述S3中对所有局部匹配事件进行排序，建立一个如图3所示的图形模型表示所有局部匹配事件之间不同的排序方式，边e有线性连接和非线性连接两种情况，若两个局部匹配事件比对至参考序列的位置相邻则属于线性连接，不相邻则属于非线性连接，例如，H1和H2中比对至参考序列的位置均为Q1，故H1和H2为非线性连接，H2和H3中比对至参考序列的位置分别为Q1和Q2，故H2和H3为线性连接，如图3所示，其中双竖线表示非线性连接，并对不同的排序形式进行打分，包括以下步骤：
S3.1、将每个组内的局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，所有的局部匹配事件按照位置从前到后依次排列；
S3.2、按照所有的局部匹配事件均可以与位置排在其后的局部匹配事件连接的条件，按照所有的可能性建立排序方式，构建有向非循环图形模型；
S3.3、设图形模型中，在前局部匹配事件和在后局部匹配事件之间的关系为边e，则每个边e的减分Length(e)＝P(e)–R(e)，一个边e的减分越少，则评分越高；
其中，P(e)代表边e的减分、R(e)代表边e的得分；
R(e)＝PID*len，PID在后的局部匹配事件的匹配质量，len是在后局部匹配事件的长度；
P(e)＝α*P_skip(e)+β*(P_ext(e)+P_sep(e))；
α和β为权重系数，可根据所分析的物种设定不同的值；
跳过中间节点的边的减分P_skip(e)＝min(PIDi*leni),其中i为所有被跳过的局部对比序列的序号，
边e有线性连接和非线性连接两种情况，若两个局部匹配事件中测序片段比对至参考序列的位置相邻则属于线性连接，不相邻则属于非线性连接，
P_ext(e)表示边e是线性连接时，测序片段和参考序列无法比对的区域，P_ext(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
P_sep(e)表示边e是非线性连接时，测序片段和参考序列无法比对的区域，P_sep(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
不同的排序方式组合成的测序序列减分为Length(p)＝Σe∈pLength(e)。
对Length(p)进行判断，选取Length(p)值最小情况下的边e对应的连接关系，确定为最优排序方式。
一种基于超长基因组的变异检测系统，包括以下模块：
局部序列比对模块，用于检测出测序片段和参考序列间存在的所有匹配，得到局部匹配事件，每一个局部匹配事件均包括参考序列上的参考片段和测序片段；
局部匹配事件分组模块，用于将所有局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，将测序片段比对至参考序列上的位置重叠或顺次相连的局部匹配事件分到一个组；
局部匹配事件排序模块，用于对每一个组内的局部匹配事件按照不同的排序方式进行评分，建立一个用来表示每一个组内的所有局部匹配事件之间不同的排序方式的图形模型，根据评分情况确定每一组内的局部匹配事件的最佳排序方式；
局部匹配事件串联模块，用于按照最佳排序方式将每一个组内的局部匹配事件中的测序片段进行串联，得到一个测序序列read；
gap添加模块，用于将每一组的测序序列read和参考序列进行对比，在测序序列read上和参考序列不匹配的区域添加gap区域，在参考序列上和测序序列read不匹配的区域添加gap区域，使测序序列read和参考序列的匹配区域的位置完全对应；
检测模块，用于对添加了gap区域后的最终测序序列进行结构变异检测。
所述局部匹配事件分组排序模块包括以下单元：
位置排序单元，用于将每个组内的局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，每一个局部匹配事件按照位置从前到后依次排列；
图形建模单元，用于按照所有的局部匹配事件均可以与位置排在其后的局部匹配事件连接的条件，按照所有的可能性建立排序方式，构建有向非循环图形模型；
评分单元，用于对不同的排序关系进行评分，设图形模型中，任一个局部匹配事件和位置排在其后的局部匹配事件之间的关系为边e，则每个边e的减分Length(e)＝P(e)–R(e)，一个边e的减分越少，则评分越高；
其中，P(e)代表边e的减分、R(e)代表边e的得分；
R(e)＝PID*len，PID位置排在后的局部匹配事件的匹配质量，len是位置排在后的局部匹配事件的长度；
P(e)＝α*P_skip(e)+β*(P_ext(e)+P_sep(e))；
α和β为权重系数，可根据所分析的物种设定不同的值；
跳过中间节点的边的减分P_skip(e)＝min(PIDi*leni),其中i为所有被跳过的局部对比序列的序号，
边e有线性连接和非线性连接两种情况，若两个局部匹配事件中测序片段比对至参考序列的位置相邻则属于线性连接，不相邻则属于非线性连接，
P_ext(e)表示边e是线性连接时，测序片段和参考序列无法比对的区域，P_ext(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
P_sep(e)表示边e是非线性连接时，测序片段和参考序列无法比对的区域，P_sep(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
不同的排序方式组合成的测序序列减分为Length(p)＝Σe∈pLength(e)，一种排序方式测序序列的减分越小，则得分越高；
决定单元，用于选取Length(p)最小的一种排序方式，确定为最优的排序方式。
本发明的VariationBlast算法的突出优势在于是它能够充分利用由不断提高(例如测序长度不断增加)的第二代以及第三代DNA测序技术所带来的可能性，利用测序结果检测各种尺度的结构变异，包括其它软件所不能够预测的大尺度插入以及所有现有软件都忽略的复合变异。VariationBlast将大大提高基因组变异和结构变异预测的灵敏度、效率和准确度，可以有效地应用于来自不同物种、不同复杂度的基因组序列，从而为各个领域的科学研究和应用提供强急需的强有力的支撑。
应用范围上，它是针对二代测序组装数据及三代测序数据开发的算法，所以对二代测序组装序列及三代测序序列特征进行了针对性优化。也就是说，它几乎能够适合目前所有完成测序或部分测序的基因组序列，包括来自真核生物与原核生物的各种复杂程度的基因组序列。基於技术限制，二代测序技术得到的reads序列较短，且对高复杂度的序列组装效果差，这也就对基因组结构变异检测软件提出了更高的要求。三代测序技术，是个新技术，有很多优点，具有很广泛的应用前景，而目前针对三代测序序列基因组结构变异检测的软件非常少，所以VariationBlast将来一定会被应用到更广泛的领域。
预测结果方面，VariationBlast可以大幅度提高基因组结构变异预测结果的准确性，可靠性，并提高输出结果可视性。预测结果的准确性，可靠性将直接影响后期实验验证，如果假阳性率太高，验证工作就会增加很多无用功，加大工作量，导致实验成本剧增。如果假阴性率太高，则可能会漏掉关键的基因组结构变异，错失解决问题的机会。更为重要的是，VariationBlast能够寻找到其他类似软件找不到的大尺度插入以及复合结构变异。
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
Claims (10)
1.一种基于超长基因组的变异检测算法，其特征在于，包括以下步骤：
S1、采用局部序列比对算法，检测出测序片段和参考序列间存在的所有匹配，得到局部匹配事件，每一个局部匹配事件均包括测序片段和参考序列上的参考片段；
S2、将所有局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，将测序片段比对至参考序列上的位置重叠或顺次相连的局部匹配事件分到一个组；
S3、对每一个组内的局部匹配事件按照不同的排序方式进行评分，建立一个用来表示每一个组内的所有局部匹配事件之间不同的排序方式的图形模型，根据评分情况确定每一组内的局部匹配事件的最佳排序方式；
S4、按照最佳排序方式将每一个组内的局部匹配事件中的测序片段进行串联，一个组得到一个测序序列read；
S5、将每一组的测序序列read和参考序列进行对比，在测序序列read上和参考序列不匹配的区域添加gap区域，在参考序列上和测序序列read不匹配的区域添加gap区域，用于使测序序列read和参考序列的匹配区域的位置完全对应；
S6、对添加了gap区域后的最终测序序列进行结构变异检测。
2.根据权利要求1所述的基于超长基因组的变异检测算法，其特征在于，所述步骤S1里采用的局部序列比对算法为BLASTn、MegaBlast、BLASTz、cross_match、BLAT和wublast算法中的一种。
3.根据权利要求1所述的基于超长基因组的变异检测算法，其特征在于，所述步骤S4中的具体实现为：
S4.1、将每一组内的局部匹配事件按照最佳排序方式进行排序，并将所有局部匹配事件中的测序片段串联成一个基础测序序列；
S4.2、对基础测序序列中的相邻的测序片段的重叠部分进行剔除，得到测序序列read。
4.根据权利要求3所述的基于超长基因组的变异检测算法，其特征在于，所述步骤S4.2中，对基础测序序列中的相邻的测序片段的重叠部分进行剔除的方法是，若相邻两个测序片段之间存在重叠部分，分别将两个测序片段的重叠部分和对应的参考片段进行对比，选取匹配度更高的一个重叠部分进行保留，并把另一个重叠部分剔除。
5.根据权利要求3所述的基于超长基因组的变异检测算法，其特征在于，步骤S4.2中对基础测序序列中的相邻的测序片段的重叠部分进行剔除后，还包括将基础测序序列和参考序列进行对比，若存在缺失片段，检测剩余未归入组内的测序片段与缺失片段的匹配度，选取匹配度最高的测序片段放入相应的基础测序序列内，再得到每一组的测序序列read的步骤。
6.根据权利要求1所述的基于超长基因组的变异检测算法，其特征在于，步骤S3中，将符合线性连接关系的局部匹配事件分到一个组时,每一个组内的局部匹配事件的重叠度不大于阈值t。
7.根据权利要求1所述的基于超长基因组的变异检测算法，其特征在于，所述步骤S3的具体实现为：
S3.1、将每个组内的局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，所有的局部匹配事件按照位置从前到后依次排列；
S3.2、按照所有的局部匹配事件均可以与位置排在其后的局部匹配事件连接的条件，按照所有的可能性建立排序方式，构建有向非循环图形模型；
S3.3、设图形模型中，在前局部匹配事件和在后局部匹配事件之间的关系为边e，则每个边e的减分Length(e)＝P(e)–R(e)，一个边e的减分越少，则评分越高；
其中，P(e)代表边e的减分、R(e)代表边e的得分；
R(e)＝PID*len，PID在后的局部匹配事件的匹配质量，len是在后局部匹配事件的长度；
P(e)＝α*P_skip(e)+β*(P_ext(e)+P_sep(e))；
α和β为权重系数，可根据所分析的物种设定不同的值；
跳过中间节点的边的减分P_skip(e)＝min(PIDi*leni),其中i为所有被跳过的局部对比序列的序号，
边e有线性连接和非线性连接两种情况，若两个局部匹配事件中测序片段比对至参考序列的位置相邻则属于线性连接，不相邻则属于非线性连接，
P_ext(e)表示边e是线性连接时，测序片段和参考序列无法比对的区域，P_ext(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
P_sep(e)表示边e是非线性连接时，测序片段和参考序列无法比对的区域，P_sep(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
不同的排序方式组合成的测序序列减分为Length(p)＝Σe∈pLength(e)，一种排序方式测序序列的减分越小，则评分越高。
8.根据权利要求7所述的基于超长基因组的变异检测算法，其特征在于，步骤S3的具体实现为：
对Length(p)进行判断，选取Length(p)值最小情况下的边e对应的连接关系，确定为最优排序方式。
9.一种基于超长基因组的变异检测系统，其特征在于，包括以下模块：
局部序列比对模块，用于检测出测序片段和参考序列间存在的所有匹配，得到局部匹配事件，每一个局部匹配事件均包括参考序列上的参考片段和测序片段；
局部匹配事件分组模块，用于将所有局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，将测序片段比对至参考序列上的位置重叠或顺次相连的局部匹配事件分到一个组；
局部匹配事件排序模块，用于对每一个组内的局部匹配事件按照不同的排序方式进行评分，建立一个用来表示每一个组内的所有局部匹配事件之间不同的排序方式的图形模型，根据评分情况确定每一组内的局部匹配事件的最佳排序方式；
局部匹配事件串联模块，用于按照最佳排序方式将每一个组内的局部匹配事件中的测序片段进行串联，得到一个测序序列read；
gap添加模块，用于将每一组的测序序列read和参考序列进行对比，在测序序列read上和参考序列不匹配的区域添加gap区域，在参考序列上和测序序列read不匹配的区域添加gap区域，使测序序列read和参考序列的匹配区域的位置完全对应；
检测模块，用于对添加了gap区域后的最终测序序列进行结构变异检测。
10.根据权利要求9所述的基于超长基因组的变异检测系统，其特征在于，所述局部匹配事件分组排序模块包括以下单元：
位置排序单元，用于将每个组内的局部匹配事件中的测序片段按照比对至参考序列的位置进行排序，每一个局部匹配事件按照位置从前到后依次排列；
图形建模单元，用于按照所有的局部匹配事件均可以与位置排在其后的局部匹配事件连接的条件，按照所有的可能性建立排序方式，构建有向非循环图形模型；
评分单元，用于对不同的排序关系进行评分，设图形模型中，任一个局部匹配事件和位置排在其后的局部匹配事件之间的关系为边e，则每个边e的减分Length(e)＝P(e)–R(e)，一个边e的减分越少，则评分越高；
其中，P(e)代表边e的减分、R(e)代表边e的得分；
R(e)＝PID*len，PID位置排在后的局部匹配事件的匹配质量，len是位置排在后的局部匹配事件的长度；
P(e)＝α*P_skip(e)+β*(P_ext(e)+P_sep(e))；
α和β为权重系数，可根据所分析的物种设定不同的值；
跳过中间节点的边的减分P_skip(e)＝min(PIDi*leni),其中i为所有被跳过的局部对比序列的序号，
边e有线性连接和非线性连接两种情况，若两个局部匹配事件中测序片段比对至参考序列的位置相邻则属于线性连接，不相邻则属于非线性连接，
P_ext(e)表示边e是线性连接时，测序片段和参考序列无法比对的区域，P_ext(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
P_sep(e)表示边e是非线性连接时，测序片段和参考序列无法比对的区域，P_sep(e)＝PIDi→j*leni→j，i和j代表没有比上参考序列的起始和终止位点；
不同的排序方式组合成的测序序列减分为Length(p)＝Σe∈pLength(e)，一种排序方式测序序列的减分越小，则得分越高；
决定单元，用于选取Length(p)最小的一种排序方式，确定为最优的排序方式。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201511005316.0A CN105483244B (zh) | 2015-12-28 | 2015-12-28 | 一种基于超长基因组的变异检测方法及检测系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201511005316.0A CN105483244B (zh) | 2015-12-28 | 2015-12-28 | 一种基于超长基因组的变异检测方法及检测系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN105483244A true CN105483244A (zh) | 2016-04-13 |
| CN105483244B CN105483244B (zh) | 2019-10-22 |
Family
ID=55670515
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201511005316.0A Active CN105483244B (zh) | 2015-12-28 | 2015-12-28 | 一种基于超长基因组的变异检测方法及检测系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN105483244B (zh) |
Cited By (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106355046A (zh) * | 2016-09-18 | 2017-01-25 | 北京百度网讯科技有限公司 | 结构体变异检测的方法及装置 |
| CN106650254A (zh) * | 2016-12-16 | 2017-05-10 | 武汉菲沙基因信息有限公司 | 一种基于转录组测序数据检测融合基因的方法 |
| CN107180166A (zh) * | 2017-04-21 | 2017-09-19 | 北京希望组生物科技有限公司 | 一种基于三代测序的全基因组结构变异分析方法和系统 |
| CN107273663A (zh) * | 2017-05-22 | 2017-10-20 | 人和未来生物科技（长沙）有限公司 | 一种dna甲基化测序数据计算解读方法 |
| CN107403076A (zh) * | 2016-05-18 | 2017-11-28 | 华为技术有限公司 | Dna序列的处理方法及设备 |
| CN108959853A (zh) * | 2018-05-18 | 2018-12-07 | 广州金域医学检验中心有限公司 | 一种拷贝数变异的分析方法、分析装置、设备及存储介质 |
| CN109266729A (zh) * | 2018-09-29 | 2019-01-25 | 中国科学院遗传与发育生物学研究所 | 一种基于基因组二代测序的大片段缺失的检测方法 |
| CN111312332A (zh) * | 2020-02-13 | 2020-06-19 | 国家卫生健康委科学技术研究所 | 基于hla基因的生物信息处理方法、装置及终端 |
| CN112669902A (zh) * | 2021-03-16 | 2021-04-16 | 北京贝瑞和康生物技术有限公司 | 检测基因组结构变异的方法、计算设备和存储介质 |
| CN113963746A (zh) * | 2021-09-29 | 2022-01-21 | 西安交通大学 | 一种模型非依赖的基因组结构变异检测系统及方法 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102460155A (zh) * | 2009-04-29 | 2012-05-16 | 考利达基因组股份有限公司 | 用于关于参考多核苷酸序列标注样本多核苷酸序列中的变异的方法和系统 |
| CN103080333A (zh) * | 2010-09-14 | 2013-05-01 | 深圳华大基因科技有限公司 | 一种基因组结构性变异检测方法和系统 |
| CN104017883A (zh) * | 2014-06-18 | 2014-09-03 | 深圳华大基因科技服务有限公司 | 组装基因组序列的方法和系统 |
| CN104160391A (zh) * | 2011-09-16 | 2014-11-19 | 考利达基因组股份有限公司 | 确定异质样本的基因组中的变异 |
| CN104603283A (zh) * | 2012-08-01 | 2015-05-06 | 深圳华大基因研究院 | 确定异常状态相关生物标志物的方法及系统 |
-
2015
- 2015-12-28 CN CN201511005316.0A patent/CN105483244B/zh active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102460155A (zh) * | 2009-04-29 | 2012-05-16 | 考利达基因组股份有限公司 | 用于关于参考多核苷酸序列标注样本多核苷酸序列中的变异的方法和系统 |
| CN103080333A (zh) * | 2010-09-14 | 2013-05-01 | 深圳华大基因科技有限公司 | 一种基因组结构性变异检测方法和系统 |
| CN104160391A (zh) * | 2011-09-16 | 2014-11-19 | 考利达基因组股份有限公司 | 确定异质样本的基因组中的变异 |
| CN104603283A (zh) * | 2012-08-01 | 2015-05-06 | 深圳华大基因研究院 | 确定异常状态相关生物标志物的方法及系统 |
| CN104017883A (zh) * | 2014-06-18 | 2014-09-03 | 深圳华大基因科技服务有限公司 | 组装基因组序列的方法和系统 |
Cited By (19)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107403076B (zh) * | 2016-05-18 | 2020-06-02 | 华为技术有限公司 | Dna序列的处理方法及设备 |
| CN107403076A (zh) * | 2016-05-18 | 2017-11-28 | 华为技术有限公司 | Dna序列的处理方法及设备 |
| CN106355046A (zh) * | 2016-09-18 | 2017-01-25 | 北京百度网讯科技有限公司 | 结构体变异检测的方法及装置 |
| CN106355046B (zh) * | 2016-09-18 | 2019-04-30 | 北京百度网讯科技有限公司 | 结构体变异检测的方法及装置 |
| CN106650254A (zh) * | 2016-12-16 | 2017-05-10 | 武汉菲沙基因信息有限公司 | 一种基于转录组测序数据检测融合基因的方法 |
| CN106650254B (zh) * | 2016-12-16 | 2018-11-20 | 武汉菲沙基因信息有限公司 | 一种基于转录组测序数据检测融合基因的方法 |
| CN107180166A (zh) * | 2017-04-21 | 2017-09-19 | 北京希望组生物科技有限公司 | 一种基于三代测序的全基因组结构变异分析方法和系统 |
| CN107273663A (zh) * | 2017-05-22 | 2017-10-20 | 人和未来生物科技（长沙）有限公司 | 一种dna甲基化测序数据计算解读方法 |
| CN107273663B (zh) * | 2017-05-22 | 2018-12-11 | 人和未来生物科技（长沙）有限公司 | 一种dna甲基化测序数据计算解读方法 |
| CN108959853A (zh) * | 2018-05-18 | 2018-12-07 | 广州金域医学检验中心有限公司 | 一种拷贝数变异的分析方法、分析装置、设备及存储介质 |
| CN108959853B (zh) * | 2018-05-18 | 2020-01-17 | 广州金域医学检验中心有限公司 | 一种拷贝数变异的分析方法、分析装置、设备及存储介质 |
| CN109266729A (zh) * | 2018-09-29 | 2019-01-25 | 中国科学院遗传与发育生物学研究所 | 一种基于基因组二代测序的大片段缺失的检测方法 |
| CN109266729B (zh) * | 2018-09-29 | 2020-11-27 | 中国科学院遗传与发育生物学研究所 | 一种基于基因组二代测序的大片段缺失的检测方法 |
| CN111312332A (zh) * | 2020-02-13 | 2020-06-19 | 国家卫生健康委科学技术研究所 | 基于hla基因的生物信息处理方法、装置及终端 |
| CN111312332B (zh) * | 2020-02-13 | 2020-10-30 | 国家卫生健康委科学技术研究所 | 基于hla基因的生物信息处理方法、装置及终端 |
| CN112669902A (zh) * | 2021-03-16 | 2021-04-16 | 北京贝瑞和康生物技术有限公司 | 检测基因组结构变异的方法、计算设备和存储介质 |
| CN112669902B (zh) * | 2021-03-16 | 2021-06-04 | 北京贝瑞和康生物技术有限公司 | 检测基因组结构变异的方法、计算设备和存储介质 |
| CN113963746A (zh) * | 2021-09-29 | 2022-01-21 | 西安交通大学 | 一种模型非依赖的基因组结构变异检测系统及方法 |
| CN113963746B (zh) * | 2021-09-29 | 2023-09-19 | 西安交通大学 | 一种基于深度学习的基因组结构变异检测系统及方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN105483244B (zh) | 2019-10-22 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN105483244A (zh) | 一种基于超长基因组的变异检测算法及检测系统 | |
| CN105886616B (zh) | 一种用于猪基因编辑的高效特异性sgRNA识别位点引导序列及其筛选方法 | |
| CN104204221B (zh) | 一种检验融合基因的方法及系统 | |
| François et al. | Spatially explicit Bayesian clustering models in population genetics | |
| CN102413029B (zh) | 基于分解的局部搜索多目标复杂动态网络社区划分方法 | |
| CN106033502A (zh) | 鉴定病毒的方法和装置 | |
| CN102206704A (zh) | 组装基因组序列的方法和装置 | |
| CN104965999B (zh) | 一种中短基因片段测序的分析拼接方法及设备 | |
| Zhang et al. | iPromoter-5mC: a novel fusion decision predictor for the identification of 5-methylcytosine sites in genome-wide DNA promoters | |
| CN102521528A (zh) | 一种基因序列数据的筛选方法 | |
| Guan et al. | Efficient iterative Hi-C scaffolder based on N-best neighbors | |
| Md Mukarram Hossain et al. | Evidence of statistical inconsistency of phylogenetic methods in the presence of multiple sequence alignment uncertainty | |
| Rachtman et al. | CONSULT: accurate contamination removal using locality-sensitive hashing | |
| Pan et al. | Novo&Stitch: accurate reconciliation of genome assemblies via optical maps | |
| Roux de Bézieux et al. | CALDERA: finding all significant de Bruijn subgraphs for bacterial GWAS | |
| Van Etten et al. | A k-mer-based approach for phylogenetic classification of taxa in environmental genomic data | |
| CN103270175B (zh) | 检测转基因外源片段插入位点的方法和系统 | |
| CN107967411A (zh) | 一种脱靶位点的检测方法、装置及终端设备 | |
| CN111767546B (zh) | 一种基于深度学习的输入结构推断方法和装置 | |
| You et al. | LTR annotator: automated identification and annotation of LTR retrotransposons in plant genomes | |
| CN102841988B (zh) | 一种对核酸序列信息进行匹配的系统和方法 | |
| CN106021998A (zh) | 单通多变体识别计算流水线 | |
| Petri et al. | isONform: reference-free transcriptome reconstruction from Oxford Nanopore data | |
| EP3663890B1 (en) | Alignment method, device and system | |
| CN106021978A (zh) | 基于光学图谱平台Irys的一种de novo测序数据组装方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| PE01 | Entry into force of the registration of the contract for pledge of patent right | ||
| PE01 | Entry into force of the registration of the contract for pledge of patent right |
Denomination of invention: A variation detection method and system based on ultra long genome Effective date of registration: 20210918 Granted publication date: 20191022 Pledgee: Wuhan area branch of Hubei pilot free trade zone of Bank of China Ltd. Pledgor: WUHAN FRASERGEN INFORMATION Co.,Ltd. Registration number: Y2021420000096 |
|
| PC01 | Cancellation of the registration of the contract for pledge of patent right | ||
| PC01 | Cancellation of the registration of the contract for pledge of patent right |
Granted publication date: 20191022 Pledgee: Wuhan area branch of Hubei pilot free trade zone of Bank of China Ltd. Pledgor: WUHAN FRASERGEN INFORMATION CO.,LTD. Registration number: Y2021420000096 |
|
| PE01 | Entry into force of the registration of the contract for pledge of patent right | ||
| PE01 | Entry into force of the registration of the contract for pledge of patent right |
Denomination of invention: A mutation detection method and detection system based on ultra long genome Granted publication date: 20191022 Pledgee: Guanggu Branch of Wuhan Rural Commercial Bank Co.,Ltd. Pledgor: WUHAN FRASERGEN INFORMATION CO.,LTD. Registration number: Y2024980021037 |