DOI： 10.3390/ijms24032595

期刊： International Journal of Molecular Sciences

中科院分区：2区生物学 

影像因子：5.6↓ 0.608

作者： Yongqing Zhang; Maocheng Wang; Zixuan Wang; Yuhang Liu; Shuwen Xiong; et al

出版日期： 2023-01-30

网址: https://www.mdpi.com/1422-0067/24/3/2595/pdf?version=1675053613

Github(数据集):https://github.com/ZhangLab312/MetaSEM 

摘要 

基因调控网络(gene regulatory network, GRNs)中的调控因子对细胞状态的识别至关重要。然而，基于 scRNA-seq 数据的 GRN 推断存在高维数和稀疏性等问题，并且需要更多的标签数据。因此， 我们提出了一个元学习 GRN 推理框架来识别调节因素。具体来说，元学习解决了高维稀疏数据特 征带来的参数优化问题。此外，采用少样本解决方案解决标签数据不足的问题。在模型中嵌入了 结构方程模型(SEM)，以识别重要的调控因子。我们将参数优化策略融入到双层优化中，提取出符 合 GRN 推理的特征。这种独特的设计使我们的模型对小规模数据具有鲁棒性。通过研究 GRN 推 断任务，我们证实了所选择的调控因子与基因表达特异性密切相关。我们进一步分析推断的 GRN，以找到细胞类型识别的重要调节因子。大量的实验结果表明，我们的模型有效地捕获了单细胞 GRN 推理中的调节器。最后，可视化结果验证了所选调控因子对细胞类型识别的重要性。

关 键词:元学习;基因调控网络推断;结构方程模型;双层的优化

1. 介绍

基因调控网络(GRNs)的推断可以更好地理解转录调控及其在细胞类型识别中的作用。基于 scRNA-seq  的 GRN 推断模型在癌症治疗[1]、细胞稳态识别[2]和单细胞多组学研究 [3]中取得了成功。然而，单细胞 RNA  测序技术存在技术噪声[4]、基因高变异性[5]、批 效应[6]等诸多局限性。那些活性在不同细胞类型之间高度可变，并预测有一小部分重要 细胞类型的必要调控因子，仍然需要更多的关注。因此，利用计算方法推断基因调控网 络来研究细胞特异性现象是生物信息学中一个具有挑战性的问题。

最近，深度学习为基于共表达的单细胞 GRN 推理带来了新的解决方案[7,8][9]。有 监督的方法有两种方法来推断 grn。一种直接设置一个确定的 ground-truth 标签作为模型 收敛目标，如 DGRN[10]和 Deep DRIM[11]。另一种是在模型中嵌入一个 ground-truth 网 络，如 GRGNN[12]和 scSGL[13]。这些模型可以预测高维 scRNA-seq 数据中潜在的基 因调控关系。然而，有监督的方法只能应用于具有标签数据的一般任务。无监督方法促  进了没有 ground-truth 标签的 GRN 推理，如 VEGA[14]和 Deep SEM[15]。SCODE[16]是 一种基于线性有序微分方程的机器学习算法。GENIE3[17]和  GRNBoost2[18]也可以在没 有标签的情况下完成 GRN 推理。但是，高度稀疏的数据特征问题仍然需要更好的解决。

元学习擅长在使用小样本[19]时解决参数初始化和数据标签不足的问题。元学习 模型由一个基础学习器和一个元学习器组成。基础学习器通过表征学习提取特征信息。元学习器通过学习基础学习器的参数， 引导基础学习器完成训练任务。元学习者综合 所有模块的训练经验， 并为新任务训练提供初始参数。通过引入元学习， 可以提高深 度学习模型的泛化能力，解决高度稀疏的数据特征。例如， 傅坤等人利用迁移学习来 缓解元学习在小样本任务[20]上的训练问题。Arkabandhu  Chowdhury  等人提出了一种 在小样本[21]上完成数据分类的元学习方法。Zitian  Chen  等通过元学习解决了图像变 形问题[22]中 one-shot 学习的低数据样本问题。然而，现有的元学习方法应该充分适应 单细胞 GRN 推理任务。

为了解决上述问题，我们提出了一个名为 Meta SEM  的元学习框架(图  1)，从 scRNA-seq 中推断 GRN。具体来说，我们采用元学习来优化每个模块的参数，以学习高维数据特征。接下来，我们采用了元解码器为编码器提供伪数据标签。然后，我们将从编码器中提取的特征向量放入训练过程中。之后，考虑到基因调控关系是一个内源变量，基因表达信息是一个外源变量，我们初始化一个结构方程模型(SEM)[23]邻接矩阵作为 GRN  层。我们将该矩阵视为调控权重矩阵， 并将其嵌入到元解码器中。最后， 使用双 层优化对所有参数进行优化。

图 1 所示。 (A)MetaSEM 概述:元解码器提取监管关系以输出伪数据标签。编码器将数据特征转换为特征向量。GRN 层是专门为嵌入 SEM 矩阵而设计的层。红色箭头表示外层循环，黄色箭头表示内层循环。MetaSEM 通过超参数优化，基于梯度实现了内外环的集成。θF 表示编码器的超参数，A θ 表示元解码器的超参数。(B)通过分析 SEM 矩阵，MetaSEM 具有三个主要功能:调节因子识别、GRN 可视化和细胞类型识别。

在这项工作中，我们从以下几个方面验证了 MetaSEM 的可靠性。我们首先将 MetaSEM 的模型性能与几种最先进的方法进行了比较。实验结果表明， MetaSEM 在 EPR、AUPR 和 AUROC 方面明显优于现有方法。接下来，为了分析 MetaSEMwe 的鲁棒性，我们探索了不同数据尺度下单细胞数据中的基本数据特征。然后， 我们在大量 RNA -seq 数据集上生成细胞类型特异性 grn，以供进一步研究。Pearson 相关分析和基因表达数据分析表明，grn 具有细胞特异性。最后，HNSCC 数据集中细胞类型特异性 grn的可视化显示了调节因子在识别细胞类型中的重要性。

2. 结果与讨论

2.1. 与现有方法的比较

为了验证 MetaSEM 的性能，我们将该模型与 BEELINE 数据集上的 DeepSEM[15]、 DGRN[10]、GENIE3[17]和 PIDC[2]四种方法进行了比较。如表 1 所示，在 1000 个基因数 据集上，MetaSEM 在三个评价指标上优于现有方法。MetaSEM 在 mHSC-L、mHSC-G 和 mHSC-E 数据集中的 EPR 分别为 1.36、1.41 和 1.21。在 mHSCs 数据集中，MetaSEM 的 EPR 平均比 DeepSEM 高 0.15。在 AUPR 和 AUROC 中，MetaSEM 比 DGRN 高 0.41。此外，由于 GENIE3 和 PIDC 是无监督的机器学习方法，因此远不如深度学习方法。结果表明，MetaSEM 能够有效记忆基因调控关系，并以此关系指导模型提取必要信息。

表 1。 四种竞争方法在 1000 个基因数据集和 500 个基因数据集上的性能比较。

图 2。 我们的模型在不同数据规模上的鲁棒性。每一列对应于一个单元格的子数据集(左 :mHSC - L，中间:mHSC-GM，右 :mHSC-E)，每一行对应于一个评价指标(上 :EPR，中 :AUPR，下:AUROC)。图中红色区域为标准差选择的结果，图中蓝色区域为随机选择的结果

2.3. MetaSEM 显示 GRN 特异性与基因表达有关

为了验证 MetaSEM  是否捕获了特定信息，我们分析了推断的 GRN  和基因表达数据。我们基于 8 个 HNSCC 子数据集生成了细胞类型特异性 GRN。图 3 显示了不同 grn之间的 Pearson  相关系数的热图。在这个矩阵中，做了两个观测。首先，平均相关系数小于 0.1，表明推断的 grn  相关性很低。第二，如图所示，相关系数最高的是成纤维细胞。这有两个原因:癌症数据集包含一些癌症相关的成纤维细胞(CAFs)[24]，成纤维细胞与癌症亚群之间的相关性最高，这与[25]的发现相对应。

图 3。不同细胞型 grn 的 Pearson 相关性。矩阵中的每个元素表示对应于两个不同细胞的 GRN 的Pearson 相关性。我们没有展示 p 值大于 0.05 的结果。

鉴于图 3 所反映的现象，我们进一步分析了数据集中基因表达的差异。如图 4 所示，成纤维细胞和内皮细胞亚群的基因表达数据与癌症亚群差异不大。然而，B 细胞亚群和肥大亚群的结果显示出显著差异。成纤维细胞亚群和内皮细胞亚群的 p  值分别为 0.17  和 0.37。相比之下，肥大细胞和 B  细胞数据集分别为 0.05  和 0.03。这些结果表明， MetaSEM 可以捕获特异性信息，这对于一般 GRN 推理任务是必不可少的。

图 4。 基因表达在不同细胞类型上的差异。红点表示正相关的基因，蓝点表示负相关的基因，黑点表示表达水平无差异的基因。灰点表示阈值以下的基因。

2.4. SEM模型中选取的调控因子具有更高的表达水平

接下来，我们分析了 MetaSEM 如何提取调节因子。在整理 GRN 层的输出时，我们发现一些基因的调控权重非常高。因此，我们收集了这些基因的调控权重。图 5  的箱线图显示了结果。这些基因在成纤维细胞、T  细胞、癌症和内皮细胞中具有较高的调节重量。如图 5 的 t-SNE 图所示，所选基因在不同样本上具有显著的重量分布。ATF4、JUN、 RPL7A 和 RPS4X 在 HNSCC 细胞上的 log2(转录本每千碱基每百万(TPM)  +  1)分别为 8.1、6.5、6.8 和 9.8。最后，我们通过交叉比较选择 ATF4、JUN、RPL7A、RPS4X 等基因作为调控因子。

图 5。八种细胞中不同基因的调控权重。提出了四种调节器:ATF4、JUN、RPL7A 和 RPS4X。箱线图显示了不同 SEM 矩阵上调控子的重量分布。t-SNE 图代表了数据集上调控因子的权重分布。

2.5. 选定的调控因子是细胞类型鉴定的主要因素

在本节中，我们将介绍选定的调控因子与细胞类型识别之间的关系。图 6 显示了癌症和成纤维细胞数据集的 grn。蓝边代表典型的调控关系。绿边和红边分别表示肿瘤 grn 和成纤维细胞 grn 的调控关系。我们标记了 grn 中的几个调控因子，如 STAT1、 JUN 和 JUNB。STAT1 编码的蛋白是 STAT 蛋白家族的一员。STAT1 介导多种基因的表达，这对于细胞在应对不同细胞刺激和病原体时的生存能力至关重要[26]。JUN 和 JUNB 属于同一基因家族，与人类恶性肿瘤有关。JUN 常发生在人类恶性肿瘤的染色体易位和缺失区域[27]。

图 6。MetaSEM 对癌症和成纤维细胞数据集的 GRN 推断可视化。节点的大小表明了调控权重。蓝色边缘为 GRN 的主体部分，表示这两个细胞之间有共同的调控关系。绿色和红色调节关系仅存在于癌症 grn 或成纤维细胞 grn 中。

为了验证调控权重对细胞类型识别的重要性。我们为每个基因的调控权重收集了细胞类型的 SEM  矩阵。然后，根据调控权重将所有基因分成等量的三部分进行细胞聚类。然后，我们使用鲁汶和 Leiden  方法对细胞类型进行聚类。如图 7 所示，聚类效果与调控权重正相关。聚类方法在归一化互信息(NMI)、v-score  和调整兰德指数(ARI)上的表现证明了这一点。

图 7。 不同监管权重的选定监管机构的可视化。每一行代表一种聚类方法( 上一行 :鲁汶，下一 行 :Leiden)。每一列通过升序排列表示所选数据的调控权重。图的降维方法为 TSNE。

3. 材料与方法

3.1. 数据准备

 BEELINE[28]数据集用于评估模型的性能。单细胞数据集包含 7  种细胞类型，包括5 种小鼠细胞和 2 种人类细胞。对于 BEELINE 数据集，我们排除了标注为低质量的细胞和少于  10%的细胞中表达的基因。然后，我们对剩余的数据进行对数归一化处理。每个细胞只保留了前 1000 个标准差基因。子数据集根据不同的细胞类型进行划分。我们以同样的方式进一步构建了 500 个基因数据集，以评估模型在标签数据不足时的表现。根据[28]中的描述对 ground-truth grn 进行预处理和归一化处理。

头颈部鳞状细胞癌(head  and  neck  squamous  cell  carcinoma,  HNSCC)数据集[25]用于研究 GRN  细胞特异性，这是一种异质性上皮肿瘤，与细胞长期暴露于酒精和烟草环境密切相关。我们根据已知的细胞类型注解(成纤维细胞、B 细胞、T 细胞、内皮细胞、树突状细胞、肥大细胞、癌症细胞、成纤维细胞、肌细胞和巨噬细胞)将数据集分为十个子集。我们丢弃了没有注释的子集，样本数量小于 50。我们还将表达在少于 30 个样本中的基因进行了定位。最后，去除成纤维细胞和肌细胞。将标准差前 1000  的基因作为训练数据集。HNSCC 数据集对应的 ground-truth 标签从 TCGA 数据库中获得[29]。

将数据集存储在以垂直轴基因表达、横轴为样本的矩阵中。一个细胞类型对应一个矩阵，矩阵中的值代表基因在样本上的表达值。与 DeepSEM 的训练过程[15]类似，我们将数据集的 64 个样本作为一个批次。每个批次都被认为是一个小样本学习任务，通过损失函数，学习目标是唯一的。在 GRN 推理中不需要划分测试集和验证集。

3.2. 模型描述

提出的MetaSEM由三部分组成:编码器、元解码器和GRN层。(i)编码器:该部分使用三层MLP将基因表达数据编码为特征向量。(ii)元解码器:该部分通过双层MLP对监管关系进行建模，并使用GRN层寻找最优伪数据标签。(iii) GRN层:GRN层通过SEM模型推断基因调控关系，并将这些关系转化为伪数据标签。

3.2.1. 我们将SEM推广为一个GRN层来模拟随机变量之间的条件依赖关系。利用元解码器迭代GRN层，提取监管信息。该模块的最终输出是一个表示GRN的邻接矩阵，矩阵中的元素描述有向边的权重。

GRN Layer的迭代公式为:A∗= θA × A + α × A(1)，其中A表示基于SEM建模得到的邻接矩阵，θA表示元解码器的模型参数，α用于控制新迭代中矩阵的学习率。

3.2.2. 我们构建了一个解码器来捕获基因表达的数据特征。编码器批量读取自然基因表达数据X。然后，使用双层全连接层获得特征向量xpi。 ypi表示来自元解码器的伪数据标签。

3.2.3. 我们将GRN层嵌入到二层MLP中来构建元解码器。

与注意机制不同[31,32]，我们使用元解码器来指导特征提取。具体来说，元解码器的目标是找到一个与标签矩阵Y一致的预测矩阵Y p i。元解码器从scRNA-seq中学习基因调控关系。然后，将调节关系存储在GRN层中。最后，元解码器输出伪数据标签来表示潜在的监管关系。该程序可以通过伪数据标签提高编码器的效率。

Y p i = JA(θA) (X p i)(4)其中JA(θA) (X p i)表示使用元编码器计算伪数据标签的过程。在每一轮训练中，从特征向量集V p计算更新的过程如式(5)所示。

4.结论

提出了一种新的基于元学习的GRN推理算法，用于分析所选调控因子在细胞类型识别中的重要性。MetaSEM从基因表达数据中学习潜在的调节关系。此外，元学习也被用来优化特征提取的过程。在不同单细胞数据集上的大量实验表明，MetaSEM在GRN推理任务中的性能优于几种先进的计算方法。最后，通过可视化推断的GRN，我们系统地分析了数据分散的重要性。我们证明了所选调节因子在细胞类型鉴定中的重要性。在未来，我们打算通过融合scATAC-seq数据和scRNA-seq数据来构建GRN，以探索GRN对单细胞的影响。