多模态单细胞测序技术从多层基因组数据中提供了丰富的细胞异质性信息。然而，在没有正确消除模态偏差的情况下去分析联合空间，往往会得到比单模态分析更差的聚类结果。如何有效利用多组学额外信息来描绘细胞状态并识别有意义的信号仍然是一个重大的挑战。因此，作者提出了一个名为SAILERX的深度学习框架，用于高效、稳健和灵活地分析多模态单细胞数据。SAILERX由一个具有不变表示学习的变分自动编码器组成，用于校正单一模态中的技术噪声，以及一个多模态数据对齐机制，用于整合来自不同模态的信息。SAILERX不是通过将两种模态投影到共享的潜在空间来进行硬对齐，而是鼓励通过成对相似性测量的两种模态的局部结构相似。这种策略对噪声的过拟合更具鲁棒性，这有助于各种下游分析，如聚类、插补和marker基因检测。此外，不变表示学习部分使SAILERX能够对多模态和单模态数据集进行综合分析，使其成为适用于更通用场景的工具。

来自：Integrated analysis of multimodal single-cell data with structural similarity

背景概述

单细胞测序提供了对单个细胞遗传信息的全基因组测量。最近的技术允许在同一细胞中同时分析多种模态，使我们能够从多个层面剖析细胞异质性。目前已经开发了多种计算方法来处理数据整合，例如在scRNA-seq空间中校正批次效应（scVI、scANVI、Scanorama、Harmony等），以及将多模态数据嵌入到联合空间中。然而，由于模态间潜在的偏差和模态内固有的批次效应，有效利用跨模态信息仍然是一个挑战。最近，已经开发了许多计算方法来分析多模态单细胞数据。许多方法使用的一种常见策略是将不同模态的数据投影到共享空间。例如现有的方法，scAI、scMM、scMVAE、BABEL和Cobolt，还有使用非负矩阵分解（NMF）或编码器-解码器类型的神经网络，将多种模态投影到共同的潜在空间。他们的基本假设是，来自不同模态的测量具有同等的信息量，并且具有共同的分布，而在现实世界中，这种假设是不成立的。（但是注意，在直觉上，考虑对比学习应该可以解决这个问题）

例如，典型的 scATAC-seq 实验通常报告在整个 32 亿碱基对基因组中，每个细胞有 1000–20000 个可定位片段，与来自同一细胞的 RNA 模态相比，这导致明显更高的丢失率。因此，一系列工作指出，模态的直接融合可能会引入严重的过拟合，导致学习的联合表征的分离度较差。

鉴于此，Signac 使用加权最近邻 WNN graph 根据来自每个细胞的两种模态的数据的可预测性生成联合嵌入。然而，信息融合是在生成单独嵌入之后完成的，而没有考虑两种模态之间的潜在交互，这可能会限制整体性能。大多数现有方法没有包含明确的机制来处理每种模态中的技术噪声，但是这在实际数据中很常见。

因此，为了解决这些问题，作者提出了一个名为SAILERX的深度学习框架，以改进单模态和多模态单细胞测序数据集的混合分析。与现有方法不同，SAILERX可以处理平行的（配对的）scRNA-seq和scATAC-seq多模态数据、单模态scATACseq数据，以及scATAC-seq与有配对关系的scRNA-seq和scATAC-seq的混合数据。为了解决模态异质性并避免过拟合，作者使用更稳健的基因表达信息作为参考模态，以规范染色质可及性模态的学习过程。具体而言，scATAC-seq数据是用变分自动编码器（VAE）建模的，scRNA-seq数据的嵌入是预先训练的，而不是在训练时明确建模。通过最小化两种模态之间嵌入空间中的成对相似性之间的距离来进一步实施正则化（见图1），这鼓励细胞的局部结构与参考模态相似，同时适应模态之间不同的技术噪声。由此产生的细胞表示隐含地包含来自两种模态的信息，并避免了过拟合风险。

SAILERX允许将具有scATAC-seq测量的数据集与具有配对scRNA-seq和scATAC-seq的数据集混合集成，有效地利用来自高质量多模态数据的信息来改进单模态数据集的分析。作者将SAILERX与现有SOTA方法进行基准测试，用于在三个具有不同测序技术和组织类型的流行单细胞数据集上进行多或单模态单细胞数据分析。

结果表明SAILERX生成的细胞表示可以提供更好的聚类和插补。作者还展示了单模态scATAC-seq数据集如何从混合数据学习中受益。对于生物学应用，这些改进有利于染色质可及性数据的下游分析。

数据集

在这项研究中，作者重点研究了具有配对scRNA-seq和scATAC-seq测量的多模态单细胞测序数据。使用了三个流行的具有不同细胞类型和测序技术的公共单细胞多组学数据集，即10x Genomics PBMC dataset，Share-seq dataset 和 SNARE-seq dataset。

PBMC dataset

10X基因组学提供多个PBMC细胞数据集，作者从10X基因组网站收集PBMC 10k Multiome和PBMC 3k。PBMC 10k 数据集主要用于对跨模态整合性能进行基准测试。对于PBMC 3k数据集，仅使用染色质可及性数据与PBMC 10k数据集进行混合分析。3k数据集的基因表达模态不用于混合分析的训练，仅使用在3k数据集中识别细胞的GT标签。对于 sc-multiome 数据集的整合，通常使用基因表达模态。对于sc-multiome的这两个数据集，通过Seurat和SeuratDisk包中的工具，使用现有的PBMC参考数据集通过标签转移对细胞类型进行注释。具体而言，使用Seurat内部的高质量数据集作为参考数据集，将细胞类型标签分别转移到PBMC 3k和PBMC 10k。

跨模态整合

对于场景一（跨模态整合），10k Multiome数据集是从10X基因组网站获取的。首先从10X基因组多组学数据集下载PBMC 10k表达矩阵和染色质可及性矩阵及其fragment文件，遵循与Signac相同的质量控制方案来过滤低质量细胞。

这保留了11,331个细胞用于进一步分析：

对于scRNA-seq，作者使用默认参数的SCTransform函数对scRNA-seq数据进行标准化。之后，使用PCA提取前50个PC，用于进一步聚类和与scATAC-seq的联合分析。
对于scATAC-seq，由于使用CellRanger识别的一组峰经常合并附近的峰，这可能会在motif富集分析等任务中造成偏差，因此，作者通过使用fragment文件以及使用MACS2软件生成独立的峰。之后，保留常染色体数据，并通过细胞矩阵获得最终的scATAC-seq峰。该矩阵进一步用于处理并与其他方法进行基准测试。例如，在Signac中，会对scATAC-seq矩阵做TF-IDF，然后对TF-IDF输出矩阵采用SVD来获得50维潜在嵌入，该嵌入进一步用于与scRNA-seq数据进行聚类和联合分析。

混合分析

关于第二种情况（混合分析），使用上述由scRNA-seq和scATAC-seq数据组成的多模态PBMC 10k数据作为参考，并使用PBMC 3k数据集的染色质可及性数据进行联合分析。作者从10X基因组中检索PBMC 3k scATAC-seq数据，并将其作为单一模态数据集。理由是，具有scATAC-seq的3k数据集比多组学的10k数据集包含的信息更少，然而，由于它们来自相同类型的细胞，我们可以使用10k多组学的数据集作为高质量reference来帮助分析3k scATAC-seq数据。

作者使用GenomicRanges包中的reduce函数合并来自scATAC-seq 10k和3k数据集的共峰，并分别为两个scATAC-seq数据重建peak-by-cell matrix（更多解释见单细胞跨模态分析综述），该矩阵进一步用于训练和评估模型，如图1B所示。

fig1

图1：SAILERX的总体设计。SAILERX结合了来自基因表达测量的信息，以改进染色质可及性的下游分析。SAILERX还可以使用一种或多种模态对多个数据集进行综合分析。
A：SAILERX将联合测量的scRNA-seq和scATAC-seq数据作为输入。scATAC-seq数据是通过VAE用不变表示学习建模的，而scRNA-seq的嵌入是在预训练期间处理的，而不是在训练过程中明确建模。通过最小化两种模态的成对余弦相似性矩阵之间的距离，施加正则化以鼓励嵌入空间中的局部结构在两种模态之间相似。潜在的scATAC-seq特征被进一步用于执行下游分析。
B：SAILERX还能够通过混合训练将单模态scATAC-seq与多模态数据集集成，这可以进一步提高对单模态数据的聚类性能。

Share-seq dataset

对于Share-seq数据集，作者从Ma等人的工作中获得Share-seq小鼠皮肤数据集，该数据集包含34,474个联合测量scRNA-seq和scATAC-seq的细胞数据。对于scRNA-seq数据，作者通过使用Signac包中默认参数的SCTransform函数，对其基因进行归一化，然后使用PCA获得前50个PC进行进一步分析。对于scATAC-seq数据，保留了Ma等人工作中使用的预处理的矩阵。

SNARE-seq dataset

对于SNARE-seq数据集，作者从 Chen 等人的工作中下载了两种模态的成人大脑皮层数据。对于scRNA-seq，遵循与之前相同的处理步骤，使用具有默认参数的 SCTransform 函数对基因进行归一化。之后，在归一化矩阵上采用 PCA，并使用前 50 个 PC 作为潜在嵌入以进行进一步分析。对于scATAC-seq，按照与BABEL相同的处理程序，在保持原始峰不变的情况下过滤掉低质量的细胞。

实验项目

去噪

通过数据重建来生成插补数据（去噪数据）。在评估插补结果时，首先使用每种方法生成插补数据。然后使用PCA将去噪数据投影到低维空间，并使用UMAP进行2D可视化。对于MAGIC的基准测试，作者使用scRNA-seq和scATAC-seq生成的数据进行公平比较。RNA图基于Seurat嵌入，ATAC图基于MAGIC自己的管道。对于scOpen的基准测试，遵循其GitHub网站上的手册生成密集插补矩阵。然后对估算的矩阵进行PCA，并用UMAP进行可视化。定量分数（ARI、NMI、Silhouette）是根据PC生成的聚类结果计算的。

NMI，Silhouette常用于评价生物异质性保护程度；

ARI（Adjusted Rand Index）是一种用于评估聚类结果的指标，它衡量了两个聚类结果之间的相似性。ARI的取值范围在 $[- 1, 1]$ 之间，ARI的值越接近1，表示聚类结果与真实标签的一致性越高；值越接近0，表示聚类结果与真实标签的独立性越高；而负值则表示聚类结果与真实标签的随机性更高。

ARI为0表示聚类结果与真实标签之间没有相关性。这并不意味着聚类结果是好或坏，而是表示聚类结果与真实标签之间的配对关系是随机的，没有明显的一致性。

fig2

图2：PBMC 10k的去噪结果。
A：SAILERX、scOpen和MAGIC生成的PBMC染色质可及性数据的UMAP可视化（MAGIC还包含了scRNA-seq的去噪），细胞由GT标签着色。
B：定量指标测量。

Marker基因分析

为了进一步评估细胞簇的质量，作者观察了PBMC数据集中标记为CD4 naive cell和B naive cell的簇以及SNARE-seq数据集中L4 cell和Pvalb cell的标记基因的表达。选择这4种细胞进行分析，是因为这些细胞与其他细胞类型相似，对它们进行聚类是具有挑战性的。在embedding空间中，CD4簇与CD8 naive和其他CD4亚型簇非常接近。L4簇靠近L2/3和L6 IT簇。特别是，单独的基因表达信息不能很好地分离出B细胞的亚型。

用于可视化的细胞类型特异性marker基因由Seurat中的FindMarker函数调用。这些基因被鉴定为标记基因，因为它们在用相应细胞类型标记的细胞中与其他细胞相比显示出显著的差异RNA表达量。细胞类型标签基于GT标签。与每种细胞类型相关的前10个选择的marker基因如表1所示。

细胞类型	marker基因
Pvalb	Erbb4, Cemip, Lrrc4c, Slit2, Cntnap4, Btbd11, Zfp536, Esrrg, Kcnc1,Cntnap5c
L4	Car10, Unc5d, Rorb, Pcdh15, Dcc, Gria4, Prkg1, Fstl4, Kcnh5, Cpne9
CD4 Naive	Bach2, Fhit, Igf1r, Ccr7, Ak5, Apba2, Lef1, Maml2, Sell,Satb1-as1
B Naive	Ighm, Ighd, Tcl1a, Bach2, Col19a1, Il4r, Skap1, Camk2D, Foxp1, Khdrbs2

表1：细胞类型与marker基因字典。

对于每种细胞类型，使用boxplot来可视化标有相应细胞类型的簇中细胞的标记基因的平均归一化表达（图3）。基因表达值通过 scTransform 归一化。 SALIERX 和其他方法之间的pairwise t-test 表明来自 SALIERX 的标记基因显示出比其他方法更高的表达。

对数据聚类到embedding空间，然后进行细胞类型注释（可以是有监督学习或KNN参考映射），也可以是人工辅助注释。

fig3

图3：对比不同方法细胞分类后的marker基因表达，ns表示不显著。

motif分析

作者对几个关键基序（key motifs）进行了基序分析，以证明在不同细胞类型之间发现细胞类型特异性基序富集的案例。

细胞类型注释是通过与上一节相同的程序通过聚类和多数投票（参考映射）来确定的。作者首先通过运行chromVAR计算每个细胞基序的活性得分，这为识别不同细胞类型之间的差异活性基序提供了方法。之后，应用MotifPlot绘制前6个基序中的4个，这4个基序代表了两种细胞类型之间差异表达最多的基序。

TF（转录因子）是一类能够结合到DNA上调控基因表达的蛋白质（可以促进也可以抑制），而motif（基序）则是指在DNA序列中存在的具有特定模式的短序列，通常与TF的结合位点相关联。研究转录因子结合位点的基序可以帮助我们理解具体基因调控机制的物理过程（可以想象为建模配体-受体的相互作用）

fig4

图4：motif富集分数。包含SNARE-seq数据集和插补数据集（由SAILERX进行插补），鉴定了Pvalb和L5-PT的细胞的基序z评分。对于每种细胞类型，都选择四个已知富集的基序（已知关键基序）。在SAILERX和所有其他方法之间进行pairwise t-test。
可以发现，SAILERX能够发现关键基序，如果在插补后的数据上，SAILERX能够发现关键基序的置信度更高。

多模态数据整合

fig5

图5：PBMC 10k多模态数据集的结果。用GT标签着色细胞。
A：通过不同方法生成的PBMC 10k Multiome embedding 的UMAP可视化。红圈显示在Seurat（仅scRNA-seq）、SAILER（仅scATAC-seq）和SAILERX（多模态）下B细胞亚簇的分离。
B：不同方法生成的ARI、NMI和Silhouette评分的聚类定量指标。误差条是通过用90%的随机子采样重复实验产生的。
C：作者主要使用Seurat的PCA作为scRNA-seq嵌入，但此处也展示了其他方法（即橙色reference：RNA模态的embedding聚类ARI分数）。蓝色为SAILERX加入对应RNA嵌入后获得的联合嵌入的聚类ARI分数。

fig6