之前介绍过使用cellphoneDB 进行细胞通讯分析scRNA分析 | 解决可能的报错,从0开始教你完成细胞通讯分析-cellphoneDB,可能会遇到一些报错。这次介绍另一款细胞通讯分析的常见方法CellChat 。CellChat是一款R包,使用更容易且可视化结果也非常不错。
一 数据输入,处理
首选安装CellChat 包,然后继续使用之前的sce2数据进行展示
#devtools::install_github("sqjin/CellChat")
library(CellChat)
library(tidyverse)
library(Seurat)
options(stringsAsFactors = FALSE)
#提取表达矩阵和细胞分类信息
#scRNA <- readRDS(url("https://zenodo.org/record/3531889/files/seuratObj.rds"))
load("sce.anno.RData")
head(sce2@meta.data)
可以使用矩阵数据、Seurat或SingleCellExperiment 对象中创建CellChat对象。如果是Seurat或SingleCellExperiment 对象,meta信息会默认使用meta.data(当然也可以指定),通过 group.by 定义分组。
这里直接使用seurat的对象sce2进行创建CellChat对象
cellchat <- createCellChat(object = sce2,
meta = sce2@meta.data,
group.by = "celltype")
cellchat
Create a CellChat object from a Seurat object
The `data` slot in the default assay is used. The default assay is RNA
Set cell identities for the new CellChat object
The cell groups used for CellChat analysis are Epi Myeloid Fibroblast T Endo un
根据分析数据的物种,可选CellChatDB.human, 或者 CellChatDB.mouse 。通过showDatabaseCategory函数可以查看该数据库的情况
CellChatDB <- CellChatDB.human
showDatabaseCategory(CellChatDB)
# use a subset of CellChatDB for cell-cell communication analysis
CellChatDB.use <- subsetDB(CellChatDB, search = "Secreted Signaling")
# set the used database in the object
cellchat@DB <- CellChatDB.use
人的数据包括61.8%的旁分泌/自分泌信号互作、21.7%的细胞外基质(ECM)受体互作和16.5%的细胞-细胞通讯互作。
注1:如果想用全部的用于cellchat分析,不进行subsetDB,直接指定cellchat@DB <- CellChatDB 即可。
注2:如果你有关心的配受体对 且 不在该数据库中,也可以自行添加上。大概步骤就是下载对应的csv(数据库),在对应的列上添加上你的配受体对信息,保存后重新读取新的csv即可,详细见https://htmlpreview.github.io/?https://github.com/sqjin/CellChat/blob/master/tutorial/Update-CellChatDB.html。
可以使用subsetData选择进行cellchat的子集,注意使用全集的话也要subsetData一下。
# This step is necessary even if using the whole database
cellchat <- subsetData(cellchat)
# do parallel ,根据配置设置
future::plan("multiprocess", workers = 1)
#识别过表达基因
cellchat <- identifyOverExpressedGenes(cellchat)
#识别过表达配体受体对
cellchat <- identifyOverExpressedInteractions(cellchat)
#project gene expression data onto PPI (Optional: when running it, USER should set `raw.use = FALSE` in the function `computeCommunProb()` in order to use the projected data)
cellchat <- projectData(cellchat, PPI.human)
cellchat@data.project[1:4,1:4]
# K16733_AAACATACTCGTTT-1 K16733_AAAGCAGAACGTTG-1 K16733_AAAGCAGACTGAGT-1 K16733_AAAGGCCTGCTCCT-1
#TNFRSF18 0.0000000 0.000000000 0.000000e+00 0.000000000
#TNFRSF4 0.0000000 0.000000000 0.000000e+00 0.177007949
#TNFRSF14 0.9381181 0.001323893 1.852638e-05 0.002543824
#TNFRSF25 0.0000000 0.000000000 0.000000e+00 0.000000000
projectData函数将配体受体对的表达值投射到PPI上为可选项,做了该步骤的话可以在data.project中查看结果。
二 推断 cell-cell communication network
前面数据和配体受体库准备好之后,就可以根据表达值推断细胞类型之间的互作了。
使用表达值推测细胞互作的概率,该步骤相对较耗时一些。
cellchat <- computeCommunProb(cellchat, raw.use = TRUE, population.size = TRUE)
# Filter out the cell-cell communication if there are only few number of cells in certain cell groups
cellchat <- filterCommunication(cellchat, min.cells = 10)
注1:raw.use = TRUE 表示使用raw数据,而不使用上一步projectData后的结果。
注2:在假设细胞数较多的群 往往比 细胞数较少的群发送更强的信号的前提下,当population.size = TRUE时候,CellChat可以在概率计算中考虑每个细胞群中细胞比例的影响。
根据需求保存结果,默认保存全部结果(推荐),设置slot.name = "netP" 保存显著的结果,指定sources.use和targets.use则能获取指定配体受体方向的结果,signaling获取指定signaling的结果。
#all the inferred cell-cell communications at the level of ligands/receptors
df.net <- subsetCommunication(cellchat)
write.csv(df.net, "cell-cell_communications.all.csv")
#access the the inferred communications at the level of signaling pathways
df.net1 <- subsetCommunication(cellchat,slot.name = "netP")
#gives the inferred cell-cell communications sending from cell groups 1 and 2 to cell groups 4 and 5.
levels(cellchat@idents)
df.net2 <- subsetCommunication(cellchat, sources.use = c("Epi"), targets.use = c("Fibroblast" ,"T"))
#gives the inferred cell-cell communications mediated by signaling WNT and TGFb.
df.net3 <- subsetCommunication(cellchat, signaling = c("CCL", "TGFb"))
使用computeCommunProbPathway计算每个信号通路的所有配体-受体相互作用的通信结果,结存存放在net 和 netP中 。
然后使用aggregateNet计算细胞类型间整合的细胞通讯结果。
#计算每个信号通路相关的所有配体-受体相互作用的通信结果
cellchat <- computeCommunProbPathway(cellchat)
#计算整合的细胞类型之间通信结果
cellchat <- aggregateNet(cellchat)
net中会有count 和 weight 两个维度,可以选择性可视化。
三 CellChat 可视化
groupSize <- as.numeric(table(cellchat@idents))
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")
根据celltype的个数,灵活调整mfrow = c(2,3) 参数 。像绘制count维度的,只需要修改下 mat <- cellchat@net$count 即可 。
mat <- cellchat@net$weight
par(mfrow = c(2,3), xpd=TRUE)
for (i in 1:nrow(mat)) {
mat2 <- matrix(0, nrow = nrow(mat), ncol = ncol(mat), dimnames = dimnames(mat))
mat2[i, ] <- mat[i, ]
netVisual_circle(mat2, vertex.weight = groupSize, weight.scale = T, edge.weight.max = max(mat), title.name = rownames(mat)[i])
}
首先根据cellchat@netP$pathways展示当前有哪些通路结果,选择感兴趣的进行展示,此处示例展示TGFb通路。levels(cellchat@idents) 查看当前的celltype顺序,然后可以通过vertex.receiver指定target 的细胞类型。
绘制层级图的话 ,需要指定layout为hierarchy ,当前版本默认下出来的是circle图。
cellchat@netP$pathways
pathways.show <- c("TGFb")
levels(cellchat@idents)
#[1] "Epi" "Myeloid" "Fibroblast" "T" "Endo" "un"
vertex.receiver = c(1,2,4,5)
netVisual_aggregate(cellchat, signaling = pathways.show,
vertex.receiver = vertex.receiver,layout = "hierarchy")
左图中间的Target是vertex.receiver选定的细胞类型,右图是除vertex.receiver选中之外的另外的细胞类型。
#Chord diagram
par(mfrow=c(1,1))
netVisual_aggregate(cellchat, signaling = pathways.show, layout = "chord")
#Heatmap
par(mfrow=c(1,1))
netVisual_heatmap(cellchat, signaling = pathways.show, color.heatmap = "Reds")
注:绘制热图的话,需要使用netVisual_heatmap函数
绘制指定受体-配体细胞类型中的全部配体受体结果的气泡图,通过sources.use 和 targets.use指定。
#指定受体-配体细胞类型
netVisual_bubble(cellchat, sources.use = c(3,5),
targets.use = c(1,2,4,6), remove.isolate = FALSE)
同时通过signaling指定展示通路
#指定TGFb和SPP1两个信号通路
netVisual_bubble(cellchat, sources.use = c(3,5), targets.use = c(1,2,4,6),
signaling = c("TGFb","SPP1"), remove.isolate = FALSE)
## Plot the signaling gene expression distribution
p1 = plotGeneExpression(cellchat, signaling = "SPP1")
p1
p2 = plotGeneExpression(cellchat, signaling = "SPP1", type = "dot")
更多的可视化展示方式详见官网:Inference and analysis of cell-cell communication using CellChat (htmlpreview.github.io) 。