单细胞空间多组学整合¶

一句话概述¶

单细胞空间多组学整合是将不同模态的单细胞数据（RNA、ATAC、蛋白、空间位置等）组合在一起分析——白话说就是：每种检测技术只能看到细胞的"一面"，整合分析就是把多面信息拼成一个完整的"3D肖像"。

核心知识点表格¶

知识点	说明
配对整合	多种模态来自同一批细胞（如CITE-seq、10x Multiome）
非配对整合	不同实验、不同细胞，需要对齐到共同空间
马赛克整合	有些细胞只有部分模态数据
WNN (Seurat)	加权近邻，用于配对多模态整合
MOFA+	多组学因子分析，可解释性强
totalVI	深度学习VAE，CITE-seq整合
MultiVI	处理缺失模态的VAE，支持RNA+ATAC+蛋白
GLUE	图链接统一嵌入，非配对整合的SOTA
SpatialGlue	空间多组学整合，利用空间信息
SIMO	2025年新方法，概率对齐的空间多组学整合

各步骤详解¶

第一步：理解多模态整合类型¶

白话解释： 多模态整合有三种场景，就像拼图的不同难度级别：

类型1 - 配对整合（简单）
  同一个细胞同时测了RNA和蛋白（如CITE-seq）
  → 数据天然配对，整合最容易
  → 工具：WNN, totalVI, MOFA+

类型2 - 非配对整合（中等）
  RNA-seq和ATAC-seq分别在不同细胞上做
  → 需要找到"锚点"来对齐两批细胞
  → 工具：Seurat v5 bridge, GLUE, Harmony

类型3 - 马赛克整合（困难）
  有些细胞测了3种模态，有些只测了1-2种
  → 需要处理大量"缺失数据"
  → 工具：MultiVI, StabMap, UINMF

第二步：配对整合——Seurat WNN¶

白话解释： 当RNA和另一种模态（蛋白/ATAC）来自同一批细胞时，WNN是最简单直接的方法。

library(Seurat)                          # 加载Seurat

# ===== 10x Multiome (RNA + ATAC) =====
# 假设已有RNA和ATAC的Seurat对象

# 各模态独立预处理
# RNA预处理
DefaultAssay(multiome) <- "RNA"          # 切换到RNA
multiome <- SCTransform(multiome)        # SCTransform归一化
multiome <- RunPCA(multiome)             # RNA的PCA

# ATAC预处理
DefaultAssay(multiome) <- "ATAC"         # 切换到ATAC
multiome <- RunTFIDF(multiome)           # TF-IDF归一化（ATAC特有）
multiome <- FindTopFeatures(multiome, min.cutoff = "q0")  # 选择特征
multiome <- RunSVD(multiome)             # LSI降维（ATAC版的PCA）

# WNN整合
multiome <- FindMultiModalNeighbors(
  multiome,                              # Seurat对象
  reduction.list = list("pca", "lsi"),   # RNA的PCA + ATAC的LSI
  dims.list = list(1:50, 2:50),          # 各用多少维（LSI跳过第1维）
  modality.weight.name = "RNA.weight"    # 权重命名
)

# 聚类和UMAP
multiome <- FindClusters(multiome, graph.name = "wsnn", resolution = 0.8)
multiome <- RunUMAP(multiome, nn.name = "weighted.nn",
  reduction.name = "wnn.umap")

# 可视化
DimPlot(multiome, reduction = "wnn.umap", label = TRUE)

# 查看RNA vs ATAC的权重分布
VlnPlot(multiome, features = "RNA.weight",
  group.by = "seurat_clusters", pt.size = 0)

第三步：非配对整合——GLUE¶

白话解释： GLUE（Graph-Linked Unified Embedding）用已知的基因调控关系（比如ATAC的peak附近有哪些基因）来把不同模态对齐到同一个坐标空间。

# 安装scglue
# pip install scglue

import anndata as ad                     # AnnData数据结构
import scanpy as sc                      # 单细胞分析
import scglue                            # GLUE整合方法
import networkx as nx                    # 网络操作

# 读取两种模态的数据
rna = sc.read_h5ad("rna_data.h5ad")      # scRNA-seq数据
atac = sc.read_h5ad("atac_data.h5ad")    # scATAC-seq数据

# 各模态独立预处理
# RNA预处理
sc.pp.normalize_total(rna)               # 归一化
sc.pp.log1p(rna)                         # log转换
sc.pp.highly_variable_genes(rna, n_top_genes=2000)  # 高变基因
sc.tl.pca(rna)                           # PCA

# ATAC预处理
sc.pp.highly_variable_genes(atac, n_top_genes=50000)  # 高变peak
scglue.data.lsi(atac, n_components=100)  # LSI降维

# 构建guidance图（基因-peak对应关系）
# 使用基因组注释信息连接RNA基因和ATAC peak
guidance = scglue.genomics.rna_atac_guidance(
    rna,                                 # RNA的AnnData
    atac,                                # ATAC的AnnData
    gtf_file="genes.gtf"                 # 基因注释文件
)

# 训练GLUE模型
scglue.models.configure_dataset(
    rna, "NB",                           # RNA用负二项分布
    use_highly_variable=True,            # 用高变基因
    use_rep="X_pca"                      # 用PCA表示
)
scglue.models.configure_dataset(
    atac, "NB",                          # ATAC也用负二项
    use_highly_variable=True,            # 用高变peak
    use_rep="X_lsi"                      # 用LSI表示
)

# 训练模型
glue = scglue.models.fit_SCGLUE(
    {"rna": rna, "atac": atac},          # 输入数据字典
    guidance,                            # guidance图
    fit_kws={"directory": "glue_model"}  # 模型保存目录
)

# 获取整合后的嵌入
rna.obsm["X_glue"] = glue.encode_data("rna", rna)   # RNA嵌入
atac.obsm["X_glue"] = glue.encode_data("atac", atac)  # ATAC嵌入

# 合并后可视化
combined = ad.concat([rna, atac])        # 合并两种模态
sc.pp.neighbors(combined, use_rep="X_glue")  # 构建近邻图
sc.tl.umap(combined)                     # UMAP
sc.pl.umap(combined, color=["modality", "cell_type"])  # 可视化

第四步：MOFA+因子分析¶

白话解释： MOFA+（Multi-Omics Factor Analysis）把多种模态分解成"隐藏因子"——每个因子可能由某种模态驱动。比如因子1主要由RNA驱动，反映细胞分化状态；因子2主要由蛋白驱动，反映细胞活化状态。

library(MOFA2)                           # 加载MOFA2

# 从Seurat对象创建MOFA对象
# 假设cite_obj是包含RNA和ADT的Seurat对象
mofa_obj <- create_mofa_from_Seurat(
  cite_obj,                              # Seurat对象
  groups = "sample",                     # 按样本分组（可选）
  assays = c("RNA", "ADT")              # 使用RNA和ADT两种模态
)

# 设置参数
data_opts <- get_default_data_options(mofa_obj)  # 数据选项
model_opts <- get_default_model_options(mofa_obj)  # 模型选项
model_opts$num_factors <- 15             # 提取15个因子

train_opts <- get_default_training_options(mofa_obj)  # 训练选项
train_opts$maxiter <- 1000               # 最大迭代次数
train_opts$seed <- 42                    # 随机种子

# 准备和训练
mofa_obj <- prepare_mofa(
  mofa_obj,                              # MOFA对象
  data_options = data_opts,              # 数据选项
  model_options = model_opts,            # 模型选项
  training_options = train_opts          # 训练选项
)
mofa_obj <- run_mofa(mofa_obj, outfile = "mofa_result.hdf5")  # 训练

# 结果可视化
# 查看各因子解释的方差比例
plot_variance_explained(mofa_obj)         # 方差解释图

# 各因子由哪种模态驱动
plot_variance_explained(mofa_obj, plot_total = TRUE)

# 因子值的UMAP可视化
plot_factor(mofa_obj, factors = 1:5)     # 显示前5个因子

第五步：空间多组学整合¶

白话解释： 当有空间信息时，整合不仅要考虑"表达相似"，还要考虑"空间位置邻近"。空间感知的整合方法比纯单细胞方法效果更好。

# ===== SIMO：空间多组学整合（2025年新方法）=====
# pip install simo-omics

# SIMO通过概率对齐整合空间转录组和单细胞多组学
import simo                              # 导入SIMO

# 读取数据
st_data = sc.read_h5ad("spatial.h5ad")   # 空间转录组
sc_rna = sc.read_h5ad("scrna.h5ad")      # 单细胞RNA
sc_atac = sc.read_h5ad("scatac.h5ad")    # 单细胞ATAC

# SIMO整合
model = simo.SIMO(
    spatial=st_data,                     # 空间数据
    single_cell=[sc_rna, sc_atac],       # 单细胞多模态列表
    modality_names=["RNA", "ATAC"]       # 模态名称
)
model.fit()                              # 训练模型

# 获取整合结果
integrated = model.get_integrated_data() # 整合后的数据
sc.pl.spatial(integrated, color="cell_type")  # 空间可视化

常见报错与解决¶

报错	原因	解决方案
`WNN: modalities have different cells`	RNA和ATAC的细胞barcode不一致	取交集：cells <- intersect(colnames(rna), colnames(atac))
`MOFA: convergence issue`	因子数太多或数据太稀疏	减少num_factors，过滤低表达基因
`GLUE: guidance graph empty`	GTF文件和peak注释不匹配	检查染色体命名是否一致（chr1 vs 1）
`MultiVI: CUDA out of memory`	GPU内存不够	减少batch_size或使用CPU
`LSI dim 1 highly correlated with depth`	ATAC的第一个LSI维度反映测序深度	FindMultiModalNeighbors时跳过dim1（dims=2:50）
`Seurat: assay not found`	切换assay名称拼写错误	用Assays(obj)查看实际assay名

速查表¶

# ========== 多模态整合方法选择指南 ==========

场景                    推荐方法         特点
配对RNA+蛋白(CITE-seq)  Seurat WNN       简单、主流、R语言
配对RNA+ATAC(Multiome)  Seurat WNN       同上
配对+批次校正           totalVI/Smmit    深度学习、处理批次效应
非配对RNA+ATAC          GLUE             利用基因组注释对齐
马赛克(部分模态缺失)    MultiVI/StabMap  可以插补缺失模态
空间+单细胞             SIMO/SpatialGlue 利用空间邻近信息
可解释性需求            MOFA+            输出可解释的因子

面试高频问题¶

Q1: 多模态整合的三种类型是什么？

(1) 配对整合(Paired/Vertical)：同一批细胞测多种模态（CITE-seq、10x Multiome），直接对应；(2) 非配对整合(Unpaired/Diagonal)：不同实验的不同细胞，需要找锚点对齐；(3) 马赛克整合(Mosaic)：部分细胞有完整模态，部分缺失。难度递增，方法各异。

Q2: WNN整合的核心思想是什么？

WNN为每个细胞的每种模态计算一个"信息量权重"。步骤：(1) 各模态独立构建KNN图；(2) 学习每个细胞在每种模态中的"局部信息量"；(3) 按权重组合多个KNN图为一个WNN图；(4) 在WNN图上做聚类和UMAP。关键优势：权重是细胞特异的——同一种细胞类型中，不同细胞的模态权重可以不同。

Q3: 空间多组学整合和普通单细胞多组学整合有什么区别？

空间多组学额外利用了"空间邻近关系"——在组织中位置相近的细胞更可能有相似的功能状态。2025年的benchmark显示，SpatialGlue、SpaMI等空间感知方法显著优于Seurat WNN和MOFA+等纯单细胞方法。因此有空间信息时，应该选用空间感知的整合方法。

Q4: MOFA+的输出怎么解读？

MOFA+输出K个"因子"，每个因子捕获数据中的一种主要变异来源。通过方差解释图，可以看出每个因子由哪种模态驱动。比如：因子1解释RNA方差30%、ADT方差5%——说明这个因子主要反映转录组层面的变异（可能是细胞分化）。因子2解释ADT方差40%——说明主要反映蛋白层面变异（可能是细胞活化状态）。

Q5: 处理多模态数据时，归一化需要注意什么？

不同模态需要不同的归一化方法：RNA用SCTransform或log-normalize；ADT/蛋白用CLR（Centered Log-Ratio）；ATAC用TF-IDF+LSI。不能把所有模态用同一种方法归一化——数据分布特征完全不同。降维时也不同：RNA用PCA，ATAC用LSI（跳过第一维），ADT直接用PCA但维度较低。