SIMBA: 单细胞特征级别的嵌入.

0. TL; DR

SIMBA (Single-cell Indexing Method with Bipartite Affinity) 是一种基于图嵌入的单细胞数据分析框架，通过将cells（细胞）与features（特征，如genes、peaks、motifs、k-mers）共同嵌入到共享的潜在空间中，实现了无需聚类的标记发现、跨批次校正和多组学整合。

该方法利用PyTorch-BigGraph进行可扩展的图嵌入，结合Softmax变换和实体类型约束，能够在统一框架下处理scRNA-seq、scATAC-seq及多模态数据，在保持生物学变异的同时消除技术批次效应，并直接推断基因调控关系。

1. 背景介绍

单细胞测序技术（single-cell sequencing）的快速发展使研究者能够在单个细胞分辨率上解析基因组（genomics）、表观基因组（epigenomics）、转录组（transcriptomics）和蛋白质组（proteomics）等多个层面的信息。然而，现有的计算方法大多针对特定任务设计，存在以下局限：

  1. 聚类中心性（cluster-centric）：传统流程（如Seurat或Scanpy）依赖降维、聚类和差异表达分析，聚类结果对分辨率参数敏感，可能导致生物学注释不一致。
  2. 任务特异性：批次校正（batch correction）、多模态分析（multimodal analysis）和多组学整合（multi-omics integration）通常需要不同的算法框架，缺乏统一解决方案。
  3. 特征关系建模不足：现有方法主要关注细胞状态学习，难以直接建模基因与调控元件（如enhancers、transcription factor (TF) motifs）之间的层次关系。

为应对这些挑战，作者提出了SIMBA。该方法将单细胞数据表示为异构图（heterogeneous graph），其中细胞和各类特征作为节点（nodes），实验测量或计算推断的关系作为边（edges）。通过多关系图嵌入技术，SIMBA将不同类型的实体映射到共同的低维空间，使得细胞与其定义性特征在嵌入空间中相邻，从而支持基于邻域查询的生物学分析，无需依赖聚类结果。

2. 方法详解

2.1 图构建（Graph Construction）

SIMBA的核心是将单细胞数据转换为多关系图 $G=(V, E)$，其中 $V$ 包含不同类型的实体，$E$ 表示实体间的关系。

根据分析任务，图中的实体可包括：

边分为两类：

  1. 实验测量边（measured edges）：如细胞-基因表达关系、细胞-peak可及性关系
  2. 计算推断边（inferred edges）：如跨批次细胞相似性、跨模态细胞对应关系、peak-motif序列包含关系

特定任务的图构建策略如下：

2.2 图嵌入（Graph Embedding）

作者采用PyTorch-BigGraph框架学习每个实体 $v \in V$ 的 $D$ 维嵌入向量 $\theta_v \in \mathbb{R}^D$（默认 $D=50$），通过优化链接预测目标函数进行训练。

对于边 $e=(u,v)$，定义分数 $s_e = \theta_u \cdot \theta_v$。优化多类log loss：

\[\mathcal{L}_e = -\log \frac{\exp(s_e)}{\sum_{e' \in \mathcal{N}} \exp(s_{e'})} w_{e'}\]

其中 $\mathcal{N}$ 为负采样（negative sampling）生成的候选边集合，$w_{e’}$ 为边权重。

由于图中只有特定实体类型间允许存在边（如细胞-基因、peak-motif），负采样时仅替换为同类型的随机实体，避免生成无效边（如基因-基因边），确保嵌入空间反映真实的生物学关系。

为处理节点度分布不均的问题，一半的负样本从同类型节点均匀采样，另一半按节点度的概率分布采样，防止嵌入过度拟合高度节点。

为防止过拟合（尤其在边-参数比率较低时），作者引入L2正则化：

\[\mathcal{L}_{\text{reg}} = \mathcal{L} + \lambda \sum_{u \in N} \sum_{d=1}^D \theta_{u,d}^2\]

其中权重衰减参数 $\lambda$ 根据训练样本量 $N_e$ 自动计算。训练过程中保留 $5\%$ 的边作为验证集，监控MRR（mean reciprocal rank）、R1、R10等指标判断收敛。

2.3 Softmax 变换与度量标准

由于不同实体类型（如细胞与peaks）具有不同的边分布，直接比较嵌入可能产生偏差。作者采用Softmax变换使特征嵌入与细胞嵌入可比：

给定细胞嵌入 ${v_{c_i}}$ 和特征嵌入 $v_{f_j}$，计算边概率：

\[p_{c_i,f_j} = \frac{\exp(v_{c_i} \cdot v_{f_j})}{\sum_{k=1}^n \exp(v_{c_k} \cdot v_{f_j})}\]

变换后的特征嵌入为细胞嵌入的加权平均：

\[\hat{v}_{f_j} = \frac{\sum_{i=1}^n p_{c_i,f_j}^{T^{-1}} v_{c_i}}{\sum_{i=1}^n p_{c_i,f_j}^{T^{-1}}}\]

其中 $T=0.5$ 为温度参数，控制分布锐度。

基于概率分布 $p_{c_i,f_j}$，作者提出四个评估特征特异性的指标：

  1. Max value：前 $k$ 个（默认 $50$）最高概率细胞的平均归一化相似度
  2. Gini index：衡量概率分布偏离均匀分布的程度，$G = \frac{\sum_{i=1}^n (2i-n-1)p_{c_i,f_j}}{n\sum_{i=1}^n p_{c_i,f_j}}$
  3. s.d.（标准差）：概率分布的离散程度
  4. Entropy（熵）：$H = -\sum p \log p$，低熵表示高特异性

通过计算null节点（保持度分布的随机重排边）的度量分布，可计算统计显著性（p-value）和假发现率（FDR）。

2.4 下游分析功能

主调控因子识别（Master Regulator Identification）： 对于转录因子，若其motif和基因在嵌入空间中距离近且均具有高细胞类型特异性（高Gini index和Max value），则判定为主调控因子。通过计算TF motif到所有基因的距离排名，结合特异性过滤，识别细胞命运决定的关键调控因子。

靶基因推断（Target Gene Inference）： 给定主调控因子，搜索其motif和基因的近邻基因（默认 $k=200$）。候选靶基因需满足：

  1. 其基因组位点附近的peaks包含该TF motif
  2. 在嵌入空间中，该基因及邻近peaks均靠近TF motif和TF基因
  3. 通过计算四种距离（基因-TF基因、基因-TF motif、peaks-TF motif、peaks-基因）的排名，筛选出符合调控逻辑的靶基因。

3. 实验分析

作者利用多个公开数据集验证SIMBA的性能，涵盖scRNA-seq、scATAC-seq、多模态（SHARE-seq、SNARE-seq、10x Multiome）、批次校正和多组学整合任务。

3.1 scRNA-seq 分析

作者在数据集 10x Genomics PBMCs（人外周血单核细胞）上评估了模型的性能。

3.2 scATAC-seq 分析

作者在一个人造血细胞数据集（含FACS分选标签）上评估模型的性能。

3.3 单细胞多模态分析

作者在SHARE-seq小鼠皮肤毛囊分化数据（同时测定RNA和ATAC）上评估模型的性能。

3.4 批次校正分析

作者在人胰腺数据集（$5$ 个批次，Baron、Muraro、Segerstolpe、Wang、Xin）和小鼠图谱数据集（$9$ 个器官系统）上评估模型的性能。

3.5 多组学整合分析

作者在10x Genomics Multiome人PBMCs数据（手动拆分为scRNA-seq和scATAC-seq以评估整合准确性）上评估模型的性能。