单细胞多样本整合之scVI和scANVI

生信菜鸟团

发布于 2023-08-23 08:50:07

4.2K0

发布于 2023-08-23 08:50:07

文章被收录于专栏：生信菜鸟团

我也有用python做过BBKNN和Harmony：

巧用Python加速单细胞分析

scVI（Single-Cell Variational Inference）是基于深度学习的整合算法，多次被评价为最佳的整合算法之一。scVI发表于2018年的Nature Methods上，后续经过多次的更新改进，文章题目是《Deep generative modeling for single-cell transcriptomics》。

而scANVI（single-cell ANnotation using Variational Inference）算法和scVI类似，是同一个团队做的算法，发表在2021年的Mol Syst Biol上，文章题目是《Probabilistic harmonization and annotation of single-cell transcriptomics data with deep generative models》。

注意，scVI和scANVI都是用于单细胞转录组数据分析的深度生成模型，它们都是scvi-tools包的一部分。scVI可以用于数据集成、降维、聚类、差异表达等任务，而scANVI可以用于半监督的细胞类型注释，即利用部分已知的细胞类型标签来推断其他细胞的状态（需要预先对每个批次的单细胞数据进行预注释）。scANVI可以看作是在scVI的基础上增加了一个分类器，因此可以从预训练的scVI模型初始化。

官网教程在https://docs.scvi-tools.org/en/stable/tutorials/notebooks/harmonization.html，
github在https://github.com/scverse/scvi-tools

需要强调的是用scVI/scANVI整合分析的话比较吃资源，如果有GPU的话速度会快5-10倍；如果没有GPU，仅用CPU运行scVI，还需谨慎。因为CPU模式的scVI运行速度实在是慢，特别是单细胞大样本数据的运行。

一. 环境部署

mamba create -n scvi python=3.9
which python
#/home/data/fuli09/miniconda3/envs/scvi2/bin/python

/home/data/fuli09/miniconda3/envs/scvi2/bin/python -m pip install scvi-tools anndata numpy scanpy scib certifi scib-metrics pymde scvi-colab
#Please be sure to install a version of PyTorch that is compatible with your GPU (if applicable).

还需要在python里运行：

from scvi_colab import install

install()

如果需要衔接Jupyter notebook使用的话，需要在scvi这个环境里安装一下几个插件，然后就能在Jupyter notebook里选择scvi这个环境（ps：我用的是JupyterLab，和Jupyter notebook是差不多的）：

mamba install -y nb_conda_kernels ipykernel
python -m ipykernel install --user --name scvi --display-name scvi

JupyterLab

二. 运行示例数据

import scanpy as sc
import scvi
from rich import print
from scib_metrics.benchmark import Benchmarker
from scvi.model.utils import mde

sc.set_figure_params(figsize=(4, 4))

%config InlineBackend.print_figure_kwargs={'facecolor' : "w"}
%config InlineBackend.figure_format='retina'

Step1. 读入数据

adata = sc.read(
    "data/lung_atlas.h5ad",
    backup_url="https://figshare.com/ndownloader/files/24539942",
)
adata

AnnData object with n_obs × n_vars = 32472 × 15148
    obs: 'dataset', 'location', 'nGene', 'nUMI', 'patientGroup', 'percent.mito', 'protocol', 'sanger_type', 'size_factors', 'sampling_method', 'batch', 'cell_type', 'donor'
    layers: 'counts'

Step2. 数据集预处理

该数据集已经按照scIB手稿中的描述进行了处理。通常，scvi工具中的模型期望以与Scanpy/Seurat相同的方式进行标准分析。

值得注意的是高度可变的基因选择问题。虽然scVI和scANVI都可以在运行时使用所有基因，但作者通常建议过滤基因以获得最佳整合性能。同时，这也有助于消除由于批次特异性基因表达而产生的批次特异性变异。

一般来说，用户可以使用Scanpy的标准流程执行高变基因选择。

adata.raw = adata  # keep full dimension safe
sc.pp.highly_variable_genes(
    adata,
    flavor="seurat_v3",
    n_top_genes=2000,
    layer="counts",
    batch_key="batch",
    subset=True,
)

Step3. 使用scVI进行整合分析

第一步，根据count数据和batch信息创建模型：

scvi.model.SCVI.setup_anndata(adata, layer="counts", batch_key="batch")
vae = scvi.model.SCVI(adata, n_layers=2, n_latent=30, gene_likelihood="nb")
#Now we train scVI. This should take a couple of minutes on a Colab session
vae.train()

将结果添加至obsm，命名为X_scVI：

adata.obsm["X_scVI"] = vae.get_latent_representation()

最后，对数据集进行聚类，并可视化：

sc.pp.neighbors(adata, use_rep="X_scVI")
sc.tl.leiden(adata)
adata.obsm["X_mde"] = mde(adata.obsm["X_scVI"])
sc.pl.embedding(
    adata,
    basis="X_mde",
    color=["batch", "leiden"],
    frameon=False,
    ncols=1,
)

image-20230723214135712

因为该数据集已做注释，所以这里可以根据注释结果来评价scVI的效果：

sc.pl.embedding(adata, basis="X_mde", color=["cell_type"], frameon=False, ncols=1)

image-20230723214212112

Step4. 使用scANVI进行整合

作者说，如果输入的数据有预注释的信息，使用scANVI可以产生更准确的整合结果。

由于上述已经训练了一个scVI模型，我们将使用它来初始化scANVI。在初始化scANVI时，我们为它提供labels_key。由于scANVI也可以用于具有部分观察到的注释的数据集，因此我们需要为它指定与未标记细胞对应的类别名称。由于我们没有未标记的细胞，我们可以给它任意一个不是现有细胞类型名称的随机名称。

lvae = scvi.model.SCANVI.from_scvi_model(
    vae,
    adata=adata,
    labels_key="cell_type",
    unlabeled_category="Unknown",
)
lvae.train(max_epochs=20, n_samples_per_label=100)

将结果添加为"X_scANVI"

adata.obsm["X_scANVI"] = lvae.get_latent_representation(adata)

结果可视化：

adata.obsm["X_mde_scanvi"] = mde(adata.obsm["X_scANVI"])
sc.pl.embedding(
    adata, basis="X_mde_scanvi", color=["cell_type"], ncols=1, frameon=False
)

image-20230723215034705

Step5. 计算整合效果

这里我们使用scIB-metrics包，它包含scIB基准测试套件中使用的度量的可伸缩实现。我们可以使用这些度量来评估集成的质量。

我们可以看到，使用标签信息和scANVI的额外训练改进了捕获生物保护的指标(cLISI, Silhouette标签)，而不会牺牲太多的批校正能力(iLISI, Silhouette批)。

bm = Benchmarker(
    adata,
    batch_key="batch",
    label_key="cell_type",
    embedding_obsm_keys=["X_pca", "X_scVI", "X_scANVI"],
    n_jobs=-1,
)
bm.benchmark()

bm.plot_results_table(min_max_scale=False)

image-20230723215710408

df = bm.get_results(min_max_scale=False)
print(df)

              Isolated labels        Leiden NMI        Leiden ARI  \
Embedding                                                           
X_pca                 0.50606          0.724169          0.474548   
X_scVI               0.540021          0.750638           0.58713   
X_scANVI             0.536124          0.771481          0.639793   
Metric Type  Bio conservation  Bio conservation  Bio conservation   

             Silhouette label             cLISI  Silhouette batch  \
Embedding                                                           
X_pca                0.557418               1.0          0.878176   
X_scVI               0.535189          0.994187          0.900951   
X_scANVI             0.588762          0.998546           0.88585   
Metric Type  Bio conservation  Bio conservation  Batch correction   

                        iLISI              KBET Graph connectivity  \
Embedding                                                            
X_pca                0.003306          0.223596           0.728504   
X_scVI               0.111625          0.339272           0.899708   
X_scANVI             0.105794          0.349628           0.922797   
Metric Type  Batch correction  Batch correction   Batch correction   

               PCR comparison Batch correction Bio conservation  \
Embedding                                                         
X_pca                     0.0         0.366716         0.652439   
X_scVI               0.891112         0.628534         0.681433   
X_scANVI             0.809051         0.614624         0.706941   
Metric Type  Batch correction  Aggregate score  Aggregate score   

                       Total  
Embedding                     
X_pca                0.53815  
X_scVI              0.660273  
X_scANVI            0.670014  
Metric Type  Aggregate score