2.基于Python的单细胞转录组数据处理-批次整合

scRNA-seq数据天然存在批次效应，批次效应是测量的基因表达水平的变化，这是测量不同组或“批次”的结果。例如，如果两个实验室从同一队列中采集样本，但这些样本的解离方式不同，则可能会出现批次效应。如果实验室 A 优化其解离方案以解离样品中的细胞，同时最大限度地减少对细胞的压力，而实验室 B 没有这样做，那么 B 组数据中的细胞很可能会表达更多与压力相关的基因。

注意，批次的划分情况很多，并不是说彻底消除批次效应就是好的，例如我们在探究肿瘤与性别的关系的时候，那么性别就应当保留，而不是当成批次效应去除。

参考：
[1] https://github.com/Starlitnightly/single_cell_tutorial
[2] https://github.com/theislab/single-cell-best-practices

整合模型的分类

批次整合模型一般包含三个步骤：1.降维，2.消除批次效应，3.生成优质嵌入。降维一般使用PCA，PCA可以降低数据的信噪比，提取细胞的差异信息。消除批次效应是模型的关键，这些方法经过长期发展被归纳为4类：全局模型，线性嵌入模型，基于图的模型，深度学习模型。

• 全局模型：把批次效应建模为加法或者乘法，ComBat是全局模型
• 线性嵌入模型：利用奇异值分解（SVD）以及其变体来获得嵌入数据，然后在嵌入中跨批次查找相似细胞的局部邻域，并以局部自适应方式（非线性）校正批次效应。常见方法包括MNN，Seurat，Harmony，Scanorama等
• 基于图的模型：这类方法通常运行速度最快。使用最近邻域图来表示每个批次的数据。通过强制连接不同批次的细胞，然后修剪细胞类型组成的边，可以校正批次效应，常见方法为BBKNN，这类方法之所以称为基于图的方法并且速度快，是因为它们不会改变细胞的embedding，而是通过修剪邻域图的边获得新的拓扑结构（新的KNN），新的KNN用于后续的聚类和可视化
• 深度学习模型：大多数深度学习模型基于VAE架构，在latent space中校正批次

下面我们比较不同的批次整合方法，我们始终要明白：不存在适用于所有场景的最佳方法。对于一些简单的批次效应校正任务来说，Seurat及Harmony的表现始终不错，而对于一些复杂的批次效应校正任务来说，scVI、scGen、scANVI和Scanorama则表现地更为出色。我们可以使用scib包来评估数据整合的效果，scib包中包含了10个指标。

查看批次效应

我们首先导入数据：

import omicverse as ov
print(f"omiverse version: {
     ov.__version__}")
import scanpy as sc
print(f"scanpy version: {
     sc.__version__}")
import scib
print(f"scib version: {
     scib.__version__}")
import scvi
print(f"scvi version: {
     scvi.__version__}")
ov.ov_plot_set()

#读取数据
adata1=ov.read('./data/neurips2021_s1d3.h5ad')
adata1.obs['batch']='s1d3'
adata2=ov.read('./data/neurips2021_s2d1.h5ad')
adata2.obs['batch']='s2d1'
adata3=ov.read('./data/neurips2021_s3d7.h5ad')
adata3.obs['batch']='s3d7'

#合并数据集
adata=sc.concat([adata1,adata2,adata3], merge='same')
adata.obs['batch'].unique()
# 删除多余的metadata
del adata.obsm
adata.obs = adata.obs[['cell_type', 'batch']]
adata.write_h5ad("./data/neurips2021_scib.h5ad", compression='gzip')

然后进行预处理：

adata = sc.read