「单细胞转录组系列」使用scCATCH进行聚类结果自动化注释

摘自：xuzhougeng

https://www.jianshu.com/p/cf7a7341b0b6

目前该软件只支持Mouse和Human，不支持其他物种，因此不是这两个物种的小伙伴可以不用看了。

scCATCH全称是single cell Cluster-based Annotation Toolkit for Cellular Heterogeneity，是一个用于实现单细胞转录组聚类结果进行注释的工具。软件核心函数是和scCATCH,findmarkergenes则是辅助用于寻找标记。

目前该软件托管在GitHub上，因此需要使用devtools::install_github()进行安装，或者直接从GitHub上下载代码使用install.packages()进行安装。我在测试这个R包发现它直接使用as.data.frame转换Seurat的稀疏矩阵，而R在转换非常大的稀疏矩阵时会报错，因此我fork了一份代码，并做了相应的修改，希望原作者能够合并我的PR。目前原作者已经修复了该问题

# 原版devtools::install_github("ZJUFanLab/scCATCH")# 修改版本#devtools::install_github("xuzhougeng/scCATCH")

我们下载10X genomics官网上的一份数据，地址为http://cf.10xgenomics.com/samples/cell-exp/3.0.2/5k_pbmc_v3/5k_pbmc_v3_filtered_feature_bc_matrix.h5

然后以https://satijalab.org/seurat/v3.1/pbmc3k_tutorial.html的步骤进行数据预处理

library(Seurat)library(scCATCH)h5_file <- "F:/5k_pbmc_v3_filtered_feature_bc_matrix.h5"# Load the PBMC dataset#pbmc.data <- Read10X(data.dir = "../data/pbmc3k/filtered_gene_bc_matrices/hg19/")pbmc.data <- Read10X_h5(h5_file)# Initialize the Seurat object with the raw (non-normalized data).pbmc <- CreateSeuratObject(counts = pbmc.data, project = "pbmc3k", min.cells = 3, min.features = 200)pbmcpbmc[["percent.mt"]] <- PercentageFeatureSet(pbmc, pattern = "^MT-")VlnPlot(pbmc, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3)pbmc <- subset(pbmc, subset = nFeature_RNA > 500 & nFeature_RNA < 5000 & percent.mt < 20)pbmc <- NormalizeData(pbmc, normalization.method = "LogNormalize", scale.factor = 10000)pbmc <- FindVariableFeatures(pbmc, selection.method = "vst", nfeatures = 2000)all.genes <- rownames(pbmc)pbmc <- ScaleData(pbmc, features = all.genes)pbmc <- RunPCA(pbmc, features = VariableFeatures(object = pbmc))ElbowPlot(pbmc)pbmc <- FindNeighbors(pbmc, dims = 1:10)pbmc <- FindClusters(pbmc, resolution = 0.2)pbmc <- RunUMAP(pbmc, dims = 1:10)DimPlot(pbmc, label = TRUE)

以UMAP对聚类结果进行可视化展示

reslution=0.2

接下来，使用findmarkergenes寻找每个cluster的差异基因。这一步的运行时间比较长，因为每个cluster都需要和其他的所有cluster按个比较，然后确定出当前cluster的特异基因。

clu_markers <- findmarkergenes(pbmc,                               species = "Human",                               cluster = 'All',                               match_CellMatch = FALSE,                               cancer = NULL,                               tissue = NULL,                               cell_min_pct = 0.25,                               logfc = 0.25,                               pvalue = 0.05)

最终得到的clu_markers是一个list，包括一个新的表达量矩阵（基于NCBI最新的Gene Symbols，并移除重复和不匹配的基因) 以及一个包括每个cluster的所有潜在标记基因。

之后使用scCATCH根据标记基因对聚类进行注释

clu_ann <- scCATCH(clu_markers$clu_markers,                   species = "Human",                   cancer = NULL,                   tissue = "Blood")

输出的clu_ann能用来修改原来的注释信息，代码如下

new.cluster.ids <- clu_ann$cell_typenames(new.cluster.ids) <- clu_ann$clusterpbmc <- RenameIdents(pbmc, new.cluster.ids)DimPlot(pbmc, reduction = "umap", label = TRUE, pt.size = 0.5) + NoLegend()

注释结果如下

注释结果

根据Seurat的教程，这个注释可能有一些问题，比如说T细胞，NK细胞没有正确注释。

注释

我觉得这可能是我选择的组织不对，所以我重新以"Blood", "Bone marrow"两个数据集进行注释

clu_ann <- scCATCH(clu_markers$clu_markers,                   species = "Human",                   cancer = NULL,                   tissue = c("Blood", "Bone marrow"))

虽然结果还是有一些问题，但是和预期结果比较接近了。

注释结果2

scCATCH文章里提到它的注释准确度其实非常高的，因此我这一次找到了文章里的测试脚本https://github.com/ZJUFanLab/scCATCH_performance_comparison/blob/master/scCATCH/scCATCH.R

根据脚本设置了参数，也就是选择'Blood','Peripheral blood','Bone marrow'作为组织类型

clu_ann <- scCATCH(clu_markers$clu_markers,                   species = "Human",                   cancer = NULL,                   tissue =  c('Blood','Peripheral blood','Bone marrow'))

这一次大部分的细胞都被正确注释

注释结果3

因此想要用好这个软件，需要非常注意组织类型的选择。如果你发现注释和预期不符，你可能没有漏了一些组织。

更多和scCATCH有关的参数阅读，请移步https://github.com/ZJUFanLab/scCATCH

引用信息:

Shao et al., scCATCH:Automatic Annotation on Cell Types of Clusters from Single-Cell RNA Sequencing Data, iScience, Volume 23, Issue 3, 27 March 2020. doi: 10.1016/j.isci.2020.100882. PMID:32062421