本发明涉及基因检测
技术领域：
，具体涉及一种基于长读长测序的多基因融合检测方法。
背景技术：
基因融合在基因组中非常普遍，也是一些类型癌症的标志。它由染色体重排而产生的，包括染色体的易位，插入，扩增，颠倒，缺失(非平衡重排)。基因融合常表现为两个不相关的基因融合形成，具有全新的功能或与两个融合前基因不同的功能。一个强启动子与一个下游功能基因(原癌基因)的融合在某些癌症中是普遍的。在生物体内发生融合基因，可导致疾病的发生。融合基因在癌症中普遍存在，与癌症的发生发展密切相关。随着近几年基于短读长测序的高通量测序技术的飞速发展和普及，高通量测序已被广泛用于基因融合检测：基于短读长测序获得数据，使用各种不断改进的算法检测基因融合。但这还是存在很大问题：1.基因组重复序列导致的多重比对使得检测结果不确定；2.无法检测大片段的多基因融合。技术实现要素：本发明的目的是提供一种基于长读长测序的多基因融合检测方法，以解决上述现有技术的问题。为达到上述目的，本发明提供了一种基于长读长测序的多基因融合检测方法，其包含以下步骤：步骤1：将通过dna长读长测序得到的读长比对到参考基因组上，得到读长的基因组坐标，并仅保留每条读长的最优比对结果；步骤2：将比对结果进行过滤，只保留可能存在基因融合突变的读长，得到候选读长数据库；步骤3：对每一条候选读长，按照其不同区域片段的比对结果，使用区域片段比对坐标序列表示读长；对全部的区域片段比对坐标序列进行聚类合并，形成包含多个融合基因读长组的候选多基因融合比对坐标序列数据库，其中每一个融合基因读长组包含多个来源于同一个基因融合突变的读长；步骤4：通过所述候选多基因融合比对坐标序列数据库中的每一个融合基因读长组确定对应的一个基因融合突变；为每一个基因融合突变确定其全部的断点坐标，形成多基因融合突变数据库。上述的基于长读长测序的多基因融合检测方法，其中，步骤1中，所述比对通过last比对算法处理。上述的基于长读长测序的多基因融合检测方法，其中，步骤1中，在进行比对前还包括步骤预处理，所述预处理为先将经过dna长读长测序得到的原始读长数据转换成fastq文件后，再通过过滤去除低质量的读长。上述的基于长读长测序的多基因融合检测方法，其中，步骤2中，所述的可能存在基因融合突变的读长是指存在2个以上的区域片段比对坐标的读长。上述的基于长读长测序的多基因融合检测方法，其中，步骤s3中，对全部的区域片段比对坐标序列进行聚类合并的具体过程为：对于任意的两条区域片段比对坐标序列a和b，且a的区域片段比对坐标数量大于b，如果对于b中的每一个区域片段比对坐标b(i)，均在a中存在域片段比对坐标a(j+i)或者a(j-i)与b(i)的左翼坐标的差值及右翼坐标的差值均小于10，则将b和a聚类合并到一组；其中，b(i)表示b中第i个区域片段比对坐标，且1≤i≤b的区域片段比对坐标总数；a(j+i)表示a中第j+i个区域片段比对坐标，且1≤j+i≤a的区域片段比对坐标总数；a(j-i)表示a中第j-i个区域片段比对坐标，且1≤j-i≤a的区域片段比对坐标总数。上述的基于长读长测序的多基因融合检测方法，其中，步骤s3中，每个融合基因读长组需要2个以上的区域片段比对坐标序列支持。上述的基于长读长测序的多基因融合检测方法，其中，该检测方法还包括步骤5：计算多基因融合突变数据库中每一个基因融合突变的概率值；如果基因融合突变的概率值大于等于基因融合突变的期望值，则将该基因融合突变按照基因融合类型进行归类，并保留在多基因融合突变数据库，否则丢弃。上述的基于长读长测序的多基因融合检测方法，其中，步骤s5中，所述基因融合类型包括易位、插入、扩增、颠倒、缺失和多基因突变。上述的基于长读长测序的多基因融合检测方法，其中，步骤s5中，所述概率值使用隐马尔夫链模型算法计算得到。上述的基于长读长测序的多基因融合检测方法，其中，步骤s5中，所述期望值是根据基因融合突变的断点数量、断点坐标的融合度能量值和测序错误率，使用隐马尔夫链模型算法计算得到的。相对于现有技术，本发明具有以下有益效果：长读长测序可以检测长片段dna/rna分子，从而直接获得融合基因全长，轻松判断融合位点。所以跟高通量测序相比，长读长测序在检测基因融合方面具有巨大优势。本发明所提供的新的基于长读长测序技术的多基因融合检测方法，具有高灵敏度和阳性预测值的优点，适合使用dna片段检测多基因融合，可以有效检测多基因融合，为临床检测疾病提供判断依据。具体实施方式以下通过具体实施例对本发明作进一步的描述，这些实施例仅用于说明本发明，并不是对本发明保护范围的限制。采用文献《nanoporesequencingdetectsstructuralvariantsincancer》(alexisl.norris等，doi：http://dx.doi.org/10.1080/15384047.2016.1139236)提供的fusion数据集作为本发明提供的一种基于长读长测序的多基因融合检测方法的检测对象。本发明提供了一种基于长读长测序的多基因融合检测方法，其包含以下步骤：步骤1：预处理和比对：先将经过dna长读长测序得到的原始读长数据转换成fastq文件后(fusion数据集的格式为fastq文件，故在此无需再进行数据转换步骤)，再通过过滤去除低质量的读长，低质量读长的定义没有固定标准，应该依据不同测序平台、该批次的测序质量、读长平均长度等因素自定义；last是一个专门用于长序列比对的软件，通过last比对算法将过滤后的读长比对到参考基因组上，得到读长的基因组坐标，并仅保留每条读长的最优比对结果；人有23条染色体，人类基因组计划通过测序得到每条染色体的碱基(actg)序列信息，将其公布出来作为参考基因组，方便后人研究；一般而言，根据公布的不同版本，主要分为grch37和grch38。last比对结果是指，将测序得到的每一条读长分别比对到参考基因组，得到它的基因组坐标，全部读长的基因组坐标即为last比对结果。步骤2：建立候选读长数据库：将last比对结果进行过滤，只保留可能存在基因融合突变的读长，得到候选读长数据库；所述的可能存在基因融合突变的读长是指存在2个以上的区域片段比对坐标(segment-alignment)的读长；如果一条读长存在一个或多个融合位点，会出现该读长上的不同区域片段(segment)比对到基因组上不同区域坐标(alignment)，且这些区域坐标不是连续的，则该读长存在2个以上的区域片段比对坐标，应该保留。如果一条读长只有一个区域片段比对坐标，则丢弃。步骤3：对候选读长进行聚类，建立候选多基因融合比对坐标序列数据库：对每一条候选读长，按照其不同区域片段的比对结果，使用区域片段比对坐标序列表示读长；对全部的区域片段比对坐标序列进行聚类合并，形成包含多个融合基因读长组的候选多基因融合比对坐标序列数据库，其中每一个融合基因读长组包含多个来源于同一个基因融合突变的读长；每个融合基因读长组需要2个以上的区域片段比对坐标序列支持；同一个基因融合突变会被测序多次，所以应该把来源于同一个基因融合突变的读长进行合并，以方便后续分析。对于融合突变，不同的基因片段连接在一起，也就是不连续的基因坐标连接在一起。举个例子，chr1:30000-50000表示1号染色体第30000到50000碱基序列。如果chr1:40000处插入一个片段chr8:7000-7800，则该突变应该表述为chr1:30000-40000；chr8:7000-7800；chr1:40001-50000。某一条读长刚好测到该区域，第1-500bp比对到chr1:38501-40000，第501-1300bp比对到chr8:7000-7800，第1301-2000bp比对到chr1:40001-40700。可以根据读长上不同区域的比对结果，转换成区域片段比对坐标序列：chr1:30000-40000；chr8:7000-7800；chr1:40001-50000。对全部的区域片段比对坐标序列进行聚类合并的具体过程为：对于任意的两条区域片段比对坐标序列a和b，且a的区域片段比对坐标数量大于b，如果对于b中的每一个区域片段比对坐标b(i)，均在a中存在域片段比对坐标a(j+i)或者a(j-i)与b(i)高度重合，如a(j+i)或者a(j-i)与b(i)的左翼坐标的差值及右翼坐标的差值均小于10，则将b和a聚类合并到一组；其中，b(i)表示b中第i个区域片段比对坐标，且1≤i≤b的区域片段比对坐标总数；a(j+i)表示a中第j+i个区域片段比对坐标，且1≤j+i≤a的区域片段比对坐标总数；a(j-i)表示a中第j-i个区域片段比对坐标，且1≤j-i≤a的区域片段比对坐标总数。所谓左翼坐标和右翼坐标是指：对于区域片段：chr1:40000-5000，其左翼坐标为chr1:40000，右翼坐标为chr1:50000。步骤4：确定断点位置，构建多基因融合突变数据库：通过所述候选多基因融合比对坐标序列数据库中的每一个融合基因读长组确定对应的一个基因融合突变；为每一个基因融合突变确定其全部的断点坐标，形成多基因融合突变数据库；本来连续的基因片段中断了，和外来基因片段连接在一起，中断又连接的位置即为断点。例如，对于区域片段比对坐标序列：chr1:30000-40000；chr8:7000-7800；chr1:40001-50000，其断点坐标为chr1:40000。在同一个融合基因读长组中，包含很多条读长。因为存在测序误差，对于每一个断点坐标，每条读长的值均不一致，一般通过取平均值或众数的方法得到确切的断点坐标。步骤5：过滤多基因融合突变数据库，降低假阳性：使用隐马尔夫链模型算法计算多基因融合突变数据库中每一个基因融合突变的概率值p；如果基因融合突变的概率值p大于等于基因融合突变的期望值e，则将该基因融合突变按照基因融合类型进行归类，并保留在多基因融合突变数据库，否则丢弃。所述期望值e是根据基因融合突变的断点数量、断点坐标的融合度能量值和测序错误率等因素，使用隐马尔夫链模型算法计算得到的，并不是一个单一阈值。断点坐标的融合度能量值是指两个基因片段断裂和重新连接融合均需要耗费的能量；两个融合片段连接处的化学键键能即为断点坐标的融合度能量值。一个融合突变含有很多特征，每个特征取值都是有一个概率分布的；根据一个突变的不同特征取值计算发生概率，通过大数据训练，形成基因融合突变期望值模型。对于某个特定融合突变，该模型给出最低阈值。低于该阈值判定该融合不可能发生。所述基因融合类型包括易位、插入、扩增、颠倒、缺失和多基因突变。使用3个融合基因检测工具(lumpy、sniffles和nanosv)分别检测fusion数据集，用于与本发明方法进行性能比较。这3个工具均为检测融合的常用工具。比较指标为：灵敏度和阳性预测值，灵敏度计算公式如下：(tp/tf)*100，阳性预测值计算公式如下：(tp/(tp+fp))*100。其中tp代表正确检测的融合突变的数量，tf代表所有的融合突变数量，fp代表假阳性融合突变的数量。fusion数据集含有18个样本，读长550–600bp，平均读长573bp，最低的融合突变浓度为1％。对每个样本分别计算阳性预测值和灵感度检测，本发明提供的方法对fusion数据集中18个样本的检测结果如表1所示：表1.本发明提供的方法对fusion数据集中18个样本的检测结果结果表明本发明提供的方法的阳性预测值均值达到100％，灵敏度均值达到100％。基于fusion数据集不同检测方法的性能(灵敏度平均值和阳性预测值平均值)比较如表2所示：表2.基于fusion数据集不同检测方法的性能比较衡量指标本发明方法lumpysnifflesnanosv灵敏度平均值/％10061.7676.4744.12阳性预测值平均值/％10043.8740.545.15由表2可知，在4种检测方法中，本发明提供的方法灵敏度平均值和阳性预测值平均值都为最优，且显著超越其余3种方法。综上所述，基于fusion数据集的检测性能比较表明，本发明提供了一种基于长读长测序的多基因融合检测方法，可有效检测多基因融合，灵敏度和阳性预测值等性能指标优于现有检测工具。尽管本发明的内容已经通过上述优选实施例作了详细介绍，但应当认识到上述的描述不应被认为是对本发明的限制。在本领域技术人员阅读了上述内容后，对于本发明的多种修改和替代都将是显而易见的。因此，本发明的保护范围应由所附的权利要求来限定。当前第1页12