本文首发于“生信补给站”公众号 https://mp.weixin.qq.com/s/csdUsCCWVO-uYwnGw29N6w
之前介绍过使用cellphoneDB 进行细胞通讯分析scRNA分析 | 解决可能的报错,从0开始教你完成细胞通讯分析-cellphoneDB,可能会遇到一些报错。这次介绍另一款细胞通讯分析的常见方法CellChat 。CellChat是一款R包,使用更容易且可视化结果也非常不错。
一 数据输入,处理
1,载入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)
2,创建CellChat对象
可以使用矩阵数据、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
3,设置参考数据库
根据分析数据的物种,可选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。
4,表达数据的预处理
可以使用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
前面数据和配体受体库准备好之后,就可以根据表达值推断细胞类型之间的互作了。
1,推断细胞通讯网络
使用表达值推测细胞互作的概率,该步骤相对较耗时一些。
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可以在概率计算中考虑每个细胞群中细胞比例的影响。
2,提取 保存结果
根据需求保存结果,默认保存全部结果(推荐),设置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"))
3,计算cell-cell communication
使用computeCommunProbPathway计算每个信号通路的所有配体-受体相互作用的通信结果,结存存放在net 和 netP中 。
然后使用aggregateNet计算细胞类型间整合的细胞通讯结果。
#计算每个信号通路相关的所有配体-受体相互作用的通信结果 cellchat <- computeCommunProbPathway(cellchat) #计算整合的细胞类型之间通信结果 cellchat <- aggregateNet(cellchat)
net中会有count 和 weight 两个维度,可以选择性可视化。
三 CellChat 可视化
文章题目就提到了可视化出众,哪能做哪些常见的图呢?
1,celltype之间通讯结果
1)根据使用netVisual_circle显示任意两个celltype之间的通讯次数(左)或总通讯强度(右)
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")
2)分别展示
根据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]) }
2,单个信号通路可视化
首先根据cellchat@netP$pathways展示当前有哪些通路结果,选择感兴趣的进行展示,此处示例展示TGFb通路。levels(cellchat@idents) 查看当前的celltype顺序,然后可以通过vertex.receiver指定target 的细胞类型。
1)层级图
绘制层级图的话 ,需要指定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选中之外的另外的细胞类型。
2)和弦图 和 热图
#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函数
3,绘制配体受体气泡图
1)指定受体-配体细胞类型
绘制指定受体-配体细胞类型中的全部配体受体结果的气泡图,通过sources.use 和 targets.use指定。
#指定受体-配体细胞类型 netVisual_bubble(cellchat, sources.use = c(3,5), targets.use = c(1,2,4,6), remove.isolate = FALSE)
2)指定受体-配体细胞类型 且 指定通路
同时通过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)
3)某条信号通路(如TGFb)的所有基因在细胞群中的表达情况展示
## Plot the signaling gene expression distribution p1 = plotGeneExpression(cellchat, signaling = "SPP1") p1 p2 = plotGeneExpression(cellchat, signaling = "SPP1", type = "dot")