本文将介绍由北京大学马剑竹课题组，清华大学交叉信息研究院曾坚阳课题组联合索尔克生物研究所，加利福尼亚大学，上海交通大学于7月22日发表在Nature Computational Science的研究成果：本文作者提出了一个深度生成模型DeepSEM，它可以推断基因调控网络(GRNs)和单细胞RNA测序数据的生物学表示。DeepSEM与最先进的方法相比，在各种单细胞计算任务上取得了优越的性能。此外，DeepSEM在小鼠皮层数据上进行验证，进一步证明了该模型的准确性和效率。因此，DeepSEM是分析细胞的scRNA-seq数据和推断GRNs的强大工具。
单细胞测序技术的快速发展为生物学家研究细胞状态提供了前所未有的机会。而在批量测序数据中未发现的实验噪声会显着降低下游生物信息学分析结果的准确性。为了解决这个问题，人们利用深度学习通过模拟基因之间复杂的相互作用模式来过滤单细胞转录组数据中的噪声。先前的基于深度学习的方法虽然能更清晰的揭示细胞异质性，捕捉转录组学相似性和细胞之间的差异，但很难评估基因调控网络(GRN)结构或数据其他内部结构。而许多计算模型是将GRN推断和单细胞数据分析模型进行结合。然而，通常需要设计复杂的实验，并且还可能引入额外的噪声。此外，基于单细胞RNA测序(scRNA-seq)数据的方法也有明显的局限性。
针对上述问题，本文作者提出基于β-VAE框架的DeepSEM模型，可以联合嵌入基因表达数据，同时构建GRNs反映单细胞内基因相互作用的内部结构，而无需依赖任何额外信息，例如TF 结合motif或单细胞ATAC测序(scATAC-seq)数据。作者在多个基准数据集上评估DeepSEM在各种单细胞任务中的性能都获得了不错的结果，并且由于VAE模型本身可以对潜在向量空间进行扰动来生成新的数据，因此在训练样本数据有限的情况下，仍能保证细胞类型分类的准确性。
给定scRNA-seq数据集作为输入，DeepSEM通过使用beta-VAE框架生成SEM来联合建模GRN和转录组(图 1)。DeepSEM有两个神经网络层，命名为GRN层和逆GRN层，以明确地对GRN结构进行建模。与传统深度学习模型将所有基因的表达一起嵌入到一个潜在空间中不同，DeepSEM的编码器函数仅将一个基因的表达作为神经网络的输入特征。不同基因的神经网络之间共享权重，或者可以将其视为使用一个神经网络来扫描所有基因。在这一步，模型中不同基因之间没有相互作用。然后另外两个全连接神经网络将这些小型神经网络的输出转换为多元高斯分布的后验均值和标准差。解耦非线性运算和基因相互作用是DeepSEM同时实现更稳健和可解释的隐含表示的关键。接下来，配备逆GRN层的解码器函数将隐含表示转换回基因表达值，这使整个框架成为自动编码器(图 2)。通过对GRN和单细胞转录组数据进行联合建模，DeepSEM作为一种多功能工具，可以通过分析不同的模块来完成单细胞数据分析中的各种任务。
图1 DeepSEM概述
左：DeepSEM两个主要模块，编码器(左下)和解码器(左上)。
右：DeepSEM通过利用不同的模块执行三个主要功能：
(1)GRN 预测(右下)，(2)scRNA-seq 数据嵌入和可视化(右中)，以及(3)scRNA-seq模拟(右上)。
DeepSEM的VAE包含四个模块：编码器、GRN层、逆GRN层和解码器。编码器和解码器都是以一个基因为输入的MLP，编码器和解码器的权重在不同基因之间共享。GRN层和逆GRN层都是基因相互作用矩阵，它们显式地对GRN网络进行建模并引导神经网络的信息流。g1-g6：基因表达数据中基因的名称；c1-c4：基因表达数据中细胞的名称；NN：神经网络。
为了评估 DeepSEM在GRN推断上的性能，作者遵循了BEELINE框架工作，该框架收集了四种不同类型的真实网络和七个scRNA-seq数据集，其中包括五种来自小鼠的细胞系和两种来自人类的细胞系。对于每个数据集，只考虑高度可变的TF和前N个(N=500和1,000)个变化最大的基因。通过BEELINE框架中使用的早期精确率(EPR)(图 3)和精确率-召回率曲线下的面积(AUPRC比率)评估性能。作者将DeepSEM与六种基线算法进行了比较，包括GENIE3、PIDC、GRNBoost2、SCODE、ppcor和SINCERITIES，这些算法已被证明在基于BEELINE评估的基准数据集上获得最先进的性能。为了从深度学习模型中获得稳定的预测，作者使用集成策略来生成最终预测。
总体而言，就EPR和AUPRC比率指标而言，DeepSEM优于scRNA-seq数据集上的所有其他基线方法(图 3)。当使用EPR进行评估时，DeepSEM实现了最佳预测性能，并且与次佳方法(GENIE3)相比至少提高了10%。在考虑AUPRC比率指标时，DeepSEM实现了最佳预测性能，与次佳方法(PIDC)相比，性能至少提高了10%。此外，DeepSEM在大多数基准测试中的表现明显优于 SCODE、ppcor和SINCERITIES。
DNA甲基化和染色质开放性会影响转录因子与顺式调控元件的结合，从而影响下游靶基因的表达。为了进一步探索DeepSEM识别的GRN的生物学意义，作者还检查了DeepSEM预测的基因调控与从细胞类型特异性表观遗传数据推断的基因调控之间的一致性。先前的研究报告称，将TF结合motif信息与表观遗传数据相结合可以细胞类型特异性方式准确预测TF结合位点。因此，作者假设，如果一个TF在给定的细胞类型中调控基因，它应该更有可能与有低CG甲基化和相应细胞类型中靶基因侧翼区域的开放染色质的TF motifs相关联。
为了验证这一假设，作者将DeepSEM框架应用于来自小鼠皮层的scRNA-seq数据集，并将结果与单核甲基胞嘧啶测序(snmC-seq)和scATAC-seq数据进行比较。为了寻找支持标记基因调控的表观遗传证据，作者使用scATAC-seq峰和差异甲基化区域(DMR)作为每个细胞类型特异性基因的潜在调控元件，然后对于每个基因，计算其预测的调节因子的比例。与作者的假设一致，作者发现不同类型的调控区域大量富集预测的TF的motifs(图 4a)，表明 DeepSEM从表观遗传方面识别的GRN具有很高的准确性。
更具体地，作者检查了与Rorb基因的相关预测，发现这些TF结合motifs位于Rorb的第一个内含子中，其中motif内的胞嘧啶在L4中特异性低甲基化，表明这些TF与Rorb之间存在细胞类型特异性调节关系(图 4b)。作为另一个例子，作者预测Syt6可能受L6 CT中的 Nf1a、Stat1和Sp3调节。因此，通过与scATAC-seq数据进行比较，作者观察到与沿着Syt6基因体的这些TF的结合motif相关的区域在L6 CT细胞中是特异性开放的(图 4c)。值得注意的是，所有这六个上游TF在目标和非目标细胞类型中都一致表达(图 4d)。总之，这些研究支持作者预测的GRN，并表明DeepSEM在研究细胞类型特异性基因相互作用网络方面的潜在效用。
先前的研究表明，通过考虑不同基因之间的相互作用，例如蛋白质-蛋白质相互作用网络、GRN、来自体细胞RNA-seq数据的共表达网络和注释通路，可以为scRNA-seq提供更好的生物学解释。特别是，将调控关系与基因表达联系起来已被证明能够有效地克服单细胞和批量测序实验中的dropout和其他技术问题。由于DeepSEM的细胞表示是从表达式到GRN活动的非线性映射，作者假设隐含表示也可以通过显式建模GRN结构来有效定义细胞状态和细胞类型。为了评估细胞表示的质量，作者应用DeepSEM来识别9个scRNA-seq数据集上的不同细胞类型。为了对DeepSEM进行基准测试，作者还将其与其他四种方法进行了比较：scVI、DCA、ZIFA 和FA。总的来说，DeepSEM 在9个基准数据集中的5个数据集上的表现优于所有四种基线方法，并且在其他数据集上也取得了较好的性能(图 5a)。
在本研究中，作者在PBMC数据集上比较了DeepSEM与其他两种基于GAN的方法 cscGAN和scGAN的模拟性能。作者发现，当使用计算方法和专家给定细胞类型标签进行聚类时，所有三种方法都可以生成与原始数据几乎无法区分的低维嵌入(图 6a)。作者训练了一个随机森林(RF)分类器以区分模拟数据和真实数据。作者的假设是，如果模拟数据生成是“真实的”，RF模型的预测性能应该接近随机。作者观察到，在区分DeepSEM生成的模拟数据和真实测试数据上，RF的分类性能与cscGAN和scGAN相比分别平均下降了1.21%和19.77%(图 6b)。DeepSEM实现了逼真的模拟，表明整合GRN可能是对scRNA-seq数据建模的关键步骤。
作者还提出了另一个概念，称为GRN一致性，用于衡量模拟单细胞的质量。GRN一致性衡量真实和模拟scRNA-seq数据之间预测GRN的差异。也就是说，如果scRNA-seq模拟是真实的，则使用模拟数据集获得的GRN预测应该与来自真实数据集的预测相匹配。作者发现，与cscGAN和scGAN相比，DeepSEM能够实现更高的GRN一致性(图 6c)。该结果表明，DeepSEM能够生成满足生物学约束的更真实的scRNA-seq数据。
图6 DeepSEM与cscGAN和scGAN模拟性能的比较
在本文中，作者介绍了一个通用的计算框架可以联合建模GRN和单细胞转录组数据。GRN的结构被明确建模为神经网络的特殊层，作为生物约束来限制参数空间。DeepSEM与最先进的方法相比，在单细胞计算任务上取得了较好的性能。未来，DeepSEM可以通过利用GRN作为“桥梁”构建共享隐空间来整合不同的单细胞模式，以及使用DeepSEM框架整合其他分子相互作用网络。与此同时，DeepSEM也存在一些局限性，例如，DeepSEM相对其他VAE模型要慢，运行时间随着所涉及的基因数量而增加。因此作者建议用户选择高度可变的基因，而不是使用整个转录组作为输入特征。
Shu, H., Zhou, J., Lian, Q. et al. Modeling gene regulatory networks using neural network architectures. Nat Comput Sci 1, 491–501 (2021).
https://doi.org/10.1038/s43588-021-00099-8
代码链接：https://github.com/HantaoShu/DeepSEM
