CN113628683A - 一种高通量测序突变检测方法、设备、装置及可读存储介质 - Google Patents
一种高通量测序突变检测方法、设备、装置及可读存储介质 Download PDFInfo
- Publication number
- CN113628683A CN113628683A CN202110972081.1A CN202110972081A CN113628683A CN 113628683 A CN113628683 A CN 113628683A CN 202110972081 A CN202110972081 A CN 202110972081A CN 113628683 A CN113628683 A CN 113628683A
- Authority
- CN
- China
- Prior art keywords
- candidate
- candidate mutation
- mutation
- mutation site
- site
- 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
- 230000035772 mutation Effects 0.000 title claims abstract description 286
- 238000001514 detection method Methods 0.000 title claims abstract description 106
- 238000012165 high-throughput sequencing Methods 0.000 title claims abstract description 53
- 238000003860 storage Methods 0.000 title claims abstract description 8
- 108090000623 proteins and genes Proteins 0.000 claims abstract description 39
- 238000007417 hierarchical cluster analysis Methods 0.000 claims abstract description 31
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 12
- 238000000034 method Methods 0.000 claims description 19
- 239000012634 fragment Substances 0.000 claims description 16
- 238000004364 calculation method Methods 0.000 claims description 12
- 238000004590 computer program Methods 0.000 claims description 9
- 108091092878 Microsatellite Proteins 0.000 claims description 6
- 238000013215 result calculation Methods 0.000 claims description 3
- 238000004458 analytical method Methods 0.000 abstract description 7
- 238000005457 optimization Methods 0.000 abstract description 2
- 238000012163 sequencing technique Methods 0.000 description 7
- 206010028980 Neoplasm Diseases 0.000 description 5
- 239000011159 matrix material Substances 0.000 description 5
- 238000002474 experimental method Methods 0.000 description 4
- 230000008569 process Effects 0.000 description 4
- 102100021569 Apoptosis regulator Bcl-2 Human genes 0.000 description 3
- 108091012583 BCL2 Proteins 0.000 description 3
- 108020004414 DNA Proteins 0.000 description 3
- 230000003321 amplification Effects 0.000 description 3
- 238000003199 nucleic acid amplification method Methods 0.000 description 3
- 108090000790 Enzymes Proteins 0.000 description 2
- 102000004190 Enzymes Human genes 0.000 description 2
- 238000012408 PCR amplification Methods 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 230000006870 function Effects 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 238000012216 screening Methods 0.000 description 2
- 238000007400 DNA extraction Methods 0.000 description 1
- 238000003705 background correction Methods 0.000 description 1
- 239000003990 capacitor Substances 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 238000007405 data analysis Methods 0.000 description 1
- 238000012217 deletion Methods 0.000 description 1
- 230000037430 deletion Effects 0.000 description 1
- 230000004069 differentiation Effects 0.000 description 1
- 238000013467 fragmentation Methods 0.000 description 1
- 238000006062 fragmentation reaction Methods 0.000 description 1
- 238000003780 insertion Methods 0.000 description 1
- 230000037431 insertion Effects 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000007481 next generation sequencing Methods 0.000 description 1
- 108020004707 nucleic acids Proteins 0.000 description 1
- 150000007523 nucleic acids Chemical class 0.000 description 1
- 102000039446 nucleic acids Human genes 0.000 description 1
- 238000002360 preparation method Methods 0.000 description 1
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
- G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
- G16B20/50—Mutagenesis
-
- 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
- G16B40/30—Unsupervised data analysis
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Engineering & Computer Science (AREA)
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Physics & Mathematics (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Medical Informatics (AREA)
- General Health & Medical Sciences (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Biophysics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Biotechnology (AREA)
- Evolutionary Biology (AREA)
- Bioethics (AREA)
- Analytical Chemistry (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Molecular Biology (AREA)
- Genetics & Genomics (AREA)
- Artificial Intelligence (AREA)
- Chemical & Material Sciences (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Evolutionary Computation (AREA)
- Epidemiology (AREA)
- Databases & Information Systems (AREA)
- Public Health (AREA)
- Software Systems (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明涉及生物技术领域，特别是涉及一种高通量测序突变检测方法、设备、装置及可读存储介质。本发明提供一种高通量测序突变检测方法，包括：获取基因样本中各候选突变位点的特征信息，所述特征信息来源于基因样本的高通量测序数据；将各候选突变位点进行层次聚类分析；获取各类候选突变位点的背景值pbg；获取各类候选突变位点的最低检测下限；确定各候选突变位点的检测阈值。本申请所提供的高通量测序突变检测方法，可以通过自适应优化的算法，针对每个样本特有的数据特点智能的设定分析的阈值，在保证阳性位点检出的前提下尽可能的去除噪音造成的假阳性，可以更加精准的从样本中检测突变，从而具有良好的产业化前景。
Description
技术领域
本发明涉及生物技术领域，特别是涉及一种高通量测序突变检测方法、设备、装置及可读存储介质。
背景技术
通过高通量测序进行肿瘤的突变检测被广泛的应用于肿瘤的基础与临床研究中。但由于在样本制备、保存、实验和分析中会引入大量的非自然来源的干扰噪音，因此，突变检测的一个关键步骤是准确区分真实的突变和不同来源的噪音信号。
目前去除噪音的实验和数据分析手段主要有以下几类：
1、针对测序过程中随机产生的噪音信号，这种类型的噪音随机低频次出现，但是可以利用测序中产生的重复，结合分子标签(CN106834275A)，虚拟分子标签(CN107944225B)等技术方式结合高深度测序来校正；
2、针对在实验过程中产生的可重复的非随机噪音，例如DNA提取、打断、捕获过程中引入的噪音，因为这种类型的噪音出现频次较高，符合一定的统计规律，因此可以利用大量阴性样本建立的背景库来建立背景校正模型来进行校正和区分(CN105574365B)；
3、针对一些特定的突变类型来设定不同的分析阈值，不同类型的变异噪音背景值不同，分类进行阈值设定，例如分别为点突变和插入缺失突变设定不同的检测阈值，可以提高分析的准确性(CN108690871A)。
上文中提到的几种技术方案，分别解决了低频次随机噪音，高频次固有噪音和不同类型突变的固有噪音的问题，但是，实际中还存在的另一种噪音类型，非随机低频次的样本特异性噪音，有多种因素可能导致这种类型的噪音出现，例如，1、样本自身DNA的损伤，在常用的福尔马林固定的肿瘤样本中常有出现，此类样本会经常出现(C>T|G>A)类型的噪音变异；2、样本片段化长度不一致(过长或过短的片段)导致的噪音，因为每个样本的实验条件不同，导致无法通过理想的背景噪音模型来复现和消除；3、样本PCR扩增轮数、扩增酶保真率不同导致的PCR错误，每个样本因为初始量不同，PCR扩增轮数不同，而每次扩增引入错误的比例又与扩增酶实验所处状态相关，样本与样本间会存在差异。这些因素的共同特点都是样本本身特异，无法在样本间有效复现，但是样本内部非随机，在同一次实验中重复出现，因此无法通过前一部分提到的技术方案来去除，影响检测结果的准确性。
发明内容
鉴于以上所述现有技术的缺点，本发明的目的在于提供一种高通量测序突变检测方法，用于解决现有技术中的问题。
为实现上述目的及其他相关目的，本发明一方面提供一种高通量测序突变检测方法，包括：
S1)获取基因样本中各候选突变位点的特征信息，所述特征信息来源于基因样本的高通量测序数据；
S2)基于各候选突变位点的特征信息和目标值s，将各候选突变位点进行层次聚类分析，其中，目标值s由各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值a和各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值b计算获得；
S3)根据层次聚类分析结果，基于候选突变位点的突变丰度，获取各类候选突变位点的背景值pbg；
S4)根据层次聚类分析结果，基于候选突变位点的深度及其所在的分类的背景值pbg，获取各类候选突变位点的最低检测下限；
S5)根据各候选突变位点所在类别的背景值pbg及各候选突变位点的最低检测下限，确定各候选突变位点的检测阈值。
在本发明一些实施方式中，还包括：将基因样本的高通量测序数据与人类参考基因组数据进行比对，以鉴别基因样本中的候选突变位点。
在本发明一些实施方式中，通过BWA算法将基因样本的高通量测序数据与人类参考基因组数据进行比对；
和/或，通过VarDict算法鉴别基因样本中的候选突变位点。
在本发明一些实施方式中，所述特征信息包括深度、变异深度、正链参考碱基深度、负链参考碱基深度、正链变异碱基深度、负链变异碱基深度、基因型、突变丰度、链偏差、读取片段上的位置、读取片段上的位置的标准差、平均碱基质量分数、碱基质量分数标准差、比对质量、高质量片段比例、高质量片段突变丰度、是否为微卫星位点、微卫星位点单元长度、片段上的错配总数、5`端的序列、3`端的序列、突变类型、重复比例中的一种或多种的组合。
在本发明一些实施方式中，目标值s的计算方法如下：
其中，a为各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值；
b为各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值。
在本发明一些实施方式中，两个位点之间的距离d为两个位点的丰度的差值，优选为两个位点的丰度的差值的绝对值。
和/或，各类候选突变位点的背景值pbg为该类中各候选突变位点的突变丰度的中位数。
和/或，各候选突变位点的最低检测下限的计算方法为：
f＝ln(1-p)/-n
其中，f为候选突变位点的最低检测下限；
p为候选突变位点所在的分类的背景值pbg；
n为候选突变位点的深度。
和/或，各候选突变位点的检测阈值为其所在类别的背景值pbg以及各候选突变位点的最低检测下限两者相比较大的数值。
在本发明一些实施方式中，还包括：根据各候选突变位点的检测阈值，获取各候选突变位点的突变检测结果。
本发明另一方面提供一种计算机可读存储介质，其上存储有计算机程序，所述程序被处理器执行时实现上述的高通量测序突变检测方法。
本发明另一方面提供一种设备，包括：处理器和存储器，所述存储器用于存储计算机程序，所述处理器用于执行所述存储器存储的计算机程序，以使所述设备执行上述的高通量测序突变检测方法。
本发明另一方面提供一种装置，所述装置包括：
特征信息获取模块，用于获取基因样本中各候选突变位点的特征信息，所述特征信息来源于基因样本的高通量测序数据；
层次聚类分析模块，用于基于各候选突变位点的特征信息和目标值s，将各候选突变位点进行层次聚类分析，其中，目标值s由各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值a和各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值b计算获得；
背景值pbg计算模块，用于根据层次聚类分析结果，基于候选突变位点的突变丰度，获取各类候选突变位点的背景值pbg；
最低检测下限计算模块，用于根据层次聚类分析结果，基于候选突变位点的深度及其所在的分类的背景值pbg，获取各类候选突变位点的最低检测下限；
检测阈值计算模块，根据各候选突变位点所在类别的背景值pbg及各候选突变位点的最低检测下限，确定各候选突变位点的检测阈值；
优选的，还包括候选突变位点鉴别模块，用于将基因样本的高通量测序数据与人类参考基因组数据进行比对，以鉴别基因样本中的候选突变位点；
优选的，还包括突变检测结果计算模块，用于根据各候选突变位点的检测阈值，获取各候选突变位点的突变检测结果。
附图说明
图1显示为本发明所提供的高通量测序突变检测方法的流程示意图。
具体实施方式
为了使本发明的发明目的、技术方案和有益技术效果更加清晰，以下结合实施例对本发明进行进一步详细说明，熟悉此技术的人士可由本说明书所揭露的内容容易地了解本申请发明的其他优点及功效。
本发明第一方面提供一种高通量测序突变检测方法，包括：
S1)获取基因样本中各候选突变位点的特征信息，所述特征信息来源于基因样本的高通量测序数据；
S2)基于各候选突变位点的特征信息和目标值s，将各候选突变位点进行层次聚类分析，其中，目标值s由各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值a和各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值b计算获得；
S3)根据层次聚类分析结果，基于候选突变位点的突变丰度，获取各类候选突变位点的背景值pbg；
S4)根据层次聚类分析结果，基于候选突变位点的深度及其所在的分类的背景值pbg，获取各类候选突变位点的最低检测下限；
S5)根据各候选突变位点所在类别的背景值pbg及各候选突变位点的最低检测下限，确定各候选突变位点的检测阈值。
本发明所提供的高通量测序突变检测方法，可以包括：获取基因样本中各候选突变位点的特征信息，所述特征信息来源于基因样本的高通量测序数据。通常来说，本领域技术人员可选择合适的方法，根据基因样本的高通量测序数据，确定所需针对的各候选突变位点，并进一步获取基因样本中各候选突变位点的特征信息。例如，还可以包括：将基因样本的高通量测序数据与人类参考基因组数据进行比对，以鉴别基因样本中的候选突变位点。基因样本的高通量测序数据可以是Fastq文件等，其通常可以从高通量测序的下机数据(例如，BCL格式的数据等)转换获得(例如，通过BCL2fastq等软件)，将基因样本的高通量测序数据与人类参考基因组数据进行比对(例如，通过BWA算法等进行对比)，比对结果可以进一步被转换为BAM文件(例如，通过samtools等软件进行转换)，根据比对结果即可以进一步鉴别基因样本中的候选突变位点(例如，可以通过VarDict算法等进行鉴别)。再例如，特征信息通常可以从适当的文件中提取(例如，BAM文件等)，特征信息具体可以包括深度、变异深度、正链参考碱基深度、负链参考碱基深度、正链变异碱基深度、负链变异碱基深度、基因型、突变丰度、链偏差、读取片段上的位置、读取片段上的位置的标准差、平均碱基质量分数、碱基质量分数标准差、比对质量、高质量片段比例、高质量片段突变丰度、是否为微卫星位点、微卫星位点单元长度、片段上的错配总数、5`端的序列、3`端的序列、突变类型、重复比例等中的一种或多种的组合，这些特征的计算方式可以参照Lai Z,Markovets A,Ahdesmaki M,Chapman B,Hofmann O,McEwen R,Johnson J,Dougherty B,Barrett JC,DryJR.VarDict:a novel and versatile variant caller for next-generationsequencing in cancer research.Nucleic Acids Res.2016Jun 20；44(11):e108.doi:10.1093/nar/gkw227.Epub 2016Apr 7.PMID:27060149；PMCID:PMC4914105.等。各特征信息的中英文对应名称及特征描述可以如表1所示。
表1
上述高通量测序突变检测方法中，所获得的各候选突变位点的特征信息，可以以矩阵的方式给出，例如，可以构成一个所有候选突变位点的特征矩阵M(Mi,j)，其中，Mi,j为第i个候选突变位点的第j个特征信息的具体数值。
本发明所提供的高通量测序突变检测方法，还可以包括：基于各候选突变位点的特征信息和目标值s，将各候选突变位点进行层次聚类分析(Hierarchical clustering)，其中，目标值s由各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值a和各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值b计算获得。在层次聚类分析中，聚类层高阈值通常可以为h，通常可以选取h采用自适应的方式进行优化，h的取值通常满足min(h)<＝h<＝max(h)，h的取值范围通常可以根据待优化的集(例如，上述的特征矩阵M(Mi,j))进行确定，上述算法可以来自于层次聚类分析，可以通过软件R的hclust函数得到，优化的目标值则如上所述可以为目标值s，当h在一定的取值范围内变化时，不同的层次聚类分析结果均会对应不同的目标值s。通常来说，目标值s越小则说明得到的聚类内部差异越小，说明聚类效果更好。例如，所提供的层次聚类分析结果为目标值s最小的情况下的层次聚类分析结果。再例如，目标值s的计算方法可以如下所示：
其中，a为各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值；
b为各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值；
上述公式中，两个位点之间的距离d通常为两个位点的丰度的差值，更具体可以为两个位点的丰度的差值的绝对值。
本发明所提供的高通量测序突变检测方法，还可以包括：根据层次聚类分析结果，基于候选突变位点的突变丰度，获取各类候选突变位点的背景值pbg。在获取了层次聚类分析结果以后，可以根据各候选突变位点的突变丰度以及他们的分类结果，基于同一类的候选突变位点的突变丰度，获取各类候选突变位点的背景值pbg。例如，各类候选突变位点的背景值pbg可以为该类中各候选突变位点的突变丰度的中位数。
本发明所提供的高通量测序突变检测方法，还可以包括：根据层次聚类分析结果，基于候选突变位点的深度及其所在的分类的背景值pbg，获取各类候选突变位点的最低检测下限。在获取了层次聚类分析结果以后，可以根据他们的分类结果，基于各候选突变位点的深度及其所在的分类的背景值pbg，获取各类候选突变位点的最低检测下限。例如，各候选突变位点的最低检测下限的计算方法可以为：
f＝ln(1-p)/-n
其中，f为候选突变位点的最低检测下限；
p为候选突变位点所在的分类的背景值pbg；
n为候选突变位点的深度。
本发明所提供的高通量测序突变检测方法，还可以包括：根据各候选突变位点所在类别的背景值pbg及各候选突变位点的最低检测下限，确定各候选突变位点的检测阈值。通常来说，各候选突变位点其所在类别的背景值pbg以及各候选突变位点的最低检测下限两者相比较大的数值，可以作为各候选突变位点的检测阈值，其原因在于，最低检测下限为该位点深度下能达到的理论最低值，背景值决定了背景噪音的大小，最低检测下限小于背景值的话以最低检测下限为阈值，否则最低只能达到背景值为下限。
本发明所提供的高通量测序突变检测方法，还可以包括：根据各候选突变位点的检测阈值，获取各候选突变位点的突变检测结果。通常来说，各候选突变位点的检测阈值可以对应于各候选突变位点特征信息中的突变丰度，并可以根据两者的比较结果，获取各候选突变位点的突变检测结果。例如，当候选突变位点的突变丰度≥候选突变位点的检测阈值，则可以认为基因样本中候选突变位点的突变为阳性。再例如，当候选突变位点的突变丰度＜候选突变位点的检测阈值，则可以认为基因样本中候选突变位点的突变为阴性。
本发明第二方面提供一种计算机可读存储介质，其上存储有计算机程序，所述程序被处理器执行时实现本发明第一方面所提供的高通量测序突变检测方法。
本发明第三方面提供一种设备，包括：处理器和存储器，所述存储器用于存储计算机程序，所述处理器用于执行所述存储器存储的计算机程序，以使所述设备执行本发明第一方面所提供的高通量测序突变检测方法。
本发明第四方面提供一种装置，所述装置包括：
特征信息获取模块，用于获取基因样本中各候选突变位点的特征信息，所述特征信息来源于基因样本的高通量测序数据；
层次聚类分析模块，用于基于各候选突变位点的特征信息和目标值s，将各候选突变位点进行层次聚类分析，其中，目标值s由各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值a和各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值b计算获得；
背景值pbg计算模块，用于根据层次聚类分析结果，基于候选突变位点的突变丰度，获取各类候选突变位点的背景值pbg；
最低检测下限计算模块，用于根据层次聚类分析结果，基于候选突变位点的深度及其所在的分类的背景值pbg，获取各类候选突变位点的最低检测下限；
检测阈值计算模块，用于根据各候选突变位点所在类别的背景值pbg及各候选突变位点的最低检测下限，确定各候选突变位点的检测阈值；
优选的，还包括候选突变位点鉴别模块，用于将基因样本的高通量测序数据与人类参考基因组数据进行比对，以鉴别基因样本中的候选突变位点；
优选的，还包括突变检测结果计算模块，用于根据各候选突变位点的检测阈值，获取各候选突变位点的突变检测结果。
本发明中，上述装置中各模块的运行原理可以参照本发明第一方面所提供的高通量测序突变检测方法，在此不做赘述。
本申请所提供的高通量测序突变检测方法，可以通过自适应优化的算法，针对每个样本特有的数据特点智能的设定分析的阈值，在保证阳性位点检出的前提下尽可能的去除噪音造成的假阳性，可以更加精准的从样本(例如，肿瘤样本)中检测突变，从而具有良好的产业化前景。
下面通过实施例对本申请予以进一步说明，但并不因此而限制本申请的范围。
实施例1
以肿瘤样本检测数据下机开始的整个分析流程为例，具体如下：
1)测序下机数据通过BCL2fastq从测序BCL数据中分离样本数据并转换成fastq文件，使用参数bcl2fastq--barcode-mismatches 1-o./demultiplex--ignore-missing-bcls--no-lane-splitting；
2)Fastq文件通过BWA算法比对到人类参考基因组并利用samtools软件转换为BAM文件，使用参数bwa mem-t 16-R"@RG\tID:DNA\tLB:DNA\tSM:S2100019497-Plasma\tPL:ILLUMINA"-Mhuman_g1k_v37_decoy.fasta；
3)使用VarDict算法鉴别出样本的所有候选突变，使用参数VarDict-b bam-p-GREF-c1-S 2-E 3-g 5；
4)对于该样本中任意的一个候选突变，使用VarDict提取表一中的所有特征，构建特征矩阵；
5)使用R,hclust函数基于特征矩阵进行层次聚类，以0.01作为h变化的步长，计算每个h下的s值，并根据s，确定最佳的分类方式；
6)定义每个分类内的变异的丰度值中位数为该分类的背景值，并结合突变深度定义AF阈值，阈值定义为所在类别的背景值pbg以及各候选突变位点的最低检测下限两者相比较大的数值；
7)筛选候选突变列表，标记高于该阈值的突变为真实突变。
8)使用已知突变位点的标准样品(例如Horizon HD780)进行标准建库实验和测序，得到的测序数据经过上述处理，比对结果如表1所示：
表1
| 处理方式 | 真阳性位点 | 假阳性位点 |
| 本专利方法 | 8/8 | 5 |
| 标准分析流程(1-3步) | 8/8 | 61 |
可见，本发明提供的筛选算法可以在检出真阳性位点的前提下明显降低假阳性位点的检出。
综上所述，本发明有效克服了现有技术中的种种缺点而具高度产业利用价值。
上述实施例仅例示性说明本发明的原理及其功效，而非用于限制本发明。任何熟悉此技术的人士皆可在不违背本发明的精神及范畴下，对上述实施例进行修饰或改变。因此，举凡所属技术领域中具有通常知识者在未脱离本发明所揭示的精神与技术思想下所完成的一切等效修饰或改变，仍应由本发明的权利要求所涵盖。
Claims (10)
1.一种高通量测序突变检测方法，包括：
S1)获取基因样本中各候选突变位点的特征信息，所述特征信息来源于基因样本的高通量测序数据；
S2)基于各候选突变位点的特征信息和目标值s，将各候选突变位点进行层次聚类分析，其中，目标值s由各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值a和各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值b计算获得；
S3)根据层次聚类分析结果，基于候选突变位点的突变丰度，获取各类候选突变位点的背景值pbg；
S4)根据层次聚类分析结果，基于候选突变位点的深度及其所在的分类的背景值pbg，获取各类候选突变位点的最低检测下限；
S5)根据各候选突变位点所在类别的背景值pbg及各候选突变位点的最低检测下限，确定各候选突变位点的检测阈值。
2.如权利要求1所述的高通量测序突变检测方法，其特征在于，还包括：将基因样本的高通量测序数据与人类参考基因组数据进行比对，以鉴别基因样本中的候选突变位点。
3.如权利要求2所述的高通量测序突变检测方法，其特征在于，通过BWA算法将基因样本的高通量测序数据与人类参考基因组数据进行比对；
和/或，通过VarDict算法鉴别基因样本中的候选突变位点。
4.如权利要求1所述的高通量测序突变检测方法，其特征在于，所述特征信息包括深度、变异深度、正链参考碱基深度、负链参考碱基深度、正链变异碱基深度、负链变异碱基深度、基因型、突变丰度、链偏差、读取片段上的位置、读取片段上的位置的标准差、平均碱基质量分数、碱基质量分数标准差、比对质量、高质量片段比例、高质量片段突变丰度、是否为微卫星位点、微卫星位点单元长度、片段上的错配总数、5`端的序列、3`端的序列、突变类型、重复比例中的一种或多种的组合。
5.如权利要求1所述的高通量测序突变检测方法，其特征在于，目标值s的计算方法如下：
其中，a为各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值；b为各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值。
6.如权利要求1所述的高通量测序突变检测方法，其特征在于，两个位点之间的距离d为两个位点的丰度的差值，优选为两个位点的丰度的差值的绝对值。
和/或，各类候选突变位点的背景值pbg为该类中各候选突变位点的突变丰度的中位数。
和/或，各候选突变位点的最低检测下限的计算方法为：
f＝ln(1–p)/–n
其中，f为候选突变位点的最低检测下限；
p为候选突变位点所在的分类的背景值pbg；
n为候选突变位点的深度。
和/或，各候选突变位点的检测阈值为其所在类别的背景值pbg以及各候选突变位点的最低检测下限两者相比较大的数值。
7.如权利要求1所述的高通量测序突变检测方法，其特征在于，还包括：根据各候选突变位点的检测阈值，获取各候选突变位点的突变检测结果。
8.一种计算机可读存储介质，其上存储有计算机程序，其特征在于，所述程序被处理器执行时实现如权利要求1-7任一权利要求所述的高通量测序突变检测方法。
9.一种设备，包括：处理器和存储器，所述存储器用于存储计算机程序，所述处理器用于执行所述存储器存储的计算机程序，以使所述设备执行如权利要求1-7任一权利要求所述的高通量测序突变检测方法。
10.一种装置，所述装置包括：
特征信息获取模块，用于获取基因样本中各候选突变位点的特征信息，所述特征信息来源于基因样本的高通量测序数据；
层次聚类分析模块，用于基于各候选突变位点的特征信息和目标值s，将各候选突变位点进行层次聚类分析，其中，目标值s由各候选突变位点与其所在的分类内部其他候选突变位点的距离d的平均值a和各候选突变位点与最近的一个分类中的候选突变位点的距离d的平均值b计算获得；
背景值pbg计算模块，用于根据层次聚类分析结果，基于候选突变位点的突变丰度，获取各类候选突变位点的背景值pbg；
最低检测下限计算模块，用于根据层次聚类分析结果，基于候选突变位点的深度及其所在的分类的背景值pbg，获取各类候选突变位点的最低检测下限；
检测阈值计算模块，根据各候选突变位点所在类别的背景值pbg及各候选突变位点的最低检测下限，确定各候选突变位点的检测阈值；
优选的，还包括候选突变位点鉴别模块，用于将基因样本的高通量测序数据与人类参考基因组数据进行比对，以鉴别基因样本中的候选突变位点；
优选的，还包括突变检测结果计算模块，用于根据各候选突变位点的检测阈值，获取各候选突变位点的突变检测结果。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110972081.1A CN113628683B (zh) | 2021-08-24 | 2021-08-24 | 一种高通量测序突变检测方法、设备、装置及可读存储介质 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110972081.1A CN113628683B (zh) | 2021-08-24 | 2021-08-24 | 一种高通量测序突变检测方法、设备、装置及可读存储介质 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113628683A true CN113628683A (zh) | 2021-11-09 |
| CN113628683B CN113628683B (zh) | 2024-04-09 |
Family
ID=78387385
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110972081.1A Active CN113628683B (zh) | 2021-08-24 | 2021-08-24 | 一种高通量测序突变检测方法、设备、装置及可读存储介质 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113628683B (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115312122A (zh) * | 2022-10-12 | 2022-11-08 | 之江实验室 | 一种CRISPR-Cas酶可突变位点推荐方法和装置 |
Citations (19)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2008005559A2 (en) * | 2006-07-07 | 2008-01-10 | Johns Hopkins University | A strategy for detecting low abundance mutations |
| JP2016077227A (ja) * | 2014-10-17 | 2016-05-16 | 国立大学法人大阪大学 | ゲノム解析装置、ゲノム解析方法及びゲノム解析プログラム |
| WO2016127944A1 (en) * | 2015-02-10 | 2016-08-18 | The Chinese University Of Hong Kong | Detecting mutations for cancer screening and fetal analysis |
| CN107944225A (zh) * | 2017-11-28 | 2018-04-20 | 慧算医疗科技（上海）有限公司 | 基因高通量测序数据突变检测方法 |
| CN107944223A (zh) * | 2017-11-10 | 2018-04-20 | 深圳裕策生物科技有限公司 | 基于二代测序的点突变检测过滤方法、装置和存储介质 |
| US20180148778A1 (en) * | 2014-12-16 | 2018-05-31 | Garvan Institute Of Medical Research | Sequencing controls |
| KR20180060764A (ko) * | 2016-11-29 | 2018-06-07 | 연세대학교 산학협력단 | 염기서열의 변이 검출방법 및 이를 이용한 염기서열의 변이 검출 디바이스 |
| KR20180060759A (ko) * | 2016-11-29 | 2018-06-07 | 연세대학교 산학협력단 | 염기서열의 변이 검출방법 및 이를 이용한 염기서열의 변이 검출 디바이스 |
| CN108603229A (zh) * | 2016-01-22 | 2018-09-28 | 格瑞尔公司 | 用于高保真测序的方法和系统 |
| CN108690871A (zh) * | 2018-03-29 | 2018-10-23 | 深圳裕策生物科技有限公司 | 基于二代测序的插入缺失突变检测方法、装置和存储介质 |
| WO2019031785A2 (ko) * | 2017-08-07 | 2019-02-14 | 연세대학교 산학협력단 | 유전자 패널에 기초한 염기서열의 변이 검출방법 및 이를 이용한 염기서열의 변이 검출 디바이스 |
| CN109439729A (zh) * | 2018-12-27 | 2019-03-08 | 上海鲸舟基因科技有限公司 | 检测低频变异用的接头、接头混合物及相应方法 |
| CN109637587A (zh) * | 2019-01-18 | 2019-04-16 | 臻悦生物科技江苏有限公司 | 检测基因融合突变的方法、装置、存储介质、处理器及转录组数据表达量标准化的方法 |
| CN109887548A (zh) * | 2019-01-18 | 2019-06-14 | 臻悦生物科技江苏有限公司 | 基于捕获测序的ctDNA占比的检测方法及检测装置 |
| CN111139291A (zh) * | 2020-01-14 | 2020-05-12 | 首都医科大学附属北京安贞医院 | 一种单基因遗传性疾病高通量测序分析方法 |
| WO2020181111A1 (en) * | 2019-03-07 | 2020-09-10 | Nvidia Corporation | Genetic mutation detection using deep learning |
| CN111718982A (zh) * | 2020-07-09 | 2020-09-29 | 深圳市海普洛斯生物科技有限公司 | 一种肿瘤组织单样本体细胞突变检测方法及装置 |
| WO2020199336A1 (zh) * | 2019-03-29 | 2020-10-08 | 北京市商汤科技开发有限公司 | 一种基因变异识别方法、装置和存储介质 |
| CN114708908A (zh) * | 2022-03-04 | 2022-07-05 | 至本医疗科技（上海）有限公司 | 检测实体瘤微小残留病灶的方法、计算设备和存储介质 |
-
2021
- 2021-08-24 CN CN202110972081.1A patent/CN113628683B/zh active Active
Patent Citations (20)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2008005559A2 (en) * | 2006-07-07 | 2008-01-10 | Johns Hopkins University | A strategy for detecting low abundance mutations |
| JP2016077227A (ja) * | 2014-10-17 | 2016-05-16 | 国立大学法人大阪大学 | ゲノム解析装置、ゲノム解析方法及びゲノム解析プログラム |
| US20180148778A1 (en) * | 2014-12-16 | 2018-05-31 | Garvan Institute Of Medical Research | Sequencing controls |
| WO2016127944A1 (en) * | 2015-02-10 | 2016-08-18 | The Chinese University Of Hong Kong | Detecting mutations for cancer screening and fetal analysis |
| CN108603229A (zh) * | 2016-01-22 | 2018-09-28 | 格瑞尔公司 | 用于高保真测序的方法和系统 |
| US20190338349A1 (en) * | 2016-01-22 | 2019-11-07 | Grail, Inc. | Methods and systems for high fidelity sequencing |
| KR20180060764A (ko) * | 2016-11-29 | 2018-06-07 | 연세대학교 산학협력단 | 염기서열의 변이 검출방법 및 이를 이용한 염기서열의 변이 검출 디바이스 |
| KR20180060759A (ko) * | 2016-11-29 | 2018-06-07 | 연세대학교 산학협력단 | 염기서열의 변이 검출방법 및 이를 이용한 염기서열의 변이 검출 디바이스 |
| WO2019031785A2 (ko) * | 2017-08-07 | 2019-02-14 | 연세대학교 산학협력단 | 유전자 패널에 기초한 염기서열의 변이 검출방법 및 이를 이용한 염기서열의 변이 검출 디바이스 |
| CN107944223A (zh) * | 2017-11-10 | 2018-04-20 | 深圳裕策生物科技有限公司 | 基于二代测序的点突变检测过滤方法、装置和存储介质 |
| CN107944225A (zh) * | 2017-11-28 | 2018-04-20 | 慧算医疗科技（上海）有限公司 | 基因高通量测序数据突变检测方法 |
| CN108690871A (zh) * | 2018-03-29 | 2018-10-23 | 深圳裕策生物科技有限公司 | 基于二代测序的插入缺失突变检测方法、装置和存储介质 |
| CN109439729A (zh) * | 2018-12-27 | 2019-03-08 | 上海鲸舟基因科技有限公司 | 检测低频变异用的接头、接头混合物及相应方法 |
| CN109637587A (zh) * | 2019-01-18 | 2019-04-16 | 臻悦生物科技江苏有限公司 | 检测基因融合突变的方法、装置、存储介质、处理器及转录组数据表达量标准化的方法 |
| CN109887548A (zh) * | 2019-01-18 | 2019-06-14 | 臻悦生物科技江苏有限公司 | 基于捕获测序的ctDNA占比的检测方法及检测装置 |
| WO2020181111A1 (en) * | 2019-03-07 | 2020-09-10 | Nvidia Corporation | Genetic mutation detection using deep learning |
| WO2020199336A1 (zh) * | 2019-03-29 | 2020-10-08 | 北京市商汤科技开发有限公司 | 一种基因变异识别方法、装置和存储介质 |
| CN111139291A (zh) * | 2020-01-14 | 2020-05-12 | 首都医科大学附属北京安贞医院 | 一种单基因遗传性疾病高通量测序分析方法 |
| CN111718982A (zh) * | 2020-07-09 | 2020-09-29 | 深圳市海普洛斯生物科技有限公司 | 一种肿瘤组织单样本体细胞突变检测方法及装置 |
| CN114708908A (zh) * | 2022-03-04 | 2022-07-05 | 至本医疗科技（上海）有限公司 | 检测实体瘤微小残留病灶的方法、计算设备和存储介质 |
Non-Patent Citations (1)
| Title |
|---|
| 邢文昊;刘永壮;王亚东;: "基于高通量测序数据的插入/删除新突变检测方法", 智能计算机与应用, no. 01, 1 January 2020 (2020-01-01), pages 265 - 269 * |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115312122A (zh) * | 2022-10-12 | 2022-11-08 | 之江实验室 | 一种CRISPR-Cas酶可突变位点推荐方法和装置 |
| CN115312122B (zh) * | 2022-10-12 | 2022-12-16 | 之江实验室 | 一种CRISPR-Cas酶可突变位点推荐方法和装置 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113628683B (zh) | 2024-04-09 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107918931B (zh) | 图像处理方法及系统及计算机可读存储介质 | |
| EP3306566B1 (en) | Method and system for processing image | |
| Rougemont et al. | Probabilistic base calling of Solexa sequencing data | |
| CN107229839B (zh) | 一种基于新一代测序数据的Indel检测方法 | |
| CN107944225B (zh) | 基因高通量测序数据突变检测方法 | |
| CN108154010B (zh) | 一种ctDNA低频突变测序数据分析方法和装置 | |
| KR102273257B1 (ko) | 리드 깊이 기반한 유전자 복제수 변이 검출 기법 및 분석장치 | |
| CN111868832A (zh) | 识别拷贝数异常的方法 | |
| JP5946030B2 (ja) | デジタル病理画像における蛍光ドット計数法 | |
| CN113628683B (zh) | 一种高通量测序突变检测方法、设备、装置及可读存储介质 | |
| CN111292807B (zh) | 一种单细胞转录组数据中分析双细胞的方法 | |
| CN117059173A (zh) | 一种拷贝数变异精确断裂点识别的方法及其应用 | |
| CN111696622A (zh) | 一种校正和评估变异检测软件检测结果的方法 | |
| JP7057973B2 (ja) | 微生物識別装置および微生物識別方法 | |
| Galinsky | Automatic registration of microarray images. I. Rectangular grid | |
| Roayaei Ardakany et al. | Efficient and accurate detection of topologically associating domains from contact maps | |
| US10443090B2 (en) | Method and apparatus for detecting translocation | |
| CN112614542B (zh) | 一种微生物鉴定方法、装置、设备及存储介质 | |
| CN114242164B (zh) | 一种全基因组复制的分析方法、装置和存储介质 | |
| Ardakany | Efficient Algorithms for the Analysis of Hi-C Contact Maps | |
| Roayaei Ardakany | Efficient Algorithms for the Analysis of Hi-C Contact Maps | |
| CN112734777B (zh) | 一种基于簇形状边界闭包聚类的图像分割方法及系统 | |
| Zinovyev | Visualizing the spatial structure of triplet distributions in genetic texts | |
| CN114664383A (zh) | 一种结合参考库先验知识的宏基因组学分箱方法及系统 | |
| CN117894367A (zh) | 微生物特异性序列保守性的筛选评估方法 |
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 | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |