细胞亚群的生物学命名

课程笔记

粉丝：有单细胞线上课程吗？

小编：什么？ 我们的单细胞转录组分析线上课程已经上线好久了，你们竟然都不知道吗，每篇推文后面的课程推荐没人看的吗，小编已哭晕在厕所

好了，戏演完了，下面郑重介绍下我们的单细胞线上课程：（详情戳下方链接） 

这个课程笔记栏目记录了学员们学习单细胞转录组课程的学习笔记

希望大家能有所收获！

作者 | 单细胞天地小编 刘小泽

课程链接在：http://jm.grazy.cn/index/mulitcourse/detail.html?cid=55
第二单元第8讲：细胞亚群的生物学命名

上一次我们好不容易得到了这个1.9G的RData，也就是作者自己做出来的Seurat对象，那么我们要怎么利用它去进一步探索呢？

加载上次运行的结果

首先还是三步走

1rm(list = ls()) 
2options(warn=-1)
3suppressMessages(library(Seurat))

然后加载进来之前保存的5.1G PBMC对象

1start_time <- Sys.time()
2load('./patient1.PBMC.output.Rdata')
3end_time <- Sys.time()
4end_time - start_time
5# Time difference of 12.6741 secs

重温上次的聚类结果

 1# 使用Seurat 2.3.4版本
 2colP<-c('green4', 
 3        'pink', 
 4        '#FF7F00', 
 5        'orchid', 
 6        '#99c9fb', 
 7        'dodgerblue2', 
 8        'grey30', 
 9        'yellow', 
10        'grey60', 
11        'grey', 
12        'red', 
13        '#FB9A99', 
14        'black'
15)
16TSNEPlot(PBMC, 
17         colors.use =  colP,
18         do.label = T)

开始新的探索

看看作者整合的数据对批次的处理

也就是把四个时间点映射到上面的tsne坐标中，并且理论上应该是：每群细胞都覆盖到四个时间点

1TSNEPlot(PBMC,group.by = "TimePoints")

再用table对比一下

1> table([e-mail]$TimePoints,PBMC@ident)
2
3                  0   1   2   3   4   5   6   7   8   9  10  11  12
4  PBMC_ARD614   665 726 572 559 420 302 457 313 283 123  17  11  68
5  PBMC_EarlyD27  43 173 245  85 120 110 543  59  91  29   7   3  84
6  PBMC_Pre      369 527 197  93 146 393   4  76  17  48  25 187   0
7  PBMC_RespD376 800 433 555 677 636 516 119 324 204 200 170  11  39

可视化一些marker基因

这些marker基因也是来源于文章的Supp Fig.7，基于他们对免疫知识的了解

 1allGenes = row.names([e-mail])
 2markerGenes <- c(
 3  "CD3D",
 4  "CD3E",
 5  "TRAC",
 6  "IL7R",
 7  "GZMA",
 8  "FCGR3A",
 9  "CD14",
10  "MS4A1",
11  "FCER1A" 
12)
13# 判断这些marker是不是存在表达矩阵中
14> markerGenes %in% allGenes
15[1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
16# 选择性保存pdf格式
17# pdf('patient1_pBMC_marker_FeaturePlot.pdf', width=10, height=15)
18FeaturePlot(object = PBMC, 
19            features.plot =markerGenes, 
20            cols.use = c("grey", "blue"), 
21            reduction.use = "tsne")

利用上面marker基因在不同细胞群的特殊表达，其实就能得到每个群的细胞命名，不过这些marker基因的获得以及和细胞名称的对应，是需要一段时间的研究才能得到的。我们这里只是展示如何操作

赋予每个cluster细胞类型

根据我们得到的cluster和原文的命名

就可以对应得到一个细胞名称列表：

 1cat >celltype-patient1-PBMC.txt
 20 "B cells"
 31 "CD4+ T cells"
 42 "Naive memory T cells"
 53 "Classical monicytes"
 64 "CD8+ effector T cells"
 75 "NK cells"
 86 "rm1"
 97 "Non-classical monocytes"
108 "Dendritic cells"
119 "rm2"
1210  "CD8+ cytotoxic T cells"
1311  "Myeloid cells"
1412  "rm3"

假如我们是先有了这个列表，核心就是要将分群的clsuter数字与细胞名称和颜色对应起来

利用Seurat V3

版本3的优势是可以直接提取出seurat对象的cluster数字，并且提供了RenameIdents函数，直接进行转换

1a=read.table('celltype-patient1-PBMC.txt')
2new.cluster.ids <- as.character(a[,2])
3names(new.cluster.ids) <- levels(PBMC_V3)
4PBMC_V3 <- RenameIdents(PBMC_V3, new.cluster.ids)
5
6DimPlot(PBMC_V3, reduction = "tsne", label = TRUE, pt.size = 0.5, cols = colP) + NoLegend()

利用Seurat V2

版本2就需要自己去对应：首先要了解它的分群信息存储在PBMC@ident中，然后要将全部12874个细胞的分群编号与celltype-patient1-PBMC.txt表中的第一列编号对应（只有这样，才能和表中的第二列对应上）

用到对应关系时，首先思考match能不能做到；如果要做，需要准备什么；对应关系搞清楚

 1# 先与表中第一列对应
 2match(as.numeric(as.character(PBMC@ident)),a[,1])
 3
 4# 第一点需要注意的是：为什么先用as.character后用as.numeric，而不是直接用as.numeric？
 5# 原因就是PBMC@ident存储的是因子型变量，直接取只会得到它们的位置信息，而不是真实的分群信息
 6> head(as.numeric(PBMC@ident))
 7[1] 2 1 2 2 2 6
 8> head(as.character(PBMC@ident))
 9[1] "1" "0" "1" "1" "1" "5"
10> head(as.numeric(as.character(PBMC@ident)))
11[1] 1 0 1 1 1 5
12
13# 第二点需要注意的是：match函数的规则是，A要在B中找到对应位置，那么就是 match(A,B)

接下来就可以得到对应的第二列，也就是细胞名称

 1labels=a[match(as.numeric(as.character(PBMC@ident)),a[,1]),2]
 2
 3# 检查一下，主要看数量的对应
 4> table(labels)
 5labels
 6                B cells            CD4+ T cells  CD8+ cytotoxic T cells 
 7                   1877                    1859                     219 
 8  CD8+ effector T cells     Classical monicytes         Dendritic cells 
 9                   1322                    1414                     595 
10          Myeloid cells                NK cells    Naive memory T cells 
11                    212                    1321                    1569 
12Non-classical monocytes                     rm1                     rm2 
13                    772                    1123                     400 
14                    rm3 
15                    191 
16> table(PBMC@ident)
17
18   0    1    2    3    4    5    6    7    8    9   10   11   12 
191877 1859 1569 1414 1322 1321 1123  772  595  400  219  212  191 

添加到metadata中，方便后面使用

[e-mail]$labels=labels
2
3TSNEPlot(PBMC, group.by = 'labels',
4         colors.use =  colP,
5         do.label = T)

但很明显，颜色标记和原文不同。因为之前只是对应了分群的编号和细胞名称，此外还需要修改颜色的顺序

1# 修改颜色顺序，思想就是：将现在的细胞名与表中的细胞名对应一下，然后这个顺序就是颜色出现的顺序
2colP=colP[match(levels(as.factor(labels)),a[,2])]
3TSNEPlot(PBMC, group.by = 'labels',
4         colors.use =  colP,
5         do.label = T)

再按时间拆分分群结果

做出文章的这张图

首先得到我们的四个时间点

1> TimePoints = [e-mail]$TimePoints
2> table(TimePoints)
3TimePoints
4  PBMC_ARD614 PBMC_EarlyD27      PBMC_Pre PBMC_RespD376 
5         4516          1592          2082          4684 

然后用一个函数`SubsetData`

举一个例子，绘制Pre时期的分群图

1# SubsetData函数其实也是利用TRUE/FALSE取子集
2PBMC_Pre  = SubsetData(PBMC,TimePoints =='PBMC_Pre')
3TSNEPlot(PBMC_Pre, 
4         colors.use = c('green4', 'pink', '#FF7F00', 'orchid', '#99c9fb', 'dodgerblue2', 'grey30', 'yellow', 'grey60', 'grey', 'red', '#FB9A99', 'black'),
5         do.label = F)
6# ggsave('PBMC_Pre_tSNE.pdf')

“

如果你对单细胞转录组研究感兴趣，但又不知道如何入门，也许你可以关注一下下面的课程

”

---