CN113035273A

CN113035273A - 一种快速、超高灵敏度的dna融合基因检测方法

Info

Publication number: CN113035273A
Application number: CN202110266339.6A
Authority: CN
Inventors: 寻雪颖; 叶雷; 邓望龙; 任用; 李诗濛; 卜范峰; 丁然; 陆光华
Original assignee: Jiangsu Xiansheng Diagnostic Technology Co ltd; Nanjing Xiansheng Diagnostic Technology Co ltd; Nanjing Xiansheng Medical Testing Co ltd
Current assignee: Jiangsu Xiansheng Diagnostic Technology Co ltd; Jiangsu Xiansheng Medical Devices Co ltd; Nanjing Xiansheng Diagnostic Technology Co ltd; Nanjing Xiansheng Medical Laboratory Co ltd
Priority date: 2021-03-11
Filing date: 2021-03-11
Publication date: 2021-06-25
Anticipated expiration: 2041-03-11
Also published as: CN113035273B

Abstract

本发明提供一种DNA融合基因检测方法，该方法在不影响高灵敏性检测前提下，达到运行时间短并具有高检出特异性，同时能够准确地表征DNA测序数据中的融合基因断点、频率及拼接信息。

Description

一种快速、超高灵敏度的DNA融合基因检测方法

技术领域

本发明涉及生物信息学分析领域，特别是涉及一种快速、超高灵敏度的DNA融合基因检测的生信分析方法。

背景技术

基因组结构变异通常指的是长度大于1kbp的基因组结构改变，主要包括大片段的缺失(Deletion)，插入(Insertion)，倒置(Inversion)及易位(Translocation)^[1]。由基因组结构变异导致的不同基因间的序列拼接，通常称为融合基因。

基因融合作为一种重要的生物标记物在肿瘤的诊断、预后及治疗中提供重要的信息。例如，ALK、ROS1、RET基因融合在非小细胞肺癌中通常是独立的致癌因子，并且是靶向药作用的靶点^[2,3]。以往，肿瘤生物标记物的检测主要是对肿瘤组织样本进行检测。液体活检是一项革命性的技术，它打开了以前意想不到的前景，主要包括检测和分离循环肿瘤细胞(CTCs)、循环肿瘤DNA(cfDNA)和外泌体^[4]。因其具有微创、可重复的检测方式，以及包含癌症患者基因组和蛋白质组学信息，在临床诊断中具有相当大的意义。基于液体活检和NGS测序技术检测基因融合事件对于检测软件的灵敏性及特异性要求更高。Guardant360 CDx伴随诊断产品是FDA批准的首个基于NGS测序技术检测cfDNA中基因突变的体外诊断试剂盒^[5]。此产品的融合检出性能可以达到投入量为5ng时最低检出限为1％，投入量为30ng时最低检出限为0.1％。

NGS测序技术可以一次同时检测多个基因融合，在临床检验上具有很大的优势。目前可以对DNA测序数据集(WGS、WES、区域捕获等)和RNA测序数据集的挖掘来识别。现有从DNA测序数据集挖掘基因融合信息的主要方法分为基于序列组装，Read pair(成对read分别比对到染色体不同位置)及Split read(同一条read比对到染色体不同位置)的方法^[5]。基于序列组装的方法通过对短序列进行从头组装或局部组装，拼接成较长的基因序列，再与参考序列比较发现基因结构变异。此方法可以检测结构变异的类型最多，但检出结果的准确性更依赖于序列组装的效果。基于Read pair方法主要通过比较discordant reads之间的距离与插入片段大小的差异来确定基因组结构变异，但此方法的灵敏度受到插入片段长度标准差的影响，并且不能给出结构变异的准确位置。与Read pair方法相比，基于splitread方法通过softclip read比对信息可以直接获得精确的断点的位置。此外，有些融合检测方法会同时运用上述两种方法，如FACTERA^[6]等

然而上述融合检测方法都存在各自的利弊。基于序列组装的方法，检测效果依赖于序列组装的质量，容易生成较多的假阳信息，并且序列组装的方法存在消耗大量计算资源的问题。基于Read pair的方法需要估计测序数据插入片段大小，及其与Read pair之间距离的差异，容易造成假阴和过多假阳信息的问题，并且此方法只能给出融合的大致区域。利用split reads确定融合断点是较为准确的方法，但不同检测方法具体实施过程中，也存在灵敏性不足，假阳检出偏高的问题。

此外，现有的融合发现方法往往在模拟数据中表现良好，但在临床样本中却高估真实肿瘤基因组中的断点，几乎都存在假阳性率高的问题。同时，在低肿瘤细胞占比情况下目前的检测方法对于低融合频率的检出敏感性存在缺陷，很难满足目前临检生产尤其是液体活检的敏感性要求。简化分析步骤，缩短运行时间，特别是低肿瘤基因组占比样本中融合事件的准确检测成为本领域面临的重大难题。

有鉴于此，提出本发明。

部分参考文献如下：[1]Alkan C,Coe B P,Eichler E E.Genome structuralvariation discovery and genotyping.[J].Nature Reviews Genetics,2011,12(5):363-76。[2]Takeuchi K,Soda M,Togashi Y,et al.RET,ROS1 and ALK fusions in lungcancer[J].Nature medicine,2012,18(3):378-381。[3]Gainor J F,Shaw A T.Noveltargets in non-small cell lung cancer:ROS1 and RET fusions[J].The oncologist,2013,18(7):865。[4]Palmirotta R,Lovero D,Cafforio P,et al.Liquid biopsy ofcancer:a multimodal diagnostic tool in clinical oncology[J].Therapeuticadvances in medical oncology,2018,10:1758835918794630。[5]FDA Approves FirstLiquid Biopsy Next-Generation Sequencing Companion Diagnostic Test。

发明内容

本发明首要目的是提供一种快速、超高灵敏度的DNA融合基因检测方法，这种融合基因的检测可以针对任意种类中的融合基因，因此该方法可应用于疾病诊断类基因和非疾病的诊断类基因的融合。

为实现上述目的，本发明首先提供了一种DNA融合基因检测的生信分析方法，所述方法包括：

1)断点查找及初步筛选；

2)融合可信度判断；

3)假阳性融合过滤；

4)融合频率计算。

在一些优选的实施方式中，还包括：

5)融合方向判断及过滤。

在一些实施方式中，所述步骤1)包括如下步骤：

a.成对断点查找：在测序获得的BAM文件中通过识别split reads的主比对信息及次比对信息(SA tag)直接搜寻融合成对断点；

b.断点初步过滤：设置成对断点间距离过滤参数为大于1k，以及过滤掉仅有1条read支持的融合断点对；

在一些优选的实施方式中，所述成对断点查找步骤如下：

搜寻BAM文件中带有soft clip区域的reads，获取全部成对断点信息；通过softclipped read的主比对位置及cigar值确定断点1的位置及CN区域；根据soft clippedread的次比对信息(SA tag)的比对位置及cigar值确定断点2的位置及CN区域；统计具有相同成对断点及CN区域的reads数量，即为成对断点的Supplyment_Support支持数；

在一些更优选的实施方式中，还包括融合Supplyment_Support支持数矫正：矫正由于PCR duplicate造成支持数偏高而导致的假阳融合的问题，记录为dupcount支持数。

在一些实施方式中，所述步骤2)包括如下步骤：

a.融合断点基因区域注释；

b.候选融合序列矫正及拼接：通过比较同一条read在两处比对位置的mapping长度，进行融合序列的矫正及拼接，并记录两个断点处read信息；

c.融合断点验证；

优选的，使用BWA软件将断点处的reads重新回比到拼接的融合序列上，若reads可以跨过融合序列拼接点长度达到给定阈值则记录为支持融合事件的reads，记录此类reads的数量，即Fusion_VD。

在一些实施方式中，所述步骤3)包括如下步骤：

a.低复杂区域标记；

b.序列相似性比较；

c.softclip序列回比；

在一些优选的实施方式中：

所述a.低复杂区域标记为：通过计算拼接序列的串联重复长度及最大单碱基占比来描述序列低复杂区域，过滤掉由于测序仪导致的假阳信息；

所述b.序列相似性比较为：截取成对断点处附近序列进行相似性比较；若序列相似，该融合是由比对算法造成的假阳信息；

所述c.softclip序列回比：截取reads的softclip区域序列回比到断点附近的参考序列上，若能回比成功该融合则是SSARs(strand-split artifact reads)嵌合序列引入的假阳信息。

在一些实施方式中，所述步骤4)的融合频率为融合支持数(Fusion_VD)与断点测序深度的比值；

在一些具体方式中，包括如下步骤：

计算融合基因两个断点处测序深度；若两断点都在有效捕获区域内，断点测序深度为两断点处测序深度的平均值；若只有一个断点在有效捕获区域内，断点测序深度为两断点处测序深度的最大值；

在一些优选的实施方式中，还通过与融合阳性标准品的融合频率比较，获得矫正系数，对融合频率进行矫正。

在一些实施方式中，所述步骤5)包括如下步骤：

a.融合方向判断：上述步骤已经获得断点处的CN区域、基因、转录本及转录方向，根据这些信息判断5’->3’的基因排列顺序；

在一些优选的实施方式中，所述判断方法为：断点1处的基因为正向转录，比对区域标记为CN，那么根据断点1处判断的基因排列顺序为gene2-gene1；若断点1处的基因为正向转录，比对区域标记为NC，则基因排列顺序为gene1-gene2；基因为负向转录时，基因排列顺序与之相反；

b.过滤并输出结果；

在一些优选的实施方式中，所述过滤并输出结果为：设置低复杂区域，序列相似性，softclip序列回比reads占比，Fusion_VD，Frequency的阈值，将满足阈值条件的信息输出到结果文件中。

本发明还提供一种DNA融合基因检测的生信分析系统，其特征在于，包括如下模块：

1)断点查找及初步筛选模块；

2)融合可信度判断模块；

3)假阳性融合过滤模块；

4)融合频率计算模块；

5)融合方向判断及过滤模块。

上述模块能够实现上述各种实施方式中的行为。

在一些实施方式中，上述述系统通过识别split reads的主比对信息及次比对信息(SA tag)直接搜寻融合成对断点，并通过断点处reads回比融合序列确认融合信息。

本发明还提供一种计算机可读介质，其存储有计算机程序，所述计算机程序被处理器执行时，实现权上述任一项方法。

本发明还提供一种电子设备，包括处理器以及存储器，所述存储器上存储一条或多条可读指令，所述一条或多条可读指令被所述处理器执行时，实现上述任一项方法。

本发明还提供上述任一所述的检测系统、计算机可读介质以及电子设备在DNA融合基因检测中的应用，尤其在涉及肿瘤相关融合基因检测中的应用。

本发明有益的技术效果：

1)本方法通过识别split reads的主比对信息及次比对信息(SA tag)直接搜寻融合成对断点并通过断点处reads回比融合序列确认融合信息。简化融合断点搜寻及确认的方法，一次性获得全部成对断点信息，保证检测方法的高灵敏性。

2)多层级过滤策略(初步过滤、判断断点处基因序列的相似性、序列复杂度以及SSAR reads占比等)，保证分析方法高灵敏性和高特异性的前提下，降低程序计算资源消耗，缩短运行时间。

3)频率计算策略，以往融合分析方法只给出融合的支持数，没有对融合频率进行过多描述，本方法提供一种融合频率计算方法，并通过融合阳性标准品进行矫正，得到相对准确的融合频率信息。

4)本方法通过对断点处reads的softclip区域(CN)及断点所在基因正负链信息，充分解析基因融合形成的具体方式，明确融合基因形成的方向及区域拼接信息。

5)相比于已有方法，本专利在不影响高灵敏性的检测前提下，达到了运行时间短并具有高检出特异性，同时能够准确地表征DNA测序数据中的融合基因断点、频率及拼接信息。

说明书附图

图1：本发明基础流程图；

图2：融合频率计算结果图；

图3：a.discordant reads、split reads示意图。b.成对断点信息查找示意图。c.融合序列拼接及矫正示意图；

图4:Sample1 RNA测序数据检出结果可视化图；

图5：Sample1 EWSR1及FLI1基因断点处DNA测序数据softclip reads可视化；

图6：断点处比对信息可视化图。

具体实施方式

下面将结合实施例对本发明的实施方案进行详细描述，但是本领域技术人员将会理解，下列实施例仅用于说明本发明，而不应视为限制本发明的范围，并且所述实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

部分术语定义

除非在下文中另有定义，本发明具体实施方式中所用的所有技术术语和科学术语的含义意图与本领域技术人员通常所理解的相同。虽然相信以下术语对于本领域技术人员很好理解，但仍然阐述以下定义以更好地解释本发明。

如本发明中所使用，术语“包括”、“包含”、“具有”、“含有”或“涉及”为包含性的(inclusive)或开放式的，且不排除其它未列举的元素或方法步骤。术语“由…组成”被认为是术语“包含”的优选实施方案。如果在下文中某一组被定义为包含至少一定数目的实施方案，这也应被理解为揭示了一个优选地仅由这些实施方案组成的组。

本发明中的术语“大约”、“大体”表示本领域技术人员能够理解的仍可保证论及特征的技术效果的准确度区间。该术语通常表示偏离指示数值的±10％，优选±5％。

在提及单数形式名词时使用的不定冠词或定冠词例如“一个”或“一种”，“所述”，包括该名词的复数形式。

此外，说明书和权利要求书中的术语第一、第二、第三、(a)、(b)、(c)以及诸如此类，是用于区分相似的元素，不是描述顺序或时间次序必须的。应理解，如此应用的术语在适当的环境下可互换，并且本发明描述的实施方案能以不同于本发明描述或举例说明的其它顺序实施。

以下术语或定义仅仅是为了帮助理解本发明而提供。这些定义不应被理解为具有小于本领域技术人员所理解的范围。

本发明中的术语“测序读出序列”：英文为“read”或“reads”，指测序平台读出的一条或一组核酸序列。

本发明中的术语“比对结果”：英文为“alignment”，指一条测序读出序列与一条参考序列之间的对应结果，一条测序读出序列可以同时具有多个比对结果。

如下实施例进一步描述，所述实施例只是为了例证本发明的特定实施方案，不应理解为以任何方式限制本发明范围之意。除非另外说明，本发明中所公开的实验方法均采用本技术领域常规技术，实施例中所用的试剂和原材料均可由市场购得。

实施例1方法建立

(1)断点查找及初步筛选

a.成对断点查找：直接搜寻BAM文件中带有soft clip区域的reads，一次性获取全部成对断点信息：I.通过soft clipped read的主比对位置及cigar值确定断点1的位置及CN区域；II.根据soft clipped read的次比对信息(SA tag)的比对位置及cigar值确定断点2的位置及CN区域；III.统计具有相同成对断点及CN区域的reads数量，即为成对断点的Supplyment_Support支持数。

b.融合Supplyment_Support支持数矫正：矫正由于PCR duplicate造成支持数偏高而导致的假阳融合的问题，记录为dupcount支持数。

c.断点初步过滤：通常认为长度大于1kbp的变异为基因结构变异，设置成对断点间距离过滤参数为大于1k，以及过滤掉仅有1条read支持的融合断点对，可以有效降低假阳检出及后续步骤计算量。

(2)融合可信度判断

a.融合断点基因区域注释：对保留的成对断点信息进行基因区域的注释，记录断点所在的基因，转录本及外显子/内含子区域。

b.候选融合序列矫正及拼接：通过比较同一条read在两处比对位置的mapping长度，如附图3c所示，进行融合序列的矫正及拼接，并记录两个断点处read信息。

c.融合断点验证：使用软件(比如BWA软件)将断点处的reads重新回比到拼接的融合序列上，若reads可以跨过融合序列拼接点长度达到给定阈值则记录为支持融合事件的reads，记录此类reads的数量，即Fusion_VD。

(3)假阳性融合过滤

a.低复杂区域标记。NGS测序仪在遇到低复杂序列时往往容易产生测序错误。通过计算拼接序列的串联重复长度及最大单碱基占比来描述序列低复杂区域，过滤掉由于测序仪导致的假阳信息。

b.序列相似性比较。截取成对断点处附近序列进行相似性比较。若序列相似，该融合是由比对算法造成的假阳信息。

c.softclip序列回比。截取reads的softclip区域序列回比到断点附近的参考序列上，若能回比成功该融合则是SSARs(strand-split artifact reads)嵌合序列引入的假阳信息。

(4)融合频率计算

I.计算融合基因两个断点处测序深度，II.若两断点都在有效捕获区域内，断点测序深度为两断点处测序深度的平均值，III.若只有一个断点在有效捕获区域内，断点测序深度为两断点处测序深度的最大值，IV.融合频率为融合支持数(Fusion_VD)与断点测序深度的比值，V.通过与融合阳性标准品的融合频率比较，获得矫正系数，对融合频率进行矫正。

(5)融合方向判断及过滤

a.融合方向判断。上述步骤已经获得断点处的CN区域、基因、转录本及转录方向。根据这些信息判断5’->3’的基因排列顺序。具体判断方法：断点1处的基因为正向转录，比对区域标记为CN，那么根据断点1处判断的基因排列顺序为gene2-gene1；若断点1处的基因为正向转录，比对区域标记为NC，则基因排列顺序为gene1-gene2。基因为负向转录时，基因排列顺序与之相反。

b.过滤并输出结果。设置低复杂区域，序列相似性，softclip序列回比reads占比，Fusion_VD，Frequency的阈值，将满足阈值条件的信息输出到结果文件中。结果信息如下表1所示。

表1结果文件内容示例

a.SampleID:表示样本名称

b.GeneFusion:表示融合基因对

c.HotFusion:YES表示在热点融合基因库中有匹配

d.Gene1/Gene2:表示融合基因1/2的名称

e.Transcript1/Transcript2:表示融合基因1/2的转录本ID

f.Strand1/Strand2:表示融合基因1/2的转录本转录方向

g.CHR1:POS1/CHR2:POS2:表示融合断点1/2的位置

h.Order1/Order2:表示融合断点1/2处reads比对情况。CN:C表示read的softclip区域，N表示read的mapping区域。

i.BreakRegion1/BreakRegion2:表示融合断点1/2在此转录本的具体区域

j.Overlap:融合断点对间的重叠或间距序列长度，负数表示重叠，正数表示插入

k.Supplyment_Support:softclipped reads支持数。

l.Fusion_VD:重比对后的reads支持数

m.Fusion_Depth:断点测序深度

n.Frequency:融合频率

o.FusionDirection:融合方向及拼接区域

p.SV_type:发生融合的结构变异类型。

实施例2组织样本模拟试验

本实施例中阳性标准品具体融合信息如下表2所示。本实施例共10例HD753为gDNA标准品，模拟新鲜组织样本，来检测本方法的检测效果。所有样本均使用靶向捕获的方式获得DNA测序数据，且靶向捕获范围覆盖标准品中的融合断点位置。测序深度2500X。

表2标准品融合基因信息

10例样本测序数据，均使用BWA(v0.7.17)MEM算法进行比对分析，设置“-Y”参数保留softclip区域序列，用于后续融合序列矫正及融合事件验证。其余参数均使用BWA MEM算法默认参数。参考序列使用人参考基因组hg19版本。经BWA比对分析得到BAM格式文件作为本检测方法的输入文件。10例样本融合检出的灵敏性和特异性结果见表3。10例融合阳性标准品中，已知融合信息共20个，频率约为5％，本检测方法全部检出，且无假阳融合报出。如表4所示，检出频率与标准品给定频率接近。

表3 5％融合阳性标准品本方法检出的结果

TP

FP

FN

Ignore

TPR

PPV

Average Runtime

20

0

30

100.00％

102.3(s)

注：Ignore信息为ROS1-SLC34A2(chr6:117658325:CN-chr4:25666629:CN)融合，在所有阳性标准品中均存在，通常认为SLC34A2-ROS1(chr4:25666628:NC-chr6:117658310:N C)有临床意义，所以将此融合信息放入Ignore中,不计入灵敏度和特异性结果计算。

表4

样本信息	融合基因	预期频率	平均频率	标准差(SD)	变异系数(CV)
						HD753	CCDC6-RET	5％	5.87％	1.2064	20.55％
HD753	SLC34A2-ROS1	5.6％	5.12％	1.0642	20.76％

实施例3液体活检样本模拟试验

液体活检对融合检测的灵敏度要求越来越高，所以本实施例使用HD786 cfDNA标准品进行测试。HD786标准品中存在频率5％左右的SLC34A2-ROS1,CCDC6-RET融合。为得到更低频率的融合阳性标准品，将HD786标准品与融合阴性样本按照不同比例混合，模拟得到预期频率1％，0.5％左右的融合阳性cfDNA样本。本实施例设置三个组别，每组10个重复样本，信息如表5。实验测序深度为10000X。具体分析过程与实施例1相同，使用本方法进行融合分析。融合检测结果如表6所示，本方法全部检出。且在相同线程(20线程)下运行时间明显缩短，如表7所示。经本检测方法得到的融合频率与标准品预期频率相近，如图2所示，相关性大于0.85。

表5融合阳性cfDNA样本信息

组别	融合预期频率	重复次数
			Group1	5％	10
Group2	1％	10
			Group3	0.5％	10

表6融合阳性cfDNA标准品本方法检出的结果

Sample	TP	FP	FN	Ignore	TPR	PPV
							Group1	20	0	0	10	100.00％	100.00％
Group2	20	0	0	10	100.00％	100.00％
							Group3	20	0	0	10	100.00％	100.00％

表7 20线程运行时间比较

检测方法	平均运行时间
		本检测方法	244.4(s)
FACTERA软件	505.1(s)

实施例4对比试验

对比1

FACTERA(v1.4.4)软件是本领域常用的检测基因组结构变异的装置。FACTERA软件首先通过discordant read pair确定融合基因对，并记录discordant reads支持数。筛选满足条件(对比例中“-m”参数设置为1，即discordant reads支持数大于1)的基因对作为候选融合基因。根据候选融合基因区域的softclipped reads信息确定断点及softclippedreads支持数。接下来，根据候选基因softclipped reads支持数进行降序排列，通过k-mer比对的方式确定排名前5的断点对是否发生融合。之后，通过blastn比对验证融合信息并计算breakpoint-spanning reads支持数，输出满足条件(对比例中“-r”参数设置为3，即breakpoint-spanning reads支持数需大于等于3)的融合。

为了比较FACTERA软件与本发明的检测方法的差别，设置参数“-s 2-m 1-r 3-p20-b 300”，其余参数使用FACTERA软件默认参数。设置上述参数，理论上融合基因检出的灵敏度与本检测方法相当。接下来对上述实施例1中10个融合阳性标准品测序数据进行分析，FACTERA软件分析结果的灵敏度及特异性如表8所示。20个融合信息检出15个，FACTERA软件的灵敏性只有75％，而本发明装置全部检出，灵敏度到达100％。FACTERA软件的特异性只有46.88％，而本发明的特异性为100.00％，远低于本方法。

表8 5％融合阳性标准品FACTERA检出的结果

TP

FP

FN

Ignore

TPR

PPV

Average Runtime

15

17

5

10

75.00％

46.88％

62.7(s)

对比2

将实施例2中样本的测序数据，使用FACTERA分析软件进行分析，检出结果如下表9，HD786 cfDNA样本未经稀释样本，融合检出敏感性只有60％，特异性只有30％。稀释至1％频率时，检出十个已知融合信息，敏感性为50％，特异性只有40％。当稀释至0.5％频率时，只检出三个已知融合信息。FACTERA软件检出的灵敏性和特异性远低于本方法。

表9融合阳性cfDNA标准品FACTERA检出的结果

Sample	TP	FP	FN	Ignore	TPR	PPV
							Group1	12	28	8	10	60.00％	30.00％
Group2	10	15	10	10	50.00％	40.00％
							Group3	3	5	17	10	15.00％	37.50％

实施例5临床试验

取4例临床样本，通过对样本的RNA进行测序分析，确认存在融合信息，样本的具体融合如下图4和表10。

表10 RNA测序分析验证融合信息表

4例样本的基因组测序数据经本发明装置分析后，全部检出，结果如图5和下表11。

表11本方法融合检出结果表

样本编号	已知融合信息	检出融合
			Sample1	EWSR1--FLI1	EWSR1—FLI1
Sample2	EWSR1--FLI1	EWSR1—FLI1
			Sample3	EWSR1--FLI1	EWSR1—FLI1
Sample4	NAB2--STAT6	NAB2-STAT6

融合序列矫正案例：

1例EZR-ROS1融合阳性临床案例，此案例EZR-ROS1融合序列中插入43bp序列，断点处比对信息可视化见附图6。此样本测序数据经本发明装置分析后可以准确矫正融合序列并且能够明确融合5’->3’的拼接顺序，如表12。而FACTERA软件并没有准确矫正，表13为FACTERA软件结果部分截图，其中Break_offset表示EZR-ROS1融合断点处序列有1bp重叠，基因组可视化截图显示插入一段序列，与实际不符。

表12本方法检出结果部分列内容

表13 FACTERA分析结果部分列内容

Region1	Region2	Break1	Break2	Break_offset
					ROS1	EZR	Chr6:117646997	Chr6:159191451	-1

注：因融合输出结果列数较多，只展示部分关键信息。

综上所述，本分析方法通过对低复杂度融合序列，融合基因断点处序列高度相似，及softclip序列回比进行过滤，可以大量减少假阳融合信息检出。实施例和对比例的检测结果证明相较于FACTERA软件本方法具有更高的灵敏性，特异性，在临床样本中表现同样优异。

最后应说明的是：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，但本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。

Claims

1.一种DNA融合基因检测的生信分析方法，其特征在于，所述方法包括如下步骤：

步骤1)断点查找及初步筛选；

步骤2)融合可信度判断；

步骤3)假阳性融合过滤；

步骤4)融合频率计算。

2.权利要求1所述的DNA融合基因检测的生信分析方法，其特征在于，所述步骤1)包括如下步骤：

b.断点初步过滤：设置成对断点间距离过滤参数为大于1k，过滤掉仅有1条read支持的融合断点对；

优选的，所述成对断点查找步骤如下：

更优选的，还包括融合Supplyment_Support支持数矫正：矫正由于PCR duplicate造成支持数偏高而导致的假阳融合的问题，记录为dupcount支持数。

3.权利要求1-2任一所述的DNA融合基因检测的生信分析方法，其特征在于，所述步骤2)包括如下步骤：

a.融合断点基因区域注释；

c.融合断点验证；优选的，将断点处的reads重新回比到拼接的融合序列上，若reads可以跨过融合序列拼接点长度达到给定阈值则记录为支持融合事件的reads，记录此类reads的数量，即Fusion_VD。

4.权利要求1-3任一所述的DNA融合基因检测的生信分析方法，其特征在于，所述步骤3)包括如下步骤：

a.低复杂区域标记；

b.序列相似性比较；

c.softclip序列回比；

优选的：

5.权利要求1-4任一所述的DNA融合基因检测的生信分析方法，其特征在于，所述步骤4)的融合频率为融合支持数(Fusion_VD)与断点测序深度的比值；具体包括如下步骤：

优选的，进一步通过与融合阳性标准品的融合频率比较，获得矫正系数，对融合频率进行矫正。

6.权利要求1-5任一所述的DNA融合基因检测方法，其特征在于，还包括：

5)融合方向判断及过滤；

优选的，所述步骤5)包括如下步骤：

a.融合方向判断：对已获得断点处的CN区域、基因、转录本及转录方向，根据这些信息判断5’->3’的基因排列顺序；

更优选的，所述判断方法为：断点1处的基因为正向转录，比对区域标记为CN，那么根据断点1处判断的基因排列顺序为gene2-gene1；若断点1处的基因为正向转录，比对区域标记为NC，则基因排列顺序为gene1-gene2；基因为负向转录时，基因排列顺序与之相反；

b.过滤并输出结果；

更优选的，所述过滤并输出结果为：设置低复杂区域，序列相似性，softclip序列回比reads占比，Fusion_VD，Frequency的阈值，将满足阈值条件的信息输出到结果文件中。

7.一种DNA融合基因检测的生信分析系统，其特征在于，包括如下模块：

1)断点查找及初步筛选模块；

2)融合可信度判断模块；

3)假阳性融合过滤模块；

4)融合频率计算模块；

优选的，还包括：

5)融合方向判断及过滤模块。

8.权利要求7所述的DNA融合基因检测的生信分析系统，其特征在于，所述系统通过识别split reads的主比对信息及次比对信息(SA tag)直接搜寻融合成对断点，并通过断点处reads回比融合序列确认融合信息。

9.一种计算机可读介质，其存储有计算机程序，所述计算机程序被处理器执行时，实现权利要求1-6中任一项所述方法。

10.一种电子设备，其特征在于，包括处理器以及存储器，所述存储器上存储一条或多条可读指令，所述一条或多条可读指令被所述处理器执行时，实现权利要求1-6中任一项所述方法。