CN119763665B - 一种基于图表示学习的基因调控网络推断方法及系统 - Google Patents

Abstract

本发明涉及一种基于图表示学习的基因调控网络推断方法及系统，该方法包括：通过单细胞RNA测序数据构建基因调控网络数据集，将已知相互作用的基因对作为正样本，未知相互作用的基因对作为负样本，并随机划分为训练集、验证集和测试集；将训练集中的先验基因调控网络和基因表达谱输入到预构建的基因调控网络推断模型中进行迭代训练，得到基因对预测值；构建包含图对比学习正则项的损失函数，计算预测值与真实值之间的损失值来更新模型参数，当训练轮次达到阈值时，输出模型训练参数，并加载到模型中，对未知相互作用的基因对进行调控关系预测。本发明解决了传统生物实验成本高和现有计算模型精度低的问题，提高了基因调控网络的推断精度。

Description

一种基于图表示学习的基因调控网络推断方法及系统

技术领域

本发明涉及基因对调控关系推断技术领域，尤其涉及一种基于图表示学习的基因调控网络推断方法及系统。

背景技术

在基因表达调控过程中，上游调控基因通过转录、翻译等步骤表达产生的转录因子与下游靶基因的特定区域结合，从而激活或抑制该基因的表达强度，实现对靶基因表达量的调节。基因表达调控对于生物细胞生长发育起到直观重要的作用，基因可通过转录和翻译等过程来调控蛋白质的合成，进而控制细胞生长发育等各项生命活动。基因调控网络通过图结构来描述细胞中基因之间的调控作用，其中节点和边分别表示基因和调控关系。基因调控网络有助于研究细胞各项生命活动的调节机制，解析复杂疾病的发病原理，为复杂疾病的药物研发提供理论依据。另外，基因调控网络可应用于代谢系统的建模与优化、微生物细胞工厂的设计与构建等。因此，合理构建基因调控网络具有重要的理论与应用价值。

然而，传统的分子相互作用实验方法，如微阵列和染色质免疫沉淀测序等，虽可用于研究基因间的相互作用，但往往成本高昂、耗时费力且高度依赖实验条件。为了克服这些局限，学者们尝试运用微分方程、布尔网络等理论对基因调控网络进行机理性建模，以实现网络的构建。然而，在面对大规模基因调控网络时，这些方法往往面临计算效率低下和精度不足的问题。

近年来，随着单细胞高通量测序技术的飞速发展，研究人员在细胞层面获取的基因表达谱数据量呈现爆炸式增长，这为机器学习、深度学习等方法在基因调控网络重构中的应用提供了重要契机。机器学习和深度学习在生物信息学领域备受瞩目，并在基因调控元件的预测与优化、蛋白质功能挖掘与设计、代谢网络分析与设计等方面取得了显著成就。

当前，已有研究利用传统机器学习方法预测基因调控网络中基因间的相互作用。例如，GENIE3模型通过随机森林算法计算调控因子与靶基因的相互作用分数，进而推断基因调控网络的结构；GRNBoost2模型则通过结合GENIE3模型与Arboreto计算框架，进一步提升了运算效率。然而，传统机器学习方法需要手工设计特征，这在一定程度上限制了基因调控网络推断性能的提升。

在此基础上，深度学习方法在基因调控网络有监督的特征表示学习方面得到了深入研究与探讨。这些方法根据先验网络与基因表达谱数据来推断基因对之间的调控关系，实现基因调控网络的重构。例如，CNNC模型利用卷积神经网络提取基因对表达值的特征，推断基因间的相互作用；DeepDRIM模型则通过多层卷积神经网络实现基因对表达值主图像及邻居图像的深层特征提取，进一步推断基因调控网络。STGRNS模型借鉴大语言模型BERT的结构，运用迁移学习思想有效提取基因对表达值的长范围特征，推断基因间潜在的相互作用。

然而，上述方法在特征提取时主要关注基因对之间的互作关系，而忽视了基因调控网络的全局结构。为此，GCNG模型引入图卷积神经网络，提取基因的低维嵌入特征，更有效地推断基因调控网络中潜在的相互作用。进而，GENELINK模型利用图注意力神经网络，通过引入注意力机制优化基因低维嵌入提取过程，进一步提升基因调控网络推断的准确性。然而，由于基因调控网络的稀疏性特点，仅利用全局连接信息进行特征提取仍面临挑战。因此，未来研究需要进一步挖掘基因调控网络中的隐式连接信息，以更有效地提取基因特征和基因间的潜在相关性，从而推动基因调控网络研究的深入发展。

发明内容

为了解决传统生物实验方法成本高且现有计算模型精度低的问题，本发明实施例提供了一种基于图表示学习的基因调控网络推断方法，该方法包括如下步骤：

S1：获取单细胞RNA测序数据，基于单细胞RNA测序数据中基因之间的调控关系，构建基因调控网络数据集，在基因调控网络数据集中，将有相互作用的基因对视为正样本，将相互作用未知的基因对视为负样本，将正样本和负样本随机划分为训练集、验证集和测试集；

S2：将训练集中由正样本产生的先验基因调控网络和基因表达谱输入到提前构建的基因调控网络推断模型中进行迭代训练，得到每次迭代计算得到的基因对预测值；

S3：基于二元交叉熵损失项和图对比学习正则项构建损失函数，利用损失函数计算基因对预测值与真实值之间的损失值来更新模型参数，判断训练轮次是否达到阈值：

若否，返回步骤S2，继续进行迭代计算；

若是，输出模型训练参数，将模型训练参数加载在基因调控网络推断模型中，对相互作用未知的基因对的调控关系进行预测，得到预测结果。

优选地，S2中，得到每次迭代计算得到的基因对预测值的方法如下：

S21：基于所述训练集中由正样本产生的先验基因调控网络的有向图，构建多个子图，每个子图对应一个邻接矩阵，对所有邻接矩阵进行数据处理，得到包括多个隐式连接矩阵的隐式连接张量；

S22：基于归一化后的基因表达谱，将多个隐式连接矩阵依次输入至图卷积神经网络中，对于每个隐式连接矩阵，输出一个基因特征矩阵，对基因特征矩阵进行层归一化处理，将归一化后的基因特征矩阵拼接起来，得到图嵌入张量；

S23：基于图嵌入张量，得到经过强化后的基因特征，计算经过强化后的基因特征中任意两个基因所对应的特征向量的哈达玛积，得到特征融合向量，将特征融合向量通过多层感知机处理，输出基因对的预测值。

优选地，S21中，得到包括多个隐式连接矩阵的隐式连接张量的方法如下：

S211：将所有邻接矩阵进行拼接后得到张量A_s∈{0,1}^5×N×N，将A_s同时输入至两个参数化层，分别得到张量和在每个参数化层中，满足如下关系：

其中，N为基因调控网络中的基因数量，Q⁽ⁱ⁾(j,:,:)为第i个参数化层输出Q⁽ⁱ⁾中第j个矩阵，B为张量Q⁽ⁱ⁾(j,:,:)中的矩阵总数量，为基因调控网络中的基因连接类型A_s(k,:,:)为输入张量A_s中第k个连接类型的子邻接矩阵，为在第i层中归一化后的训练参数；

S212：将两个参数化层输出的张量Q⁽¹⁾和Q⁽²⁾在对应维度上作内积，得到隐式连接张量该过程可表示为：

A_L(j,:,:)＝Q⁽¹⁾(j,:,:)·Q⁽²⁾(j,:,:),j＝1,2,...,B

其中，表示输出隐式连接张量A_L中第j个隐式连接矩阵。

优选地，S22中，图嵌入张量中每个张量X_IE(j,:,:)的表达式如下：

其中，H为图嵌入后的基因特征维度，表示输出隐式连接张量A_L中第j个隐式连接矩阵，N为基因调控网络中的基因数量，I为单位矩阵，为的度矩阵，X_F为归一化后的基因表达谱，W为图卷积神经网络的训练参数，LayerNorm(·)表示层归一化函数，B为图嵌入张量的张量数量。

优选地，S23中，得到经过强化后的基因特征的方法如下：

S231：对图嵌入张量进行维度变换，得到

其中，M为H的平方根，由B个特征张量组成，i＝1,2,...,B，并且每个张量可表示为N个通道的特征矩阵i＝1,2,...,B,j＝1,2,...,N；

S232：对的每一个特征张量中所有N个通道上的特征矩阵分别进行全局最大池化和全局平均池化，得到两个通道特征，将这两个通道特征输入到参数共享中进行卷积操作后，进行逐元素求和，并经过Sigmoid激活函数得到通道注意力

S233：将特征张量通过广播机制与通道注意力特征N_c,i作哈达玛积后，得到输出特征该过程表示为：

S234：沿的通道维度N分别进行全局最大池化及全局平均池化得到两个空间特征，并在通道维度N上进行拼接，进行卷积操作后，得到空间注意力特征该过程表示为：

S235：将通过广播机制与空间注意力特征N_s,i作哈达玛积后，得到强化后的特征图

S236：将所有经过强化后的特征图和原始特征图各自进行残差连接，得到M×M维的输入矩阵，并将输入矩阵展平成向量在维度B上求均值，得到经过强化后的基因特征

其中，Reshape(·)表示维度变换操作。

优选地，基因对的预测值的计算方法如下：

其中，y为基因对的标签值，和为基因调控网络推断模型输出的预测值，n为训练批次中样本数量。

优选地，S3中，所述损失函数的表达式如下：

Loss＝αLoss_bc+βN_gc

其中，Loss表示总损失函数值，Loss_bc表示二元交叉熵损失项，N_gc表示图对比学习正则项，α和β均为加权系数。

优选地，所述二元交叉熵损失项Loss_bc的获取方法如下：

其中，y为基因对的标签值，和为基因调控网络推断模型输出的预测值，n为训练批次中样本数量，y_s为第s个基因对的标签值，为第s个基因对的预测值。

优选地，所述图对比学习正则项N_gc的获取方法如下：

S31：随机丢弃归一化后的基因表达谱X_F与隐式连接矩阵A中部分值后，将其输入到图卷积神经网络中，得到在显式连接下的基因嵌入

其中，A^d为A丢弃部分值后的连接矩阵，I为单位矩阵，为的度矩阵，为X_F丢弃部分值的基因表达谱；

S32：对于经过强化后的基因特征X_GE和基因嵌入X_EE中第m个基因的特征向量u_m＝X_GE(m,:)和v_m＝X_EE(m,:)，计算u_m与v_m之间的第一基因对函数P(u_m,v_m)和v_m与u_m之间的第二基因对函数P(v_m,u_m)，如下：

其中，θ为神经网络拟合参数，τ代表温度系数，1_[m≠n]代表指示函数，当m≠n时为1，当m＝n时为0；

S33：基于第一基因对函数P(u_m,v_m)和第二基因对函数P(v_m,u_m)，得到图对比学习正则项N_gc：

其中，N代表基因的数量。

基于同一发明构思，本发明实施例还提供了一种基于图表示学习的基因调控网络推断系统，该系统用于实现上述的基于图表示学习的基因调控网络推断方法，具体包括：

基因调控网络数据集构建模块，用于获取单细胞RNA测序数据，基于单细胞RNA测序数据中基因之间的调控关系，构建基因调控网络数据集，在基因调控网络数据集中，将有相互作用的基因对视为正样本，将相互作用未知的基因对视为负样本，将正样本和负样本随机划分为训练集、验证集和测试集；

基因调控网络推断模型构建模块，用于将训练集中由正样本产生的先验基因调控网络和基因表达谱输入到提前构建的基因调控网络推断模型中进行迭代训练，得到每次迭代计算得到的基因对预测值；

基因对调控关系预测模块，用于基于二元交叉熵损失项和图对比学习正则项构建损失函数，利用损失函数计算基因对预测值与真实值之间的损失值来更新模型参数，若训练轮次达到阈值时，输出模型训练参数，将模型训练参数加载在基因调控网络推断模型中，对相互作用未知的基因对的调控关系进行预测，得到预测结果。

从以上技术方案可以看出，本发明申请具有以下有益效果：

第一，解决传统生物实验高的问题：本发明主要依赖于单细胞RNA测序数据，而非传统的生物实验方法，如微阵列和染色质免疫沉淀测序等。这样可以显著降低实验成本，单细胞RNA测序数据可以通过公开数据集或者高通量测序技术获得，相对于传统的实验方法更为经济；通过构建基因调控网络数据集，充分利用了单细胞RNA测序数据中的基因调控关系，将有相互作用的基因对视为正样本，将相互作用未知的基因对视为负样本，这种方法提高了数据的利用率，使得在有限的数据下能够获得更多的信息。

第二，解决现有基因对关系计算模型精度低下的问题：

(1)在构建基因调控网络推断模型时，引入了隐式连接的概念，通过构建多个子图和隐式连接矩阵，能够捕捉到基因调控网络中更为复杂的相互作用关系，从而提高模型的推断精度；

(2)将图卷积神经网络(GCN)应用于基因调控网络推断中，能够学习到基因的低维嵌入表示，这些表示不仅保留了基因的原始特征信息，还融入了其在网络中的上下文信息，从而增强了模型的表达能力；

(3)使用通道注意力和空间注意力机制来强化基因特征，以及通过全局最大池化和全局平均池化等技术来提取关键特征，这些步骤有助于去除冗余信息，保留最具代表性的特征，从而提升模型的推断精度；

(4)综合考虑基因对预测值与真实值之间的差异以及基因在网络中的位置关系，构建包含图对比学习正则项的损失函数，引导模型学习到更为准确的基因调控关系。

附图说明

为了更清楚地说明本发明实施案例或现有技术中的技术方案，下边将对实施例中所需要使用的附图做简单说明，通过参考附图会更清楚的理解本发明的特征和优点，附图是示意性的而不应该理解为对本发明进行任何限制，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，可以根据这些附图获得其他的附图。其中：

图1为本发明的实施例一中提供的一种基于图表示学习的基因调控网络推断方法的流程图；

图2为本发明的实施例一中提供的基于图表示学习的基因调控网络推断方法的具体流程图；

图3为本发明的实施例一中提供的构建基因调控网络数据集以及样本集划分示意图；

图4中A表示提取基因调控网络结构中的隐式连接，B表示提取基因调控网络在隐式连接下的基因特征，C表示强化基因调控网络中的基因特征提取，D表示输出基因对调控关系的预测结果；

图5是利用注意力机制模块得到强化后的特征图的流程图；

图6是获取图对比学习正则项的流程图；

图7中A表示在不同规模的金标准网络框架下，本发明所提出的基因调控网络推断模型(简称为“GRLGRN模型”)在七种细胞类型中的AUROC指标，B表示GRLGRN与其他对比模型在数据集上的AUROC指标的分布情况；

图8中A表示在不同规模的金标准网络框架下，本发明所提出的基因调控网络推断模型(简称为“GRLGRN模型”)在七种细胞类型中的AUPRC指标，B表示GRLGRN与其他对比模型在数据集上的AUPRC指标的分布情况；

图9为实施例二中提供的一种基于图表示学习的基因调控网络推断系统的框图；

说明书附图标记说明：

100、基因调控网络数据集构建模块；200、基因调控网络推断模型构建模块；300、基因对调控关系预测模块。

具体实施方式

为使本发明实施例的目的、技术方案与优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

实施例一

如图1和图2所示，本发明实施例提出一种基于图表示学习的基因调控网络推断方法，该方法包括：

S1：获取单细胞RNA测序数据，基于单细胞RNA测序数据中基因之间的调控关系，构建基因调控网络数据集，在基因调控网络数据集中，将有相互作用的基因对视为正样本，将相互作用未知的基因对视为负样本，将正样本和负样本随机划分为训练集、验证集和测试集；

S2：将训练集中由正样本产生的先验基因调控网络和基因表达谱输入到提前构建的基因调控网络推断模型中进行迭代训练，得到每次迭代计算得到的基因对预测值；

S3：基于二元交叉熵损失项和图对比学习正则项构建损失函数，利用损失函数计算基因对预测值与真实值之间的损失值来更新模型参数，判断训练轮次是否达到阈值：

若否，返回步骤S2，继续进行迭代计算；

若是，输出模型训练参数，将模型训练参数加载在基因调控网络推断模型中，对相互作用未知的基因对的调控关系进行预测，得到预测结果。

从上述技术方案可知，本发明通过利用单细胞RNA测序数据，构建了基因调控网络数据集，并在此基础上进行模型训练和预测，显著减少了直接进行生物实验的需求，从而降低了研究成本；通过引入图表示学习，将基因调控网络转化为图结构数据，并利用图对比学习正则项优化损失函数，来增强模型的泛化能力。同时，通过迭代训练和损失函数的优化，使得模型能够更准确地预测基因对之间的调控关系。这种方法不仅考虑了基因表达谱的信息，还利用了基因调控网络的拓扑结构信息，从而提高了模型的推断精度。

具体地，在步骤S1中，在公开数据中获取人胚胎干细胞(hESC)、人成熟肝细胞(hHEP)、小鼠树突状细胞(mDC)、小鼠胚胎干细胞(mESC)、小鼠红系造血干细胞(mHSC-E)、粒细胞-单核细胞谱系小鼠造血干细胞(mHSC-GM)和淋巴样谱系的小鼠造血干细胞(mHSCL)总计7种细胞类型的单细胞RNA测序(scRNA-seq)数据。其中，每种细胞类型均对应cell-type-specific ChIP-seq、non specific ChIP-seq和STRING数据库提供的3个不同类型的标准网络，7种细胞共计21个标准网络。针对每个标准网络，根据基因表达值的变化幅度，筛选出了含转录因子在内的表达值变化最为显著的前500个和前1000个基因。基于这些关键基因及其间的调控关系，构成两种不同规模的基因调控网络(共计42个)，以及描述每个基因在时间序列上表达强度的基因表达谱矩阵。

进而，如图3所示，将有相互作用的基因对视为正样本，标签记为1，将相互作用未知的基因对视为负样本，标签记为0。将正样本和负样本均按照3:1:1的比例随机划分为训练集、验证集和测试集，通过这样划分出来的正负样本被称为困难样本，目的是为了以加速训练过程、降低计算复杂性和更好的提取基因特征。

优选地，S2中，得到每次迭代计算得到的基因对预测值的方法如下：

S21：如图4中A所示，基于所述训练集中由正样本产生的先验基因调控网络的有向图，构建包括“调控基因-目标基因”的有向子图及其反向图、“调控基因-调控基因”的有向子图及其反向图、基因自连接子图这5个子图，每个子图对应一个邻接矩阵，对所有邻接矩阵进行数据处理，得到包括多个隐式连接矩阵的隐式连接张量；

S22：基于归一化后的基因表达谱，将多个隐式连接矩阵依次输入至图卷积神经网络中，对于每个隐式连接矩阵，输出一个基因特征矩阵，对基因特征矩阵进行层归一化处理，将归一化后的基因特征矩阵拼接起来，得到图嵌入张量；

S23：基于图嵌入张量，得到经过强化后的基因特征，计算经过强化后的基因特征中任意两个基因所对应的特征向量的哈达玛积，得到特征融合向量，将特征融合向量通过多层感知机处理，输出基因对的预测值。

优选地，S21中，得到包括多个隐式连接矩阵的隐式连接张量的方法如下：

S211：将所有邻接矩阵进行拼接后得到张量为同时考虑不同类型的隐式连接信息，将A_s同时输入至两个独立的参数化层，分别得到张量和在每个参数化层中，满足如下关系：

其中，N为基因调控网络中的基因数量，Q⁽ⁱ⁾(j,:,:)为第i个参数化层输出Q⁽ⁱ⁾中第j个矩阵，B为张量Q⁽ⁱ⁾(j,:,:)中的矩阵总数量，为基因调控网络中的基因连接类型A_s(k,:,:)为输入张量A_s中第k个连接类型的子邻接矩阵，为在第i层中归一化后的训练参数；的归一化过程满足：

S212：将两个参数化层输出的张量Q⁽¹⁾和Q⁽²⁾在对应维度上作内积，得到隐式连接张量该过程表示为：

A_L(j,:,:)＝Q⁽¹⁾(j,:,:)·Q⁽²⁾(j,:,:),j＝1,2,...,B

其中，表示输出隐式连接张量A_L中第j个隐式连接矩阵。

进一步地，在步骤S22中，如图4中B所示，将归一化后的基因表达谱记为其中N为基因调控网络中的基因数量，D为基因表达谱的维度，基于归一化后的基因表达谱，将B个隐式连接矩阵依次输入至图卷积神经网络中，对于每个隐式连接矩阵，输出一个基因特征矩阵，对基因特征矩阵进行层归一化处理，将归一化后的基因特征矩阵拼接起来，得到图嵌入张量其中的每个张量X_IE(j,:,:)的表达式如下：

其中，H为图嵌入后的基因特征维度，表示输出隐式连接张量A_L中第j个隐式连接矩阵，N为基因调控网络中的基因数量，I为单位矩阵，为的度矩阵，X_F为归一化后的基因表达谱，W为图卷积神经网络的训练参数，LayerNorm(·)表示层归一化函数，B为图嵌入张量的张量数量。

进一步地，在步骤S23中，如图4中C所示，基于图嵌入张量得到经过强化后的基因特征的方法如下：

S231：对图嵌入张量进行维度变换，得到

其中，M为H的平方根，由B个特征张量组成，i＝1,2,...,B，并且每个张量可表示为N个通道的特征矩阵

i＝1,2,...,B,j＝1,2,...,N；

S232：如图5所示，将输入至由通道注意力和空间注意力构成的卷积块注意力模块(CBAM)，对的每一个特征张量中所有N个通道上的特征矩阵分别进行全局最大池化和全局平均池化，得到两个通道特征，将这两个通道特征输入到参数共享中进行卷积操作。其中，该神经网络由卷积操作Conv₁、非线性激活函数(ReLU)以及卷积操作Conv₂组成。两个通道特征经过该神经网络后，进行逐元素求和，并经过Sigmoid激活函数得到通道注意力

其中，将这两个通道特征输入到参数共享中进行卷积操作满足如下：

Conv_c(·)＝Conv₂(ReLU(Conv₁(·)))

S233：将特征张量通过广播机制与通道注意力特征N_c,i作哈达玛积后，得到输出特征该过程表示为：

S234：沿的通道维度N分别进行全局最大池化及全局平均池化得到两个空间特征，并在通道维度N上进行拼接，得到一个大小为2×M×M的空间分数。进行卷积操作Conv₃后，利用Sigmoid激活函数得到空间注意力特征该过程表示为：

S235：将通过广播机制与空间注意力特征N_s,i作哈达玛积后，得到强化后的特征图

S236：将所有经过强化后的特征图和原始特征图各自进行残差连接，得到M×M维的输入矩阵，并将输入矩阵展平成向量在维度B上求均值，得到经过强化后的基因特征

其中，Reshape(·)表示维度变换操作。卷积操作Conv₁、Conv₂和Conv₃的参数设置如表1所示。

表1

卷积操作名称	参数设置
		Conv1	卷积核参数数量＝1×N×N/r，步长＝1
Conv2	卷积核参数数量＝1×N/r×N，步长＝1
		Conv3	卷积核参数数量＝1×2×M×M，步长＝1

进一步地，如图4中D所示，计算经过强化后的基因特征中任意两个基因i和j所对应的特征向量X_GE(i,:)和X_GE(j,:)的哈达玛积，得到特征融合向量，将特征融合向量输入至非线性激活函数ReLU，通过多层感知机处理，得到基因对的预测值如下：

其中，y为基因对的标签值，和为基因调控网络推断模型输出的预测值，n为训练批次中样本数量。上述的多层感知机由两个全连接层构成，分别含有ReLU函数和Sigmoid函数。

进一步地，在步骤S3中，为避免由Non-Specific CHIP-seq提供的基准网络极度稀疏导致提取信息不充分的情况出现，采用如图6所示的图对比学习框架来最大化X_GE和X_EE中同一基因特征向量的一致性，得到图对比学习正则项N_gc，具体方法如下：

S31：随机丢弃归一化后的基因表达谱X_F与隐式连接矩阵A中部分值后，将其输入到图卷积神经网络中，得到在显式连接下的基因嵌入

其中，A^d为A丢弃部分值后的连接矩阵，I为单位矩阵，为的度矩阵，为X_F丢弃部分值的基因表达谱；

S32：对于经过强化后的基因特征X_GE和基因嵌入X_EE中第m个基因的特征向量u_m＝X_GE(m,:)和v_m＝X_EE(m,:)，计算u_m与v_m之间的第一基因对函数P(u_m,v_m)和v_m与u_m之间的第二基因对函数P(v_m,u_m)，如下：

其中，θ为神经网络拟合参数，τ代表温度系数，1_[m≠n]代表指示函数，当m≠n时为1，当m＝n时为0；模型超参数见表2。

表2

S33：基于第一基因对函数P(u_m,v_m)和第二基因对函数P(v_m,u_m)，得到图对比学习正则项N_gc：

其中，N代表基因的数量。

基于上述的图对比学习正则项N_gc，定义所述基因调控网络推断模型在训练阶段中采用的损失函数Loss的表达式如下：

Loss＝αLoss_bc+βN_gc

其中，Loss表示总损失函数值，Loss_bc表示二元交叉熵损失项，N_gc表示图对比学习正则项，α和β均为加权系数。其中，所述二元交叉熵损失项Loss_bc的获取方法如下：

其中，y为基因对的标签值，和为基因调控网络推断模型输出的预测值，n为训练批次中样本数量，y_s为第s个基因对的标签值，为第s个基因对的预测值。

采用Adam优化策略对模型参数实施迭代更新，期间密切监测模型在验证集上的性能表现，据此精细调整超参数配置，以期达到最优的模型效能。随后，统一在相同的数据集上对所有模型进行训练与测试，并依据AUROC与AUPRC两项评估指标进行综合比较。

其中，AUROC指标通过计算在不同阈值下的真正例率(TPR，也称为召回率)与假正例率(FPR)来绘制ROC曲线，并计算该曲线下的面积来得出。AUROC主要用于衡量模型区分正负样本的能力，其取值范围在0到1之间，AUROC值越接近1，模型的性能越好。

AUPRC指标通过计算在不同召回率下的平均精度来绘制PR曲线，并计算该曲线下的面积来得出。AUPRC在处理不平衡数据集的二分类问题时，能够更好地反映模型在少数类上的性能表现，其取值范围在0到1之间，数值越高表示模型的性能越好。

在现有模型体系中，诸如GNNLINK，其通过多层图卷积神经网络高效整合基因调控网络中的高阶连接信息；GENELINK，则运用图注意力神经网络计算基因连接间的注意力权重，进而推断基因调控网络；STGRNS，创造性地将基因表达谱进行分割与拼接，并引入BERT的Encoder架构以推断基因调控网络；GENIE3，以随机森林算法为核心，对调控因子与靶基因的相互作用分数进行排序；GRNBoost2，巧妙结合GENIE3与Arboreto计算框架，显著加速了基因调控网络的推理进程；以及GNE模型，利用多层感知机与基因调控网络的全局结构信息，对基因进行深入的聚类分析。这些模型各具特色，构成了当前基因调控网络推断领域的代表性方法。

针对AUROC指标，本发明所提出的GRLGRN模型与其他模型的性能对比，通过可视化手段在图7与图8中得以直观展现。具体而言，图7中A与图8中A分别展示了在不同规模的金标准网络框架下，GRLGRN模型在七种细胞类型中的AUROC与AUPRC指标。数据显示，GRLGRN模型在绝大多数数据集上均展现出了超越其他模型的卓越性能。在对应的三种不同类型的金标准网络数据集中，GRLGRN模型分别在85.7％(12/14)、57.1％(8/14)和92.8％(13/14)的数据集上优于排名第二的模型，性能提升幅度分别达到约6.93％、4.3％和6.43％。就AUPRC指标而言，GRLGRN模型在73.8％(31/42)的数据集上表现优于次优模型，且在金标准网络下的平均提升幅度至少达到约7.84％。

图7中B与图8中B则分别揭示了GRLGRN模型与其他对比模型在数据集上的AUROC和AUPRC指标的分布情况，基准数据集涵盖了TFs+500(上方)和TFs+1000(下方)的真值网络。最终，表3汇总了GRLGRN模型及其他六个对比模型在42个数据集上的AUROC和AUPRC得分的平均值，进一步验证了GRLGRN模型在基因调控网络推断领域的显著优势。

表3

实施例二

基于和实施例一方法的同一发明构思，本发明实施例还提供了一种基于图表示学习的基因调控网络推断系统，该系统用于实现上述实施例一中的基于图表示学习的基因调控网络推断方法，如图9所示，系统具体包括：

基因调控网络数据集构建模块100，用于获取单细胞RNA测序数据，基于单细胞RNA测序数据中基因之间的调控关系，构建基因调控网络数据集，在基因调控网络数据集中，将有相互作用的基因对视为正样本，将相互作用未知的基因对视为负样本，将正样本和负样本随机划分为训练集、验证集和测试集；

基因调控网络推断模型构建模块200，用于将训练集中由正样本产生的先验基因调控网络和基因表达谱输入到提前构建的基因调控网络推断模型中进行迭代训练，得到每次迭代计算得到的基因对预测值；

基因对调控关系预测模块300，用于基于二元交叉熵损失项和图对比学习正则项构建损失函数，利用损失函数计算基因对预测值与真实值之间的损失值来更新模型参数，若训练轮次达到阈值时，输出模型训练参数，将模型训练参数加载在基因调控网络推断模型中，对相互作用未知的基因对的调控关系进行预测，得到预测结果。

本实施例的一种基于图表示学习的基因调控网络推断系统，用于实现前述的基于图表示学习的基因调控网络推断方法，因此基于图表示学习的基因调控网络推断系统中的具体实施方式可见前文基于图表示学习的基因调控网络推断方法的实施例部分，例如，基因调控网络数据集构建模块100、基因调控网络推断模型构建模块200和基因对调控关系预测模块300，分别用于实现上述基于图表示学习的基因调控网络推断方法中步骤S1、S2和S3，所以，其具体实施方式可以参照相应的各个部分实施例的描述，为了避免冗余，在此不再赘述。

本领域内的技术人员应明白，本申请的实施例可提供为方法、系统、或计算机程序产品。因此，本申请可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。

本申请是参照根据本申请实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。

这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。

显然，上述实施例仅仅是为清楚地说明所作的举例，并非对实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式变化或变动。这里无需也无法对所有的实施方式予以穷举。而由此所引申出的显而易见的变化或变动仍处于本发明创造的保护范围之中。

Claims (9)

1.一种基于图表示学习的基因调控网络推断方法，其特征在于，包括以下步骤：

S1：获取单细胞RNA测序数据，基于单细胞RNA测序数据中基因之间的调控关系，构建基因调控网络数据集，在基因调控网络数据集中，将有相互作用的基因对视为正样本，将相互作用未知的基因对视为负样本，将正样本和负样本随机划分为训练集、验证集和测试集；

S2：将训练集中由正样本产生的先验基因调控网络和基因表达谱输入到提前构建的基因调控网络推断模型中进行迭代训练，得到每次迭代计算得到的基因对预测值；

S3：基于二元交叉熵损失项和图对比学习正则项构建损失函数，利用损失函数计算基因对预测值与真实值之间的损失值来更新模型参数，判断训练轮次是否达到阈值：

若否，返回步骤S2，继续进行迭代计算；

若是，输出模型训练参数，将模型训练参数加载在基因调控网络推断模型中，对相互作用未知的基因对的调控关系进行预测，得到预测结果；

其中，S2中，得到每次迭代计算得到的基因对预测值的方法如下：

S21：基于所述训练集中由正样本产生的先验基因调控网络的有向图，构建包括调控基因-目标基因的有向子图及其反向图、调控基因-调控基因的有向子图及其反向图、基因自连接子图这些子图，每个子图对应一个邻接矩阵，对所有邻接矩阵进行数据处理，得到包括多个隐式连接矩阵的隐式连接张量；

S22：基于归一化后的基因表达谱，将多个隐式连接矩阵依次输入至图卷积神经网络中，对于每个隐式连接矩阵，输出一个基因特征矩阵，对基因特征矩阵进行层归一化处理，将归一化后的基因特征矩阵拼接起来，得到图嵌入张量；

S23：基于图嵌入张量，使用通道注意力和空间注意力机制来强化基因特征，得到经过强化后的基因特征，计算经过强化后的基因特征中任意两个基因所对应的特征向量的哈达玛积，得到特征融合向量，将特征融合向量通过多层感知机处理，输出基因对的预测值。

2.根据权利要求1的基于图表示学习的基因调控网络推断方法，其特征在于，S21中，得到包括多个隐式连接矩阵的隐式连接张量的方法如下：

S211：将所有邻接矩阵进行拼接后得到张量A_s∈{0,1}^5×N×N，将A_s同时输入至两个参数化层，分别得到张量和在每个参数化层中，满足如下关系：

其中，N为基因调控网络中的基因数量，Q⁽ⁱ⁾(j,:,:)为第i个参数化层输出Q⁽ⁱ⁾中第j个矩阵，B为张量Q⁽ⁱ⁾(j,:,:)中的矩阵总数量，为基因调控网络中的基因连接类型A_s(k,:,:)为输入张量A_s中第k个连接类型的子邻接矩阵，为在第i层中归一化后的训练参数；

S212：将两个参数化层输出的张量Q⁽¹⁾和Q⁽²⁾在对应维度上作内积，得到隐式连接张量该过程可表示为：

A_L(j,:,:)＝Q⁽¹⁾(j,:,:)·Q⁽²⁾(j,:,:),j＝1,2,...,B

其中，表示输出隐式连接张量A_L中第j个隐式连接矩阵。

3.根据权利要求1的基于图表示学习的基因调控网络推断方法，其特征在于，S22中，图嵌入张量中每个张量X_IE(j,:,:)的表达式如下：

其中，H为图嵌入后的基因特征维度，表示输出隐式连接张量A_L中第j个隐式连接矩阵，N为基因调控网络中的基因数量，I为单位矩阵，为的度矩阵，X_F为归一化后的基因表达谱，W为图卷积神经网络的训练参数，LayerNorm(·)表示层归一化函数，B为图嵌入张量的张量数量。

4.根据权利要求1的基于图表示学习的基因调控网络推断方法，其特征在于，S23中，得到经过强化后的基因特征的方法如下：

S231：对图嵌入张量进行维度变换，得到

其中，M为H的平方根，由B个特征张量组成，i＝1,2,...,B，并且每个张量可表示为N个通道的特征矩阵

S232：对的每一个特征张量中所有N个通道上的特征矩阵分别进行全局最大池化和全局平均池化，得到两个通道特征，将这两个通道特征输入到参数共享中进行卷积操作后，进行逐元素求和，并经过Sigmoid激活函数得到通道注意力

S233：将特征张量通过广播机制与通道注意力特征N_c,i作哈达玛积后，得到输出特征该过程表示为：

S234：沿的通道维度N分别进行全局最大池化及全局平均池化得到两个空间特征，并在通道维度N上进行拼接，进行卷积操作后，得到空间注意力特征该过程表示为：

S235：将通过广播机制与空间注意力特征N_s,i作哈达玛积后，得到强化后的特征图

S236：将所有经过强化后的特征图和原始特征图各自进行残差连接，得到M×M维的输入矩阵，并将输入矩阵展平成向量在维度B上求均值，得到经过强化后的基因特征

其中，Reshape(·)表示维度变换操作。

5.根据权利要求1的基于图表示学习的基因调控网络推断方法，其特征在于，基因对的预测值的计算方法如下：

其中，y为基因对的标签值，和为基因调控网络推断模型输出的预测值，n为训练批次中样本数量。

6.根据权利要求1的基于图表示学习的基因调控网络推断方法，其特征在于，S3中，所述损失函数的表达式如下：

Loss＝αLoss_bc+βN_gc

其中，Loss表示总损失函数值，Loss_bc表示二元交叉熵损失项，N_gc表示图对比学习正则项，α和β均为加权系数。

7.根据权利要求6的基于图表示学习的基因调控网络推断方法，其特征在于，所述二元交叉熵损失项Loss_bc的获取方法如下：

其中，y为基因对的标签值，和为基因调控网络推断模型输出的预测值，n为训练批次中样本数量，y_s为第s个基因对的标签值，为第s个基因对的预测值。

8.根据权利要求6的基于图表示学习的基因调控网络推断方法，其特征在于，所述图对比学习正则项N_gc的获取方法如下：

S31：随机丢弃归一化后的基因表达谱X_F与隐式连接矩阵A中部分值后，将其输入到图卷积神经网络中，得到在显式连接下的基因嵌入

其中，A^d为A丢弃部分值后的连接矩阵，I为单位矩阵，D为的度矩阵，为X_F丢弃部分值的基因表达谱；

S32：对于经过强化后的基因特征X_GE和基因嵌入X_EE中第m个基因的特征向量u_m＝X_GE(m,:)和v_m＝X_EE(m,:)，计算u_m与v_m之间的第一基因对函数P(u_m,v_m)和v_m与u_m之间的第二基因对函数P(v_m,u_m)，如下：

其中，θ为神经网络拟合参数，τ代表温度系数，1_[m≠n]代表指示函数，当m≠n时为1，当m＝n时为0；

S33：基于第一基因对函数P(u_m,v_m)和第二基因对函数P(v_m,u_m)，得到图对比学习正则项N_gc：

其中，N代表基因的数量。

9.一种基于图表示学习的基因调控网络推断系统，其特征在于，系统用于实现权利要求1至8任意一项所述的基于图表示学习的基因调控网络推断方法，具体包括：

基因调控网络数据集构建模块，用于获取单细胞RNA测序数据，基于单细胞RNA测序数据中基因之间的调控关系，构建基因调控网络数据集，在基因调控网络数据集中，将有相互作用的基因对视为正样本，将相互作用未知的基因对视为负样本，将正样本和负样本随机划分为训练集、验证集和测试集；

基因调控网络推断模型构建模块，用于将训练集中由正样本产生的先验基因调控网络和基因表达谱输入到提前构建的基因调控网络推断模型中进行迭代训练，得到每次迭代计算得到的基因对预测值；

基因对调控关系预测模块，用于基于二元交叉熵损失项和图对比学习正则项构建损失函数，利用损失函数计算基因对预测值与真实值之间的损失值来更新模型参数，若训练轮次达到阈值时，输出模型训练参数，将模型训练参数加载在基因调控网络推断模型中，对相互作用未知的基因对的调控关系进行预测，得到预测结果。