【DeepMaps】文献精讲：Single-cell biological network inference using a heterogeneous graph transformer
标题： Single-cell biological network inference using a heterogeneous graph transformer
（使用异构图Transformer的单细胞生物网络推理）
作者团队： 山东大学刘丙强、美国俄亥俄州立大学Qin Ma、美国密苏里大学Dong Xu
摘要： 单细胞多组学（scMulti-omics）允许同时对多种模式进行量化，以捕获复杂分子机制和细胞异质性的复杂性。现有的工具不能有效地推断不同细胞类型中活跃的生物网络以及这些网络对外部刺激的反应。文章提出了用于从scMulti-omics推断生物网络的DeepMAPS。它在异构图中建模scMulti-omics，并使用多头图转换器以鲁棒的方式在局部和全局上下文中学习细胞和基因之间的关系。
基准测试结果表明，DeepMAPS在细胞聚类和生物网络构建方面优于现有工具。它还展示了在肺肿瘤白细胞CITE-seq数据和匹配的弥漫性小淋巴细胞淋巴瘤scRNA-seq和scATAC-seq数据中获得细胞类型特异性生物网络的竞争能力。此外，我们还部署了一个具有多种功能和可视化的DeepMAPS web服务器，以提高scMulti-omics数据分析的可用性和可重复性。
文献链接： https://www.nature.com/articles/s41467-023-36559-0
代码链接： https://github.com/OSU-BMBL/deepmaps
1. 背景
单个单细胞模态仅反映了遗传特征的快照，并部分描述了细胞的特性，导致复杂生物系统中的表征偏差。单细胞多组学（scMultiomics）允许同时定量多种模式，以充分捕捉复杂分子机制和细胞异质性的复杂性。当这些分析与强大的计算分析方法相结合时，可以促进各种生物学研究。现有的scMulti-omics数据综合分析工具都没有明确考虑细胞和模态之间的拓扑信息共享。因此，它们不能有效地推断出与细胞聚集同时存在的多种细胞类型的活性生物网络，也无法阐明这些复杂网络在特定细胞类型对外部刺激的反应。
因此，团队开发了DeepMAPS（基于深度学习的单细胞数据多组学分析平台），这是一个异构图转换器框架，用于从scMulti-omics数据中推断细胞类型特异性生物网络。
该框架采用了一种先进的GNN模型，即异构图转换器（HGT），具有以下优点：
- (1) 构建了以细胞和基因为节点，以细胞和基因之间的关系为边的一体化异构图。
- (2) 该模型同时捕获细胞和基因之间的邻居和全局拓扑特征，以同时构建细胞-细胞关系和基因-基因关系。
- (3) 该HGT模型中的注意机制能够估计基因对特定细胞的重要性，从而可用于区分基因贡献并增强生物学可解释性。
- (4) 该模型是无假设的，不依赖于基因共表达的约束，因此有可能推断出其他工具通常无法找到的基因调控关系。
2. 方法
总的来说，DeepMAPS是一个端到端的、无假设的框架，可以从scMulti-omics数据中推断出细胞类型特异性的生物网络。在DeepMAPS框架中有五个主要步骤，如下图所示。
2.1 Preprocessing 数据预处理
对数据预处理部分，去除低质量细胞和低表达基因，然后根据具体数据类型采用不同的归一化方法。生成一个集成的细胞-基因矩阵来表示每个细胞中每个基因的组合活性。针对不同的scMulti-omics数据类型，采用了不同的数据集成方法。
DeepMAPS的分析以 多个scRNA-seq（多基因表达矩阵）、CITE-seq（基因和表面蛋白表达矩阵）和scRNA-ATAC-seq（基因表达和染色质可及性矩阵） 数据的原始计数矩阵为输入。对于每个数据矩阵，将模态表示定义为行（基因、蛋白质或峰区），将细胞定义为列，除非有例外情况。在每个数据矩阵中，如果一行或一列包含的非零值少于0.1%，则删除该行或一列。数据质量控制由Seurat v3进行，包括但不限于总读取次数、线粒体基因比例、黑名单比例。
2.1.1 Multiple scRNA-seq data 多基因表达矩阵
对多基因表达矩阵，使用对数归一化。
为了学习细胞和基因的联合表示，首先生成一个整合了输入scMultiomics数据信息的细胞-基因矩阵。使用Seurat v3从每个矩阵中选择前2000个高度可变的基因。如果基质中少于2000个基因，则全部被选择进行整合。然后在Seurat中应用广泛使用的典型相关分析（CCA）来对齐这些矩阵并协调scRNA-seq数据，从而得到矩阵
X
=
x
i
j
∣
i
=
1
,
2
,
.
.
.
,
I
;
j
=
1
,
2
,
.
.
.
,
J
\mathcal{X} = {x_{ij} | i=1,2,...,I; j=1,2,...,J}
X=xij∣i=1,2,...,I;j=1,2,...,J，
x
i
j
x_{ij}
xij表示J细胞中I基因的归一化表达值。
2.1.2 CITE-seq data 基因和表面蛋白表达矩阵
对基因和表面蛋白表达矩阵，使用对数归一化。
选择前
I
1
I_1
I1个高度可变的基因（
I
1
=
2000
I_1 = 2000
I1=2000）和
I
2
I_2
I2个蛋白质（
I
2
I_2
I2为基质中蛋白质的总数）。然后将这两个矩阵垂直连接，得到一个矩阵
X
=
x
i
j
∣
i
=
1
,
2
,
.
.
.
,
I
1
+
I
2
;
j
=
1
,
2
,
.
.
.
,
J
\mathcal{X} = {x_{ij} | i=1,2,...,I_1+I_2; j=1,2,...,J}
X=xij∣i=1,2,...,I1+I2;j=1,2,...,J，
x
i
j
x_{ij}
xij表示J细胞中
I
1
I_1
I1基因和
I
2
I_2
I2蛋白的归一化表达值。（可视为J细胞中的
I
1
+
I
2
=
I
I_1 + I_2 = I
I1+I2=I基因）。对
X
\mathcal{X}
X进行中心对数比（CLR）变换如下：
CLR
(
x
i
j
)
=
log
(
1
+
(
x
i
j
exp
(
∑
k
j
j
log
(
1
+
x
i
j
)
∣
ξ
j
∣
)
)
)
\operatorname{CLR}\left(x_{i j}\right)=\log \left(1+\left(\frac{x_{i j}}{\exp \left(\frac{\sum_{k j_{j}}^{\log \left(1+x_{i j}\right)}}{\left|\xi_{j}\right|}\right)}\right)\right)
CLR(xij)=log
1+
exp(∣ξj∣∑kjjlog(1+xij))xij
其中，
Z
j
Z_j
Zj表示细胞J中非零基因的索引集，
∣
⋅
∣
|·|
∣⋅∣表示集合中元素的个数。
2.1.3 scRNA-ATAC-seq 基因表达和染色质可及性矩阵
基因表达矩阵
X
=
x
i
j
∣
i
=
1
,
2
,
.
.
.
,
I
;
j
=
1
,
2
,
.
.
.
,
J
\mathcal{X} = {x_{ij} | i=1,2,...,I; j=1,2,...,J}
X=xij∣i=1,2,...,I;j=1,2,...,J，
x
i
j
x_{ij}
xij表示J细胞中I基因的归一化表达值。
然后，通过对细胞群体中潜在调控信号如何控制基因表达的建模，使用左截断混合高斯（LTMG）模型来提供所有细胞中每个基因的定性表示。具体来说，如果基因i可以用所有J细胞的
G
i
G_i
Gi高斯分布来表示，这意味着可能存在
G
i
G_i
Gi调节信号来调节该基因。可以生成与
X
R
X^R
XR具有相同维数的矩阵
X
R
=
x
i
j
R
′
X^R = {x^{R'}_{ij}}
XR=xijR′，其中基因表达用
x
i
j
R
′
=
1
,
2
,
.
.
,
G
i
x^{R'}_{ij}=1,2,..,G_i
xijR′=1,2,..,Gi的离散值来标记。
2.2 Data intergration 数据整合构建异构图
从集成矩阵中构建异质图，以细胞和基因为节点，以细胞中存在的基因为边。
构建具有细胞节点和基因节点的异构图，其中未加权的细胞-基因边缘表示细胞中基因活性的存在，并且通过双层GNN图自编码器从基因-细胞集成矩阵中学习每个节点的初始嵌入（方法）。这种异构图为清晰地表示和有机地集成scMulti-omics数据提供了机会，从而可以协同学习具有生物学意义的特征。
首先通过多组学数据整合得到矩阵
X
=
x
i
j
∣
i
=
1
,
2
,
.
.
.
,
I
;
j
=
1
,
2
,
.
.
.
,
J
\mathcal{X} = {x_{ij} | i=1,2,...,I; j=1,2,...,J}
X=xij∣i=1,2,...,I;j=1,2,...,J，
x
i
j
x_{ij}
xij 表示J细胞中I基因的归一化表达值（用于多个scRNA-seq和CITE-seq）或GAS（用于scRNA-ATAC-seq）。随后通过两个自编码器计算基因和细胞的初始嵌入。细胞自编码器将每个细胞的基因维数从I维降至512维，再降至256维；基因编码器将每个基因的细胞维度从J维减少到512维，再降至256维。由此，每个细胞和基因都有256维的相同初始嵌入。输出层是与
X
X
X 具有相同维数的重构矩阵
X
^
\hat{X}
X^，细胞自编码器的损失函数为
X
X
X 与
X
^
\hat{X}
X^的均方误差MSE：
l
o
s
s
=
M
S
E
(
X
,
X
^
=
∑
i
(
X
−
X
^
)
2
)
loss = MSE(X, \hat{X} = \sum_i(X - \hat{X})^2)
loss=MSE(X,X^=i∑(X−X^)2)
基因自编码器从所有细胞中学习基因的低维特征，具有类似于细胞自编码器的编码器、隐空间和解码器，输入
X
T
X^T
XT为
X
X
X的转置矩阵，基因自编码器的损失函数为：
l
o
s
s
=
M
S
E
(
X
T
,
X
^
T
=
∑
i
(
X
T
−
X
^
T
)
2
)
loss = MSE(X^T, \hat{X}^T = \sum_i(X^T - \hat{X}^T)^2)
loss=MSE(XT,X^T=i∑(XT−X^T)2)
其中
X
^
T
\hat{X}^T
X^T是输出层的重构矩阵，与
X
T
X^T
XT具有相同维数。给定整合矩阵
X
X
X，作者构建具有两种节点类型（细胞节点和基因节点）和一种边类型（基因-细胞）的异构图
G
=
(
V
,
E
)
G=(V,E)
G=(V,E)用于下一节的计算。其中
V
=
V
C
∪
V
G
V = V^C \cup V^G
V=VC∪VG,
V
G
=
v
i
G
∣
i
=
1
,
2
,
.
.
.
,
I
V^G = {v_i^G | i=1,2,...,I}
VG=viG∣i=1,2,...,I表示所有的基因，
V
C
=
v
j
C
∣
j
=
1
,
2
,
.
.
.
,
J
V^C = {v_j^C | j=1,2,...,J}
VC=vjC∣j=1,2,...,J表示所有的细胞。
E
=
e
i
,
j
E = {e_{i, j}}
E=ei,j表示
v
i
G
v_i^G
viG和
v
j
C
v_j^C
vjC之间的边。
2.3 The HGT framework HGT模型工作框架
建立HGT模型，共同学习细胞和基因的低维嵌入，并生成一个关注分数来表示基因对细胞的重要性。
这里运用一个无监督的HGT框架来学习所有节点的图嵌入并挖掘基因和细胞之间的关系。HGT的输入是整合矩阵
X
X
X，输出的是细胞核基因的嵌入以及表示基因对细胞重要性的注意力分数。
2.3.1 多头注意力机制与向量线性映射
在这里，DeepMAPS采用异构多头注意力机制对异构图上的整体拓扑信息（全局关系）和邻居消息传递（局部关系）进行建模。
令
H
l
[
v
t
]
\mathcal{H}^l[v_t]
Hl[vt]和
H
l
−
1
[
v
s
]
\mathcal{H}^{l-1}[v_s]
Hl−1[vs]表示
v
t
v_t
vt和
v
s
v_s
vs第l层HGT的嵌入。采用多头注意力机制将
H
l
[
v
t
]
\mathcal{H}^l[v_t]
Hl[vt]和
H
l
−
1
[
v
s
]
\mathcal{H}^{l-1}[v_s]
Hl−1[vs]等分入
H
H
H个头。多头注意力允许模型共同注意来自不同嵌入子空间的信息，并且每个头可以并行运行一个注意机制以减少计算时间。
对于第
l
l
l层HGT的第
h
h
h个头，
H
l
[
v
t
]
\mathcal{H}^l[v_t]
Hl[vt]更新自
H
l
−
1
[
v
t
]
\mathcal{H}^{l-1}[v_t]
Hl−1[vt]和
H
l
−
1
[
v
s
]
\mathcal{H}^{l-1}[v_s]
Hl−1[vs]。
H
0
[
v
t
]
\mathcal{H}^0[v_t]
H0[vt]和
H
0
[
v
s
]
\mathcal{H}^0[v_s]
H0[vs]是
v
t
v_t
vt和
v
s
v_s
vs的初始嵌入。具体来说，
Q
l
i
n
e
a
r
τ
(
v
t
)
h
Q_{{linear}_{\tau(v_t)}}^h
Qlinearτ(vt)h函数将
v
t
v_t
vt线性映射到第
h
h
h个查询向量
Q
h
(
v
t
)
Q^h(v_t)
Qh(vt)，维度为 $\mathbb{R} ^{d} \to \mathbb{R} ^{\frac{d}{H} } $
其中
d
d
d是
H
l
−
1
[
v
t
]
\mathcal{H}^{l-1}[v_t]
Hl−1[vt]的维数，
d
H
\frac{d}{H}
Hd是每个头的向量维数。
类似的，
K
l
i
n
e
a
r
τ
(
v
s
)
h
K_{{linear}_{\tau(v_s)}}^h
Klinearτ(vs)h和
V
l
i
n
e
a
r
τ
(
v
s
)
h
V_{{linear}_{\tau(v_s)}}^h
Vlinearτ(vs)h函数将源节点
v
s
v_s
vs映射到第
h
h
h个键向量
K
h
(
v
s
)
K^h(v_s)
Kh(vs)和第
h
h
h个值向量
V
h
(
v
s
)
V^h(v_s)
Vh(vs)。其中计算公式为：
Q
h
(
v
t
)
=
Q
linear
τ
(
v
t
)
h
(
H
l
−
1
[
v
t
]
)
,
K
h
(
v
s
)
=
K
linear
τ
(
v
s
)
h
(
H
l
−
1
[
v
s
]
)
,
V
h
(
v
s
)
=
V
linear
τ
(
v
v
s
h
(
H
l
−
1
[
v
s
]
)
\begin{array}{l} Q^{h}\left(v_{t}\right)=Q_{\text {linear }_{\tau\left(v_{t}\right)}}^{h}\left(\mathcal{H}^{l-1}\left[v_{t}\right]\right), \\ K^{h}\left(v_{s}\right)=K_{\text {linear }_{\tau\left(v_{s}\right)}}^{h}\left(\mathcal{H}^{l-1}\left[v_{s}\right]\right), \\ V^{h}\left(v_{s}\right)=V_{\text {linear }_{\tau_{\left(v v_{s}\right.}}}^{h}\left(\mathcal{H}^{l-1}\left[v_{s}\right]\right) \end{array}
Qh(vt)=Qlinear τ(vt)h(Hl−1[vt]),Kh(vs)=Klinear τ(vs)h(Hl−1[vs]),Vh(vs)=Vlinear τ(vvsh(Hl−1[vs])
每种类型的节点都有一个唯一的线性投影，以最大限度地模拟分布差异。
2.3.2 HGT框架的构成
见HGT文献精讲
2.3.3 子图上的HGT训练
为了提高HGT模型在大型异构图（数万个节点和数百万条边）上的效率和能力，作者部署了一种改进的HGT采样方法，用于子图选择和HGT训练。
对于含有I个基因和J个细胞的图G，子图的并集应覆盖基因和细胞节点的
a
a%
a（设为30%）节点，以保证训练效率。因此，采样器从给定的异构图G构建了许多小的子图（DeepMAPS中为50个），并使用多个gpu将这些子图以不同的批次馈送到HGT模型中。每个图应该包括
a
%
×
I
/
50
a\% \times I/50
a%×I/50个基因和
a
%
×
J
/
50
a\% \times J/50
a%×J/50个细胞。
以细胞j为目标节点，对应基因i作为源节点，计算其边
e
s
,
t
e_{s,t}
es,t的概率为：
Prob
(
e
s
,
t
)
=
x
(
v
s
,
v
t
)
∑
v
∈
N
(
v
t
)
x
(
v
,
v
t
)
\operatorname{Prob}\left(e_{s, t}\right)=\frac{x\left(v_{s}, v_{t}\right)}{\sum_{v \in \mathcal{N}\left(v_{t}\right)} x\left(v, v_{t}\right)}
Prob(es,t)=∑v∈N(vt)x(v,vt)x(vs,vt)
式中，
x
(
v
s
,
v
t
)
=
x
i
j
x(v_s, v_t)=x_{ij}
x(vs,vt)=xij表示基因i在聚合矩阵x中在细胞j中的表达或GAS值。因此，对于每个目标节点
v
t
v_t
vt，基于采样概率
P
r
o
b
(
e
s
,
t
)
{Prob}(e_{s, t})
Prob(es,t)随机选择
a
%
×
I
/
50
a
%
×
J
/
50
\frac{a\% \times I/50}{a\% \times J/50}
a%×J/50a%×I/50个邻居基因。
HGT中的超参数，例如
W
ϕ
(
e
s
,
t
)
A
T
T
,
W
ϕ
(
e
s
,
t
)
M
S
G
W^{ATT}_{\phi (e_{s,t})}, W^{MSG}_{\phi (e_{s,t})}
Wϕ(es,t)ATT,Wϕ(es,t)MSG和
θ
\theta
θ，将在一个周期内从子图1到50依次训练和更新。子图训练以无监督的方式使用图自编码器(GAE)进行。HGT是编码器层，嵌入的内积是解码器层。将GAE的损失函数计算为重构矩阵
X
^
\hat{X}
X^和聚合矩阵
X
X
X的Kullback-Leibler散度（KL）：
l
o
s
s
=
K
L
(
s
o
f
t
m
a
x
(
X
^
,
s
o
f
t
m
a
x
(
X
)
)
)
loss = KL(softmax(\hat{X}, softmax(X)))
loss=KL(softmax(X^,softmax(X)))
损失被控制或达到100次epoch，则子图训练完成.
初始图决定了消息传递的路径以及如何在DeepMAPS中计算注意力分数。
在每个HGT层中，每个节点（细胞或基因）被视为目标节点，其邻居被视为源节点。DeepMAPS根据节点嵌入的协同作用（即注意力分数）评估其邻居节点的重要性和可以传递给目标的信息量。因此，具有高度正相关嵌入的细胞和基因更有可能在彼此内部传递信息，从而最大化嵌入的相似性和差异性。
作为一个重要的训练结果，给出了一个注意力分数来表示基因对细胞的重要性。一个基因对一个细胞的关注度越高，表明该基因对确定细胞身份和表征细胞异质性具有相对重要的意义。
2.4 Clustering 预测细胞聚类和功能基因模块
作者比较了DeepMAPS中不同的聚类方法（即Leiden， Louvain和SLM），并比较了聚类分辨率（即0.4,0.8,1.2和1.6）对细胞聚类结果的影响。发现这些聚类方法之间没有显著差异，Louvain的性能略好于其他两种。最后，在选择相同的聚类分辨率时，DeepMAPS的得分高于其他工具。我们还发现，在大多数情况下，分辨率越高，细胞聚类预测得分越低；因此，我们在DeepMAPS中选择分辨率为0.4作为默认参数。
2.5 Network inference 识别细胞聚类中基因关联网络
在每种细胞类型中推断出多种生物网络，例如基因调控网络（GRN）和基因关联网络。作者使用SFP模型来选择对细胞聚类特征具有重要贡献的基因，并构建细胞聚类中火星基因关联网络。
定义一个新的异构图
G
~
=
(
V
,
E
~
)
\tilde{G} = (V, \tilde{E} )
G~=(V,E~)，其中
V
∈
V
G
∪
V
C
V \in V^G \cup V^C
V∈VG∪VC, $\tilde{E} \in \tilde{E}^1 \cup \tilde{E}^2
，有
，有
，有\tilde{E}^1
表示基因
−
基因关系，
表示基因-基因关系，
表示基因−基因关系，\tilde{E}^2
表示基因
−
细胞关系。设
表示基因-细胞关系。 设
表示基因−细胞关系。设 Z
为
L
o
u
v
a
i
n
聚类预测的聚类个数，
为Louvain聚类预测的聚类个数，
为Louvain聚类预测的聚类个数，V^{C[z]} = {v^{C[z]}}
表示与
表示与
表示与z=1,2,…,Z$的聚类标签中的单元集对应的节点集。
3. 结果
3.1 DeepMAPS在scMulti-omics数据的细胞聚类和生物网络推断方面性能卓越
文章在10个sc多组学数据集上对DeepMAPS的细胞聚类性能进行了基准测试，包括3个多个scRNA-seq数据集、3个ite -seq数据集以及4个匹配的scRNA-seq和scATAC-seq （scRNA-ATAC-seq）数据集，这些数据集来自同一个细胞。这些数据集涵盖的细胞数量从3,009到32,029不等；平均读取深度（仅考虑scRNA-seq数据）范围为2,885至11,127；零表达率（仅考虑scRNA-seq数据）为82%至96%。
作者将DeepMAPS与四种基准测试工具（Seurat v3和v45、18、MOFA + 6、TotalVI8、Harmony7和GLUE19 (Methods)）在平均轮廓宽度（ASW）、Calinski-Harabasz （CH）、Davies-Bouldin指数（DBI）和调整后的Rand指数（ARI）方面进行了比较，以评估细胞聚类性能。对于每个数据集，在36个参数组合上训练DeepMAPS，包括头的数量、学习率和训练epoch的数量。为了确保公平性，还使用不同的参数组合对每个基准测试工具进行了调优（方法）。在ARI和ASW上，比较所有基准工具在所有测试数据集中的最佳性能。还注意到，Seurat是性能第二好的工具，在所有基准数据集中，不同参数选择的差异很小。
每个框显示了使用不同参数设置的工具的最小、第一四分位数、中位数、第三四分位数和最大ARI或AWS结果（DeepMAPS: n = 96， RNA-RNA和CITEseq的Seurat: n = 16， RNA-ATAC的36，Harmony: n = 36, MOFA +: n = 36, TotalVI: n = 48, GLUE: n = 72）。点表示异常值。
文章根据参数组合在网格搜索基准测试中的性能选择每个数据类型的默认参数。所有基准数据集中平均ARI/ASW分数中位数最高的参数组合被视为相应数据类型的默认参数。
文章进一步在五个独立的数据集上独立测试了我们的默认参数选择，分别是R-test， C-test, A-test-1， - 2和- 3，通过将结果与使用其默认参数的相同基准测试工具进行比较。对于三个带有基准单元标签的测试数据集，DeepMAPS在ARI得分方面表现最好，而对于两个没有单元标签的scRNA-ATAC-seq数据集，比较中的基准测试工具实现了相似的性能。
为了评估DeepMAPS的鲁棒性，对三个带有基准标签的独立测试数据集进行了留一测试。首先基于基准标签去除细胞簇中的所有细胞，然后对剩余的细胞应用DeepMAPS和其他工具。对于每个数据集，DeepMAPS的leave-one结果都优于ARI得分较高的其他工具，这表明DeepMAPS使用的消息传递和注意机制以稳健的方式维持了细胞-细胞关系。
在三个具有基准标记的独立数据集上的细胞聚类UMAP表明，在DeepMAPS中获得的潜在表征可以更好地保留scRNA-seq数据的异质性（下图）。对于r测试数据集，所有工具都显示出分离间充质细胞、白细胞和内皮细胞的能力，但未能分离尿路上皮基底细胞和膀胱细胞。然而，DeepMAPS UMAP上的细胞更紧凑，膀胱细胞（红点）的分组比MOFA +和Seurat（图d）更好。对于C-test数据集，同一簇中的细胞更加有序和紧凑（例如，B细胞簇和NK细胞簇），而来自不同簇的细胞在DeepMAPS UMAP上更加彼此分开（例如，CD8细胞簇和CD4细胞簇）。(图e)。对于A-test-1数据集，DeepMAPS是唯一能够准确分离每种细胞类型的工具。相反，Seurat和MOFA +错误地将PDX1或PDX2群体分为两个簇，并包含更多的不匹配（图f）。
3.2 DeepMAPS可以从scMulti-omics数据中推断出具有统计学和生物学意义的基因关联网络
作者在中心性评分和功能富集方面评估了DeepMAPS可以推断的两种生物网络（基因关联网络和GRN）。对于R-test数据集（下图a）和C-test数据集（下图b），作者使用了两种中心性评分来比较所有方法识别的基因关联网络，即接近中心性（CC）和特征向量中心性（EC）。CC反映了网络中一个节点与所有其他节点的平均连通性，EC则根据其连接的节点反映了节点的重要性。节点中心度越高的基因关联网络表明检测到的基因越有可能参与关键和功能性的生物系统。作者还计算了细胞簇中基因表达的Pearson相关系数，使用数据集中的所有基因构建基因共表达网络。实验表明，DeepMAPS中产生的基因关联网络不仅是共表达的，而且对细胞的注意力影响很大。因此，网络中的基因往往对细胞类型更重要。
为了评估 DeepMAPS 是否能够在特定细胞类型中识别具有生物学意义的基因调控网络，作者使用三个公共的功能数据库 Reactome、DoRothEA 和 TRRUST v2 对基因调控模块进行了富集分析。
为避免比较中出现偏差，作者将 DeepMAPS 推断出的细胞类型特异性 GRN 与（1）IRIS3 和 SCENIC（在 scRNA-seq 矩阵上）（2）基于 DeepMAPS 的基因活性分数（GAS）计算的 IRIS3 和 SCENIC（在基因-细胞矩阵上）（3）MAESTRO（在 scATAC-seq 矩阵上）（4）MAESTRO（在原始 scRNA-seq 和 scATAC-seq 矩阵上）进行比较。作者使用六个从人类组织收集的数据集，结果显示 DeepMAPS 鉴定的 GRN 包含的独特转录因子（TF）调控比其他方法更多（下图c）。
作者还比较了不同工具富集到一个功能/通路的细胞类型特异性调节子（ CTSR ）的数量。在大多数scRNA-ATAC-seq数据集上，DeepMAPS 在仅富集一个功能/通路的调节子数量和富集评分方面优于其他方法（下图d，e）。
3.3 DeepMAPS准确地识别PBMC和肺肿瘤免疫CITE-seq数据中细胞类型并推断细胞间通讯
文章将DeepMAPS应用于已发表的混合外周血单个核细胞（PBMC）和肺肿瘤白细胞CITE-seq数据集以证明在表征细胞身份方面建模scMulti-omics的能力。该数据集包括在3485个细胞上测量的rna和蛋白质。DeepMAPS确定了13个细胞簇。同时通过可视化标记的maker基因和蛋白质的表达水平来注释每个簇。与仅使用蛋白质或RNA鉴定的细胞类型相比，模型分离或准确地注释了无法使用个体模态分析表征的细胞群。例如，只有使用整合的蛋白质和RNA才能成功识别DC簇。通过结合从RNA和蛋白质中捕获的信号，DeepMAPS成功地在CITE-seq数据中识别出生物学上合理且有意义的细胞类型。
文章还比较了两种细胞类型之间的模态相关性。利用记忆B细胞和浆细胞之间的顶级差异表达基因和蛋白，对相关矩阵进行分层聚类。结果清楚地将这些特征分为两个不相关的模块：一个与记忆B细胞相关，另一个与浆细胞相关。（图c）此外，还发现两个模块中的特征与我们的HGT包埋捕获的成熟轴显著相关例如，一个HGT包埋（第51个）在浆细胞和记忆B细胞之间显示出显著差异（图d， e）。在比较EM CD8+ T细胞和TRM CD8+ T细胞时也观察到类似的发现（图f）。然而，有可能确定一个具有代表性的HGT嵌入（第56位），该嵌入维持了两组明确分离的嵌入信号（图h）。这些结果指向任何两个由多个基因和蛋白质的协调激活和抑制组成的细胞簇，导致细胞状态的逐渐转变，可以通过DeepMAPS潜在HGT空间的特定维度捕获。
3.4 DeepMAPS在弥漫性小淋巴细胞淋巴瘤scRNA-seq和scATAC-seq数据中识别特异性GRN
为了进一步扩展DeepMAPS对GRN推断的能力，作者使用了10 × Genomics网站上提供的单细胞Multiome ATAC+基因表达数据集。同时，通过RNA速度平衡细胞中基因每种模态的权重来整合基因表达和染色质可及性。此外，为了建立TF-基因链接，作者考虑了基因表达、基因可及性、TF-基序结合亲和力、峰-基因距离和TF编码基因表达。在一个细胞簇中发现受相同TF调节的基因被归类为调节子，而具有较高中心性评分的调节子对细胞簇的表征具有更显著的影响。
DeepMAPS在DSLL数据中识别出11个细胞簇，所有的簇都是基于基因标记进行人工注释的（图b）。对三个B细胞群进行的基于RNA速度的伪时间分析假设两个DSLL状态来自正常B细胞，并且状态1比状态2更早出现，尽管这两个状态似乎部分混合（图c）。作者进一步选择了三个细胞群中具有最高正则子中心性得分的前20个转录因子（图d）。实验可知，这些转录因子显示出了正常状态和两种DSLL状态之间的差异，并推断出两种DSLL状态中的变异调控模式。
作者构建了一个由四种细胞类型特异性调节子组成的GRN（图e），在DSLL状态-1中，RAS明显高于正常B细胞和DSLL状态-2（图f）。当放大单个调节子时，可以观察到不同的调节模式（图g）。JUN是DSLL状态-1中最活跃的调节子，能够调节5个独特的下游基因和12个与DSLL状态-2共有的基因。
3.5 DeepMAPS为分析scMulti-omics数据提供了一个多功能且用户友好的门户网站
由于单细胞测序数据的复杂性，过去三年开发了许多的网络服务器和对接器。然而，这些工具中的大多数仅提供最小的功能，例如细胞聚类和差异基因分析。它们不支持scMulti-omics数据的联合分析，尤其缺乏对生物网络推断的充分支持。另一方面，作者还记录了DeepMAPS和基准方法在不同数据集上的运行时间，这些数据集的细胞数量范围从1000 到160000。由于深度学习模型（DeepMAPS和TotalVI）的运行时间比Seurat和MOFA+更长。为此，作者提供了一个无代码、交互式和非编程的接口，以减轻scMulti-omics数据的编程负担（图6a）。网络服务器支持使用DeepMAPS分析多个RNA-seq数据、CITE-seq数据和scRNA-ATAC-seq数据（图6b）。该服务器主要包括三个步骤：数据预处理、细胞聚类和注释以及网络构建。此外，DeepMAPS服务器支持实时计算和交互式图形表示。除了上述进展，DeepMAPS网络服务器还强调了一个额外的功能，即阐明特定细胞类型对外部刺激的反应中的复杂网络。
4.总结
DeepMAPS是一个深度学习框架，实现了异构图表示学习和图转换器，用于从scMulti-omics数据中研究生物网络。通过构建包含细胞和基因的异构图，DeepMAPS可以同时识别它们的联合嵌入，并能够在完整的框架中推断细胞类型特异性生物网络以及细胞类型。此外，应用异构图转换器将细胞-基因关系建模为可解释的统一多关系。通过这种方式，可以大大缩短图中的训练和学习过程，从而从更远的距离考虑细胞的影响。
扫一扫
1046
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
