作者，追风少年i

开头先放一张marker表，供大家参考

这一篇内容很简单

接上一篇，上一篇文章单细胞 & 空间整合去批次方法比较介绍了以下几种方法整合去批次的代码

• CCA
• merge
• SCT
• merge & SCT
• merge加harmony
• SCT && harmony

关于上述方法我只强调一点，就是ScaleData的时候vars.to.regress的作用，这个要引起大家的重视。

这一篇我们要进行补充，因为上述的方法仍然存在一定的局限性。面对大数据集，几十上百万的细胞量，上述方法无能为力，R语言的原因，处理起来很慢。

• rpca的方法：

ifnb.list <- SplitObject(ifnb, split.by = "stim") # normalize and identify variable features for each dataset independently ifnb.list <- lapply(X = ifnb.list, FUN = function(x) { x <- NormalizeData(x) x <- FindVariableFeatures(x, selection.method = "vst", nfeatures = 2000) }) # select features that are repeatedly variable across datasets for integration run PCA on each # dataset using these features features <- SelectIntegrationFeatures(object.list = ifnb.list) ifnb.list <- lapply(X = ifnb.list, FUN = function(x) { x <- ScaleData(x, features = features, verbose = FALSE) x <- RunPCA(x, features = features, verbose = FALSE) }) ####rpca immune.anchors <- FindIntegrationAnchors(object.list = ifnb.list, anchor.features = features, reduction = "rpca" # this command creates an 'integrated' data assay immune.combined <- IntegrateData(anchorset = immune.anchors)

这里大家要注意rpca的方法，关于rpca，大家可以参考文章10X单细胞数据整合分析Seurat之rpca（large data，细胞量超过20万），原理我都已经说的很清楚了。

• rpca + reference

bm280k.list <- SplitObject(bm280k, split.by = "orig.ident") bm280k.list <- lapply(X = bm280k.list, FUN = function(x) { x <- NormalizeData(x, verbose = FALSE) x <- FindVariableFeatures(x, verbose = FALSE) }) features <- SelectIntegrationFeatures(object.list = bm280k.list) bm280k.list <- lapply(X = bm280k.list, FUN = function(x) { x <- ScaleData(x, features = features, verbose = FALSE) x <- RunPCA(x, features = features, verbose = FALSE) }) anchors <- FindIntegrationAnchors(object.list = bm280k.list, reference = c(1, 2), reduction = "rpca", dims = 1:50) bm280k.integrated <- IntegrateData(anchorset = anchors, dims = 1:50)

注意这里的reference在指定的时候一般指定为control样本，或者预后好的样本。

• bbknn（scanpy，python），这是很多高分文章采用的大细胞量整合方法

import scanpy as sc import pandas as pd import seaborn as sns adata_ref = sc.datasets.pbmc3k_processed() # this is an earlier version of the dataset from the pbmc3k tutorial adata = sc.datasets.pbmc68k_reduced() ####注意这里的不同，高变基因取交集 var_names = adata_ref.var_names.intersection(adata.var_names) adata_ref = adata_ref[:, var_names] adata = adata[:, var_names] ####ref数据预处理 sc.pp.pca(adata_ref) sc.pp.neighbors(adata_ref) sc.tl.umap(adata_ref) sc.tl.ingest(adata, adata_ref, obs='louvain') adata.uns['louvain_colors'] = adata_ref.uns['louvain_colors'] # fix colors adata_concat = adata_ref.concatenate(adata, batch_categories=['ref', 'new']) adata_concat.obs.louvain = adata_concat.obs.louvain.astype('category') adata_concat.obs.louvain.cat.reorder_categories(adata_ref.obs.louvain.cat.categories, inplace=True) # fix category ordering adata_concat.uns['louvain_colors'] = adata_ref.uns['louvain_colors'] # fix category colors ####采用bbknn矫正批次 sc.tl.pca(adata_concat) sc.external.pp.bbknn(adata_concat, batch_key='batch') ####然后进行标准过程即可，多样本就是不停的迭代这个过程。

python版本大细胞量harmony矫正。

import scanpy as sc import scanpy.external as sce ####大家读取自己的数据 adata = sc.datasets.pbmc3k() sc.pp.recipe_zheng17(adata) sc.tl.pca(adata) ####矫正批次 sce.pp.harmony_integrate(adata, 'sample') sc.pp.neighbors(adata) sc.tl.umap(adata)

其中sce.pp.harmony_integrate(adata, 'sample')这句其实就是下面这个调用下面的语句

import harmonypy harmony_out = harmonypy.run_harmony(adata.obsm["X_pca"], adata.obs, 'batch') adata.obsm[adjusted_basis] = harmony_out.Z_corr.T ###obsm就是降维的数据

大样本量还是推荐scanpy做分析，bbknn或者harmony矫正。

生活很好，有你更好