Dictys：单细胞多组学分析发育连续性的动态基因调控网络

编译|程昭龙 审稿|林荣鑫，王静

本文介绍由哈佛医学院的Luca Pinello通讯预印在bioRxiv的研究成果：基因调控网络(GRN)是细胞功能和特性的关键决定因素，并且会在发育和疾病期间动态重组。尽管经过了几十年的发展，GRN推理仍然面临诸多挑战，如动态重组、因果推理、反馈回路建模和上下文特异性。为了解决这些问题，作者开发了一种动态GRN推断和分析方法Dictys，该方法利用了染色质可及性、基因表达的多组学单细胞分析、上下文特异性转录因子(TF)足迹、随机过程网络和scRNA-seq读取计数的高效概率模型。Dictys提高了GRN重建的准确性和再现性，并能够跨发育环境对特定上下文和动态GRN进行推断和比较分析。Dictys通过细胞类型特异性和动态GRN进行网络分析，恢复了人类血液和小鼠皮肤发育的独特见解。其动态网络可视化可以对发育驱动因子TF及其调控目标进行时间分辨的发现和研究。同时，Dictys是一个免费、开源和用户友好的Python包。

1

简介

生物分子网络包含一系列依赖于上下文的物理相互作用，化学反应，以及 mRNA、蛋白质、DNA 调节元件等之间的其他因果依赖性。它们共同定义了细胞特性、基因功能和最终的表型特征。在发育过程中，高度有序、动态进化的GRN由一系列细胞类型特异性TF组成。它们通过与其近端（即启动子）和远端（即增强子）顺式调控DNA元件的结合来调控细胞命运并决定基因表达，从而建立调控程序。因此，GRN已被用于模拟TF活动，确定关键调控因子的优先级，并了解它们如何动态地重塑调控方案。

在过去的几十年中，已经提出了许多方法来从bulk转录组数据推断GRN，但这些方法在区分因果关系与相关性、直接影响与间接影响方面面临挑战。同时，这些基于bulk的方法需要相当大的成本来获得足够的转录组样本和统计能力，这严重限制了它们实际的负担能力，更不用说诸如上下文特异性或动态GRN的比较研究了。

单细胞RNA测序(scRNA-seq) 通过经济高效地并行分析数千个细胞的转录程序，克服了其中的一些局限性，并促进了GRN推断的新计算发展。然而，单细胞GRN推断仍然面临着一些挑战，其中一些已经存在数十年之久，且尚未得到解决。实际上，一些基于稳态回归模型的方法，不仅无法考虑反馈循环或检测噪声，反而会因单细胞稀疏性加剧和混淆反馈循环或检测噪声。单细胞多组学分析可以联合分析基因表达和染色质可及性，但它们仍有待用于GRN推断的机制研究。

因此，作者提出了Dictys来应对这些挑战，并利用最新的多组学技术，使用上下文特定的TF足迹、随机过程模型、概率规划和动态网络。这些进步提高了GRN推断的质量，并通过专用的网络可视化，实现了对特定上下文的GRN和沿差异化路径的动态GRN的比较分析。同时，作者通过重新发现已知的关键调控因子以及发现和优先考虑新的调控候选因子及其靶基因，证明了它们在血液和皮肤系统中作为发育环境的效用。

2

结果

Dictys概述

Dictys是一个可用于网络推理、分析和可视化的集成Python包，它采用scRNA-seq和bulk或scATAC-seq数据的联合或单独配置文件，来发现与发育等连续过程相关的主调控子、目标基因及其重组过程。Dictys可以在标注的细胞群或生物环境(例如集群或分类种群)和沿连续过程(例如伪时间、RNA速度或基于克隆的轨迹)的动态网络中推断上下文特定网络。Dictys还提供了基于网络的分析和可视化，用于发现和检查每个TF的调控活动变化，其中TF的表达水平只是一个代表，不能完全捕获。

为了给每组细胞重建上下文特定的 GRN，Dictys 首先从pseudo-bulk或bulk染色质可及性数据的 TF 足迹推断调节区（即启动子和增强子）中的 TF 结合位点（图 1a，b）。与染色质可及性峰相比，TF 足迹区域要短很多，并且可以减轻假阳性结合位点。该步骤是根据推断的结合位点和接近度，将上下文特异性调控 TF-靶基因链作为 TF 结合网络优先考虑。

Dictys使用单细胞转录组数据细化初始的TF结合网络(图1c)。作者使用Ornstein-Uhlenbeck (OU) 过程对单细胞转录动态进行建模以考虑反馈循环，其稳态分布表征了单细胞表达的生物学变化。相反，单细胞技术变异/噪声采用稀疏二项式抽样建模。该生成过程在scRNA-seq读取计数上进行训练，以使用概率编程框架Pyro来推断所有的动态和随机参数，包括GRN。得到的GRN进一步归一化，以解决由于单细胞稀疏性造成的方差低估偏见。这种动态GRN还可以将每个基因的扰动模拟为基础转录率的变化，并分析推导出总效应(直接+间接)GRN作为其他基因稳态表达的相应变化。

图1 Dictys结构图

Dictys包含一套用于理解和比较特定上下文网络的功能。它通过识别每个TF(调节子)的靶基因组，可以根据恢复的网络将TF调控活性量化为靶基因数，这与仅基于TF表达水平定义的基因水平分析形成对比(CPM，图1d)。通过在全局水平上比较上下文特定网络，Dictys可以基于特定上下文的目标计数对调控标记TF进行优先级排序，并在点图中可视化。在两个特定上下文网络之间，Dictys可以根据恢复的差异调节和差异表达揭示TF活性变化的不同模式。这种关系可以在散点图中可视化，以发现调节活性发生强烈变化但表达不发生变化的独特TF，如果仅根据表达信息将会遗漏这些变化。通过进一步整合这两个差分轴， Dictys提供了一个可视化的条形图来对TF排名。在单网络级别中，Dictys还以网络图或热图格式可视化每个调控子，以便进行深入研究。

Dictys可以推断和分析(伪)时间分辨的动态GRN，并通过实验剖析连续过程中的基因调控变化。沿着提供的轨迹，Dictys首先定义一个移动窗口，将细胞子集划分为重叠的小亚群，然后为每个亚群重建一个静态GRN，进而用高斯核平滑重建动态GRN(图1e)。通过动态GRN, Dictys将每个TF的调节活性曲线定义为它们随时间的调节活性变化。然后，Dictys在单调或瞬态模式下发现具有高度可变调节活性曲线的TF，并使用综合网络可视化(INV)对单个基因和调控进行调查分析 (图1d)。

总体而言，Dictys为从单细胞转录组和染色质可及性获得的特定上下文和动态GRN提供了一个推理、可视化和分析框架，解决了现有方法在上下文特异性、时间分辨率、反馈循环和单细胞噪声检测方面的几个局限性。

Dictys推断的细胞类型特异性GRN可用于发现人类造血过程中的TF

为了证明Dictys在发育环境中的效用，作者利用一个包含骨髓单核细胞scRNA-seq和scATAC-seq数据的人类血液数据集进行实验(图2a)。首先确定细胞类型的特异性TF足迹，并通过质量评估确定足够的足迹，从而可以在人类造血过程中重建12种细胞类型的细胞类型特异性GRN。

图2 Dictys识别人类造血过程中的TF

通过比较所有这些细胞类型的GRN，Dictys识别出具有明确细胞类型特异性调控活性的调控标记TF（图2b）。此外，Dictys根据这些调控标记TF的顶部激活靶点，揭示了它们的细胞类型特异性调控程序（图2c）。这表明细胞类型特异性GRN可以识别调控标记TF及其调控程序，而不是基于均值表达的分析。作者接下来对比了造血早期到晚期的 GRN，以分析基于差异调节和差异表达（CPM）的每个 TF 的调控变化。尽管红细胞与祖细胞的差异分析恢复了许多已建立的 TF，但一些具有红细胞特异性功能的 TF 在差异调节中表现出更强的 logFC（图 2d）。因此，采用差异表达和差异调节的均值 logFC 综合 TF 排序可以提高排序质量，特别是对于已知单轴较弱的 TF（图 2e）。综上所述，差异调节可以识别在差异表达中无法检测到的具有细胞类型特异性功能的 TF，并为比较 TF 排序提供独立信息。

Dictys 通过利用小鼠皮肤中的多组学数据和转录组-染色质可及性关联来改进 GRN 推断

为证明 Dictys 如何利用最近的转录组-染色质可及性联合数据来改进 GRN 推断。作者重新分析了关于小鼠皮肤发育的 SHARE-seq 数据集（图 3a）。利用这种多模态数据，作者将初始 TF 结合网络限制在染色质峰值可及性和靶基因表达之间具有群体水平相关性的TF结合网络中。Dictys 鉴定出了与之前的表达标记不同的调控标记 TF（图 3b）。例如，基底表皮细胞和毛囊 TAC-1 之间的差异调节分析揭示了已知但排名较低或在差异表达中具有相反效果的 TF（图 3c）。通过综合 TF 排名，这些 TF比单独的差异表达更具细胞类型特异性（图 3d）。因此，Dictys 可以与最近的多组学数据一起使用，并利用转录组-染色质可及性关联来改进 GRN 推断。

图3 Dictys 利用多组学和转录组-染色质可及性联合数据推断GRN

Dictys在定量基准测试中优于现有方法

由于金标准完整性和正确性的限制，以及假设和问题表述的差异，GRN 推断基准测试仍然具有挑战性。因此，作者基于或独立于金标准，在上述血液和皮肤数据集上建立了五个基准来全面和定量地评估推断的GRN 质量。主要比较了 SCENIC、CellOracle 和 Dictys 这几种方法，因为这些方法可与scRNA-seq 数据一起明确地模拟 TF 结合。

首先，在 TF 结合评估中，作者从 Cistrome 数据库中收集了 512 个人类血液（114 个 TF）和 33 个小鼠皮肤（16个 TF）的 ChIP-seq 实验，以测试每种方法是否可以恢复这些实验支持的 TF-靶基因连接。在 TF 结合 + 染色质环评估中，作者通过将多个 Erythroid 特异性 ChIP-seq 实验与染色质构象数据相交，进一步整理了更可能具有调控作用的 TF 结合位点。实验结果表明，在所有评估指标中，Dictys 在两个数据集上均优于现有方法（图 4a-d）。

图4 Dictys基准测试

TF 结合并不总是调控附近的基因，因此作者接下来利用敲除(KO)或敲低(KD)实验从KnockTF数据库中提取造血TF进行扰动评估。为了量化再现差异基因表达 (DGE) 统计数据的能力，作者计算了这些统计数据与在重建的连续网络上传播的计算机模拟 TF 扰动预测的下游效应之间的 Pearson 相关性。为了捕捉间接效应，作者使用总效应 GRN 进行比较，或者作为 Dictys 的每个 TF 转录率变化的稳态效应，或者使用 CellOracle 描述的三步传播（图 4e）。通过考虑反馈循环和随机过程网络的无限步传播，总效应 GRN 表现出明显的优势（图 4f）。总之，通过系统基准测试表明，Dictys 在五种不同的评估中优于现有方法，同时不同参数选择模拟了广泛的生物应用，提供了稳健和特定细胞类型的网络推断。

动态GRN在时间分辨下研究TF及其在人类造血中的调控

原则上，单细胞技术能够研究和确定诸如发育或疾病进展等连续过程的特征。然而，这只能通过专门设计利用这种可能性的计算方法来实现。尽管现有的 GRN 推断方法已应用于离散细胞组，但 Dictys 可以使用任何连续的细胞顺序（例如时间、伪时间、RNA 速度或谱系数据）来重建动态 GRN 并揭示网络的连续重组。为了在造血中展示这一分析，作者首先利用 STREAM 推断了三个发育谱系（红细胞、B 细胞和单核细胞）的伪时间轨迹，并从相匹配的 scATAC-seq 与 ArchR 中整合了染色质可及性信息（图 5a）。通过使用伪时间作为发育顺序的替代，作者为每个谱系推断出一个动态 GRN，用于在(伪)时间分辨下发现和研究与发育相关的TF和调控因素。

图5 Dictys重构动态GRN和发现TF及调控

对于每个 TF，动态 GRN 用曲线直接量化其时间依赖性的调控活性。该调控活性曲线包含了给定 TF 潜在作用机制的重要信息，这些信息可以通过曲线的几何特征（例如距离和面积）进行总结（图 5a，b）。已知的 TF 证明了调控活性曲线在离散细胞类型之外的效用，以及无偏分析是如何优先考虑其他不太了解或新的与发育相关的 TF。总体而言，Dictys的动态GRN重建对造血功能提供了精细的时间分辨率和调控见解，这些信息是无法从粗粒度细胞簇或基于平均表达的分析中获得的。

3

总结

在这项研究中，作者提出了一种 GRN 推断方法 Dictys，该方法可用于从单细胞转录组和染色质可及性数据中重建、分析和可视化上下文特异性和动态 GRN。通过综合基准测试表明，Dictys 与现有的静态 GRN 推断方法相比具有卓越的性能，尤其是在细胞类型特异性和数据驱动的可重复性方面。Dictys 使得从发现到调查、从基因到调控、从离散组到连续过程的 GRN 比较分析成为可能。此外，Dictys 还揭示了当前基于平均表达变化和细胞簇的单细胞分析无法捕捉到的血液和皮肤系统的生物学见解。

GRN 推断是生物学中一个长期存在的基本问题，有许多挑战。Dictys 将 GRN 的稳态线性模型扩展到具有稳态观测的随机过程模型。然而，该方法没有考虑非线性、非稳态或机械模型，也没有包含其他可能改进推断的数据模式，例如未剪接的 mRNA、蛋白质、DNA 甲基化、单细胞扰动或细胞谱系测量。Dictys 假定处于稳态，但仍揭示了非稳态发育系统中的丰富生物学。同时，当细胞数量有限时，GRN 推断具有很高的方差，可以通过核平滑和基于目标计数的分析来稳定这一点。TF 结合预测的准确性可能受到细胞计数和缺乏生化/染色质结合分析的直接单细胞测量限制。此外，由于金标准的质量和数量有限，对 GRN 推断方法的评估仍然存在一定困难。

参考资料

Wang L, Trasanidis N, Wu T, et al. Dictys: dynamic gene regulatory network dissects developmental continuum with single-cell multi-omics[J]. bioRxiv, 2022.

数据链接：

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE139369

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE140203

代码链接：

https://github.com/pinellolab/dictys