基于有向图卷积的上皮细胞基因调控关系预测方法

技术领域

本发明涉及智能信息处理和生物信息领域，具体来说是一种基于有向图卷积的上皮细胞基因调控关系的预测方法。

背景技术

对于细胞来说，转录调节是非常重要的，它允许细胞适应不断变化的环境，并且对外界刺激做出反应来动态的调节基因的表达。转录因子(TF)通过调控自身及其对应靶基因(Target)发挥作用，转录因子与靶基因以及他们之间的相互作用构成了基因调控网络(GRN)。GRN旨在捕捉转录因子和靶基因这些实体分子之间的依赖关系。基因调控网络通过控制基因表达在生物体的生长发育中有着至关重要的作用，具体表现为：GRN是维持生命过程、控制生化反应、调节化合物水平等的核心元素，在各种生物体和系统中发挥着重要作用。因此，准确的构建GRN可以提高我们对基因的调控机制以及生物体分子机制的理解。

为了更好的了解基因间的调控机制，近年来开发了很多算法用于GRN的构建：基于信息论的方法以基因的相关性为基础来量化基因之间的统计依赖性，从得到的基因间的相关矩阵中来预测GRN的结构。基于布尔网络的方法是将使用布尔变量来表示基因的状态，然后使用布尔函数来计算基因之间的关系，可以捕捉到基因之间的动态相互作用。特征选择也逐渐应用于推断GRN，例如GENIE3、dynGENIE3、Jump3、BiXGBoost等方法。随着深度学习的发展，越来越多的深度学习方法应用于GRN推断，例如CNN、TDL、DeepDRIM等方法。

上述方法都是基于传统方法例如信息论和传统机器学习，或者是使用基础的卷积神经网络或与循环神经网络相结合。对于GRN这样的有向图数据来说，图神经网络是可以更高效便捷，但是因为在构建数据的过程中更为繁琐所以很少有人使用。其次仅使用基因表达数据的方法中最终的预测效果都会因为无法挖掘更深层的信息而较差。同时GRN存在的大量低度基因节点会导致以往模型的无法得到完全的表达。

发明内容

本发明是为了解决上述现有技术存在的不足之处，提出一种基于有向图卷积的上皮细胞基因调控关系预测方法，以期能更高效准确的预测出基因之间的调控关系，从而有助于生物研究人员更有效地研究基因之间的调控关系。

本发明为达到上述发明目的，采用如下技术方案：

本发明一种基于有向图卷积的上皮细胞基因调控关系预测方法的特点在于，包含如下步骤：

步骤1：获取上皮细胞的基因集合X并进行预处理，形成n×c维的基因表达特征矩阵F

步骤2：构建条件变分自动编码器，并用于对上皮细胞的基因控制图G中的基因节点的特征及其邻居基因节点的特征进行处理，得到n×h维邻域特征矩阵F

步骤3：获取上皮细胞的n×t维基因序列矩阵S，并输入到双向门控循环单元中进行处理，得到n×t维序列特征矩阵F

通过生物学公式对n×t维基因序列矩阵S进行计算，得到n×e维生物特征矩阵F

步骤4：将F

步骤5：构建有向图卷积神经网络并用于对F

步骤6：利用随机梯度下降法对所述有向图卷积神经网络进行训练，并计算交叉熵损失函数以更新网络参数，直到交叉熵损失函数收敛为止，从而得到训练好的上皮细胞基因调控模型；

步骤7：将上皮细胞的基因表达特征矩阵和基因序列矩阵输入上皮细胞基因调控模型中进行预测，并输出基因之间存在调控关系的预测分数，并通过设定的阈值最终判定是否存在调控关系。

本发明所述的基于有向图卷积的上皮细胞基因调控关系的预测方法的特点也在于，所述条件变分自动编码器包括两部分：编码器网络和解码器网络，其中，编码器网络和解码器网络均由多层感知机构成；

以上皮细胞的基因控制图G中任意一个基因节点q的特征i

y

式(1)中，f是激活函数，W

z

式(2)中，NN

所述步骤5包括：

步骤5.1：通过式(3)、式(4)和式(5)分别得到基因控制图G中基因节点q的度信息，包括节点q的一阶度信息A

A

式(3)中，A(q，a)表示基因控制图G中与基因节点q之间存在有控制关系的基因节点a的邻接矩阵；

式(4)中，P

式(5)中，P

步骤5.2：将基因控制图G中所有基因节点的度信息构成基因控制图G的一阶入度矩阵A

步骤5.3：利用式(6)，式(7)和式(8)得到各个基因节点之间调控的潜在关系，包括：一阶度信息的潜在关系Z

式(6)中，

式(7)中，

式(8)中，

步骤5.4：利用式(9)得到基因控制图G中所有基因节点之间存在调控关系的预测分数矩阵R：

式(9)中，α和β是两个权重参数，Concat是级联操作，Softmax是激活函数。

本发明一种电子设备，包括存储器以及处理器的特点在于，所述存储器用于存储支持处理器执行所述基于有向图卷积的上皮细胞基因调控关系预测方法的程序，所述处理器被配置为用于执行所述存储器中存储的程序。

本发明一种计算机可读存储介质，计算机可读存储介质上存储有计算机程序的特点在于，所述计算机程序被处理器运行时执行所述基于有向图卷积的上皮细胞基因调控关系预测方法的步骤。

与现有技术相比，本发明的有益效果在于：

1、GRN中存在大量的低度基因，这些基因的存在会影响现有的神经网络模型表达效果，使模型无法高效的预测基因之间的调控关系，因此本发明使用条件变分自动编码器(CVAE)进行图神经网络的局部增强，通过增强图神经网络中节点的局部特征，以增强图神经网络的表达能力，从而解决难以挖掘低度基因节点在基因调控网络的关联信息的问题；

2、目前先用的神经网络模型仅仅使用细胞基因的表达值作为神经网络模型输入的特征矩阵，但是现有的深度学习计算方法很难再从基因的表达数据中挖掘出更多有用的信息。所以本发明使用循环神经网络提取基因的序列特征值，并通过具有生物学意义的序列特征描述符，以此从不同角度提取基因的特征值，从而提高了上皮细胞基因调控模型的预测性能；

3、基于神经网络的方法中，大部分使用的都是传统的卷积神经网络或者与循环神经网络相结合使用，其次使用神经网络的方法中事先规定了基因的调控方向是从调控因子到靶基因，这会丢失一些环状的调控关系，因此本发明使用有向图卷积神经网络，在高效处理有向图数据的同时更准确的预测基因调控关系。

附图说明

图1为本发明方法流程图；

图2为本发明条件变分自动编码器结构图；

图3为本发明在DREAM4模拟数据中的五个不同的网络中运行十次五折交叉验证下每种方法的受试者工作特征曲线下面积和精确召回曲线下面积图与其他现有方法的比较结果；

图4为本发明在大肠杆菌真实数据中的三个不同的网络中运行十次五折交叉验证下每种方法的受试者工作特征曲线下面积和精确召回曲线下面积图与其他现有方法的比较结果；

图5为本发明在膀胱尿路上皮异常细胞中的受试者工作特征曲线下面积和精确召回曲线下面积图与其他现有方法的比较结果。

具体实施方式

本实施例中、一种基于有向图卷积的基因调控关系预测方法，如图1所示，是按照如下步骤进行：

步骤1：获取上皮细胞的基因集合X并进行预处理，形成n×c维的基因表达特征矩阵F

步骤2：构建条件变分自动编码器，包括两部分：编码器网络和解码器网络，如图2所示。编码器网络和解码器网络都由多层感知机构成；

步骤2.1：以上皮细胞的基因控制图G中任意一个基因节点q的特征i

y

式(1)中，f是激活函数，W

z

式(2)中，NN

如图2所示，条件变分自动编码器的编码器结构和解码器结构都是由两层的多层感知机组成。其中编码器以基因节点q的特征i

步骤3：获取上皮细胞的n×t维基因序列矩阵S，并输入到双向门控循环单元中进行处理，得到n×t维序列特征矩阵F

其中基因的序列矩阵S为基因的碱基序列(如ATCGGCT...)，每一个基因的碱基序列长度并不相同，为了能够通过神经网络(双向门控循环单元)对基因的碱基序列进行处理，我们在较短的序列中进行补0操作，使其与最长的基因序列长度相同，从而使各基因的序列可以拼接成一个矩阵(即序列矩阵S)输入进神经网络中进行处理。

通过生物学公式对n×t维基因序列矩阵S进行计算，得到n×e维生物特征矩阵F

其中生物学公式可以通过收集得到，也就是说目前模型输入的生物特征并不唯一，可以再添加任意数量的生物特征，并且收集得到的该类公式越多生物意义越有效，最终得到的预测结果也会更好。

步骤4：将F

步骤5：构建有向图卷积神经网络；

步骤5.1：通过式(3)、式(4)和式(5)分别得到基因控制图G中基因节点q的度信息，包括节点q的一阶度信息A

A

式(3)中，A为基因控制图G的邻接矩阵，a表示基因控制图G中与基因q之间存在有控制关系的基因；A(q，a)表示基因控制图G中与基因节点q之间存在有控制关系的基因节点a的邻接矩阵；

式(4)中，P

式(5)中，P

其中二阶入度信息的计算公式中，分子部分表示同时指向基因q和基因b的基因k，分母部分表示基因k指向的所有基因节点，公式最终计算得到的就是某一同时指向基因q和基因b的节点占它指向的全部节点数量的比重。

步骤5.2：将基因控制图G中所有基因的节点度信息整理成矩阵，得到基因控制图G节点的一阶入度矩阵A

步骤5.3：利用式(6)，式(7)和式(8)得到各个基因节点之间调控的潜在关系，包括：通过一阶度信息的潜在关系Z

式(6)中，

式(7)中，

式(8)中，

其中使用自环矩阵是为了在聚合邻居节点信息的时候消除无向图结构边的影响，如果不加入自环，那么会因为无向图的边不具有方向的特点，使得节点将无法聚合到自身的特征信息。

步骤5.4：利用式(9)得到预测分数矩阵R，表示基因控制图G中基因之间存在关联的可能性大小：

式(9)中，α和β是两个权重参数，Concat是级联操作，Softmax是激活函数；

预测分数矩阵R中，如果表示任意两个节点之间存在关联的预测分数大于所设定的阈值，则表示相应两个基因之间存在调控关系，且调控方向为调控因子到靶基因，否则，表示相应两个基因之间不存在调控关系；

步骤6：利用随机梯度下降法对所述有向图卷积神经网络进行训练，并计算交叉熵损失函数以更新模型的超参数，直到交叉熵损失函数收敛为止。同时不断修改设定阈值，从而得到训练好的上皮细胞基因调控模型；

步骤7：将上皮细胞的基因与序列特征信息输入上皮细胞基因调控模型中进行预测，并输出基因之间存在调控关系的预测分数，并通过设定的阈值最终判定是否存在确定的调控关系。

本实施例中，一种电子设备，包括存储器以及处理器，该存储器用于存储支持处理器执行上述方法的程序，该处理器被配置为用于执行该存储器中存储的程序。

本实施例中，一种计算机可读存储介质，是在计算机可读存储介质上存储有计算机程序，该计算机程序被处理器运行时执行上述方法的步骤。

为了评估模型的性能，将结果与单个分类器和其他方法进行比较，评价指标为受试者工作特征曲线下面积(AUROC)。AUROC为受试者工作特征曲线下面积，这个面积的数值介于0到1之间，能够直观的评价出分类器的好坏，AUC的值越大，分类器效果越好。对实验采用十次五折交叉验证，并计算平均值作为模型结果。

图3展示了本申请的模型在DREAM4模拟数据集中的结果与其他预测模型Jump3、dynamical GEne Network Inference with Ensemble of trees(dynGENIE3)、Bidirectional eXtreme Gradient Boosting(BiXGBoost)和Multi-source Multi-modelFusion method(MMFGRN)进行比较。实验结果表明，本发明模型在DREAM4的网络1到网络5的五个网络数据中的结果都优于其他预测模型，验证了本模型的有效性和优越性。

图4展示了本申请的模型在大肠杆菌真实数据集中的结果与其他预测模型Jump3、dynGENIE3和BiXGBoost进行比较。在冷压力、热压力和氧气应激压力三个网络数据中，本发明模型在真实数据集上的预测结果也优于其他预测模型，这表明了模型具有优越性，也表示模型能够更准确的预测基因之间的调控关系。

图5展示了本申请模型在膀胱尿路上皮异常细胞数据中实施案例结果，与GENIE-Random Forest(GENIE-RF)、GENIE-Extrem Tree(GENIE-ET)、Network Inference usingMultiple Ensemble Feature Importance Algorithms(NIMEFI)、GEne NetworkInference method based on a Multi-level Strategy(GENIMS)、Partial LeastSquares NET(PLSNET)、NARROMI和Gene Networks Inference based on Projection andLagged Regression(GNIPLR)方法进行比较。实验结果表明，本发明模型在实际案例的膀胱尿路上皮异常细胞数据中的预测结果优于其他预测模型，验证本发明模型在实际应用中的可行性以及有效性。