CN114999573B - 一种基因组变异检测方法及检测系统 - Google Patents

一种基因组变异检测方法及检测系统 Download PDF

Info

Publication number
CN114999573B
CN114999573B CN202210391526.1A CN202210391526A CN114999573B CN 114999573 B CN114999573 B CN 114999573B CN 202210391526 A CN202210391526 A CN 202210391526A CN 114999573 B CN114999573 B CN 114999573B
Authority
CN
China
Prior art keywords
sequence
copy number
sequencing
sequencing sequence
reference sequence
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
CN202210391526.1A
Other languages
English (en)
Other versions
CN114999573A (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.)
Harbin Yinji Technology Co ltd
Original Assignee
Harbin Yinji Technology 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 Harbin Yinji Technology Co ltd filed Critical Harbin Yinji Technology Co ltd
Priority to CN202210391526.1A priority Critical patent/CN114999573B/zh
Publication of CN114999573A publication Critical patent/CN114999573A/zh
Application granted granted Critical
Publication of CN114999573B publication Critical patent/CN114999573B/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
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G16B30/10—Sequence alignment; Homology search
    • 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/20—Allele or variant detection, e.g. single nucleotide polymorphism [SNP] detection
    • 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/30—Detection of binding sites or motifs
    • 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
    • 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

  • Life Sciences & Earth Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Health & Medical Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Medical Informatics (AREA)
  • General Health & Medical Sciences (AREA)
  • Theoretical Computer Science (AREA)
  • Biophysics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Biotechnology (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Analytical Chemistry (AREA)
  • Chemical & Material Sciences (AREA)
  • Molecular Biology (AREA)
  • Genetics & Genomics (AREA)
  • Bioethics (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Data Mining & Analysis (AREA)
  • Databases & Information Systems (AREA)
  • Epidemiology (AREA)
  • Evolutionary Computation (AREA)
  • Public Health (AREA)
  • Software Systems (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

本发明提出了一种基因组变异检测方法和检测系统，将样本的基因组划分成连续非重叠窗口，将每个窗口参考序列分别和测序序列进行序列比对，并进行Reads深度信息的标准化处理，根据所述序列比对结果以及深度信息，确定所述基因组的潜在变异区域；应用De‑Bruijn图和哈希表，将测序序列与参考序列进行匹配，得到最佳单倍体；根据中间比对结果与最终比对结果之间的映射关系，获得每个测序序列与参考序列的校正对齐结果；根据每个测序序列与所述参考序列的校正对齐结果，确定断点；对断点进行筛选，确定拷贝数变异区域，利用深度信息计算拷贝数变异区域内的拷贝数变异比例，得到准确的拷贝数变异检测结果。

Description

一种基因组变异检测方法及检测系统
技术领域
本发明涉及生物信息学领域，更具体地涉及基因组变异检测方法及检测系统。
背景技术
基因研究近年来已经看到快速的进展。已经测序包括一些个体人类的特定有机体的完整基因组，并且其作为参考变得可用，即用作研究相同物种的成员的参考的基因序列。在基因研究、基因测试、个体化医药和多个其它应用中，通常有用的是获取基因材料的样本、确定该样本的序列、以及关于一个或多个现有参考分析该样本以鉴定序列变异或获取关于样本的其它有用信息。
基因测试的现有方法通常为定位或映射长、邻接样本序列到参考中的位置。然而，一些用于获得样本序列的技术产生包括具有预定空间关系的多个较短序列的多核苷酸序列，在一些情况下包括具有可变距离的空间关系的多个较短序列的多核苷酸序列。在后者的情况下，在这些较短序列中碱基的相对基因组位置仅是近似已知，并且通常是以具有可变但受限数量的间隔或重叠(称为缺口距离)的短邻接读出(read)的形式。
CNV是基因组结构变异的一种形式。CNV的狭义定义通常是指染色体中DNA片段的拷贝数变化。这种形式的基因组结构变异的类型和原因可以包括：缺失(末端缺失、间质缺失)；易位(相互易位、罗伯逊易位)；反转；环状染色体；双着丝粒染色体；CNV的更广泛的定义还包括例如结构变异，例如染色体非整倍性和部分非整倍性。
目前可用的检测拷贝数变异的方法主要包括高分辨率染色体核型分析、FISH(荧光原位杂交)、阵列CGH(阵列比较基因组杂交)、MLPA(多重连接依赖性探针扩增)、PCR(聚合酶链反应)等，其中FISH检测被认为是遗传诊断的黄金标准，其可以有效地用于检测大多数已知的染色体缺失或重复。然而，这些方法通常具有低效率，特别是当用于全基因组扫描时，这可能消耗大量资源或可能无法检测精确的CNV。因此，迫切需要一种检测基因组变异的新方法，以检测的精度。
发明内容
本发明为了解决上述技术问题，提出了一种检测拷贝数变异的方法及其系统，该基因组变异检测方法包括如下步骤：
步骤1，将样本的基因组划分成长度相等的连续非重叠窗口，将每个窗口内的测序序列分别和参考序列进行序列比对；
步骤2，进行reads深度信息的标准化处理；
步骤3，根据所述序列比对结果以及深度信息，确定所述基因组的潜在变异区域；
步骤4，对潜在变异区域的标准化后的测序序列应用De-Bruijn图和哈希表，使用隐马氏匹配算法依次将各测序序列的与参考序列进行匹配，得到N个最佳单倍体；
步骤5，分别将潜在变异区域内的测序序列与N个最佳单倍体进行比对后得到中间比对结果，N个最佳单倍体与参考序列片段进行比对后的最终比对结果；
步骤6，根据中间比对结果与最终比对结果之间的映射关系，获得每个测序序列与参考序列片段的校正对齐结果；
步骤7，根据每个测序序列与所述参考序列的校正对齐结果，确定断点；
步骤8，对断点进行筛选，确定拷贝数变异区域，利用深度信息计算拷贝数变异区域内的拷贝数变异比例，并拷贝数变异比例与阈值进行比较，得到准确的拷贝数变异检测结果。
进一步地，所述步骤3具体包括：采用基于中位数归一化的方法对GC含量和Mappability值进行校正。
进一步地，所述步骤5中：在测序序列与N个最佳单倍体进行比对时的不匹配的区域添加第一空白区域，在N个最佳单倍与参考序列进行比对时的不匹配的区域添加第二空白区域，用于使测序序列、N个最佳单倍体和参考序列的匹配区域的位置完全对应。
进一步地，所述步骤6中：采用预先训练的神经网络模型对数据映射关系进行处理,采用径向基核函数的支持向量机模型进行多次训练，正向映射之后的输出值与目标值进行比较，得到预测误差之后，将误差逆映射至各个映射参数，直到映射参数接近训练集目标。
进一步地，所述步骤8中：对断点进行筛选具体包括：将显著值最小的m个断点作为候选断点进行排列，从每轮循环和迭代中移除具有最小显著值的候选断点，并且重新计算两个相邻断点的显著值，循环迭代，直至所有候选断点的显著值都小于终止阈值。
进一步地，所述步骤8中：在已获得筛选的断点之后，确定拷贝数变异区域。设两个断点之间的窗口为，(xL，xR),待测试样品相对于正常样品的拷贝数变异比例按照下式计算：
Figure BDA0003595778510000031
其中与参考序列比对的正常样品的Read和待测样品的Read的数量分别是
aN和aT，并且落入窗口(xL，xR)内的Read的数量分别是N(xL，xR)和T(xL，xR),将所述拷贝数变异比例≤0.75和≥1.25分别作为染色体片段缺失和重复的检测阈值。
该基因组变异检测系统，用于实现基因组变异检测方法，包括：
序列比对单元，用于将每个窗口内的测序序列分别和参考序列进行序列比对；
数据处理单元，对测序序列进行reads深度信息的标准化处理；
确定单元，用于根据所述序列比对结果以及深度信息，确定所述基因组的潜在变异区域；
匹配单元，用于对潜在变异区域标准化后的测序序列应用De-Bruijn图和哈希表，使用隐马氏匹配算法依次将各测序序列的与参考序列进行匹配，得到N个最佳单倍体；
映射单元，用于将潜在变异区域内的测序序列与N个最佳单倍体进行比对后得到中间比对结果，将N个最佳单倍体与参考序列片段进行比对后的最终比对结果；根据中间比对结果与最终比对结果之间的映射关系，获得每个测序序列与参考序列片段的校正对齐结果；
断点检测单元，用于根据每个测序序列与所述参考序列的校正对齐结果，确定断点；
变异检测单元，用于确定拷贝数变异区域，利用深度信息计算拷贝数变异区域内的拷贝数变异比例，并拷贝数变异比例与阈值进行比较，得到准确的拷贝数变异检测结果。
进一步地，所述数据处理单元，具体采用基于中位数归一化的方法对GC含量和Mappability值进行校正。
进一步地，所述映射单元，具体采用预先训练的神经网络模型对数据映射关系进行处理,采用径向基核函数的支持向量机模型进行多次训练，正向映射之后的输出值与目标值进行比较，得到预测误差之后，将误差逆映射至各个映射参数，直到映射参数接近训练集目标。
进一步地，所述变异检测单元通过对断点进行筛选，确定拷贝数变异区域，所述对断点进行筛选具体包括：将显著值最小的m个断点作为候选断点进行排列，从每轮循环和迭代中移除具有最小显著值的候选断点，并且重新计算两个相邻断点的显著值，循环迭代，直至所有候选断点的显著值都小于终止阈值。
本发明与现有技术中的基因组变异检测的分析方法对比，提高了检测的准确率，可精确检测到微缺失区域；利用训练的神经网络模型最大的优势是可以和任何一种检测方法相结合，该检测过程可以从机器学习过程中自动的获得权值参数，而不需要人为的设定，所以避免了人为的错误，提高了检测的精度。
附图说明
附图1为本发明的检测拷贝数变异的方法的流程示意图；
附图2为本发明的基因组变异检测系统示意图；
附图3为本发明的为序列比对示意图；
附图4为GC含量对reads深度的影响点图；
附图5为Mappability对reads深度的影响点图；
附图6为隐藏状态转移概率示意图；
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述。详细描述本方案之前，为了便于理解本发明实施例所述的技术方案，对本发明实施例中的一些重要名词进行详细解释，需要理解的是，下述名词仅是本发明实施例技术人员为了描述方便进行的命名，并不代表或暗示所指的系统或元件必须有此命名，因此不能理解为对本发明实施例的限制：
基因组：包含人的所有遗传信息，而不单单只是一些外在和内在特征，也包含很多目前而言不明其功能性的碱基序列，基因是基因组的一个子集，基因是控制性状的遗传单元，性状为个体的各种外在和内在特征。在基因组的碱基序列中，按照碱基的排列顺序对每个碱基进行顺序编码，使每一个碱基具有一个相匹配的编码，单个编码代表基因组中的一个碱基对，连续的编码区间代表基因组中的一个碱基片段。
碱基序列：是一种由核苷酸分子构成的长链聚合物，这种长链聚合物具有方向性，核苷酸分子也可称为碱基，不同的碱基可以相互组合，使得DNA序列形成双螺旋结构。碱基序列根据核苷酸分子上携带的四种不同碱基类型可以抽象为由字符集{A，C，G，T}构成的字符串，可以结合的碱基称为互补碱基，互补规则为A与T互补，G与C互补。
测序：确定碱基序列内容的过程，如：将完整的样本碱基序列打碎，从中筛选出满足特定长度(通常为数百bp)的片段，然后在每个片段的一端或者两端各读取一段长度为数十至数百bp的序列，这些读取出的序列长度通常远远小于被测样本碱基序列的长度，使得全部短序列的总长度达到样本DNA长度的数倍至数十倍，从而使获得样本碱基序列成为可能。
样本：在测序中提供被测序碱基序列的个体。
测序序列：是指一段特定长度的DNA片段，该DNA片段可以从打碎的碱基序列上读取出来。
参考序列：“人类基因组计划”中所组装出来的人类基因组，是通过拼接方法获得的一套人类碱基序列。参考基因组并不是某个人类个体的基因组序列，它来自于多个具有代表性的人类种族个体的序列，是一套综合的序列。由于双链碱基序列呈互补关系，根据其中一条碱基序列即可以获得另一条序列内容，参考序列仅包括双链碱基序列中的一条。其他物种，如果蝇、小鼠等同样有参考序列，本发明实施例是针对人类基因组进行研究，故在本发明实施例中参考序列仅指人类基因组参考序列。其中，参考序列与碱基序列相同，都可以看作为定义在字符集{A，C，G，T}上的字符串，其中“A”、“C”、“G"、“T”分别对应核苷酸所携带的四种碱基。
基因组变异：是指基因组中碱基对组成或排列顺序的改变，主要包括：单核苷酸变异和indel(Insertion和Deletion的简称)两种；单核苷酸变异：通常称为单核苷酸多态性(英文：Single Nucleotide Polymorphism，SNP)，通俗的说法就是单个DNA碱基的不同。indel：指的是基因组的某个位置上所发生的小片段序列的插入或删除，其长度通常在100bp以下。
单倍型，是单倍体基因型的简称，在遗传学上是指在同一染色体上进行共同遗传的多个基因座上等位基因的组合；通俗的说法就是若干个决定同一性状的紧密连锁的基因构成的基因性。按照某一指定基因座上基因重组发生的数量，单倍型甚至可以指至少两个基因座或整个染色体。
本实施例中从样品的核酸分子的至少一部分获得读段。可以对样品中的核酸分子的至少一部分或整个核酸分子进行测序以获得Read。例如，将来自样品的基因组DNA分子随机片段化以获得DNA片段，然后对其进行测序以获得具有一定长度的读段。获得的读段的长度可以在一定范围内，并且可以通过截断获得具有固定长度的读段。DNA片段的长度可以为50bp/1500bp，例如50bp/150bp、150bp/350bp、350bp/500bp、500bp/700bp、700bp/1000bp或1000bp/1500bp。例如，DNA片段可以具有50bp、90bp、100bp、150bp、300bp、350bp、500bp、700bp、1000bp、1500bp的长度。在一个示例中，300bp和700bp是优选的，350bp和500bp是更优选的。
读段的长度可能由于不同的测序仪而具有大的差异，例如，诸如illumina-solexa和life technologies-solid的设备的一般序列长度在300bp的范围内，而通过Roche-454、常规Sanger测序、超现代单分子测序系统获得的序列长度可以是大约或超过1000bp。为了满足独特比对的要求，当选择独特映射的读段时，通常选择具有20bp或更长的长度的序列，优选地，所选择的序列具有26bp或更长的长度。
如图1所示，为本发明的检测拷贝数变异的方法的流程示意图，首先将样本基因组划分成长度相等的连续非重叠窗口(窗口长度默认为200bp)，将基因组的多窗口测序序列分别和参考序列进行序列比对。
对测序序列进行数据进行reads深度的标准化。在基因组的某些区域，异常的reads深度不仅可能是拷贝数变异的反映，也可能是测序过程中其它因素引起的偏差，其中GC含量和Mappability是影响最大的两个因素，如图4和图5所示，分别为GC含量对reads深度的影响点图和Mappability对reads深度的影响点图；reads深度的标准化即是指校正GC含量和Mappability等因素引起的偏差，使得基因组特定区域的异常reads深度能够真实地反映具体的拷贝数。
本实施例中采用基于中位数归一化的方法对GC含量进行校正：
Figure BDA0003595778510000071
其中ri表示校正前某一窗口的reads深度值，r′i，GC表示GC含量校正后该窗口的reads深度值，m是所有窗口的reads深度的中位数，mGC表示和当前窗口拥有有同样GC含量的所有窗口的reads深度的中位数。
采用如下方法进行Mappability校正：
Figure BDA0003595778510000072
其中ri表示校正前某一窗口的reads深度值，r′i，MAP表示Mappability校正后该窗口的reads深度值，m是所有窗口的reads深度的中位数，mMAP表示和当前窗口拥有同样Mappability的所有窗口的reads深度的中位数。
在GC含量校正和Mappability校正的同时，需要对:reads深度信息构建混合概率分布来拟合真实测序数据的reads深度。
根据所述序列比对结果以及深度信息，确定所述基因组的潜在变异区域；如图3所示为序列比对示意图。
对标准化后的测序序列应用De-Bruijn图和哈希表，依次读取各测序序列的与参考序列进行匹配，使用Pair HMM算法(隐马氏匹配算法)对序列进行重新比对，并结合其单倍型的概率信息得到N个最佳单倍体(N为大于等于1的整数)，
具体地，对于序列比对来说，隐藏状态存在碱基匹配和Indel这两种状态，M、X、Y状态也就对应着碱基匹配、碱基缺失以及碱基插入。参考图6，其中转移概率δ、ε分别表示出现Indel和维持Indel的概率，而每种隐藏状态对应单倍型的输出概率则由测序对齐的质量分数给出。利用维比特算法，我们得到如下所示的迭代公式：
Figure BDA0003595778510000081
Figure BDA0003595778510000082
Figure BDA0003595778510000083
通过上述迭代，我们可以得到一个给定Reads深度行为，以单倍型为列的概率矩阵:
Figure BDA0003595778510000091
此时，给定Reads深度行为的单个位点等位基因的似然值就由包含该等位基因的单倍型对应的最大概率给出，如果等位基因特异性拷贝数估计值的似然得分低于阈值，则从所述数据集中去除等位基因特异性拷贝数估计值，得最佳单倍体。
分别将潜在变异区域内的测序序列片段与最佳单倍体进行比对后得到中间比对结果，最佳单倍体与参考序列片段进行比对后的最终比对结果，在测序序列片段与最佳单倍体不匹配的区域添加第一空白区域，在最佳单倍体与参考序列片段不匹配的区域添加第二空白区域，用于使测序序列片段、精选单倍体和参考序列的匹配区域的位置完全对应。
根据中间比对结果与最终比对结果之间的映射关系，获得每个测序序列片段与参考序列片段的校正对齐结果。对于映射关系的设定采用预先训练的神经网络模型对数据映射关系进行处理,基于基因组的特征、测序序列的整体信息、测序序列所属物种基因组的多态发生情况、全体测序序列的长度分布情况和总体质量值分布情况，设计映射参数或训练参数，所述测序序列所属物种基因组的多态发生情况包括SNP、插入/删除发生的频率和插入/删除的长度分布情况。
训练集以更新后的样本数据将样本数据随机划分为训练集和预测集，对训练集和预测集基于采用径向基核函数的支持向量机模型进行多次训练，根据预测集的结果是否正确分别调整训练参数。正向映射之后的输出值要与目标值进行比较，得到预测误差之后，将误差逆映射至前面各个参数，以更新参数，直到参数接近训练集目标，达到收敛。
利用训练的神经网络模型最大的优势是可以和任何一种检测方法相结合，该检测过程可以从机器学习过程中自动的获得权值参数，而不需要人为的设定，所以避免了人为的错误，提高了检测的精度。
最后，根据每个测序序列片段与所述参考序列片段的校正对齐结果，确定断点。
对于测序序列片段中的每个断点，分别从断点的两侧选择n个窗口(例如100个窗口)作为两个群体用于统计检验，通过计算位点两侧之间的reads数量信息的差异获得对应于每个断点的显著值，具有最小显著值的m个断点(例如3000个位点)作为候选断点。
筛选断点：所有排列的候选断点记录为BC＝{b1，b2，...bs}，每个候选断点存在于两个连续片段之间，其中这两个片段是分别从前一个断点到所述断点和从所述断点到下一个断点的区域，移除具有最小显著值的候选断点。
从每轮循环和迭代中移除具有最小显著值的候选断点，并且重新计算两个相邻断点的显著值，循环迭代，直至所有候选断点的显著值都小于终止阈值。终止阈值通常预先设定。例如，通过对由正常样本组成对照样本集进行分析处理获得该终止阈值。
在已获得筛选的断点之后，确定拷贝数变异区域。设两个断点之间的窗口为，(xL，xR),待测试样品相对于正常样品的拷贝数变异比例按照下式计算：
Figure BDA0003595778510000101
其中与参考序列比对的正常样品的Read和待测样品的Read的数量分别是
aN和aT，并且落入窗口(xL，xR)内的Read的数量分别是N(xL，xR)和T(xL，xR),将所述拷贝数变异比例≤0.75和≥1.25分别作为染色体片段缺失和重复的检测阈值。
另一方面，如图2所示，本发明还提出了用于实现上述基因组变异检测方法的检测系统200，该检测系统200包括：
序列比对单元10，用于将每个窗口内的测序序列分别和参考序列进行序列比对；
数据处理单元20，用于对测序序列进行reads深度信息的标准化处理；具体采用基于中位数归一化的方法对GC含量和Mappability值进行校正；
确定单元30，用于根据所述序列比对结果以及深度信息，确定所述基因组的潜在变异区域；
匹配单元40，用于对潜在变异区域标准化后的测序序列应用De-Bruijn图和哈希表，使用隐马氏匹配算法依次将各测序序列的与参考序列进行匹配，得到N个最佳单倍体；
映射单元50，用于将潜在变异区域内的测序序列与N个最佳单倍体进行比对后得到中间比对结果，将N个最佳单倍体与参考序列片段进行比对后的最终比对结果；根据中间比对结果与最终比对结果之间的映射关系，获得每个测序序列与参考序列片段的校正对齐结果；具体采用预先训练的神经网络模型对数据映射关系进行处理,采用径向基核函数的支持向量机模型进行多次训练，正向映射之后的输出值与目标值进行比较，得到预测误差之后，将误差逆映射至各个映射参数，直到映射参数接近训练集目标。
断点检测单元60，用于根据每个测序序列与所述参考序列的校正对齐结果，确定断点；
变异检测单元70，用于确定拷贝数变异区域，计算拷贝数变异区域内的拷贝数变异比例，并拷贝数变异比例与阈值进行比较，得到准确的拷贝数变异检测结果。
所述变异检测单元通过对断点进行筛选，确定拷贝数变异区域，所述对断点进行筛选具体包括：将显著值最小的m个断点作为候选断点进行排列，从每轮循环和迭代中移除具有最小显著值的候选断点，并且重新计算两个相邻断点的显著值，循环迭代，直至所有候选断点的显著值都小于终止阈值。
所属领域的技术人员可以清楚地了解到，为描述的方便和简洁，上述描述的单元和系统的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。
在本申请所提供的几个实施例中，应该理解到，所揭露的系统，设备和方法，可以通过其它的方式实现。例如，以上所描述的设备实施例仅仅是示意性的，例如，所述单元的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。
所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本实施例方案的目的。
最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (7)

1.一种基因组变异检测方法，其特征在于，包括如下步骤：
步骤1，将样本的基因组划分成长度相等的连续非重叠窗口，将每个窗口内的参考序列分别和测序序列进行序列比对；
步骤2，对测序序列进行Reads深度信息的标准化处理；
步骤3，根据序列比对结果以及深度信息，确定所述基因组的潜在变异区域；
步骤4，对潜在变异区域的标准化后的测序序列应用De-Bruijn图和哈希表，使用隐马氏匹配算法依次将测序序列与参考序列进行匹配，得到N个最佳单倍体；
步骤5，分别将潜在变异区域内的测序序列与N个最佳单倍体进行比对后得到中间比对结果，N个最佳单倍体与参考序列进行比对后的最终比对结果；
在测序序列与N个最佳单倍体进行比对时的不匹配的区域添加第一空白区域，在N个最佳单倍与参考序列进行比对时的不匹配的区域添加第二空白区域，用于使测序序列、N个最佳单倍体和参考序列的匹配区域的位置完全对应；
步骤6，根据中间比对结果与最终比对结果之间的映射关系，获得每个测序序列与参考序列的校正对齐结果；
采用预先训练的神经网络模型对数据映射关系进行处理,基于基因组的特征、测序序列的整体信息、测序序列所属物种基因组的多态发生情况、全体测序序列的长度分布情况和总体质量值分布情况，设计训练参数；
将样本数据随机划分为训练集和预测集，对训练集和预测集基于采用径向基核函数的支持向量机模型进行多次训练，根据预测集的结果是否正确分别调整训练参数；
正向映射之后的输出值与目标值进行比较，得到预测误差之后，将误差逆映射至各个映射参数，以更新参数，直到参数接近训练目标，达到收敛；
步骤7，根据每个测序序列与所述参考序列的校正对齐结果，确定断点；
步骤8，对断点进行筛选，确定拷贝数变异区域，利用深度信息计算拷贝数变异区域内的拷贝数变异比例，并拷贝数变异比例与阈值进行比较，得到准确的拷贝数变异检测结果。
2.根据权利要求1所述的基因组变异检测方法，其特征在于，所述步骤2具体包括：采用基于中位数归一化的方法对GC含量和Mappability值进行校正。
3.根据权利要求1所述的基因组变异检测方法，其特征在于，所述步骤8中:对断点进行筛选具体包括：将显著值最小的m个断点作为候选断点进行排列，从每轮循环和迭代中移除具有最小显著值的候选断点，并且重新计算两个相邻断点的显著值，循环迭代，直至所有候选断点的显著值都小于终止阈值。
4.根据权利要求1所述的基因组变异检测方法，其特征在于，所述步骤8中:在已获得筛选的断点之后，确定拷贝数变异区域；设两个断点之间的窗口为
Figure QLYQS_1
,待测试样品相对于正常样品的拷贝数变异比例按照下式计算：/>
Figure QLYQS_2
其中与参考序列比对的正常样品的Read和待测样品的Read的数量分别是
Figure QLYQS_3
和/>
Figure QLYQS_4
，并且落入窗口/>
Figure QLYQS_5
内的Read的数量分别是/>
Figure QLYQS_6
和/>
Figure QLYQS_7
,将所述拷贝数变异比例≤0.75和≥1.25分别作为染色体片段缺失和重复的检测阈值。
5.一种基因组变异检测系统，用于实现如权利要求1-4中任意一项所述的基因组变异检测方法，其特征在于，该检测系统包括：
序列比对单元，用于将每个窗口内的测序序列分别和参考序列进行序列比对；
数据处理单元，用于对测序序列进行reads深度信息的标准化处理；
确定单元，用于根据序列比对结果以及深度信息，确定所述基因组的潜在变异区域；
匹配单元，用于对潜在变异区域标准化后的测序序列应用De-Bruijn图和哈希表，使用隐马氏匹配算法依次将各测序序列的与参考序列进行匹配，得到N个最佳单倍体；
映射单元，用于将潜在变异区域内的测序序列与N个最佳单倍体进行比对后得到中间比对结果，将N个最佳单倍体与参考序列片段进行比对后的最终比对结果；
在测序序列与N个最佳单倍体进行比对时的不匹配的区域添加第一空白区域，在N个最佳单倍与参考序列进行比对时的不匹配的区域添加第二空白区域，用于使测序序列、N个最佳单倍体和参考序列的匹配区域的位置完全对应；
根据中间比对结果与最终比对结果之间的映射关系，获得每个测序序列与参考序列片段的校正对齐结果；
采用预先训练的神经网络模型对数据映射关系进行处理,基于基因组的特征、测序序列的整体信息、测序序列所属物种基因组的多态发生情况、全体测序序列的长度分布情况和总体质量值分布情况，设计训练参数；
将样本数据随机划分为训练集和预测集，对训练集和预测集基于采用径向基核函数的支持向量机模型进行多次训练，根据预测集的结果是否正确分别调整训练参数；
正向映射之后的输出值与目标值进行比较，得到预测误差之后，将误差逆映射至各个映射参数，以更新参数，直到参数接近训练目标，达到收敛；
断点检测单元，用于根据每个测序序列与所述参考序列的校正对齐结果，确定断点；
变异检测单元，用于确定拷贝数变异区域，利用深度信息计算拷贝数变异区域内的拷贝数变异比例，并拷贝数变异比例与阈值进行比较，得到准确的拷贝数变异检测结果。
6.根据权利要求5所述的基因组变异检测系统，其特征在于，所述数据处理单元，具体采用基于中位数归一化的方法对GC含量和Mappability值进行校正。
7.根据权利要求5所述的基因组变异检测系统，其特征在于，所述变异检测单元通过对断点进行筛选，确定拷贝数变异区域，所述对断点进行筛选具体包括：将显著值最小的m个断点作为候选断点进行排列，从每轮循环和迭代中移除具有最小显著值的候选断点，并且重新计算两个相邻断点的显著值，循环迭代，直至所有候选断点的显著值都小于终止阈值。
CN202210391526.1A 2022-04-14 2022-04-14 一种基因组变异检测方法及检测系统 Active CN114999573B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210391526.1A CN114999573B (zh) 2022-04-14 2022-04-14 一种基因组变异检测方法及检测系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210391526.1A CN114999573B (zh) 2022-04-14 2022-04-14 一种基因组变异检测方法及检测系统

Publications (2)

Publication Number Publication Date
CN114999573A CN114999573A (zh) 2022-09-02
CN114999573B true CN114999573B (zh) 2023-07-07

Family

ID=83024054

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210391526.1A Active CN114999573B (zh) 2022-04-14 2022-04-14 一种基因组变异检测方法及检测系统

Country Status (1)

Country Link
CN (1) CN114999573B (zh)

Families Citing this family (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115641911B (zh) * 2022-10-19 2023-05-23 哈尔滨工业大学 一种针对序列间重叠检测的方法
CN115602244B (zh) * 2022-10-24 2023-04-28 哈尔滨工业大学 一种基于序列比对骨架的基因组变异检测方法
CN115375682B (zh) * 2022-10-24 2023-01-20 湖南自兴智慧医疗科技有限公司 一种染色体罗氏易位异常检测方法、系统及存储介质
CN115631789B (zh) * 2022-10-25 2023-08-15 哈尔滨工业大学 一种基于泛基因组的群体联合变异检测方法
CN115602246B (zh) * 2022-10-31 2023-06-20 哈尔滨工业大学 一种基于群体基因组的序列比对方法
CN115910199B (zh) * 2022-11-01 2023-07-14 哈尔滨工业大学 一种基于比对框架的三代测序数据结构变异检测方法
CN115762633B (zh) * 2022-11-23 2024-01-23 哈尔滨工业大学 一种基于三代测序的基因组结构变异基因型校正方法

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105392894B (zh) * 2012-01-20 2018-05-29 深圳华大基因股份有限公司 确定样本基因组中是否存在拷贝数变异的方法、系统和计算机可读介质
AU2012376134B2 (en) * 2012-04-05 2016-03-03 Bgi Genomics Co., Ltd. Method and system for detecting copy number variation
CN108121897B (zh) * 2016-11-29 2020-05-08 华为技术有限公司 一种基因组变异检测方法及检测装置
CN111916150A (zh) * 2019-05-10 2020-11-10 北京贝瑞和康生物技术有限公司 一种基因组拷贝数变异的检测方法和装置
JP2023521893A (ja) * 2020-04-15 2023-05-25 チャン ザッカーバーグ バイオハブ， インコーポレイテッド 機械学習モデルを用いた局所的祖先推論
CN113270141B (zh) * 2021-06-10 2023-02-21 哈尔滨因极科技有限公司 一种基因组拷贝数变异检测整合算法
CN113555062B (zh) * 2021-07-23 2022-07-12 哈尔滨因极科技有限公司 一种用于基因组碱基变异检测的数据分析系统及分析方法

Also Published As

Publication number Publication date
CN114999573A (zh) 2022-09-02

Similar Documents

Publication Publication Date Title
CN114999573B (zh) 一种基因组变异检测方法及检测系统
Liu et al. Interrogating the “unsequenceable” genomic trinucleotide repeat disorders by long-read sequencing
US8725422B2 (en) Methods for estimating genome-wide copy number variations
CN108256289B (zh) 一种基于目标区域捕获测序基因组拷贝数变异的方法
CN111081315B (zh) 一种同源假基因变异检测的方法
US20130316915A1 (en) Methods for determining absolute genome-wide copy number variations of complex tumors
CN115631789B (zh) 一种基于泛基因组的群体联合变异检测方法
CN113249453B (zh) 一种检测拷贝数变化的方法
CN108137642A (zh) 分子质量保证方法在测序中的应用
CN108595912B (zh) 检测染色体非整倍性的方法、装置及系统
Ahsan et al. A survey of algorithms for the detection of genomic structural variants from long-read sequencing data
CN111276189B (zh) 基于ngs的染色体平衡易位检测分析系统及应用
CN115394359B (zh) 一种通过转录组检测单细胞染色体拷贝数变异方法
CN108694304B (zh) 一种身份关系鉴定方法、装置、设备及存储介质
WO2019132010A1 (ja) 塩基配列における塩基種を推定する方法、装置及びプログラム
CN116246703A (zh) 一种核酸测序数据的质量评估方法
CN116434843A (zh) 一种碱基测序质量评估方法
Roy et al. NGS-μsat: bioinformatics framework supporting high throughput microsatellite genotyping from next generation sequencing platforms
CN114420214A (zh) 核酸测序数据的质量评估方法和筛选方法
US20220364080A1 (en) Methods for dna library generation to facilitate the detection and reporting of low frequency variants
Schwender Statistical analysis of genotype and gene expression data
JP7166638B2 (ja) 多型検出法
CN117672354B (zh) 比较哺乳动物近源物种完整基因组组装质量的方法和装置
Prodanov Read Mapping, Variant Calling, and Copy Number Variation Detection in Segmental Duplications
CN112599189B (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