CN115762646B - 一种基于gan样本平衡的泛癌公共驱动通路识别方法 - Google Patents

Abstract

本发明公开了一种基于GAN样本平衡的泛癌公共驱动通路识别方法，包括如下步骤：1)生成对应癌症的符合真实数据分布的体细胞突变数据；2)最小化各个癌症总权重之间离散度的模型CDP‑HA；3)引入单亲遗传算法对该模型CDP‑HA进行求解。这种方法成为识别泛癌公共驱动通路的有用工具，扩展性强，实用性强。

Description

一种基于GAN样本平衡的泛癌公共驱动通路识别方法

技术领域

本发明涉及生物信息学领域，用于识别癌症驱动通路，具体是一种基于GAN样本平衡的泛癌公共驱动通路识别方法。

背景技术

癌症是一种威胁人类健康且十分复杂的疾病，其病因包含多种遗传和环境因素。从分子水平了解致癌机理是一个巨大的挑战，有利于医学中癌症的诊断、治疗和药物设计。随着新一代测序技术(NGS)的快速发展，研究人员可以更好地描述癌症的分子特征。目前，几个大型的癌症基因组项目(癌症基因组图谱(TCGA)、国际癌症基因组联盟(ICGC)、癌症细胞系百科全书(CCLE))已经产生并分析了海量的数据，为进一步了解癌症的分子机制和致癌机制提供了前所未有的机会。已有的研究表明，只有功能驱动突变会促进癌症的发展，而乘客突变对癌症发展的影响很小。区分功能驱动突变和乘客突变已成为研究癌症发病机制的一项重要任务。

早期研究大多致力于设计能够有效识别出突变率显著较高的个体驱动基因。然而，由于癌症具有巨大的突变异质性，即同一癌症中的不同驱动基因突变都会导致癌症发生。因此识别单个驱动基因并不能有效了解癌症发展机制。进一步的研究表明，癌症的发生通常由部分的通路被破坏而导致，不同的驱动突变组合可能会扰乱这些通路(细胞信号或调控通路)。因此，识别驱动通路成为在通路水平上了解癌症致癌机制的关键途径。当前，驱动通路识别问题可以划分三个方向：识别单个驱动通路，识别协作驱动通路，识别泛癌的公共驱动通路和特异驱动通路。本文主要研究识别泛癌的公共驱动通路问题。

在泛癌规模上识别的公共驱动通路是为了研究不同癌症类型之间可能具有的共性，有利于增强对癌症致病机制的理解。TCGA泛癌项目已经收集了12种癌症类型的数千名癌症患者产生的多平台突变数据，为进一步研究这类问题提供了机会。最近，一类基于先验知识的研究方法被提出，这些方法通常使用基因-基因相互作用(GGI)网络、蛋白质-蛋白质相互作用(PPI)网络以及通路-通路相互作用(PaPaI)网络。虽然它们能够有较好的识别效果，但依赖于先验知识一方面会遗漏更好的突变基因组合的发现，另一方面则会限制寻找通路的范围，因为现有的先验知识不够完善且包含部分通路信息，如Leiserson等人提出基于定向热扩散模型的HotNet2方法，它试图通过结合蛋白质—蛋白质相互作用网络进行分析得到通路和蛋白复合物；Kim等人研究了多种癌症类型之间的不同类型的互斥性，并基于HumanNet网络提出了用于识别子网络/通路的MEMCover方法。Hajkarim等人研究基于大型基因-基因相互作用网络提出了DAMOKLE算法，它试图识别两种癌症中样本突变频率显著差异的子网络。另一类则是从头识别方法。Zhang等人提出了ComMDP方法，该方法利用驱动通路的两种特性：高互斥性和高覆盖度，再将适用于单一癌症的最大权重子矩阵问题模型直接扩展到是用于多种癌症类型，试图通过累加绝对权重值识别公共驱动通路。Wu等人引入了CDP-V模型，该模型使用相对比例代替绝对数字，并利用方差最小化每个比例的离散度，同时提出了使用调和均值最小化每个比例的离散度的CDP-H模型，减少了参数的使用。试图通过相对权重值识别公共驱动通路。

在上述研究方法中，ComMDP方法中使用绝对权重值，它通过把每个癌种的绝对权重值累加得出一组权重值最大的基因集，没有考虑到癌种之间样本量不均衡的问题，因此每个癌种样本量差距较大时，会导致识别结果更倾向于样本量较大的癌种，从而会遗漏某些共性更高的驱动通路。此外，PGA-V方法中使用相对权值计算的方式解决样本量不均衡的问题，但是却引进了一个人为参数，参数的确定需要大量的实验，而实验的过程会比较繁琐，不利于在实际应用中的扩展。然而，样本量不均衡的问题都未得到实际的解决。

随着深度学习的快速发展，解决数据不平衡的问题愈发急迫。在2012年,数据增强策略被提出，它是通过从已有的真实数据转换生成符合真实数据分布的附加数据项。但同一时期内并没有较好的数据生成方式。直到2014年Goodfellow等人提出了一类基于博弈论的强大生成模型——生成对抗网络(GANs)。虽然早期，诞生了大量基于深度学习的生成模型，但GANs是最成功的生成模型之一，并已成功应用于增强各个领域的数据。

发明内容

本发明的目的是针对现有技术的不足，而提供一种基于GAN样本平衡的泛癌公共驱动通识别方法，这种方法成为识别泛癌公共驱动通路的有用工具，扩展性强，实用性强。

实现本发明目的的技术方案是：

一种基于GAN样本平衡的泛癌公共驱动通路识别方法，包括如下步骤：

1)生成对应癌症的符合真实数据分布的体细胞突变数据：

1.1)设置对抗生成网络框架：

假设一个具有m_r个样本，n_r个基因的示例训练集，SNV-GANs的生成器网络定义为G(z)，生成器的输入是：z～norm(0,1)，生成器的定义如下：

1.1.1)输入层用GFC1将噪声向量z映射为维度(1,128)的张量zn；

1.1.2)隐藏层将步骤1.1.1)的张量zn放入GFC2进行映射，得到的结果放入GFC3进行映射，GFC4同理，最终映射为维度(1,1024)张量zn′；

1.1.3)输出层通过GFC5将张量zn′映射为维度(1,m_r*n_r)的张量gn，再将gn重置成维度为(m_r,n_r)的张量张量为生成器的输出；

其中输入层和隐藏层均使用了dropout函数对部分神经元进行冻结，并采用由式(1)定义的激活函数ReLU，输出层使用了由式(2)定义的激活函数Sigmod，辨别器网络定义为D(x)，辨别器的输入是真实数据x～P_real或者生成数据x代表一组体细胞突变数据样本，辨别器定义如下：

1.1.4)输入层用DFC1将x映射为维度(m_r,256)的张量xn；

1.1.5)隐藏层将步骤1.1.4)的张量xn放入DFC2进行映射，得到的结果放入DFC3进行映射，DFC4同理，最终映射为维度(m_r,16)张量xn′；

1.1.6)输出层通过DFC5将张量xn′映射为维度(m_r,1)的张量dn，即辨别器的输出；其中输入层和隐藏层均使用了由式(1)定义的激活函数ReLU，输出层使用了由式(2)定义的激活函数Sigmod；

ReLU(x)＝max(0,x) (1)，

1.2)SNV-GANs的训练过程：

1.2.1)给定一个体细胞突变矩阵A_r(m_r×n_r)和一个随机抽取样本的比例参数∝(∝<1)，按照比例参数∝从矩阵A_r随机抽取样本数量为m，m＝m_r*∝的子矩阵M_r，大小为m×n_r，一共抽取64个子矩阵M_r构造成训练集X，并输入到生成对抗网络中训练；

1.2.2)初始化辨别器D(.)的参数θ_d、生成器G(.)的参数θ_g；

1.2.3)令当前轮次epoch＝1，随机生成一个1×100的高斯分布的噪声向量z；

1.2.4)将步骤1.2.3)中的z作为生成器的输入，得到一个大小为m×n_r的向量 1.2.5)根据式(3)计算生成器损失值，然后更新生成器的参数θ_g：

其中G(z⁽ⁱ⁾)表示生成器通过噪声向量z⁽ⁱ⁾生成的生成数据，D(G(z⁽ⁱ⁾))表示辨别器生成数据判定为真实数据的概率，loss_G越小越好；

1.2.6)从训练集X中随机抽取一个样本组x；

1.2.7)随机生成一个1×100的高斯分布的噪声向量z；

1.2.8)将步骤1.2.7)中的z作为生成器的输入，得到一个大小为m×n_r的向量 1.2.9)根据式(4)计算生成器损失值，然后更新生成器的参数θ_g：

其中D(x⁽ⁱ⁾)表示辨别器将生成数据判定为真实数据的概率，1-D(G(z⁽ⁱ⁾))表示D将生成数据判定为生成数据的概率，loss_D越大越好；

1.2.10)判断当前轮次epoch是否达到设定的最大轮次：若是，则停止训练；否则，返回步骤1.2.3)最终得到训练好的生成器G(.)；

1.3)数据处理：

1.3.1)随机生成一个1×100的高斯分布的噪声向量z；

1.3.2)将步骤1.3.1)中的向量z输入到通过训练得到的生成器G(.)，得到生成数据G_data＝G(z)；

1.3.3)将G_data中大于等于0.85的值置为1，小于0.85的值置为0，得到一个新的二元矩阵A_fakedata；

1.3.4)取R种癌症体细胞突变矩阵中的样本数量最大的癌症的样本量为m_max，0<max<R；再将矩阵A_r插入扩增矩阵此时，扩增矩阵的样本数量为

1.3.5)若需要扩增的样本数量大于矩阵A_r的样本数量m_r，执行步骤1.3.6)，若需要扩增的样本数量小于矩阵A_r的样本数量m_r，执行步骤1.3.8)；

步骤1.3.6)、从步骤1.3.3)中的A_fakedata中随机抽取数量为m_r的样本分别计算矩阵和矩阵A_r中每个基因的突变率，分别得到两个对应突变概率集合V和Q；

1.3.7)将步骤1.3.6)中得到的两个集合V和Q输入由式(5)定义的JS散度公式得出一个分布值，分布值越小抽取矩阵与矩阵A_r中的突变率越相似，所以分布值小于等于0.09，即将矩阵插入扩增矩阵中，并更新反之则重复步骤1.3.6)；

1.3.8)从步骤1.3.3)中的A_fakedata中随机抽取一个样本直接加入扩增矩阵中，更新

1.3.9)若当前扩增矩阵的样本量等于最大样本量m_max，样本补充结束；反之，执行步骤1.3.5)最终得到一个与最大样本量m_max相等的扩增矩阵再使矩阵2)最小化各个癌症总权重之间离散度的模型CDP-HA：

设有R，R≥2种癌症类型，对于每种癌症类型，一个二进制体细胞突变矩阵表示为用于记录样本中的基因是否突变，其具有m_r行和n_r列，行代表样本或者患者，列代表基因，r＝1,2,3,..,R.，a_i-表示矩阵A_r中的第i个样本，a_-j表示矩阵A_r中的第j个基因，在第r种癌症的突变矩阵中的第i个样本的第j个基因发生突变时，否则给定一组大小为k的基因集S，表示对应矩阵A_r中大小为m_r×k的子矩阵，表示子矩阵中基因a发生突变的的样本，表示中覆盖的样本总数，用于衡量基因集S的覆盖度， 重叠覆盖的样本总和，衡量基因集S的互斥度；

根据上一段符号和问题的定义，构造了非线性最大化权重函数模型CDP-HA：给定R种癌症类型的m_r行，n_r列的二元体细胞突变矩阵A_r，一个参数K，令W_C(S)为最大权重和函数，确定一个m×K的子矩阵具体公式(6)如下：

其中表示基因集S在第r个癌种中的绝对权重值；

3)引入单亲遗传算法对模型CDP-HA进行求解：

3.1)设定适应度函数：

假设给定染色体E，令M_E代表与染色体对应的一个子矩阵，矩阵M_E的规模大小为m×K，个体适应度函数Fitness(E)的定义如下公式(7)所示，个体适应度函数值越大，代表可行解方案越好；

Fitness(E)＝W_C(M_E) (7)；

3.2)设定选择算子：

采用轮盘赌选择和精英策略来产生新一代种群，适应度最高的个体直接从父代遗传到子代，然后运用轮盘赌选择算子来生成其余N-1个个体；

3.3)设定重组算子：

采用一种基于贪心策略的重组算子，步骤如下：首先，给定一个父代染色体E＝{e₁,e₂,...,e_k}(e_i＝1,2,...,n)，其中e_i代表一个基因序号，因此E也称为基因集，由此确定候选基因集合其次，从基因集E中随机删除一个基因，得到基因集E^′；最后，基于贪心策略，从候选集合中选出最佳基因即并产生最终的新后代

3.4)设定参数：

输入R种癌症增强后的体细胞突变矩阵A_r，基因数g_number，参数k，种群规模N，算法执行次数t，最大演化代数maxg；

3.5)构造初始种群：

染色体用十进制的编码方式进行编码，一条染色体就代表一个个体，被用来代表问题的解向量，在单亲遗传算法中，用K个基因组成的集合作为一个问题解，即E＝{e₁,e₂,...,e_k}(e_i＝1,2,...,n)，种群中的个体初始化的方法为：随机产生1至n的自然数据集，每个数字表示突变矩阵中的一个基因，按顺序对n个基因进行分组，获得n/k个基因集S₁,S₂,...,S_n/k。令选择基因集S_max的基因构成初始染色体，通过产生N个初始染色体来产生初始种群；选取前K个数字作为一条初始染色体，生成一个初始种群pop₀，种群大小为N，计算pop₀种群染色体的适应值，将pop₀中最优的染色体相比较，保存最好的个体到变量best中，初始迭代次数step＝0；

3.6)执行迭代操作：

3.6.1)若step>maxg，转入步骤3.6.5)，得到大小为K的公共驱动通路，否则转入步骤3.6.2)；

3.6.2)对于种群pop_step，首先将pop_step中适应值最高的best染色体放入pop_step+1，然后执行轮盘赌选择算子选出其余N-1个染色体放入pop_step+1；

3.6.3)若step<0.7*maxg或Fitness(E′)>Fitness(E)，更新染色体E＝E′，否则不更新，保留E，step＝step+1；

3.6.4)取pop_step+1中适应值最高的染色体，若该染色体适应值大于best染色体的适应值，则更新best染色体，即best＝pop_step+1的最优染色体；

3.6.5)将best染色体转换为基因集，由此得到子矩阵M，并将子矩阵M输出，输出的M即为大小为K的公共驱动通路S。

本技术方案提出使用对抗生成网络(Generative Adversarial Networks，简称GAN)生成少样本癌症的体细胞突变数据的样本，使多个癌种之间的样本达到平衡，再利用调和均值累积每个癌种的绝对权重值的数学模型CDP-HA，最后使用单亲遗传算法求解该模型。对于几种类型的癌症，本技术方案提出的方法能够有效的对少样本数量癌症进行补充，解决样本数量差异问题。同时，基于提出的模型识别的基因集不仅在这些癌症的大部分样本中发生突变，在单独的癌症中的突变样本的比例很接近。此外，该方法检测到了一些具有生物学意义的基因集，这些基因集是其他方法中缺失的。因此，本技术方案提出的基于对抗生成网络缩小不同种类癌症之间的样本数据量差异成为可能，为识别癌症公共驱动通路提供了新的思路。

与现有技术相比，本技术方案具有以下优势：

(1)设计了一个非线性最小化离散度的最大化权重函数来度量多种癌症类型的相对权重。

(2)提出了一种适用于癌症体细胞突变数据的数据增强方法SNV-GANs，具有较高的使用价值，也是首次将GANs应用于癌症体细胞突变数据上。

(3)整个技术方案找出的泛癌公共驱动通路，包含更多富集在同一条重要信号通路的基因，以及识别出的基因富集在更多的重要信号通路。

这种方法不仅能成为生物癌症体细胞突变数据增强的有用工具，及识别泛癌公共驱动通路的有用工具，还能提供更多的生物信息，扩展性前，实用性强，而且能识别更多富集在重要信号通路的基因，以及识别出的基因富集在更多的重要信号通路。

附图说明

图1为实施例中公共驱动通路示例图；

图2为实施例中单亲遗传算法的示例图；

图3为实例中SNV-GANs的网络模型；

图4为实例中SNV-GANs的训练流程示意图；

图5为实例中SNV-GANs训练过程伪代码；

图6为实例中生成数据的实际效果示例图。

具体实施方式

下面结合附图及具体实施例对本发明作进一步的详细描述，但不是对本发明的限定。

实施例：

本例实验步骤1)在一台Linux服务器(Intel(R)Xeon(R)Gold 6230 2.10GHz CPU，内存为256G，显存为32G)，编译运行环境为Python 3.7.9。步骤2)和3)在一台电脑(Intel(R)Core(TM)i5-6500 3.20GHz CPU，内存为32G)上进行，操作系统为Windows 10，编译运行工具为Eclipse 4.23，编译环境Java 1.8.0。

本例针对泛癌公共驱动通路识别问题进行阐述。

一种基于GAN样本平衡的泛癌公共驱动通路识别方法，包括如下步骤：

假设有三个癌症类型的体细胞突变数据A_r：A₁(COADCORE具有95个样本，211个基因)，A₂(BLCA具有95个样本，211个基因)、A₃(LUAD具有95个样本，211个基因)。

1)生成对应癌症的符合真实数据分布的体细胞突变数据：

1.1)设置对抗生成网络框架：

假设使用BLCA癌症体细胞突变数据A₂作为示例训练集，SNV-GANs的生成器网络定义为G(z)，生成器的输入是：z～norm(0,1)，生成器的定义如下，如图3所示：

1.1.1)输入层用GFC1将噪声向量z映射为维度(1,128)的张量zn；

1.1.2)隐藏层将步骤1.1.1)的张量zn放入GFC2进行映射，得到的结果放入GFC3进行映射，GFC4同理，最终映射为维度(1,1024)张量zn′；

1.1.3)输出层通过GFC5将张量zn′映射为维度(1,95*211)的张量gn，再将gn重置成维度为(95,211)的张量张量为生成器的输出；

其中输入层和隐藏层均使用了dropout函数对部分神经元进行冻结，并采用由式(1)定义的激活函数ReLU，输出层使用了由式(2)定义的激活函数Sigmod，辨别器网络定义为D(x)，辨别器的输入是真实数据x～P_real或者生成数据x代表一组体细胞突变数据样本，辨别器定义如下，如图3所示：

1.1.4)输入层用DFC1将x映射为维度(95,256)的张量xn；

1.1.5)隐藏层将步骤1.1.4)的张量xn放入DFC2进行映射，得到的结果放入DFC3进行映射，DFC4同理，最终映射为维度(95,16)张量xn′；

1.1.6)输出层通过DFC5将张量xn′映射为维度(95,1)的张量dn，即辨别器的输出；其中输入层和隐藏层均使用了由式(1)定义的激活函数ReLU，输出层使用了由式(2)定义的激活函数Sigmod；

ReLU(x)＝max(0,x) (1)，

1.2)SNV-GANs的训练过程：

1.2.1)给定一个体细胞突变矩阵A₂(95×211)和一个随机抽取样本的比例参数∝＝0.7，按照比例参数∝从矩阵A_r随机抽取样本数量为m＝95*0.7≈67的子矩阵M₂，大小为67×211，一共抽取64个子矩阵M₂构造成训练集X，并输入到生成对抗网络中训练，如图4所示；

1.2.2)初始化辨别器D(.)的参数θ_d、生成器G(.)的参数θ_g；

1.2.3)令当前轮次epoch＝1，随机生成一个1×100的高斯分布的噪声向量z；

1.2.4)将步骤1.2.3)中的z作为生成器的输入，得到一个大小为67×211的向量

1.2.5)根据式(3)计算生成器损失值，然后更新生成器的参数θ_g：

其中G(z⁽ⁱ⁾)表示生成器通过噪声向量z⁽ⁱ⁾生成的生成数据，D(G(z⁽ⁱ⁾))表示生成器将生成数据判定为真实数据的概率，loss_G越小越好；

1.2.6)从训练集X中随机抽取一个样本组x；

1.2.7)随机生成一个1×100的高斯分布的噪声向量z；

1.2.8)将步骤1.2.7)中的z作为生成器的输入，得到一个大小为95×211的向量

1.2.9)根据式(4)计算生成器损失值，然后更新生成器的参数θ_g：

其中D(x⁽ⁱ⁾)表示辨别器将生成数据判定为真实数据的概率，1-D(G(z⁽ⁱ⁾))表示D将生成数据判定为生成数据的概率，loss_D越大越好；

1.2.10)判断当前轮次epoch是否达到10000次：若是，则停止训练；否则，返回步骤1.2.3)最终得到训练好的生成器G(.)；

1.3)数据处理：

1.3.1)随机生成一个1×100的高斯分布的噪声向量z；

1.3.2)将步骤1.3.1)中的向量z输入到通过训练得到的生成器G(.)，得到生成数据G_data＝G(z)；

1.3.3)将G_data中大于等于0.85的值置为1，小于0.85的值置为0，得到一个新的二元矩阵A_fakedata；

1.3.4)取3种癌症体细胞突变矩阵中的样本数量最大的癌症的样本量为m_max＝m₁＝489；再将矩阵A₂插入扩增矩阵此时，扩增矩阵的样本数量为

1.3.5)若需要扩增的样本数量大于矩阵A₂的样本数量m₂＝95，执行步骤1.3.6)，若需要扩增的样本数量小于矩阵A₂的样本数量m₂＝95，执行步骤1.3.8)；

1.3.6)从步骤1.3.3)中的A_fakedata中随机抽取数量为m₂的样本分别计算矩阵和矩阵A₂中每个基因的突变率，分别得到两个对应突变概率集合V和Q；1.3.7)将步骤1.3.6)中得到的两个集合V和Q输入由式(5)定义的JS散度公式得出一个分布值，分布值越小抽取矩阵与矩阵A_r中的突变率越相似，所以分布值小于等于0.09，即将矩阵插入扩增矩阵中，并更新反之则重复步骤1.3.6)；

1.3.8)从步骤1.3.3)中的A_fakedata中随机抽取一个样本直接加入矩阵A_r中，更新

1.3.9)若当前扩增矩阵的样本量等于最大样本量m_max＝489，样本补充结束；反之，执行步骤1.3.5)最终得到一个与最大样本量m_max＝489相等的矩阵得到新的扩增矩阵的突变率与原始的矩阵A₂的突变率基本一致，再使矩阵如图6所示；训练过程伪代码如图5；

2)最小化各个癌症总权重之间的离散度的模型CDP-HA：

设有R＝3种癌症类型，对于每种癌症类型，一个二进制体细胞突变矩阵表示为用于记录样本中的基因是否突变，其具有m_r行和n_r列，行代表样本或者患者，列代表基因，r＝1,2,3,..,R.，a_i-表示矩阵A_r中的第i个样本，a_-j表示矩阵A_r中的第j个基因，在第r种癌症的突变矩阵中的第i个样本的第j个基因发生突变时，否则给定一组大小为k的基因集S，表示对应矩阵A_r中大小为m_r×k的子矩阵，表示子矩阵中基因a发生突变的的样本，表示中覆盖的样本总数，用于衡量基因集S的覆盖度， 重叠覆盖的样本总和，衡量基因集S的互斥度；

根据以上定义，构造了非线性最大化权重函数模型CDP-HA：给定R种癌症类型的m_r行，n_r列的二元体细胞突变矩阵A_r，一个参数K(0<K≤10)，令W_C(S)为最大权重和函数，确定一个m×K的子矩阵具体公式(6)如下：

其中表示基因集S在第r个癌种中的绝对权重值；如图1所示是三个癌症的突变矩阵，有两个子矩阵S₁和S₂，规模大小为K为3，根据本例提出的模型，可得发现S₂的权重更高，而且S₂更符合公共驱动通路所需要的共性；

3)引入单亲遗传算法对模型CDP-HA进行求解：

3.1)设定适应度函数：

假设给定染色体E，令M_E代表与染色体对应的一个子矩阵，矩阵M_E的规模大小为m×K，个体适应度函数Fitness(E)的定义如下公式(7)所示，个体适应度函数值越大，代表可行解方案越好；

Fitness(E)＝W_C(M_E) (7)；

3.2)设定选择算子：

采用轮盘赌选择和精英策略来产生新一代种群，适应度最高的个体直接从父代遗传到子代，然后运用轮盘赌选择算子来生成其余N-1个个体；

3.3)设定重组算子：

采用一种基于贪心策略的重组算子，步骤如下：首先，给定一个父代染色体E＝{e₁,e₂,...,e_k}(e_i＝1,2,...,n)，其中代表一个基因序号，因此E也称为基因集,由此确定候选基因集合其次，从基因集E中随机删除一个基因，得到基因集E^′；最后，基于贪心策略，从候选集合中选出最佳基因即并产生最终的新后代

3.4)设定参数：

输入3种癌症增强后的体细胞突变矩阵A_r，基因数g_number＝211，参数k＝3，种群规模N＝20，算法执行次数step＝10，最大演化代数maxg＝1000；

3.5)构造初始种群：

染色体用十进制的编码方式进行编码，一条染色体就代表一个个体，被用来代表问题的解向量，在单亲遗传算法中，用K＝3个基因组成的集合作为一个问题解，即E＝{e₁,e₂,...,e_k}(e_i＝1,2,...,n)，种群中的个体初始化的方法为：随机产生1至20的自然数据集，每个数字表示突变矩阵中的一个基因，按顺序对20个基因进行分组，获得n/k＝20/3≈7个基因集S₁,S₂,...,S_n/k。令选择基因集S_max的基因构成初始染色体，通过产生N个初始染色体来产生初始种群；选取前K＝3个数字作为一条初始染色体，生成一个初始种群pop₀，种群大小为N，计算pop₀种群染色体的适应值，将pop₀中最优的染色体相比较，保存最好的个体到变量best中，初始迭代次数step＝0，如图2所示；

3.6)执行迭代操作：

3.6.1)若step>maxg，转入步骤3.6.5)，得到大小为K的公共驱动通路，否则转入步骤3.6.2)；

3.6.2)对于种群pop_step，首先将pop_step中适应值最高的best染色体放入pop_step+1，然后执行轮盘赌选择算子选出其余N-1＝20-1＝19个染色体放入pop_step+1；3.6.3)若step<700或Fitness(E′)>Fitness(E)，更新染色体E＝E′，否则不更新，保留E，step＝step+1；

3.6.4)取pop_step+1中适应值最高的染色体，若该染色体适应值大于best染色体的适应值，则更新best染色体，即best＝pop_step+1的最优染色体；

3.6.5)将best染色体转换为基因集，由此得到子矩阵M，并将子矩阵M输出，输出的M即为大小为K＝3的公共驱动通路S。

Claims (1)

1.一种基于GAN样本平衡的泛癌公共驱动通路识别方法，其特征在于，包括如下步骤：

1)生成对应癌症的符合真实数据分布的体细胞突变数据：

1.1)设置对抗生成网络框架：

假设一个具有m_r个样本，n_r个基因的示例训练集，SNV-GANs的生成器网络定义为G(z)，生成器的输入是：z～norm(0,1)，生成器的定义如下：

1.1.1)输入层用GFC1将噪声向量z映射为维度(1,128)的张量zn；

1.1.2)隐藏层将步骤1.1.1)的张量zn放入GFC2进行映射，得到的结果放入GFC3进行映射，GFC4同理，最终映射为维度(1,1024)张量zn′；

1.1.3)输出层通过GFC5将张量zn′映射为维度(1,m_r*n_r)的张量gn，再将gn重置成维度为(m_r,n_r)的张量张量为生成器的输出；

其中输入层和隐藏层均使用了dropout函数对部分神经元进行冻结，并采用由式(1)定义的激活函数ReLU，输出层使用了由式(2)定义的激活函数Sigmod，辨别器网络定义为D(x)，辨别器的输入是真实数据x～P_real或者生成数据x代表一组体细胞突变数据样本，辨别器定义如下：

1.1.4)输入层用DFC1将x映射为维度(m_r,256)的张量xn；

1.1.5)隐藏层将步骤1.1.4)的张量xn放入DFC2进行映射，得到的结果放入DFC3进行映射，DFC4同理，最终映射为维度(m_r,16)张量xn′；

1.1.6)输出层通过DFC5将张量xn′映射为维度(m_r,1)的张量dn，即辨别器的输出；

其中输入层和隐藏层均使用了由式(1)定义的激活函数ReLU，输出层使用了由式(2)定义的激活函数Sigmod；

ReLU(x)＝max(0,x) (1)，

1.2)SNV-GANs的训练过程：

1.2.1)给定一个体细胞突变矩阵A_r(m_r×n_r)和一个随机抽取样本的比例参数∝，∝<1，按照比例参数∝从矩阵A_r随机抽取样本数量为m，m＝m_r*∝的子矩阵M_r，大小为m×n_r，一共抽取64个子矩阵M_r构造成训练集X，并输入到生成对抗网络中训练；

1.2.2)初始化辨别器D(.)的参数θ_d、生成器G(.)的参数θ_g；

1.2.3)令当前轮次epoch＝1，随机生成一个1×100的高斯分布的噪声向量z；

1.2.4)将步骤1.2.3)中的z作为生成器的输入，得到一个大小为m×n_r的向量

1.2.5)根据式(3)计算生成器损失值，然后更新生成器的参数θ_g：

其中G(z⁽ⁱ⁾)表示生成器通过噪声向量z⁽ⁱ⁾生成的生成数据，D(G(z⁽ⁱ⁾))表示辨别器将生成数据判定为真实数据的概率，loss_G越小越好；

1.2.6)从训练集X中随机抽取一个样本组x；

1.2.7)随机生成一个1×100的高斯分布的噪声向量z；

1.2.8)将步骤1.2.7)中的z作为生成器的输入，得到一个大小为m×n_r的向量

1.2.9)根据式(4)计算生成器损失值，然后更新生成器的参数θ_g：

其中D(x⁽ⁱ⁾)表示辨别器将生成数据判定为真实数据的概率，1-D(G(z⁽ⁱ⁾))表示D将生成数据判定为生成数据的概率，loss_D越大越好；

1.2.10)判断当前轮次epoch是否达到设定的最大轮次：若是，则停止训练；否则，返回步骤1.2.3)最终得到训练好的生成器G(.)；

1.3)数据处理：

1.3.1)随机生成一个1×100的高斯分布的噪声向量z；

1.3.2)将步骤1.3.1)中的向量z输入到通过训练得到的生成器G(.)，得到生成数据G_data＝G(z)；

1.3.3)将G_data中大于等于0.85的值置为1，小于0.85的值置为0，得到一个新的二元矩阵A_fakedata；

1.3.4)取R种癌症体细胞突变矩阵中的样本数量最大的癌症的样本量为m_max，0<max<R；再将矩阵A_r插入扩增矩阵此时，扩增矩阵的样本数量为1.3.5)若需要扩增的样本数量大于矩阵A_r的样本数量m_r，执行步骤1.3.6)，若需要扩增的样本数量小于矩阵A_r的样本数量m_r，执行步骤1.3.8)；

1.3.6)从步骤1.3.3)中的A_fakedata中随机抽取数量为m_r的样本分别计算矩阵和矩阵A_r中每个基因的突变率，分别得到两个对应突变概率集合V和Q；

1.3.7)将步骤1.3.6)中得到的两个集合V和Q输入由式(5)定义的JS散度公式得出一个分布值，分布值越小抽取矩阵与矩阵A_r中的突变率越相似，所以分布值小于等于0.09，即将矩阵插入扩增矩阵中，并更新反之则重复步骤1.3.6)；

1.3.8)从步骤1.3.3)中的A_fakedata中随机抽取一个样本直接加入扩增矩阵中，更新

1.3.9)若当前扩增矩阵的样本量等于最大样本量m_max，样本补充结束；反之，执行步骤1.3.5)最终得到一个与最大样本量m_max相等的的扩增矩阵再使矩阵2)最小化各个癌症总权重之间离散度的模型CDP-HA：

设有R，R≥2种癌症类型，对于每种癌症类型，一个二进制体细胞突变矩阵表示为用于记录样本中的基因是否突变，其具有m_r行和n_r列，行代表样本或者患者，列代表基因，r＝1,2,3,..,R.，a_i-表示矩阵A_r中的第i个样本，a_-j表示矩阵A_r中的第j个基因，在第r种癌症的突变矩阵中的第i个样本的第j个基因发生突变时，否则给定一组大小为k的基因集S，表示对应矩阵A_r中大小为m_r×k的子矩阵，表示子矩阵中基因a发生突变的的样本，表示中覆盖的样本总数，用于衡量基因集S的覆盖度，重叠覆盖的样本总和，衡量基因集S的互斥度；

根据上一段符号和问题的定义，构造了非线性最大化权重函数模型CDP-HA：给定R种癌症类型的m_r行，n_r列的二元体细胞突变矩阵A_r，一个参数K，令W_C(S)为最大权重和函数，确定一个m×K的子矩阵具体公式(6)如下：

其中表示基因集S在第r个癌种中的绝对权重值；3)引入单亲遗传算法对该模型CDP-HA进行求解：

3.1)设定适应度函数：

假设给定染色体E，令M_E代表与染色体对应的一个子矩阵，矩阵M_E的规模大小为m×K，个体适应度函数Fitness(E)的定义如下公式(7)所示，个体适应度函数值越大，代表可行解方案越好；

Fitness(E)＝W_C(M_E) (7)；

3.2)设定选择算子：

采用轮盘赌选择和精英策略来产生新一代种群，适应度最高的个体直接从父代遗传到子代，然后运用轮盘赌选择算子来生成其余N-1个个体；

3.3)设定重组算子：

采用一种基于贪心策略的重组算子，步骤如下：首先，给定一个父代染色体E＝{e₁,e₂,...,e_k}，e_i＝1,2,...,n，其中e_i代表一个基因序号，因此E也称为基因集，由此确定候选基因集合其次，从基因集E中随机删除一个基因，得到基因集E^′；最后，基于贪心策略，从候选集合中选出最佳基因即并产生最终的新后代

3.4)设定参数：

输入R种癌症增强后的体细胞突变矩阵A_r，基因数g_number，参数k，种群规模N，算法执行次数t，最大演化代数maxg；

3.5)构造初始种群：

染色体用十进制的编码方式进行编码，一条染色体就代表一个个体，被用来代表问题的解向量，在单亲遗传算法中，用K个基因组成的集合作为一个问题解，即E＝{e₁,e₂,...,e_k}，e_i＝1,2,...,n，种群中的个体初始化的方法为：随机产生1至n的自然数据集，每个数字表示突变矩阵中的一个基因，按顺序对n个基因进行分组，获得n/k个基因集S₁,S₂,...,S_n/k，令选择基因集S_max的基因构成初始染色体，通过产生N个初始染色体来产生初始种群；选取前K个数字作为一条初始染色体，生成一个初始种群pop₀，种群大小为N，计算pop₀种群染色体的适应值，将pop₀中最优的染色体相比较，保存最好的个体到变量best中，初始迭代次数step＝0；

3.6)执行迭代操作：

3.6.1)若step>maxg，转入步骤3.6.5)，得到大小为K的公共驱动通路，否则转入步骤3.6.2)；

3.6.2)对于种群pop_step，首先将pop_step中适应值最高的best染色体放入pop_step+1，然后执行轮盘赌选择算子选出其余N-1个染色体放入pop_step+1；

3.6.3)若step<0.7*maxg或Fitness(E′)>Fitness(E)，更新染色体E＝E′，否则不更新，保留X，step＝step+1；

3.6.4)取pop_step+1中适应值最高的染色体，若该染色体适应值大于best染色体的适应值，则更新best染色体，即best＝pop_step+1的最优染色体；

3.6.5)将best染色体转换为基因集，由此得到子矩阵M，并将子矩阵M输出，输出的M即为大小为K的公共驱动通路S。