CN112687341B - 一种以断点为中心的染色体结构变异鉴定方法 - Google Patents
一种以断点为中心的染色体结构变异鉴定方法 Download PDFInfo
- Publication number
- CN112687341B CN112687341B CN202110268544.6A CN202110268544A CN112687341B CN 112687341 B CN112687341 B CN 112687341B CN 202110268544 A CN202110268544 A CN 202110268544A CN 112687341 B CN112687341 B CN 112687341B
- Authority
- CN
- China
- Prior art keywords
- breakpoint
- reads
- read
- structural variation
- positions
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 71
- 210000001726 chromosome structure Anatomy 0.000 title description 5
- 210000000349 chromosome Anatomy 0.000 claims abstract description 45
- 238000004590 computer program Methods 0.000 claims abstract description 12
- 230000002759 chromosomal effect Effects 0.000 claims abstract description 8
- 238000012163 sequencing technique Methods 0.000 claims description 43
- 108091028043 Nucleic acid sequence Proteins 0.000 claims description 18
- 150000007523 nucleic acids Chemical group 0.000 claims description 16
- 108020004414 DNA Proteins 0.000 claims description 15
- 238000001914 filtration Methods 0.000 claims description 15
- 108090000623 proteins and genes Proteins 0.000 claims description 10
- 238000012217 deletion Methods 0.000 claims description 5
- 230000037430 deletion Effects 0.000 claims description 5
- 230000005945 translocation Effects 0.000 claims description 5
- 238000003780 insertion Methods 0.000 claims description 4
- 230000037431 insertion Effects 0.000 claims description 4
- 230000035772 mutation Effects 0.000 claims description 3
- 238000012545 processing Methods 0.000 claims description 2
- 230000000875 corresponding effect Effects 0.000 description 13
- 238000001514 detection method Methods 0.000 description 13
- 230000004927 fusion Effects 0.000 description 10
- 108091032973 (ribonucleotides)n+m Proteins 0.000 description 8
- 206010028980 Neoplasm Diseases 0.000 description 8
- 230000035945 sensitivity Effects 0.000 description 8
- 201000011510 cancer Diseases 0.000 description 4
- 238000010276 construction Methods 0.000 description 3
- 201000010099 disease Diseases 0.000 description 3
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 3
- 239000000284 extract Substances 0.000 description 3
- 108091008794 FGF receptors Proteins 0.000 description 2
- 102100023600 Fibroblast growth factor receptor 2 Human genes 0.000 description 2
- 101710182389 Fibroblast growth factor receptor 2 Proteins 0.000 description 2
- 238000012896 Statistical algorithm Methods 0.000 description 2
- 230000004931 aggregating effect Effects 0.000 description 2
- 238000004422 calculation algorithm Methods 0.000 description 2
- 238000012733 comparative method Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 102000052178 fibroblast growth factor receptor activity proteins Human genes 0.000 description 2
- PCHJSUWPFVWCPO-UHFFFAOYSA-N gold Chemical compound [Au] PCHJSUWPFVWCPO-UHFFFAOYSA-N 0.000 description 2
- 108020004707 nucleic acids Proteins 0.000 description 2
- 102000039446 nucleic acids Human genes 0.000 description 2
- 239000002773 nucleotide Substances 0.000 description 2
- 125000003729 nucleotide group Chemical group 0.000 description 2
- 238000003745 diagnosis Methods 0.000 description 1
- 238000002405 diagnostic procedure Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 229940079593 drug Drugs 0.000 description 1
- 239000003814 drug Substances 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 102000052116 epidermal growth factor receptor activity proteins Human genes 0.000 description 1
- 108700015053 epidermal growth factor receptor activity proteins Proteins 0.000 description 1
- 239000012634 fragment Substances 0.000 description 1
- 230000014509 gene expression Effects 0.000 description 1
- 230000003862 health status Effects 0.000 description 1
- 108091033319 polynucleotide Proteins 0.000 description 1
- 102000040430 polynucleotide Human genes 0.000 description 1
- 239000002157 polynucleotide Substances 0.000 description 1
- 230000003252 repetitive effect Effects 0.000 description 1
- 238000007480 sanger sequencing Methods 0.000 description 1
- 241000894007 species Species 0.000 description 1
- 210000001519 tissue Anatomy 0.000 description 1
- 210000004881 tumor cell Anatomy 0.000 description 1
- 238000012049 whole transcriptome sequencing Methods 0.000 description 1
Images
Landscapes
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明涉及一种以断点为中心的染色体结构变异鉴定方法。具体而言,本发明涉及一种可由计算机实施的以断点为中心的鉴定染色体结构变异的方法。本发明还涉及用于实施所述方法的计算机系统或装置、以及存储有能够实施所述方法的计算机程序的计算机可读介质。
Description
技术领域
本发明涉及核酸测序领域。本发明涉及一种以断点为中心的染色体结构变异鉴定方法。具体而言,本发明涉及一种可由计算机实施的以断点为中心的鉴定染色体结构变异的方法。本发明还涉及用于实施所述方法的计算机系统或装置、以及存储有能够实施所述方法的计算机程序的计算机可读介质。
背景技术
染色体结构变异是指染色体改变其原有的基因组序列和结构,如发生插入、缺失、倒置、易位等变异,从而产生新的染色体结构。染色体结构变异在肿瘤细胞中经常出现,如EML4-ALK基因融合现象、EGFR激酶结构域重复等,检测染色体结构变异对于指导肿瘤治疗以及用药具有重要临床意义。
染色体结构变异检测的难点在于准确鉴定结构变异断点产生的位置以及方向,进而才能准确分析结构变异的类型、可靠性、功能性以及变异丰度等特征。
检测染色体结构变异断点的难点之一在于目前常用的二代测序技术的读长有限,而断点位置在建库打断的DNA片段中随机分布,使得大部分测序读段都不包含有断点信息或仅包含有限的跨断点的序列,从而导致难以从单一读段中获取染色体结构变异两端的断点信息。
检测染色体结构变异断点的难点之二在于染色体结构变异常发生在基因组的重复序列区域,导致测序读段在断点附近的比对质量较低或测序错误率较高,从而增加了断点检测的难度。
检测染色体结构变异断点的难点之三在于对于低频染色体结构变异,测序产生的数据中包含的支持结构变异的读段数目少,其包含断点的读段更加有限,从而增加了断点检测的难度。
当前染色体结构变异断点的鉴定方法主要有两种,其一是对测序产生的原始读段先组装成更长的序列,进而再比对到基因组上获得准确的断点位置,其二是将测序产生的原始读段直接比对到基因组上,进而从具有二级比对的读段中提取断点位置信息。前者的难点在于组装过程资源耗费多、组装结果准确度低,从而导致检测灵敏性与特异性都较低,后者的难点在于若无具有二级比对结果的读段将无法鉴定结构变异的具体断点,以及该种方法存在遗漏部分支持结构变异的读段,以及由于建库过程中不可避免的会产生大量具有二级比对结果的假阳性读段等,从而导致该种方法的检测灵敏性与特异性依然存在不足。
为解决前述方法中存在的可能遗失支持结构变异的读段、对低频结构变异检测灵敏度低、检出结果中特异性低等问题,本发明提供了一种以结构变异断点为中心的检测方法,可以显著提高基因组结构变异检测性能。
发明内容
本发明通过提供一种以断点为中心的结构变异鉴定算法解决了上述技术问题。具体而言,本发明方法首先提取了所有潜在的结构变异支持读段,避免了其它方法中存在的支持读段丢失的可能;其次,本发明独特的以断点为中心的结构变异支持读段统计算法统计了所有支持该断点的读段同时又避免了假阳性读段,并充分利用了读段的位置信息、正负链信息、软剪切的序列和方向信息,从而保证了本发明敏感性和特异性。
在一个方面,本发明涉及一种可由计算机实施的以断点为中心的鉴定染色体结构变异的方法,所述方法按顺序包括以下步骤:
(1) 将样本DNA文库经双端测序后所产生的成对读段比对到参考基因组上,获得读段比对结果;
(2) 提取读段比对结果中的比对位置在目标区域内且支持结构变异的所有成对读段,其中提取的成对读段依据其包含的两个读段比对位置是否跨越断点来划分读段类型,其包括跨越断点的读段以及未跨越断点的读段,并分别依据读段是否包含软剪切序列提取可能的断点精确位置;其中对于跨越断点的读段,要求两个成对读段或/和其二级比对结果比对在不同染色体上或者比对在相同染色体上且二者之间的距离大于一定阈值;其中对于未跨越断点的读段,要求两个成对读段比对在相同染色体上且二者之间的距离小于一定阈值且包含软剪切序列;
(3) 将跨越断点的读段利用成对读段的比对位置进行聚类,聚类时要求成对读段与其它成对读段对应的比对位置均在一定距离范围内,聚类后将所有跨越断点的读段分配到不同的子类中,每一个子类包含了所有支持基因组上对应两个区域间结构变异的所有跨越断点的读段;
(4) 以从步骤(3)中产生的每一个读段子类为单位鉴定其中可能的结构变异断点位置,鉴定结构变异断点位置的规则为:
I) 在子类中存在具有二级比对结果的读段的情况下,则该读段直接提供了一种可能的结构变异位置,即第I)类断点;
II) 在子类中不存在第I)类断点但存在两个读段比对位置都跨越断点的成对读段,且两个读段都不存在二级比对结果,则可通过验证两个读段末端是否存在重叠区域判断该断点的可靠性,以及软剪切序列是否存在于另一个读段的末端中从而判断断点间是否存在其它插入片段,即第II)类断点;
III) 在子类中不存在第I)和II)类断点但在断点两边同时存在不同成对读段的单个读段包含软剪切序列,则可将来自不同成对读段的包含软剪切序列的单个读段视为来源于同一成对读段,按上述第II)类断点处理方式鉴定断点,即第III)类断点;
IV) 在子类中不存在第I)至III)类断点但在断点单边存在单个读段包含软剪切序列,则此时将鉴定结构变异单端断点,并认为另一端断点具体位置不可知,即第IV)类断点;
V) 在子类中不存在第I)至IV)类断点但存在成对读段比对位置跨越断点,则此时认为结构变异两端断点具体位置均不可知,即第V)类断点;
(5) 获得结构变异断点位置后,以断点为中心统计支持该结构变异的读段数目以及类型,
其中在断点两端都存在精确的断点位置的情况下,用该断点过滤子类中所有成对读段,并根据读段类型统计支持该结构变异的读段数目以及类型,
其中在仅有单端断点位置确定的情况下,仅对断点已知的一端根据成对读段类型过滤子类中读段,断点未知的一端仅对正负链的信息过滤,和
其中在断点两端位置都未知的情况下,读段两端均只对正负链信息过滤;和
(6) 对子类中每个可能的结构变异断点进行读段数目过滤,并将满足过滤条件的结构变异断点比对到基因上,并推断结构变异产生的结果。
在一个方面,在步骤(1)中将样本DNA文库经二代测序后所产生的读段比对到人类参考基因组上,获得包含读段比对结果的BAM文件,其中一种可用的比对软件为BWA,和一种可用的人类参考基因组版本为hg19;和在步骤(2)中提取BAM文件中的比对位置在目标区域内且支持结构变异的所有读段。
在一个方面,读段类型包括:读段类型1):两个读段都不包含断点;读段类型2):仅一个读段包含断点,其跨越断点的长度不足以产生二级比对结果;读段类型3):两个读段都包含断点,但其跨越断点的长度均不足以产生二级比对结果;读段类型4):至少一个读段包含断点且其跨越断点的长度足以产生二级比对结果,此时对另一个读段的比对位置及是否包含断点并无限定;和读段类型5):至少一个读段包含断点,且其跨越断点的长度不足以产生二级比对结果;其中,所述5种类型中,前3种类型的两个读段比对位置都跨越断点,第4种类型的一个读段的一级比对结果与其二级比对结果的比对位置跨越断点,对其中另一个读段无限定,第5种类型的两个读段都比对在断点的同一侧;其中,两个成对读段比对位置跨越断点是指成对读段的两端的一级比对位置分布在断点两侧或单个读段的一级比对结果与其二级比对结果的比对位置分布在断点两侧;其中,读段包含断点是指断点位置在该读段的内部,读段的一部分序列跨越断点,在比对后将形成软剪切序列。
在一个方面,在步骤(2)中,要求跨越断点的成对读段所包含的两个读段比对在相同染色体上时二者间比对位置之间的距离的阈值是500bp至500Kb、特别是1Kb至100Kb、更特别是5Kb至50Kb,例如10Kb;和/或要求非跨越断点的成对读段所包含的两个读段比对在相同染色体上时二者间比对位置之间的距离的阈值是50bp至5Kb、特别是100bp至2Kb、更特别是200bp至1Kb,例如600bp。
在一个方面,在步骤(3)中,聚类时要求成对读段与其它成对读段对应的比对位置间的距离为小于50bp至10Kb的范围、特别是小于100bp至5Kb的范围、更特别是小于200bp至1Kb的范围,例如小于500bp的范围。
在一个方面,在步骤(6)中,过滤方式为要求总的特异的读段数目大于1-20、特别是大于2-10,例如大于3,且存在至少两种及以上的支持所述断点的读段类型。
在一个方面,在步骤(1)中,用比对软件将样本核酸序列经双端测序生成的读段比对到对应参考基因组产生的BAM文件,以该结果作为输入文件。
在一个方面,在步骤(6)中,依据以下方式推断结构变异产生的结果:当跨断点的成对读段的两个读段分别比对在相同染色体的正链、负链上,且正链的比对位置小于负链的比对位置,则表征该读段支持结构变异缺失类型;当跨断点的成对读段的两个读段分别比对在相同染色体的正链、负链上,且正链的比对位置大于负链的比对位置,则表征该读段支持结构变异重复类型;当跨断点的成对读段的两个读段都比对在相同染色体的正链或负链上,则表征该读段支持结构变异倒置类型;和当跨断点的成对读段的两个读段分别比对在不同染色体上,则表征该读段支持结构变异易位类型。在该方面中,推断结构变异产生的结果如图2所示。
在一个方面,步骤(1)还包括从基因测序装置获得样本核酸序列经双端测序生成的读段,并且将读段信息经过输入设备输入到计算机中。
在一个方面,在步骤(1)中计算机处理器接收并处理输入的读段信息,从而获得读段的比对结果;在步骤(2)中计算机处理器提取样本核酸序列中的目标区域支持结构变异的所有读段,划分读段类型,并提取可能的断点精确位置;在步骤(3)中计算机处理器将跨越断点的读段利用成对读段两端的比对位置进行聚类;在步骤(4)中计算机处理器以从步骤(3)中产生的每一个读段子类为单位鉴定其中可能的结构变异断点位置;在步骤(5)中计算机处理器以断点为中心统计支持结构变异的读段数目以及类型;和/或在步骤(6)中计算机处理器对子类中每个可能的结构变异断点进行读段数目过滤,并将满足过滤条件的结构变异断点比对到基因上,并推断结构变异产生的结果。
在一个方面,在步骤(5)中,还以断点为中心统计支持该结构变异的两个断点间缺失或插入的序列长度和/或断点附近的序列复杂度的信息。
在一个方面,计算机处理器将步骤(6)中的结构变异产生的结果传输并显示到输出设备中。
在一个方面,本发明涉及用于实施本发明方法的计算机系统,其包括:输入设备,用于输入样本核酸序列的测序读段信息;计算机存储器,用于存储计算机程序指令;计算机处理器,用于执行所述计算机程序指令,其中所述计算机程序指令实施步骤(1)至(6),对样本核酸序列的测序读段信息进行处理并生成结构变异产生的结果,并将结果传输到输出设备;和输出设备,用于显示所述结果。
在一个方面,本发明涉及一种计算机可读介质,其中所述计算机可读介质存储有计算机程序,其中所述计算机程序能被计算机处理器执行以实施本发明的方法。
在一个方面,本发明涉及用于实施本发明方法的装置,其包括:基因测序数据输入模块,用于输入样本核酸序列的测序读段信息;获取测序读段比对结果模块,用于实施步骤(1);提取目标区域结构变异支持读段模块,用于实施步骤(2);跨越断点读段聚类模块,用于实施步骤(3);从读段子类鉴定结构变异断点模块,用于实施步骤(4);以断点为中心统计支持结构变异读段数模块,用于实施步骤(5);结构变异过滤模块,用于实施步骤(6);和显示结果模块,用于显示从步骤(6)获得的结果。
在一个方面,本发明方法可以不是疾病诊断方法。具体而言,本发明方法的直接目的可以不是获得诊断结果或健康状况,而只是从患者或受试者获取作为中间结果的信息,或处理该信息。在一个方面,所述中间结果信息可为染色体结构变异的结果。在另一方面,所述中间结果信息可为图2中所示的染色体结构变异类型,包括缺失、重复、倒置和易位。在一个方面,根据现有技术中的医学知识和本申请公开的内容,医师在获得和知晓所述中间结果信息,不一定能获得直接得出疾病的诊断结果或健康状况。例如,医师可通过本发明方法获得受试者的染色体结构变异的结果,但是该受试者自身不一定已经罹患癌症或肿瘤。又例如,医师可通过本发明方法获得受试者的染色体结构变异的结果,但不会直接判断出患者已经罹患癌症或肿瘤。
在一个方面,患者或受试者是人。
附图说明
图1显示了支持基因组结构变异的读段类型。
图2是不同结构变异类型及跨越断点的读段方向示意图。
图3是以断点为中心检测结构变异流程图。
图4显示了从读段子类中鉴定结构变异断点位置的断点类型。
具体实施方式
术语"测序"通常是指确定多核苷酸分子的核苷酸序列的过程。通常,对基因组进行测序需要首先用常规手段将基因组序列进行片段化,以构建核酸片段文库。
术语"读段"通常是指测序过程的输出,即已通过任何测序方法产生的核酸的核苷酸或碱基序列信息。
术语"成对读段"通常是指双端测序所产生的两个配对读段。
术语“读长”通常是指测序能够准确读取的核酸序列的最大长度。
术语“比对”通常是指将目标序列(例如读段)与参考序列比较,从而确定参考序列是否含有目标序列(例如读段)的过程。
术语“软剪切序列”通常是指:在读段进行比对后,若存在部分序列比对到参考序列某位置,另一部分比对到参考序列另一位置或不能比对到参考序列,则该读段被称为软剪切序列。
术语“聚类”通常是指:将比对在染色体上的不同读段依据其比对的位置之间的距离聚集为不同的读段类的过程,即将读段两端比对位置的距离少于给定阈值的读段聚成同一类。
术语“断点”通常是指:染色体上发生结构变异的染色体位置。
术语“提取”或“提取读段”通常是指:从比对文件如BAM文件中获得相关读段的过程。
从基因组结构变异产生的生物学原理以及二代测序建库测序的原理来看,如不考虑其它复杂情况,支持结构变异的成对读段在测序比对结果中可以分为以下几种类型(图1):1)两个读段都不包含断点(图1中读段类型1);2)仅一个读段包含断点,其跨越断点的长度不足以产生二级比对结果(图1中读段类型2);3)两个读段都包含断点,但其跨越断点的长度均不足以产生二级比对结果(图1中读段类型3);4)至少一个读段包含断点且其跨越断点的长度足以产生二级比对结果,此时对另一个读段的比对位置及是否包含断点并无限定(图1中读段类型4);5)至少一个读段包含断点,且其跨越断点的长度不足以产生二级比对结果(图1中读段类型5)。上述5种类型中,前3种的两个读段比对位置都跨越断点,第4种的一个读段的一级比对结果与其二级比对结果的比对位置跨越断点,对其中另一个读段无限定,第5种的两个读段都比对在断点的同一侧。
从基因组结构变异产生的生物学原理以及二代测序建库测序的原理来看,如不考虑其它复杂情况,支持结构变异的成对读段比对位置如跨越断点,即两个读段分别比对在断点的两端,且比对位置间的距离远大于建库时DNA模板打断长度(如大于10Kb)或比对在不同染色体上,此时单个读段比对的正负链信息可表征结构变异中断点的方向,结合成对读段的比对位置,可表征结构变异的类型。
如图2所示,当跨断点的成对读段的两个读段分别比对在相同染色体的正链、负链上,且正链的比对位置小于负链的比对位置,则表征该读段支持结构变异缺失类型;当跨断点的成对读段的两个读段分别比对在相同染色体的正链、负链上,且正链的比对位置大于负链的比对位置,则表征该读段支持结构变异重复类型;当跨断点的成对读段的两个读段都比对在相同染色体的正链或负链上,则表征该读段支持结构变异倒置类型;当跨断点的成对读段的两个读段分别比对在不同染色体上,则表征该读段支持结构变异易位类型。
本发明方法的详细流程和步骤如图3所示。
S01 获取样本测序读段的比对结果,一般为双端测序生成的读段经由比对软件(如bwa)比对到对应参考基因组产生的BAM文件,以该结果作为本发明的输入文件。
S02 提取目标区域支持结构变异的所有读段,读段的类型为图1中描述的5种类型,其中,对于比对结果中两个读段比对位置跨越结构变异断点的读段,要求两个读段比对在不同染色体上或者比对在相同染色体上且二者比对位置之间的距离大于一定阈值(如10Kb);对于比对结果中两个读段比对位置非跨越断点的成对读段,要求两个读段比对在相同染色体上且二者比对位置之间的距离小于一定阈值(如600bp)。提取的读段依据是否包含软剪切序列提取可能的断点精确位置。
S03 对于比对结果中两个读段比对位置跨越结构变异断点的读段,利用两个读段的比对位置进行聚类,聚类时要求两个读段与其它成对读段对应的比对位置均在一定距离范围内(如都小于500bp)。聚类后,所有比对位置跨越断点的读段将分配到不同的子类中,每一个子类包含了所有支持该结构变异的所有跨越断点的读段。
S04 以S03产生的每一个读段子类为单位鉴定其中可能的结构变异断点位置,具体步骤如下(图4):1)如果子类中存在具有二级比对结果的读段,则该读段直接提供了一种可能的结构变异位点(对应图1中读段类型4,图4中断点类型1(或第I)类断点));2)如果子类中不存在第I类断点但存在两个读段比对位置都跨越断点的读段,且两个读段都不存在二级比对结果,则可通过验证两个读段末端是否存在重叠区域判断该断点的可靠性,以及软剪切序列是否存在于另一个读段的末端中从而判断断点间是否存在其它插入片段(对应图1中读段类型3,图4中断点类型2(或第II)类断点));3)如果子类中不存在第I、II类断点但存在不同成对读段的单个读段包含软剪切序列,则可将来自不同成对读段的包含软剪切序列的单个读段视为同一成对读段,按上述断点类型2处理鉴定断点(对应图1中读段类型2,图4中断点类型3(或第III)类断点));4)如果子类中不存在第I-III类断点但存在单个读段包含软剪切序列且不同成对读段中包含软剪切序列的单个读段无法组成新的成对读段,则此时将鉴定结构变异单端断点,并认为另一端断点具体位置不可知(对应图1中读段类型2,图4中断点类型4、5(或第IV)类断点));5)如果子类中不存在第I-IV类断点但存在成对读段比对位置跨越断点,则此时认为结构变异两端断点具体位置均不可知(对应图1中读段类型1,图4中断点类型6(或第V)类断点))。
S05 获得结构变异断点位置后,以断点为中心统计支持该结构变异的读段数目以及类型和其它特征。如结构变异的断点两端都存在精确的断点位置,则用该断点过滤子类中所有读段,并按图1中读段类型统计支持该结构变异的读段数目以及类型;如仅有单端断点位置确定,则仅对断点已知的一端按图1中读段类型过滤子类中读段,断点未知的一端仅对正负链的信息过滤;如断点两端位置都未知,则读段两端均只对正负链信息过滤。
S06 对子类中每个可能的结构变异断点进行读段数目过滤,一种可行的过滤方式为要求总的特异的读段数目大于3,且存在至少两种及以上的图1中支持该断点的读段类型。将满足过滤条件的结构变异断点比对到基因上,并依据图2推断结构变异产生的结果。
本发明通过提供一种以断点为中心的结构变异鉴定算法实现了提高的敏感性和特异性。具体而言,本发明方法首先提取了所有潜在的结构变异支持读段,避免了其它方法中存在的支持读段丢失的可能;其次,本发明独特的以断点为中心的结构变异支持读段统计算法统计了所有支持该断点的读段同时又避免了假阳性读段,并充分利用了读段的位置信息、正负链信息、软剪切的序列和方向信息,从而保证了本发明敏感性和特异性。本申请下面的实施例充分验证了本发明方法所能达到的技术效果。然而,该实施例不应理解成限制本发明实施和本发明效果的限制条件。技术人员理解,可以使用各种类型的癌症样品或者其他疾病的样品实施本发明的方法并可实现相应的本发明效果。
实施例
取实体瘤患者手术切除组织所制备的FFPE样本604例,分别提取DNA及RNA建库,DNA文库经捕获后进行二代测序,RNA文库进行全转录组测序。测序得到各样本的DNA序列及RNA序列,其中DNA序列经BWA比对到人类参考基因组版本hg19上,RNA序列经STAR比对到人类参考基因组版本hg19对应的转录本上。样本DNA序列经比对获得的BAM文件作为本发明的输入文件。任意公开发表的DNA融合检测软件可作为本发明的DNA检测结果对比方法。一种具体的可用的DNA融合检测软件如GeneFuse可作为本发明的对比方法。作为本发明检测结果的第二种对比方法,样本RNA序列经比对获得的BAM文件利用STAR-fusion检测RNA融合结果。
在所实施的604例样本中,利用本发明所述融合鉴定方法检测了与FGFR2相关的融合,利用GeneFuse作为DNA检测FGFR2相关融合的对比方法,利用RNA作为检测FGFR2相关融合的金标准,结果表明(表1),本发明所述方法相比于RNA金标准的检出灵敏度可达89%,对比方法检出灵敏度仅46%。进一步的,对于本发明所报告的3例DNA阳性融合而RNA阴性的样本,DNA融合序列经Sanger测序全部获得验证,其DNA与RNA检测结果不一致可能受基因表达水平影响。本发明所报告的FGFR2在该癌种中融合阳性率更符合公开数据报导的结果,相比于对比方法具有显著的临床意义。
以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明,不能认定本发明的具体实施方式仅限于此,对于本发明所属技术领域的普通技术人员来说,在不脱离本发明构思的前提下,还可以做出若干简单的推演或替换,都应当视为属于本发明由所提交的权利要求书确定专利保护范围。
Claims (19)
1.一种可由计算机实施的以断点为中心的鉴定染色体结构变异的方法,所述方法按顺序包括以下步骤:
(1) 将样本DNA文库经双端测序后所产生的成对读段比对到参考基因组上,获得读段比对结果;
(2) 提取读段比对结果中的比对位置在目标区域内且支持结构变异的所有成对读段,其中提取的成对读段依据其包含的两个读段比对位置是否跨越断点来划分读段类型,其包括跨越断点的读段以及未跨越断点的读段,并分别依据读段是否包含软剪切序列提取可能的断点精确位置;其中对于跨越断点的读段,要求两个成对读段或/和其二级比对结果比对在不同染色体上或者比对在相同染色体上且二者之间的距离大于一定阈值;其中对于非跨越断点的读段,要求两个成对读段比对在相同染色体上且二者之间的距离小于一定阈值且包含软剪切序列;
(3) 将跨越断点的读段利用成对读段的比对位置进行聚类,聚类时要求成对读段与其它成对读段对应的比对位置均在一定距离范围内,聚类后将所有跨越断点的读段分配到不同的子类中,每一个子类包含了所有支持基因组上对应两个区域间结构变异的所有跨越断点的读段;
(4) 以从步骤(3)中产生的每一个读段子类为单位鉴定其中可能的结构变异断点位置,鉴定结构变异断点位置的规则为:
I) 在子类中存在具有二级比对结果的读段的情况下,则该读段直接提供了一种可能的结构变异位置,即第I)类断点;
II) 在子类中不存在第I)类断点但存在两个读段比对位置都跨越断点的成对读段,且两个读段都不存在二级比对结果,则可通过验证两个读段末端是否存在重叠区域判断该断点的可靠性,以及软剪切序列是否存在于另一个读段的末端中从而判断断点间是否存在其它插入片段,即第II)类断点;
III) 在子类中不存在第I)和II)类断点但在断点两边同时存在不同成对读段的单个读段包含软剪切序列,则可将来自不同成对读段的包含软剪切序列的单个读段视为来源于同一成对读段,按上述第II)类断点处理方式鉴定断点,即第III)类断点;
IV) 在子类中不存在第I)至III)类断点但在断点单边存在单个读段包含软剪切序列,则此时将鉴定结构变异单端断点,并认为另一端断点具体位置不可知,即第IV)类断点;
V) 在子类中不存在第I)至IV)类断点但存在成对读段比对位置跨越断点,则此时认为结构变异两端断点具体位置均不可知,即第V)类断点;
(5) 获得结构变异断点位置后,以断点为中心统计支持该结构变异的读段数目以及类型,
其中在断点两端都存在精确的断点位置的情况下,用该断点过滤子类中所有成对读段,并根据读段类型统计支持该结构变异的读段数目以及类型,
其中在仅有单端断点位置确定的情况下,仅对断点已知的一端根据成对读段类型过滤子类中读段,断点未知的一端仅对正负链的信息过滤,和
其中在断点两端位置都未知的情况下,读段两端均只对正负链信息过滤;和
(6) 对子类中每个可能的结构变异断点进行读段数目过滤,并将满足过滤条件的结构变异断点比对到基因上,并推断结构变异产生的结果。
2.如权利要求1所述的方法,其中在步骤(1)中将样本DNA文库经二代测序后所产生的读段比对到人类参考基因组上,获得包含读段比对结果的BAM文件,其中一种可用的比对软件为BWA,和一种可用的人类参考基因组版本为hg19;和在步骤(2)中提取BAM文件中的比对位置在目标区域内且支持结构变异的所有读段。
3.如权利要求1或2所述的方法,其中所述读段类型包括:
读段类型1):两个读段都不包含断点;
读段类型2):仅一个读段包含断点,其跨越断点的长度不足以产生二级比对结果;
读段类型3):两个读段都包含断点,但其跨越断点的长度均不足以产生二级比对结果;
读段类型4):至少一个读段包含断点且其跨越断点的长度足以产生二级比对结果,此时对另一个读段的比对位置及是否包含断点并无限定;
和读段类型5):至少一个读段包含断点,且其跨越断点的长度不足以产生二级比对结果;
其中,所述5种类型中,前3种类型的两个读段比对位置都跨越断点,第4种类型的一个读段的一级比对结果与其二级比对结果的比对位置跨越断点,对其中另一个读段无限定,第5种类型的两个读段都比对在断点的同一侧;
其中,两个成对读段比对位置跨越断点是指成对读段的两端的一级比对位置分布在断点两侧或单个读段的一级比对结果与其二级比对结果的比对位置分布在断点两侧;和
其中,读段包含断点是指断点位置在该读段的内部,读段的一部分序列跨越断点,在比对后将形成软剪切序列。
4.如权利要求1或2所述的方法,其中在步骤(2)中,要求跨越断点的成对读段比对在相同染色体上时二者间比对位置之间的距离的阈值是500bp至500Kb;和/或要求非跨越断点的成对读段比对在相同染色体上时二者间比对位置之间的距离的阈值是50bp至5Kb。
5.如权利要求1或2所述的方法,其中在步骤(2)中,要求跨越断点的成对读段比对在相同染色体上时二者间比对位置之间的距离的阈值是1Kb至100Kb;和/或要求未跨越断点的成对读段比对在相同染色体上时二者间比对位置之间的距离的阈值是100bp至2Kb。
6.如权利要求1或2所述的方法,其中在步骤(2)中,要求跨越断点的成对读段比对在相同染色体上时二者间比对位置之间的距离的阈值是5Kb至50Kb;和/或要求未跨越断点的成对读段比对在相同染色体上时二者间比对位置之间的距离的阈值是200bp至1Kb。
7.如权利要求1或2所述的方法,其中在步骤(3)中,聚类时要求成对读段与其它成对读段对应的比对位置间的距离为小于某个设定的数值,该设定的数值选自50bp至10Kb的范围中的某个值。
8.如权利要求1或2所述的方法,其中在步骤(3)中,聚类时要求成对读段与其它成对读段对应的比对位置间的距离为小于某个设定的数值,该设定的数值选自100bp至5Kb的范围中的某个值。
9.如权利要求1或2所述的方法,其中在步骤(3)中,聚类时要求成对读段与其它成对读段对应的比对位置间的距离为小于某个设定的数值,该设定的数值选自200bp至1Kb的范围中的某个值。
10.如权利要求1或2所述的方法,其中在步骤(6)中,过滤方式为要求总的特异的读段数目大于某个设定的数值,该设定的数值选自1-20中的某个值,且存在至少两种及以上的支持所述断点的读段类型。
11.如权利要求1或2所述的方法,其中在步骤(1)中,用比对软件将样本核酸序列经双端测序生成的读段比对到对应参考基因组产生的BAM文件,以该结果作为输入文件。
12.如权利要求1或2所述的方法,其中在步骤(6)中,依据以下方式推断结构变异产生的结果:
当跨断点的成对读段的两个读段分别比对在相同染色体的正链、负链上,且正链的比对位置小于负链的比对位置,则表征该读段支持结构变异缺失类型;
当跨断点的成对读段的两个读段分别比对在相同染色体的正链、负链上,且正链的比对位置大于负链的比对位置,则表征该读段支持结构变异重复类型;
当跨断点的成对读段的两个读段都比对在相同染色体的正链或负链上,则表征该读段支持结构变异倒置类型;和
当跨断点的成对读段的两个读段分别比对在不同染色体上,则表征该读段支持结构变异易位类型。
13.如权利要求1或2所述的方法,其中步骤(1)还包括从基因测序装置获得样本核酸序列经双端测序生成的读段,并且将读段信息经过输入设备输入到计算机中。
14.如权利要求1或2所述的方法,其中
在步骤(1)中计算机处理器接收并处理输入的读段信息,从而获得读段的比对结果;
在步骤(2)中计算机处理器提取样本核酸序列中的目标区域支持结构变异的所有读段,划分读段类型,并提取可能的断点精确位置;
在步骤(3)中计算机处理器将跨越断点的成对读段利用读段的比对位置进行聚类;
在步骤(4)中计算机处理器以从步骤(3)中产生的每一个读段子类为单位鉴定其中可能的结构变异断点位置;
在步骤(5)中计算机处理器以断点为中心统计支持结构变异的读段数目以及类型;和/或
在步骤(6)中计算机处理器对子类中每个可能的结构变异断点进行读段数目过滤,并将满足过滤条件的结构变异断点比对到基因上,并推断结构变异产生的结果。
15.如权利要求1或2所述的方法,在步骤(5)中,还以断点为中心统计支持该结构变异的两个断点间缺失或插入的序列长度和/或断点附近的序列复杂度的信息。
16.如权利要求1或2所述的方法,其中计算机处理器将步骤(6)中的结构变异产生的结果传输并显示到输出设备中。
17.用于实施如权利要求1-16中任一项所述的方法的计算机系统,其包括:
输入设备,用于输入样本核酸序列的测序读段信息;
计算机存储器,用于存储计算机程序指令;
计算机处理器,用于执行所述计算机程序指令,其中所述计算机程序指令实施步骤(1)至(6),对样本核酸序列的测序读段信息进行处理并生成结构变异产生的结果,并将结果传输到输出设备;和
输出设备,用于显示所述结果。
18.一种计算机可读介质,其中
所述计算机可读介质存储有计算机程序,
其中所述计算机程序能被计算机处理器执行以实施权利要求1-16中任一项的方法。
19.用于实施如权利要求1-16中任一项所述的方法的装置,其包括:
基因测序数据输入模块,用于输入样本核酸序列的测序读段信息;
获取测序读段比对结果模块,用于实施步骤(1);
提取目标区域结构变异支持读段模块,用于实施步骤(2);
跨越断点读段聚类模块,用于实施步骤(3);
从读段子类鉴定结构变异断点模块,用于实施步骤(4);
以断点为中心统计支持结构变异读段数模块,用于实施步骤(5);
结构变异过滤模块,用于实施步骤(6);和
显示结果模块,用于显示从步骤(6)获得的结果。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110268544.6A CN112687341B (zh) | 2021-03-12 | 2021-03-12 | 一种以断点为中心的染色体结构变异鉴定方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110268544.6A CN112687341B (zh) | 2021-03-12 | 2021-03-12 | 一种以断点为中心的染色体结构变异鉴定方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112687341A CN112687341A (zh) | 2021-04-20 |
CN112687341B true CN112687341B (zh) | 2021-06-04 |
Family
ID=75455450
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110268544.6A Active CN112687341B (zh) | 2021-03-12 | 2021-03-12 | 一种以断点为中心的染色体结构变异鉴定方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112687341B (zh) |
Families Citing this family (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN113628680B (zh) * | 2021-09-06 | 2022-06-21 | 哈尔滨师范大学 | 一种基于基准集的基因组结构变异性能检测方法 |
CN114005490B (zh) * | 2021-12-30 | 2022-04-22 | 北京优迅医疗器械有限公司 | 基于二代测序技术的循环肿瘤dna融合检测方法 |
Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2014055920A1 (en) * | 2012-10-05 | 2014-04-10 | The Regents Of The University Of California | Targeted detection of recurrent genomic rearrangements |
CN104302781A (zh) * | 2013-05-15 | 2015-01-21 | 深圳华大基因科技有限公司 | 一种检测染色体结构异常的方法及装置 |
WO2016070230A1 (en) * | 2014-11-05 | 2016-05-12 | University Of South Australia | Detecting sequence mutations in leukaemic fusion genes |
CN107992721A (zh) * | 2017-11-10 | 2018-05-04 | 深圳裕策生物科技有限公司 | 用于检测目标区域基因融合的方法、装置和存储介质 |
WO2019200328A1 (en) * | 2018-04-13 | 2019-10-17 | Guardant Health, Inc. | Methods for detecting and suppressing alignment errors caused by fusion events |
CN111081318A (zh) * | 2019-12-06 | 2020-04-28 | 人和未来生物科技(长沙)有限公司 | 一种融合基因检测方法、系统和介质 |
CN111243663A (zh) * | 2020-02-26 | 2020-06-05 | 西安交通大学 | 一种基于模式增长算法的基因变异检测方法 |
Family Cites Families (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108073790B (zh) * | 2016-11-10 | 2022-03-01 | 安诺优达基因科技(北京)有限公司 | 一种染色体变异检测装置 |
CN106611106B (zh) * | 2016-12-06 | 2019-05-03 | 北京荣之联科技股份有限公司 | 基因变异检测方法及装置 |
CN109280702A (zh) * | 2017-07-21 | 2019-01-29 | 深圳华大基因研究院 | 确定个体染色体结构异常的方法和系统 |
CN110010193B (zh) * | 2019-05-06 | 2021-09-03 | 西安交通大学 | 一种基于混合策略的复杂结构变异检测方法 |
-
2021
- 2021-03-12 CN CN202110268544.6A patent/CN112687341B/zh active Active
Patent Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2014055920A1 (en) * | 2012-10-05 | 2014-04-10 | The Regents Of The University Of California | Targeted detection of recurrent genomic rearrangements |
CN104302781A (zh) * | 2013-05-15 | 2015-01-21 | 深圳华大基因科技有限公司 | 一种检测染色体结构异常的方法及装置 |
WO2016070230A1 (en) * | 2014-11-05 | 2016-05-12 | University Of South Australia | Detecting sequence mutations in leukaemic fusion genes |
CN107992721A (zh) * | 2017-11-10 | 2018-05-04 | 深圳裕策生物科技有限公司 | 用于检测目标区域基因融合的方法、装置和存储介质 |
WO2019200328A1 (en) * | 2018-04-13 | 2019-10-17 | Guardant Health, Inc. | Methods for detecting and suppressing alignment errors caused by fusion events |
CN111081318A (zh) * | 2019-12-06 | 2020-04-28 | 人和未来生物科技(长沙)有限公司 | 一种融合基因检测方法、系统和介质 |
CN111243663A (zh) * | 2020-02-26 | 2020-06-05 | 西安交通大学 | 一种基于模式增长算法的基因变异检测方法 |
Non-Patent Citations (3)
Title |
---|
《Robust and exact structural variation detection with paired-end and soft-clipped alignments: SoftSV compared with eight algorithms》;Bartenhagen C et al;《ResearchGate》;20160131;第17卷(第1期);全文 * |
《基因组高通量测序数据结构变异识别算法》;王春宇等;《智能计算机与应用》;20150514;第5卷(第1期);全文 * |
《癌症基因组和转录组测序数据分析方法研究》;李倩倩;《中国优秀硕士学位论文全文数据库 医药卫生科技辑》;20150715(第2015年第07期);全文 * |
Also Published As
Publication number | Publication date |
---|---|
CN112687341A (zh) | 2021-04-20 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
JP7385686B2 (ja) | 無細胞核酸の多重解像度分析のための方法 | |
CN108368546B (zh) | 无细胞dna分析中基因融合检测的方法和应用 | |
CN107077537B (zh) | 用短读测序数据检测重复扩增 | |
EP2926288B1 (en) | Accurate and fast mapping of targeted sequencing reads | |
CN109767810B (zh) | 高通量测序数据分析方法及装置 | |
DE202013012824U1 (de) | Systeme zum Erfassen von seltenen Mutationen und einer Kopienzahlvariation | |
CN105408496A (zh) | 检测稀有突变和拷贝数变异的系统和方法 | |
CN112349346B (zh) | 检测基因组区域中的结构变异的方法 | |
CN108690871A (zh) | 基于二代测序的插入缺失突变检测方法、装置和存储介质 | |
CN112687341B (zh) | 一种以断点为中心的染色体结构变异鉴定方法 | |
CN112599188B (zh) | 一种融合驱动基因单端锚定的dna融合断点注释方法 | |
CN116064755B (zh) | 一种基于连锁基因突变检测mrd标志物的装置 | |
CN116469462B (zh) | 一种基于双重测序的超低频dna突变识别方法和装置 | |
CN110622250A (zh) | 用于检测插入和缺失的方法和系统 | |
CN115989544A (zh) | 用于在基因组的重复区域中可视化短读段的方法和系统 | |
CN116356001B (zh) | 一种基于血液循环肿瘤dna的双重背景噪声突变去除方法 | |
CN109461473B (zh) | 胎儿游离dna浓度获取方法和装置 | |
CN109920480B (zh) | 一种校正高通量测序数据的方法和装置 | |
CN114093424B (zh) | 病变特异性数据筛选及处理方法、装置、设备及存储介质 | |
CN105528532B (zh) | 一种rna编辑位点的特征分析方法 | |
CN104428423A (zh) | 确定外源基因在人类基因组中整合方式的方法和系统 | |
CN116153417B (zh) | 甲基化特征筛选方法及装置 | |
CN117577178B (zh) | 一种结构变异精确断裂信息的检测方法、系统及其应用 | |
EP3409788B1 (en) | Method and system for nucleic acid sequencing | |
CN119091965A (zh) | 一种基于DNA-based测序数据检测融合基因的方法 |
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 |