本期内容同步发布于
单细胞天地微信公众号。
回顾——细胞通讯网络构建
前面的帖子中我们已经成功地进行了细胞间通讯网络的构建，总的来看借助下面简易的分析流程即可完成：
library(CellChat)
library(Seurat)
library(SeuratData)
library(future)
#load data
data("pbmc3k.final")
#create CellChat object
cellchat <- createCellChat(object = GetAssayData(object = pbmc3k.final, slot = 'data'),
meta = pbmc3k.final@meta.data,
group.by = 'seurat_annotations')
## The cell groups used for CellChat analysis are Naive CD4 T Memory CD4 T CD14+ Mono B CD8 T FCGR3A+ Mono NK DC Platelet
#add CellChat database
cellchat@DB <- CellChatDB.human #human
#cell communication
plan(strategy = 'multiprocess', workers = 4)
cellchat <- subsetData(cellchat) #save time and memory
cellchat <- identifyOverExpressedGenes(cellchat)
cellchat <- identifyOverExpressedInteractions(cellchat)
cellchat <- computeCommunProb(cellchat)
cellchat <- filterCommunication(cellchat, min.cells = 10)
#pathway level
cellchat <- computeCommunProbPathway(cellchat)
cellchat <- aggregateNet(cellchat)
今天继续介绍：
（1）我们如何来理解这些分析结果；
（2）我们如何将这些结果有效地展示出来，也就是可视化。
CellChat分析结果可视化
目前CellChat提供了多种图形来对结果进行可视化，包括但不限于：
Hierarchical plot （层次聚类图）
Chord diagram （和弦图）
Circle plot （环状图）
Bubble plot （气泡图）
此外结果的可视化也分为在单个受体-配体层次与代谢通路（多个受体-配体整合）层次，上游分析的相关结果分别存储在cellchat@net和cellchat@netP中。
单个受体-配体层次
这个部分的可视化主要是统计不同细胞类群之间的受体-配体对数和通讯强度（strength）。
library(dplyr)
groupSize <- table(cellchat@idents) %>% as.numeric()
par(mfrow = c(1, 2), xpd = TRUE)
netVisual_circle(cellchat@net$count,
vertex.weight = groupSize,
weight.scale = T,
label.edge = F,
title.name = "Number of interactions")
netVisual_circle(cellchat@net$weight,
vertex.weight = groupSize,
weight.scale = T,
label.edge = F,
title.name = "Interaction weights/strength")
首先简单介绍一下par()函数，这是R基础绘图系统中一个比较常见的函数，在CellChat中主要是为了排列多个图形，例如这里的mfrow参数的形式就为c(nr, nc)，表示将后续绘制的图形按照nr行nc列排布，而mfrow表示先横向后纵向，对应的也有mfcol表示先纵向后横向。
再来详细介绍一下netVisual_circle()函数，这个环状图中最关键的两个概念就是edge和vertex，前者指的是不同细胞类群之间相互关联的“连线”，后者指的是位于环状图边界上的各个细胞群的节点的大小，显然在上面的例子中我们分别使用了不同细胞类群之间的受体-配体对数和通讯强度作为“连线”强度的度量，而不同细胞类群中的细胞数量作为节点大小的度量。我们会对这些信息scale一下，防止差异过大影响可视化的结果。
为了让信息不显得那么冗余，我们也可以单独可视化以每个细胞类群为信号起点的环状图，实现的思路很简单，以cellchat@net$weight为例，这个矩阵中的每一行代表信号发出的细胞，每一列代表接收信号的细胞，所以我们可以通过分别可视化每一行的数据来研究以每个细胞群为信号起点的细胞间通讯：
groupSize <- table(cellchat@idents) %>% as.numeric()
par(mfrow = c(3, 3), xpd = TRUE)
for(i in 1:nrow(cellchat@net$weight)){
mat <- matrix(0,
nrow = nrow(cellchat@net$weight),
ncol = ncol(cellchat@net$weight),
dimnames = dimnames(cellchat@net$weight))
mat[i, ] <- cellchat@net$weight[i, ]
netVisual_circle(mat,
vertex.weight = groupSize,
weight.scale = T,
edge.weight.max = max(cellchat@net$weight),
title.name = rownames(mat)[i])
}
值得注意的是，这里多了一个参数：edge.weight.max = max(cellchat@net$weight)，这主要是为了让不同的图之间具有可比性，因为都是按照总体的最大值进行scale的。
信号通路层次可视化
单个受体-配体层次的分析只能让我们知道哪些细胞类群之间可能存在通讯，至于是何种类型的通讯、会有什么生物学意义却比较难知道，这个时候如果能够将这种通讯放到信号通路的层次上就能较好地解决这个问题。
首先来查看一下我们分析的结果中有哪些信号通路：
cellchat@netP$pathways
## [1] "MHC-I" "CD99" "MIF" "MHC-II" "GALECTIN" "ITGB2"
## [7] "APP" "ANNEXIN" "SELPLG" "LCK" "ICAM" "CD40"
## [13] "ESAM" "IL16"
Hierarchical plot
Hierarchical plot一个比较突出的优点在于能够比较清晰有层次地展示细胞类群之间的通讯，我们可以通过指定vertex.receiver来自定义信号接收细胞群，这个参数可以是细胞类群的名称，例如CD8 T，也可以是索引，例如c(1, 2, 3)。
cellchat@meta$seurat_annotations %>% head()
## [1] Memory CD4 T B Memory CD4 T CD14+ Mono NK
## [6] Memory CD4 T
## 9 Levels: Naive CD4 T Memory CD4 T CD14+ Mono B CD8 T FCGR3A+ Mono NK ... Platelet
vertex.receiver <- 1:3 #Naive CD4 T, Memory CD4 T, CD14+ Mono
cellchat@netP$pathways
## [1] "MHC-I" "CD99" "MIF" "MHC-II" "GALECTIN" "ITGB2"
## [7] "APP" "ANNEXIN" "SELPLG" "LCK" "ICAM" "CD40"
## [13] "ESAM" "IL16"
pathways.show <- 'MHC-I'
netVisual_aggregate(cellchat,
signaling = pathways.show,
vertex.receiver = vertex.receiver,
layout = 'hierarchy')
这里我们可视化了Naive CD4 T, Memory CD4 T, CD14+ Mono作为信号接收细胞群的MHC-I信号通路，详细来看Hierarchical plot：
左侧部分展示的是以Naive CD4 T, Memory CD4 T, CD14+ Mono作为信号接收细胞群的MHC-I信号通路，其中左半部分是自分泌相关信号，这里的所谓自分泌信号就是指这几类细胞类群自己释放的信号作用于自己这几类细胞，相应的右半部分就是展示的旁分泌信号，也就是其它类的细胞所释放的信号作用于这几类细胞。
右侧部分展示是以除了Naive CD4 T, Memory CD4 T, CD14+ Mono以外的其它细胞群作为信号接收细胞群的MHC-I信号通路，只不过此时左右半部分的信息反过来了：左半部分是旁分泌信号，右半部分是自分泌信号。
Chord diagram
netVisual_aggregate(cellchat,
signaling = pathways.show,
vertex.receiver = vertex.receiver,
layout = 'chord')
改一下layout就可以实现Chord diagram，不过这个图的信号很杂，将所有和MHC-I的信号糅合在了一起，从这一点上来看效果不如Hierarchical plot。不过好在它还是可以和环状图一样，接受sources.use和targets.use参数，大大简化图形信息，帮助我们更精准地定位到我们想要的内容：
netVisual_aggregate(cellchat,
signaling = pathways.show,
sources.use = 1:2, #naive CD4 T, Memory CD4 T
targets.use = 5:6, #CD8 T, FCGR3A+ Mono
layout = 'chord')
还有一个玩法在于，有的时候我们会在细胞类型注释的时候对某个细胞类群进行细分，但是我们在可视化细胞通讯的时候可能会想将这一大类细胞单独“聚”在一起，这个时候Chord diagram就来了：
#combine Naive CD4 T, Memory CD4 T and CD8 T
group.celltype = c(rep('T', 2), 'CD14+ Mono', 'B', 'T', 'FCGR3A+ Mono', 'NK', 'DC', 'Platelet')
names(group.celltype) = levels(cellchat@meta$seurat_annotations)
netVisual_aggregate(cellchat,
signaling = pathways.show,
group = group.celltype,
layout = 'chord')
可以看到这三类细胞已经聚在一起了。
Circle plot
这个部分和前面一样，只不过现在是在信号通路的层次上进行的罢了：
netVisual_aggregate(cellchat,
signaling = pathways.show,
layout = 'circle')
Bubble plot
netVisual_bubble(cellchat,
sources.use = 1:3,
targets.use = 4:7,
signaling = 'MHC-I')
## Comparing communications on a single object
气泡图的横轴标明了信号起始和接收的细胞类群（以 -> 标明）。同时纵轴标上了属于这个MHC-I信号通路的具体细胞间相互作用（受体-配体对），你可以通过pairLR.use传入字符串向量来选择性可视化哪些具体的相互作用。还有很多其它实用的参数，例如使用remove.isolate = TRUE去掉空白行和空白列，让整个图形仅仅保留有信号的部分。
除了这几类主要的可视化以外，我们还可以使用热图：
netVisual_heatmap(cellchat, signaling = pathways.show, color.heatmap = 'Reds')
## Do heatmap based on a single object
左侧已经标明了信号发出的细胞群，上方和右方的柱状图分别表示每一列和每一行数值的总和。
无论怎么变，最终总还是在对cellchat@netP$prob这个array在可视化。此外这些函数的帮助文档都比较通俗易懂，大家可以多用?来获取帮助信息。
基因表达量可视化
我们在Seurat对象中很容易进行基因表达的可视化，CellChat同样为我们提供了可能性，我们只需要使用plotGeneExpression()函数就可以了，我们可以指定我们想可视化的信号通路（此时会自动可视化该信号通路内富集的重要基因），也可以指定一些具体的基因，还可以指定可视化图形的类型："violin"，"dot"，"bar"。
plotGeneExpression(cellchat, signaling = 'MHC-I')
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
可视化的简单分享就到这里啦~大家赶紧行动起来吧！
