肿瘤纯度和平均倍体信息的预测方法、系统、存储介质

摘要

本发明属于信息处理技术领域，公开了一种肿瘤纯度和平均倍体信息的推测方法、系统、存储介质。从fasta文件中读取标准序列，从bam文件中读入测序读段；利用BWA软件对read段比对对齐，使用SAMtools提取readcount的值；定义Bin的长度，除去无用的位置和去除噪音，求出每个不包含缺失位置的bin的readdepth值；使用Free‑C软件检测出拷贝数变异CNV区域；使用高斯混合模型确定loss区域的类型；确定迭代范围和步长；采用偏差最小策略来筛选最优的结果，确定此时的肿瘤纯度和平均倍体为最终的结果。本发明防止噪音在模型中的传播，得到较为可靠的结果；缩短了需要比对的数据量。

说明书

技术领域

本发明属于信息处理技术领域，尤其涉及一种肿瘤纯度和平均倍体信息的预测方法、系统、存储介质。

背景技术

高斯混合模型是一种对数据进行聚类的算法，它的基本假设是总体的分布中含有K个子高斯分布。整个算法最主要的是以下几个概念：高斯分布：又叫正态分布，其数学期望μ决定了分布的位置，其标准差δ决定了分布的幅度。期望最大算法：在含有隐变量概率参数模型中，求解似然方程。该算法用于迭代求解高斯混合模型。极大似然估计：一种常见的统计方法，用于推测与一组样本相关的概率密度函数的参数。

目前，针对新一代测序技术的肿瘤纯度和平均倍体的检测方法，主要有以下几种：1.基于拷贝数的方法；通过计算肿瘤和正常基因组的拷贝数来估计肿瘤纯度；2.基于SNV杂合位点的方法：考虑肿瘤样本中带有体细胞突变的杂合位点的等位基因的突变部分，来估计肿瘤的纯度，或者同时估计出肿瘤的倍性；3.结合拷贝数变化和等位基因频率的方法：结合拷贝数变化和正常基因组中杂合位点上等位基因频率，来推断肿瘤纯度和倍性；

上述方法中，结合拷贝数改变的方法，在实践中往往能获得更好的准确度，同时方法本身也表现出更强的鲁棒性。这是因为，在NGS数据中，总的读数是非常大的，这就使得这种方法具有一定的统计学基础。现阶段使用该方法来检测肿瘤纯度和倍性的主要思想可分为以下几种：(1)通过确定肿瘤区域中映射到各个部分的reads服从的分布，建立概率模型来推断肿瘤的纯度和倍性。(2)通过结合全外显子组测序(WES)数据来推断肿瘤的纯度和倍性。(3)利用读段计数，并结合一定的先验知识，建立模型来推断肿瘤纯度。

通过上述分析，现有技术存在的问题及缺陷为：现有技术在肿瘤纯度和平均倍体检测中，对异常信息敏感，时间复杂度过大，计算时间长。

解决以上问题及缺陷的难度为：

1、如何在减少算法复杂度的基础上，确保算法结果能满足一定的准确度是困难的。2、数据带有噪音是普遍现象，当噪音比较大的时候，保证算法依然能给出一个让人接受的结果是困难的。

解决以上问题及缺陷的意义为：算法比较快的执行速度和较为可靠的结果，可以满足在某些特定应用场景下的需求，这些应用场景一般要求算法尽快给出结果。

发明内容

针对现有技术存在的问题，本发明提供了一种肿瘤纯度和平均倍体信息的预测方法、系统、存储介质。

本发明是这样实现的，一种肿瘤纯度和平均倍体信息的预测方法，所述肿瘤纯度和平均倍体信息的预测方法包括：

从fasta文件中读取标准序列，从bam文件中读入测序读段；

利用BWA软件对read段比对对齐，使用SAMtools提取readcount的值；定义Bin的长度，去除无用的位置和去除噪音，求出每个不包含缺失位置的bin的readdepth值；使用Free-C软件检测出拷贝数变异CNV区域；

使用高斯混合模型确定loss区域的类型；

确定迭代范围和步长；

采用偏差最小策略来筛选最优的结果，即当某一步计算出来的结果与实际的观测值的偏差小于阈值时，确定此时的肿瘤纯度和平均倍体为最终的结果。

进一步，所述使用高斯混合模型确定loss区域的类型包括：从确定的CNV区域中，取出loss类型的区域[l

进一步，所述确定迭代范围和步长包括：得到一组分类标签label∈[0，1]，对于标签为0的纯和缺失区域，利用其标签消去目标方程中的平均倍体p，迭代肿瘤的纯度；确定肿瘤纯度的范围是[0.05,0.95]，步长为0.01；对于标签为1的杂合缺失区域，代入计算得到的肿瘤纯度，利用目标方程迭代平均倍体，确定平均倍体的范围是[1.7,2.7]，步长为0.01。

本发明的另一目的在于提供一种计算机设备，所述计算机设备包括存储器和处理器，所述存储器存储有计算机程序，所述计算机程序被所述处理器执行时，使得所述处理器执行如下步骤：

从fasta文件中读取标准序列，从bam文件中读入测序读段；

利用BWA软件对read段比对对齐，使用SAMtools提取readcount的值；定义Bin的长度，去除无用的位置和去除噪音，求出每个不包含缺失位置的bin的readdepth值；使用Free-C软件检测出拷贝数变异CNV区域；

使用高斯混合模型确定loss区域的类型；

确定迭代范围和步长；

采用偏差最小策略来筛选最优的结果，即当某一步计算出来的结果与实际的观测值的偏差小于阈值时，确定此时的肿瘤纯度和平均倍体为最终的结果。

本发明的另一目的在于提供一种计算机可读存储介质，存储有计算机程序，所述计算机程序被处理器执行时，使得所述处理器执行如下步骤：

从fasta文件中读取标准序列，从bam文件中读入测序读段；

利用BWA软件对read段比对对齐，使用SAMtools提取readcount的值；定义Bin的长度，去除无用的位置和去除噪音，求出每个不包含缺失位置的bin的readdepth值；使用Free-C软件检测出拷贝数变异CNV区域；

使用高斯混合模型确定loss区域的类型；

确定迭代范围和步长；

采用偏差最小策略来筛选最优的结果，即当某一步计算出来的结果与实际的观测值的偏差小于阈值时，确定此时的肿瘤纯度和平均倍体为最终的结果。

本发明的另一目的在于提供一种实施所述肿瘤纯度和平均倍体信息的预测方法的肿瘤纯度和平均倍体信息的预测系统，所述肿瘤纯度和平均倍体信息的预测系统包括：

读取数据模块，从fasta文件中读取标准序列，从bam文件中读入测序读段；

确定CNV区域模块，用于检测出拷贝数变异CNV区域；

确定loss区域的类型模块，用于使用高斯混合模型确定loss区域的类型；

确定迭代范围和步长模块，用于确定算法的搜索范围；

确定肿瘤纯度和平均倍体模块，用于采用偏差最小策略来筛选最优的结果，即当某一步计算出来的结果与实际的观测值的偏差小于阈值时，确定此时的肿瘤纯度和平均倍体为最终的结果。

本发明的另一目的在于提供一种终端，所述终端搭载所述的肿瘤纯度和平均倍体信息的预测系统。

结合上述的所有技术方案，本发明所具备的优点及积极效果为：本发明通过减少推测变量的个数，控制变量之间的相互依赖，并引入校正值，这样可以防止噪音在模型中的传播，进而得到较为可靠的结果；聚焦于拷贝数变异区域的loss部分，大大缩短了需要比对的数据量，可以大幅提高算法的速度。

本发明可以控制变量之间的相互影响，降低模型对异常值的敏感性。同时，通过减少比对的区域，使得模型在准确率有保证的前提下，获得更高的效率。

附图说明

为了更清楚地说明本申请实施例的技术方案，下面将对本申请实施例中所需要使用的附图做简单的介绍，显而易见地，下面所描述的附图仅仅是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下还可以根据这些附图获得其他的附图。

图1是本发明实施例提供的肿瘤纯度和平均倍体信息的预测方法流程图。

图2是本发明实施例提供的肿瘤纯度和平均倍体信息的预测系统的结构示意图；

图2中：1、读取数据模块；2、确定CNV区域模块；3、确定loss区域的类型模块；4、确定迭代范围和步长模块；5、确定肿瘤纯度和平均倍体模块。

图3是本发明实施例提供的肿瘤纯度和平均倍体信息的预测方法实现流程图。

图4是本发明实施例提供对比了本发明与已知最好的算法的结果示意图；

图4(a)是各算法在仿真样本上输出的肿瘤纯度的值；(b)代表的是各算法在仿真样本上输出的平均倍体的值。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

针对现有技术存在的问题，本发明提供了一种肿瘤纯度和平均倍体信息的预测方法、系统、存储介质，下面结合附图对本发明作详细的描述。

如图1所示，本发明提供的肿瘤纯度和平均倍体信息的预测方法包括以下步骤：

S101：从fasta文件中读取标准序列，从bam文件读取测序读段；

S102：利用BWA软件对read段比对对齐，使用SAMtools提取readcount的值，；定义Bin的长度，除去无用的位置和去除噪音，求出每个不包含缺失位置的bin的readdepth值；使用Free-C软件检测出拷贝数变异CNV区域；

S103：使用高斯混合模型确定loss区域的类型；

S104：确定迭代范围和步长；

S105：采用偏差最小策略来筛选最优的结果，即当某一步计算出来的结果与实际的观测值的偏差小于阈值时，确定此时的肿瘤纯度和平均倍体为最终的结果。

本发明提供的肿瘤纯度和平均倍体信息的预测方法业内的普通技术人员还可以采用其他的步骤实施，图1的本发明提供的肿瘤纯度和平均倍体信息的预测方法仅仅是一个具体实施例而已。

如图2所示，本发明提供的肿瘤纯度和平均倍体信息的预测系统包括：

读取数据模块1，从fasta文件中读取标准序列，从bam文件中读入测序读段。

确定CNV区域模块2，用于检测出拷贝数变异CNV区域。

确定loss区域的类型模块3，用于使用高斯混合模型确定loss区域的类型。

确定迭代范围和步长模块4，用于确定算法的搜索范围。

确定肿瘤纯度和平均倍体模块5，用于采用偏差最小策略来筛选最优的结果，即当某一步计算出来的结果与实际的观测值的偏差小于阈值时，确定此时的肿瘤纯度和平均倍体为最终的结果。

下面结合附图对本发明的技术方案作进一步的描述。

如图3所示，本发明提供的肿瘤纯度和平均倍体信息的预测方法包括以下步骤：

(1)数据读入

从fasta文件中读取标准序列，从bam文件中读取测序读。

(2)确定CNV区域

利用BWA软件对read段比对对齐，使用SAMtools提取readcount的值，；定义Bin的长度，除去无用的位置和去除噪音，求出每个不包含缺失位置的bin的readdepth值；使用Free-C软件检测出拷贝数变异CNV区域；(3)

(3)使用高斯混合模型确定loss区域的类型

从(2)中确定的CNV区域中，取出loss类型的区域[l

(4)确定迭代范围和步长

通过(3)的计算，本发明得到了一组分类标签label∈[0，1]，对于标签为0的纯和缺失区域，本发明利用其标签消去目标方程中的平均倍体p，从而迭代肿瘤的纯度。在本发明中，本发明确定肿瘤纯度的范围是[0.05,0.95]，步长为0.01。对于标签为1的杂合缺失区域，本发明代入上一步计算得到的肿瘤纯度，此时利用目标方程迭代平均倍体。本发明确定平均倍体的范围是[1.7,2.7]，步长为0.01。

(5)确定肿瘤纯度和平均倍体

采用偏差最小策略来筛选最优的结果，即当某一步计算出来的结果与实际的观测值的偏差小于阈值时，确定此时的肿瘤纯度和平均倍体为最终的结果。

本发明研究拷贝数变异区域的loss部分，通过减少比对的区域，减少模型的计算时间；本发明在技术上使用高斯混合模型，对loss区域进行分类，通过结合不同缺失区域的特点，可以迭代求出肿瘤的纯度和平均倍体；

下面结合仿真对本发明的技术效果作详细的描述。

为了验证本发明的有效性，本发明在仿真数据上对比了本发明与已知最好的算法的结果。如图4所示：

图4(a)是各算法在仿真样本上输出的肿瘤纯度的值，其中蓝色代表的是从groundtruth中估计出来的值，可作为一个参考标准。棕色代表的是本发明的输出值，深绿代表的是现存算法ABSOLUTE的输出值。通过对比本发明可以发现:1、随着覆盖度的增加，本发明的输出值与估计值有接近的趋势；2、本发明的输出结果与ABSOLUTE的输出结果相差不大，这说明本发明通过合理减少比对数据可以实现与现存算法相当的效果；图4(b)代表的是各算法在仿真样本上输出的平均倍体的值。本发明可以发现：1、除了覆盖度为15x时，本发明的输出值与估计值有较小的偏差之外，在其他的覆盖度下，本发明的方法都很好的近似了平均倍体；2、总体来看，ABSOLUTE用于推测平均倍体的时候，产生了较大的偏差。

应当注意，本发明的实施方式可以通过硬件、软件或者软件和硬件的结合来实现。硬件部分可以利用专用逻辑来实现；软件部分可以存储在存储器中，由适当的指令执行系统，例如微处理器或者专用设计硬件来执行。本领域的普通技术人员可以理解上述的设备和方法可以使用计算机可执行指令和/或包含在处理器控制代码中来实现，例如在诸如磁盘、CD或DVD-ROM的载体介质、诸如只读存储器(固件)的可编程的存储器或者诸如光学或电子信号载体的数据载体上提供了这样的代码。本发明的设备及其模块可以由诸如超大规模集成电路或门阵列、诸如逻辑芯片、晶体管等的半导体、或者诸如现场可编程门阵列、可编程逻辑设备等的可编程硬件设备的硬件电路实现，也可以用由各种类型的处理器执行的软件实现，也可以由上述硬件电路和软件的结合例如固件来实现。

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，都应涵盖在本发明的保护范围之内。