跟着Theis Lab学单细胞: 格式转换，你需要知道底层逻辑（下）
跟着Theis Lab学单细胞: 格式转换，你需要知道底层逻辑（下）
单细胞最佳实践
基于内存的互操作性
本章导读 SECTION PREVIEW
1. 使用Seurat内置函数在SingleCellExperiment和Seurat对象间转换
2. 使用rpy2和reticulate包实现R和Python间的相互调用
3. 了解多模态数据的转换
本节对应原书 5.4 In-memory interoperability
上一篇文章中我们提到，格式转换分为磁盘上、内存中转化两种策略。在这一篇中，我们将讨论基于内存的单细胞格式转化的方式和原理（即“内存互操作” in-memory interoperability）。这类方法常常需要同时运行两个编程语言环境，并让它们访问同一个对象，或在需要时进行转换。通常情况下，一种语言（如R）作为主要环境，另一种语言（如Python）则通过接口来访问数据。在进行需要实时交互的分析时，这种策略使得研究者能在两种语言中同时进行工作。需要创建多语言的文档时（比如此教程）内存互操作也经常使用。然而，缺点也很明显：研究者需要对两种语言都相当熟悉，而且某些过于复杂的对象可能并不被两种编程语言支持。同时，内存互操作对硬件要求和负担也更大，因为数据被复制时会占用更多空间（因此对大数据不太友好）。
01
R生态系统（sce、seurat）内部的互操作
在讨论R和Python间的转换前，我们先来看看R内部的两个单细胞生态系统（Seurat、SingleCellExperiment）之间的转化。这个操作比较简单，因为Seurat包就提供了这种转换函数。
Seurat --> SingleCellExperiment:
sce_from_seurat <- Seurat::as.SingleCellExperiment(seurat)
sce_from_seurat
class: SingleCellExperiment
dim: 2000 100
metadata(0):
assays(2): X logcounts
rownames(2000): Gene-0 Gene-1 ... Gene-1998 Gene-1999
rowData names(0):
colnames(100): Cell_0 Cell_1 ... Cell_98 Cell_99
colData names(9): n_genes_by_counts log1p_n_genes_by_counts ...
pct_counts_in_top_500_genes ident
reducedDimNames(2): PCA UMAP
mainExpName: NULL
altExpNames(0):
SingleCellExperiment --> Seurat:
seurat_from_sce <- Seurat::as.Seurat(sce_from_seurat)
seurat_from_sce
An object of class Seurat
2000 features across 100 samples within 1 assay
Active assay: RNA (2000 features, 0 variable features)
2 dimensional reductions calculated: PCA, UMAP
这里困难的部分在于两种对象结构的差异，我们需要确保设置正确的参数，这样转换函数才能知道该转换哪些信息以及将它们放在哪里。不过，很多时候可能不需要完全转换Seurat对象到SingleCellExperiment，因为许多核心的Bioconductor单细胞分析包也接受矩阵作为输入：
# 使用Bioconductor的scuttle包计算CPM，输入是Seurat对象中的矩阵
cpm <- scuttle::calculateCPM(Seurat::GetAssayData(seurat, slot = "counts"))
02
在Python中访问R
rpy2这个包允许我们在Python中访问R的函数和对象。比如：
counts_mat = adata.layers["counts"].T
rpy2.robjects.globalenv["counts_mat"] = counts_mat
cpm = rpy2.robjects.r("scuttle::calculateCPM(counts_mat)")
cpm
<2000x100 sparse matrix of type '<class 'numpy.float64'>'
with 126146 stored elements in Compressed Sparse Column format>
如果你使用Jupyter notebook，还可以利用IPython"接口"创建原生R代码的单元格：
%%R -i counts_mat -o magic_cpm
# 这里是R代码
magic_cpm <- scuttle::calculateCPM(counts_mat)
# 回到Python代码
magic_cpm # 现在可以在Python中使用R计算的结果
对于单细胞数据，anndata2ri包特别有用。它是rpy2的扩展，允许R将AnnData对象视为SingleCellExperiment对象，避免了不必要的转换，使R代码可以直接在Python对象上运行：
%%R -i adata
qc <- scuttle::perCellQCMetrics(adata)
head(qc)
sum detected total
Cell_0 2005 1274 2005
Cell_1 1941 1233 1941
Cell_2 2011 1270 2011
Cell_3 1947 1268 1947
Cell_4 1933 1265 1933
Cell_5 2031 1289 2031
不过，如果对象包含不支持的数据类型，可能需要先修改对象才能正确读取。
03
在R中访问Python
访问Python从R会话类似于从Python访问R，但这里的接口是由{reticulate}包提供的：
# R
# 将R对象转换为Python对象
reticulate_list <- reticulate::r_to_py(LETTERS)
print(reticulate_list)
# 使用Python内置函数
py_builtins <- reticulate::import_builtins()
py_builtins$zip(letters, LETTERS)
如果在RMarkdown或Quarto文档中工作，也可以使用reticulate的Python引擎编写原生Python代码块。这样我们就可以使用魔术变量r和py访问另一种语言中的对象：
# R
# 一个访问Python对象的R代码块
print(py$py_object)
# Python
# 一个访问R对象的Python代码块
print(r$r_object)
与anndata2ri不同，目前没有R包提供直接接口让Python将SingleCellExperiment或Seurat对象视为AnnData对象。但我们仍然可以通过reticulate访问AnnData的大部分内容：
# R
# 打印Python环境中的AnnData对象
py$adata
# 访问obs槽，pandas DataFrames会自动转换为R的data.frames
head(adata$obs)
前文提到，R包{anndata}也提供了在R中访问
AnnData的接口，但多数分析工具仍不支持这种方式。因此，在进行需要对整个对象进行复杂分析时，仍有需要将对象从R完全转换到Python（或相反）的必要。由于创建了整个数据的副本，这种方式需要耗费大量内存，但是相对的对各种分析工具的支持性也大大提升了，比如之前提到的{zellkonverter}包：
# 将AnnData转换为SingleCellExperiment
sce <- zellkonverter::AnnData2SCE(adata, verbose = TRUE)
# 反向转换
adata2 <- zellkonverter::SCE2AnnData(sce, verbose = TRUE)
04
多模态数据的互操作性
多模态数据的复杂性为互操作性带来了额外的挑战。不同的数据结构对多模态数据的支持也不同：
- • SingleCellExperiment：通过"alternative experiments"存储多个模态，要求在细胞维度（列）上匹配
- • Seurat：使用不同的"assays"存储多个模态
- • AnnData：仅限于单模态数据，无法直接存储多个模态
为了解决这个限制，MuData对象应运而生。
开发者在设计中考虑了互操作性。虽然MuData的主要平台是Python，但作者提供了：
{MuDataSeurat} R包：将磁盘上的H5MU格式读取为Seurat对象{MuData} R包：将H5MU格式读取为Bioconductor MultiAssayExperiment对象
虽然这种官方支持非常实用，但由于对象结构的差异，仍然存在一些不一致之处。
Python中使用MuData
在Python中，MuData的使用非常直观：
# 读取文件
mdata = mudata.read_h5mu("datasets/original.h5mu")
print(mdata)
# 保存新文件
python_h5mu_file = Path(temp_dir.name) / "python.h5mu"
mdata.write_h5mu(python_h5mu_file)
MuData object with n_obs × n_vars = 411 × 56
obs: 'louvain', 'leiden', 'leiden_wnn', 'celltype'
var: 'gene_ids', 'feature_types', 'highly_variable'
obsm: 'X_mofa', 'X_mofa_umap', 'X_umap', 'X_wnn_umap'
...
2 modalities
prot: 411 x 29
...
rna: 411 x 27
...
R中使用MuData
在R中，可以通过两种主要方式使用MuData：
- 1. Bioconductor
使用{MuData}包将MuData读入为MultiAssayExperiment对象：
# 读取为MultiAssayExperiment对象
mae <- MuData::readH5MU("datasets/original.h5mu")
print(mae)
# 保存为H5MU文件
bioc_h5mu_file <- tempfile(fileext = ".h5mu")
MuData::writeH5MU(mae, bioc_h5mu_file)
A MultiAssayExperiment object of 2 listed
experiments with user-defined names and respective classes.
Containing an ExperimentList class object of length 2:
[1] prot: SingleCellExperiment with 29 rows and 411 columns
[2] rna: SingleCellExperiment with 27 rows and 411 columns
05
其他语言的互操作性
除了R和Python，单细胞数据还可以与其他编程语言进行互操作：
1. Julia
- • Muon.jl：实现AnnData/MuData及H5AD/H5MU读写
- • scVI.jl：支持H5AD格式交互
JavaScript
- • Vitessce：支持Zarr格式的AnnData可视化
- • kana工具集：可解析H5AD及RDS格式
Rust
- • anndata-rs：高性能H5AD读写实现
06
后记
本节是“导言”部分的最后一章，希望读者朋友们通过这部分的学习，对单细胞分析的流程和工具有了基础的了解。下一届开始，我们将进入正式的下游分析流程，从讨论质量控制开始，了解如何选择合适的质控标准，去除低质量细胞、双细胞/多细胞和环境RNA的影响。本章英文原文及完整代码，可以在后台回复scbp05获取。
----------
单细胞领域大牛编写的最强硬核单细胞教程来啦！想要追更学习笔记，请给我们加一个小✨，这样最新的推送才会显示在主页：