学习笔记总结于『生信技能树』马拉松课程

下文使用代码读取的文件“obj.Rdata”是正式开始注释的上一步——降维聚类分群的结果，详见 文章复现学习 | ROS（2）单细胞流程 一文

1.整理数据

library(Seurat)
library(celldex)
library(SingleR)
library(BiocParallel)
library(ggplot2)
load("obj.Rdata")
p1 <- UMAPPlot(sce.all);p1

得到umap图如图1所示

2.SingleR不同参考数据集的注释结果比较

singleR 里面，人类总共5个数据，放入循环

fs = dir("../supp/single_ref/",pattern = "ata")[-c(4,6)]# 4和6是小鼠的，所以将其去掉
fs
## [1] "ref_BlueprintEncode.RData"                 
## [2] "ref_DatabaseImmuneCellExpressionData.RData"
## [3] "ref_Human_all.RData"                       
## [4] "ref_MonacoImmuneData.Rdata"                
## [5] "ref_NovershternHematopoieticData.Rdata"

发现有多余的前缀和后缀，只需要留下中间的关键词如BlueprintEncode作为标题即可，其余的去掉

\\. 中的“.”号单独出现是正则表达式里面的语法，我们不想用它的语法，而是只需要它作为“.”放在这里，那前面就要加两个斜杠

library(stringr)
fs2 = str_remove_all(fs,"ref_|\\.Rdata|\\.RData")

用5个参考数据集循环做了一遍自动注释

scRNA = list()# 用于存储注释的结果
new.cluster.ids = list()# 每次重命名的变量名称
p = list()# 用于存储图片
for(i in 1:length(fs)){
  ref <- get(load(paste0("../supp/single_ref/",fs[[i]]))) # 虽然fs是向量，但此处使用两个中括号是for循环的一种保险措施
  scRNA[[i]] = sce.all
  pred.scRNA <- SingleR(test = scRNA[[i]] @assays$RNA@data, 
                        ref = ref,
                        labels = ref$label.main, 
                        clusters = scRNA[[i]] @active.ident)
  new.cluster.ids[[i]] <- pred.scRNA$pruned.labels
  names(new.cluster.ids[[i]]) <- levels(scRNA[[i]] )
  scRNA[[i]]  <- RenameIdents(scRNA[[i]] ,new.cluster.ids[[i]])
  p[[i]] = UMAPPlot(scRNA[[i]],label = T)+
    labs(title = fs2[[i]])
}
library(patchwork)
wrap_plots(c(p,list(p1)),nrow = 2)

R包 patchwork 的 wrap_plots() 函数，适用于拼多张图，操作对象可以是列表，那样就不用写p1+p2+p3+p4+p5+p6了

从图2也能看出DatabaseImmuneCellExpressionData的NA最多，而这部分也是最具争议的簇，五次注释得到了五个不同的结果