使用异构图Transformer进行单细胞生物网络推断.

paper：Single-cell biological network inference using a heterogeneous graph transformer

0. TL; DR

DeepMAPS 是一个用于从scMulti-omics数据中推断生物学网络的全新框架。该框架的核心思想是将scMulti-omics数据建模为一个异构图（heterogeneous graph），并利用一个多头图Transformer（multi-head graph transformer），以一种鲁棒的方式，同时学习细胞和基因在局部和全局上下文中的关系。

基准测试表明，DeepMAPS在细胞聚类和生物学网络构建方面均优于现有工具。作者在肺肿瘤白细胞CITE-seq数据和淋巴瘤scRNA-seq/scATAC-seq配对数据上的案例研究，充分展示了其推断细胞类型特异性生物学网络的强大能力。

1. 背景介绍

单细胞测序技术，特别是scRNA-seq和scATAC-seq，已经彻底改变了我们研究细胞异质性的方式。然而，任何单一的组学模态都只能提供细胞状态的一个快照。为了真正理解细胞内部复杂的分子机器是如何运转的，需要同时观察多个层面——这正是单细胞多组学（scMulti-omics）技术的用武之地。

近年来，涌现出许多用于整合scMulti-omics数据的计算工具，如Seurat、MOFA+、Harmony和totalVI。这些工具在细胞类型和状态的预测、批次效应的去除以及不同模态间的对齐等方面取得了巨大成功。

然而，这些方法大多存在一个共同的局限性：它们主要关注于细胞层面的整合，即如何将不同的细胞正确地分群和对齐。它们在很大程度上忽略了数据中蕴含的、更深层次的拓扑信息，即细胞与细胞之间、基因与基因之间、以及细胞与基因之间的复杂关系。

因此，现有方法虽然能告诉我们有哪些细胞类型，但很难回答更进一步的问题：在特定的细胞类型中，哪些生物学网络是活跃的？当细胞受到外界刺激时（如药物处理或疾病状态），这些网络是如何响应的？

为了回答这些问题，需要一个能够同时实现细胞聚类和网络推断的全新框架。

2. DeepMAPS 方法

DeepMAPS的整体框架是一个端到端的、无需先验假设的流程，其核心是构建一个异构图并利用异构图Transformer（HGT）进行表征学习和网络推断。

2.1 异构图的构建与初始化

DeepMAPS的第一步是进行数据预处理和模态整合。针对不同类型的scMulti-omics数据（如多个scRNA-seq、CITE-seq、scRNA-ATAC-seq），它采用了不同的、领域内成熟的整合策略（如CCA、CLR变换、基于RNA velocity的加权等），最终生成一个统一的细胞-基因活性矩阵 $X$。

基于这个整合后的矩阵 $X$，作者构建了一个异构二分图（bipartite heterogeneous graph） $G$。图中有两种类型的节点：细胞节点 $V_C$ 和基因节点 $V_G$。如果基因i在细胞j中的活性值 $x_{ij} > 0$，则在它们之间连接一条边。这是一个无向、无权重的图，只表示存在性。

为了获得每个节点的初始特征表示，作者使用了两个独立的图自编码器（Graph Autoencoders），分别学习细胞和基因的初始嵌入。细胞自编码器输入是矩阵 $X$，将每个细胞（行）从基因数降维到256维。基因自编码器输入是转置矩阵 $X^T$，将每个基因（行）从细胞数降维到256维。这样，每个细胞和基因节点都有了一个相同维度的初始嵌入向量。

2.2 异构图Transformer (HGT)

DeepMAPS利用HGT来学习图中所有节点（细胞和基因）的最终联合嵌入，并同时计算出基因-细胞对的注意力分数。HGT的每一层都包含三个关键步骤：异构互注意力、异构信息传递和目标节点聚合。

⚪ 异构互注意力 (Heterogeneous Mutual Attention)

HGT的核心是其多头注意力机制。对于一个目标节点 $v_t$（可以是细胞或基因）和它的一个邻居源节点 $v_s$，模型需要计算 $v_s$ 对 $v_t$ 的重要性，即注意力权重。

首先，通过三个独立的线性变换，将源节点和目标节点的上一层嵌入 $H^{l-1}$，分别投影为Query（$Q$）、Key（$K$）和Value（$V$）向量。然后，第h个注意力头的权重计算如下：

\[\text{ATT\_head}^h(v_s, e_{s,t}, v_t) = (K^h(v_s) W_{ATT}^{\phi(e_{s,t})}) (Q^h(v_t))^T \cdot \frac{\mu_{\tau(v_s), \phi(e_{s,t}), \tau(v_t)}}{\sqrt{d}}\]

这里的权重矩阵 $W_{ATT}$ 和先验重要性张量 $\mu$ 都是与边的类型（meta-relation）相关的。这意味着HGT能够为不同类型的节点间交互（如细胞-基因）学习不同的注意力模式。

⚪ 异构信息传递 (Heterogeneous Message Passing)

源节点 $v_s$ 传递给目标节点 $v_t$ 的“信息”，是通过对Value向量进行线性变换得到的：

\[\text{MSG\_head}^h(v_s, e_{s,t}, v_t) = V^h(v_s) W_{MSG}^{\phi(e_{s,t})}\]

同样，权重矩阵 $W_{MSG}$ 也是边类型相关的。

⚪ 目标节点聚合 (Target-specific Aggregation)

目标节点 $v_t$ 的新嵌入，是通过聚合其所有邻居传递来的、经过注意力加权的信息，并与自身的旧嵌入进行残差连接得到的：

\[\tilde{H}^l(v_t) = \text{Aggregate}_{v_s \in N(v_t)} (\text{Attention}(v_s, e_{s,t}, v_t) \cdot \text{Message}(v_s, e_{s,t}, v_t)) \\ H^l(v_t) = \sigma(\text{ReLU}(\tilde{H}^l(v_t))) \oplus H^{l-1}(v_t)\]

通过堆叠多层HGT，每个节点的嵌入就能融合来自其多跳邻居的信息，从而捕捉到图的全局拓扑结构。

2.3 模型训练与输出

为了高效地训练大规模异构图，DeepMAPS采用了子图采样策略。它从全图中随机采样出多个小规模的子图（例如50个），并在这些子图上以mini-batch的方式进行训练。

HGT被置于一个图自编码器（GAE）的框架中进行无监督训练。其目标是最小化重构损失，即要求解码器（通过节点嵌入的内积计算）重构出的邻接矩阵，与原始图的邻接矩阵尽可能相似。

训练完成后，模型输出两个关键信息：

节点嵌入：所有细胞和基因的最终低维嵌入向量。细胞嵌入可用于聚类，基因嵌入可用于功能分析。
注意力分数：在最后一层HGT中计算出的、每个基因-细胞对的注意力分数。这个分数 $a_{i,j}$ 代表了基因i对于定义细胞j身份的重要性。

2.4 细胞类型特异性生物学网络推断

基于HGT的输出，DeepMAPS通过一个基于斯坦纳森林问题（Steiner Forest Problem, SFP）的模型来推断每个细胞簇的细胞类型活性基因关联网络。

SFP模型的目标是在由基因和细胞构成的图中，找到一个连接该簇所有细胞的、总权重最小的子图（一个森林）。图中的边权重由基因-基因嵌入相似度和基因-细胞注意力分数共同决定。

最终，SFP求解出的最优子图中所包含的基因，就被认为是该细胞类型中最重要的活性基因，它们之间的连接构成了该细胞类型特异性的基因关联网络。

3. 实验分析

作者在十个不同来源、不同类型的scMulti-omics数据集上，对DeepMAPS的性能进行了全面的基准测试，并与Seurat, MOFA+, Harmony, TotalVI, GLUE等SOTA工具进行了比较。

3.1 实验一：细胞聚类性能

在所有十个基准数据集上，无论是对于多来源的scRNA-seq数据、CITE-seq数据，还是配对的scRNA-ATAC-seq数据，DeepMAPS在细胞聚类性能上（通过ARI或ASW指标评估）均取得了最佳表现 (图a)。

Seurat是表现第二好的工具，但其性能在不同参数选择下波动较大。通过留一法交叉验证，作者进一步证明了DeepMAPS的鲁棒性。即使在移除某个细胞簇后，其在剩余数据上的聚类性能依然显著优于其他工具（图c）。

UMAP可视化也直观地显示，DeepMAPS生成的细胞嵌入，其类内更紧凑，类间边界更清晰（图Ad-f）。

3.2 实验二：生物学网络推断的评估

作者进一步评估了DeepMAPS推断出的生物学网络的质量。

作者比较了DeepMAPS与IRIS3生成的基因关联网络。结果显示，DeepMAPS生成的网络具有显著更高的中心性得分（包括紧密中心性CC和特征向量中心性EC），这意味着其识别出的基因在网络中扮演着更核心、更重要的角色（图a, b）。

作者比较了DeepMAPS与SCENIC, MAESTRO等工具推断的GRN。结果显示DeepMAPS推断出的GRN包含更多独特的、在功能数据库（如Reactome, DoRothEA）中得到富集的转录因子（TF）（图c）。

更重要的是，DeepMAPS能够识别出更多高度细胞类型特异性的调控子（regulon）——即那些只在一个生物学功能/通路中显著富集的调控子。其F1分数也显著高于其他方法（图d, e）。

这些结果表明，DeepMAPS不仅能进行高质量的细胞聚类，更能推断出具有更高统计显著性和生物学意义的细胞类型特异性网络。

3.3 案例研究一：肺肿瘤免疫微环境CITE-seq数据

作者将DeepMAPS应用于一个包含PBMC和肺肿瘤白细胞的CITE-seq数据集。

通过整合RNA和蛋白质信息，DeepMAPS成功地识别出13个细胞类型，包括先前仅用单一模态难以区分的树突状细胞（DC）亚群（图a）。

作者发现，HGT学习到的某些潜在维度能够清晰地捕捉到细胞状态的连续转变。例如，第51个嵌入维度清晰地分开了浆细胞和记忆B细胞，第56个维度则分开了效应记忆T细胞和组织驻留记忆T细胞。沿着这些维度的基因/蛋白表达变化，揭示了细胞分化成熟过程中的分子动态（图4 d-h）。

基于DeepMAPS的细胞分型，作者利用CellChat推断了细胞间的通讯网络，并识别出了如DC与T细胞之间的CD6-ALCAM信号通路，以及肿瘤相关巨噬细胞（TAM）与T细胞之间的NECTIN-TIGIT免疫抑制通路等已知的、关键的互作关系（图i）。

3.4 案例研究二：淋巴瘤scRNA-ATAC-seq数据

作者进一步将DeepMAPS应用于一个弥漫性小淋巴细胞淋巴瘤（DSLL）的scRNA-ATAC-seq数据集，以展示其在疾病研究和GRN推断中的威力。

DeepMAPS成功地将DSLL B细胞分为了两个不同的恶性状态（state-1和state-2）。RNA velocity分析进一步揭示了从正常B细胞到这两个恶性状态的潜在发育轨迹（图c）。

作者识别出了在正常B细胞和两个DSLL状态中差异活跃的TF和调控子。例如，JUN, KLF16等与癌症发生相关的调控子在DSLL state-1中特异性高活性（图e, f）。

通过整合细胞间通讯分析，作者发现巨噬细胞通过BAFF信号通路，能够激活DSLL细胞中的AP-1（其亚基包括JUN）转录复合体，进而驱动CDK6等下游癌基因的表达，从而揭示了一个从细胞间互作到下游基因调控的、驱动DSLL发展的完整分子机制（图h）。