单细胞亚群鉴定

生信喵

1 细胞亚群鉴定
1.1 细胞亚群鉴定原理
      细胞亚群鉴定是进行单细胞转录组分析的最基础一步，是赋予细胞数据以生物学意义的关键过程。细胞亚群鉴定主要借助 marker 基因在各个细胞亚群的表达情况来判断细胞亚群所属细胞类型。
      利用 Seurat 分析之后，我们可以将细胞分成很多 cluster 类，也称为细胞亚群，我们可以找到每一亚群包含的 Marker 基因。

      
      t-SNE 绘图

      接下来就可以根据据亚群上调表达的基因进行细胞亚型鉴定。细胞类型的鉴定原理并不复杂，主要是根据细胞中表达基因的差别进行分类。例如在 Seurat 案例中，表 5 和下表中所示，如果是 MS4A1 基因高表达，则是 B 细胞。

      每个 cluster 标记基因与细胞类型

Cluster ID	Markers	Cell Type
0	IL7R, CCR7	Naive CD4+ T
1	CD14, LYZ	CD14+ Mono
2	IL7R, S100A4	Memory CD4+
3	MS4A1	B
4	CD8A	CD8+ T
5	FCGR3A, MS4A7	FCGR3A+ Mono
6	GNLY, NKG7	NK
7	FCER1A, CST3	DC
8	PPBP	Platelet

      
      marker 基因与细胞类型

1.2 人工鉴定与自动化鉴定
      亚细胞群细胞类型鉴定可以分为人工更坚定方法与自动化鉴定两种。人工方法就是我们找到每个亚群中的 Marker 基因，然后与已知的数据库进行比对，对细胞亚群进行鉴定，这个过程类似于物种鉴定。但人工的方法有几个问题：
      1、首先，完全依赖数据库信息，如果数据库中没有，则无法鉴定，这就会忽略掉一些中间态细胞；
      2、第二，鉴定准确性也与数据库高度相关，如果数据库有错误，则错误一直累积；
      3、第三、人工方法比较耗时；
      虽然人工方法有缺点，但人工方法结合背景信息，比单纯软件数据计算更准确。
自动化细胞注释的局限性
      自动化鉴定则是软件根据一定的聚类结果自动与 Marker 数据库进行比对，根据一定的算法规则进行判断。这一类软件替代的是细胞注释过程的机械工作，也就是通过 marker 基因分布判断细胞亚群所属细胞类型。所以，这一类软件本身也是围绕 marker 基因来展开分析的。而不同的软件在 marker 基因的决定上会有不同的策略：一类是基于基因在不同细胞类型的表达量差异去判定 marker 基因；一类是基于人工提供的 marker 基因信息做后续分析。虽然自动化软件比较方便，但软件也有一定局限性；
      第一：自动化鉴定同样以来数据库；
      第二：对于相似度很高的细胞类型鉴定准确性会降低；
      正是因为自动化细胞注释软件的局限性，它们更多的还是作为细胞注释中的辅助工具，人工注释的介入才能更好地提高细胞注释的准确性。


1.3 常用人工鉴定数据库
      目前已经有很多数据库可用于细胞类型的鉴定。例如

      1、Cell Marker（人、小鼠，）：http://bio-bigdata.hrbmu.edu.cn/CellMarker/
      2、panglaoDB（人、小鼠）：https://panglaodb.se/
      3、MCA（小鼠）：http://bis.zju.edu.cn/MCA/
      4、CancerSEA：http://biocc.hrbmu.edu.cn/CancerSEA/
      5、LnCeCell：http://www.bio-bigdata.net/LnCeCell/
      6、ColorCells：http://rna.sysu.edu.cn/colorcells/
      7、SC2disease：http://easybioai.com/sc2disease/
      8、HCA(Human Cell Atlas)：https://data.humancellatlas.org/
      9、SCP(Single Cell Portal)：https://singlecell.broadinstitute.org/single_cell
      10、SCEA(Single Cell Expression Atlas)：https://www.ebi.ac.uk/gxa/sc/home
      11、HCL(Human Cell Landscape)：https://db.cngb.org/HCL/
      12、CDCP：https://db.cngb.org/cdcp/
      13、scRNASeqDB：https://bioinfo.uth.edu/scrnaseqdb/
      14、Cell BLAST：https://cblast.gao-lab.org/
      15、signatureDB：https://lymphochip.nih.gov/signaturedb/
      16、VascularSingleCells：http://betsholtzlab.org/VascularSingleCells/database.html
      17、GRNdb：http://www.grndb.com/
      18、SpatialDB：https://www.spatialomics.org/SpatialDB/


3 自动化鉴定原理
      自动化细胞类型鉴定主要就是减少人工的大量操作，可以自动化完成工作。但无论人工还是自动化鉴定，本质上都是根据已知信息来推测未知，这都属于机器学习范畴。属于一种分类工作。目前自动化细胞鉴定原理不同，可以分为基于基因表达量矩阵方法与基于 marker 基因的方法。
3.1 基于基因表达量矩阵
      这一类软件的一大特点就是低人工占比，直接输入细胞基因表达矩阵，然后等待结果即可。
      这是因为软件知道数据库，自动完成比对判断。
      我们需要做的就是输入待鉴定细胞群，然后看输出结果，最后根据输出结果做些验证和矫正就可以完成细胞注释。
      首先，根据基因表达量矩阵寻找可靠的 marker 基因；
      第二，根据 marker 基因计算未鉴定细胞/细胞亚群与所有细胞类型的相关性；

      第三： 根据相关性进行细胞类型注释
      不过不同的软件在第二步计算时有所不同，例如singleR以基因差异倍数来寻找marker基因，以 spearman 相关系数来计算相关性。

      
      singleR 原理图

      celaref 以差异基因来寻找 marker 基因，通过 marker 基因在待鉴定细胞亚群中上调倍数排名来体现相关性。

      
      celaref 原理图

      CHETAH 首先以细胞类型构建层级聚类树，将分支特异性基因作为 marker 基因，以细胞为某一条分支的可能性来体现相关性。

      
      CHETAH 原理图

3.2 基于已知 marker 基因
      该方法利用机器学习方法，对 marker 基因进行机器学习分类器训练。使用最广泛的分类器训练方式微人工神经网络模型，例如软件 SuperCT 和 scANVI，它们将 marker 基因和细胞类型为节点，以 marker 基因在细胞类型在不同细胞类型的表达量和基因是否为相应细胞类型的 marker 基因两个依据作为节点间的权重评判，以此构建人工神经网络模型；最后，以该模型去判断待鉴定细胞的细胞类型。

      
      SuperCT 原理图


4 利用 singleR 进行细胞鉴定
4.1 SingleR 简介
      Single R 可以自动完成单细胞类型的鉴定。软件主要分为三个步骤：
      首先、输入没有注释的单细胞转录组数据；
      第二、基于参考库中的 marker 基因对输入的细胞进行鉴定。
      第三、输出鉴定结果。
      文档：http://bioconductor.org/books/release/SingleRBook/
      https://www.bioconductor.org/pac ... st/doc/SingleR.html
      celldex 文档：https://www.bioconductor.org/pac ... /doc/userguide.html

1. library(SingleR)
   
2. # BiocManager::install("SingleR",lib = "/opt/R/4.2.1/lib/R/library",destdir = '/home/xhs/Rpack/download')

复制代码

数据库下载：
      singleR 进行细胞分类需要依赖数据库，数据库可以从网上下载，也可以自行构建，目前有七个参考数据库。

      目前 singler 已经将这些数据库打包到 celldex 包中，可以通过安装 celldex 包进行安装和下载。下载的数据为 SummarizedExperiment 类对象，可以通过添加Ensemble=TRUE，下载 Ensemble 类型基因 ID，通过添加 cell.ont = c("all", "nonna", "none")添加 Cell Ontology terms 信息。

      中国大陆地区的部分网络其实是没办法联通singleR数据库中心，所以文件我尝试了好几次都是下载失败的，于是我切去了香港服务器下载了数据。我将数据上传到百度云吧， 获取以下7个数据库请在公众号"生信喵实验柴"后台回复“20220728”。

生信喵

1. library(celldex)
   
2. ref <- HumanPrimaryCellAtlasData()
   
3. save(ref,file = 'HumanPrimaryCellAtlasData.Rdata')
   
4. ref <- BlueprintEncodeData()
   
5. save(ref,file = 'BlueprintEncodeData.Rdata')
   
6. ref <- MouseRNAseqData()
   
7. save(ref,file = 'MouseRNAseqData.Rdata')
   
8. ref <- ImmGenData()
   
9. save(ref,file = 'ImmGenData.Rdata')
   
10. ref <- DatabaseImmuneCellExpressionData()
    
11. save(ref,file = 'DatabaseImmuneCellExpressionData.Rdata')
    
12. ref <- NovershternHematopoieticData()
    
13. save(ref,file = 'NovershternHematopoieticData.Rdata')
    
14. ref <- MonacoImmuneData()
    
15. save(ref,file = 'MonacoImmuneData.Rdata')

复制代码

可以通过常用函数探索 SummarizedExperiment 类。

1. class(ref)
   
2. dim(ref)
   
3. str(ref)
   
4. assayNames(ref)
   
5. assays(ref)
   
6. rowData(ref)
   
7. colData(ref)
   
8. rownames(ref@assays@data@listData$logcounts)colData(x, ...)

复制代码

4.2 SingleR 使用

1. #读入Seurat数据
   
2. pbmc <- readRDS("pbmc3k_final.rds")
   
3. 
   
4. #获得表达矩阵
   
5. query <- pbmc@assays$RNA@data
   
6. class(query)
   
7. 
   
8. #SingleR注释
   
9. pred.pbmc <- SingleR(test=query, assay.type.test=1, ref=ref,clusters = pbmc@active.ident,labels=ref$label.main)
   
10. head(pred.pbmc)
    
11. # ref$label.fine
    
12. 
    
13. #查看注释结果
    
14. View(as.data.frame(pred.pbmc))
    
15. pred.pbmc$pruned.labels
    
16. #统计
    
17. table(pred.pbmc$pruned.labels)
    
18. 
    
19. #注释结果验证
    
20. library(pheatmap)
    
21. plotScoreHeatmap(pred.pbmc)
    
22. plotDeltaDistribution(pred.pbmc, ncol = 3)
    
23. 
    
24. 
    
25. new.cluster.ids <- pred.pbmc$labels
    
26. names(new.cluster.ids)
    
27. names(new.cluster.ids) <- levels(pbmc)
    
28. pbmc <- RenameIdents(pbmc,new.cluster.ids)
    
29. levels(pbmc)
    
30. #重新绘图
    
31. DimPlot(pbmc, reduction = "umap", label = TRUE, pt.size = 0.5) + NoLegend()
    
32. DimPlot(pbmc, reduction = "tsne", label = TRUE, pt.size = 0.5) + NoLegend()

复制代码

singleR 注释完成之后细胞进行重命名

生信喵

帖子中的资源，移步关注公众号："音云商城"，回复关键字20220728即可获取。