这个教程概述了使用CellChat对单个数据集进行细胞间通信网络的推理、分析和可视化的步骤。我们通过将CellChat应用于来自患者的病变(LS,患病)人类皮肤细胞的scRNA-seq数据来展示CellChat的多种功能。
CellChat需要细胞的基因表达数据作为用户输入,并通过将基因表达与信号配体(ligands)、受体(receptors)及其辅助因子(cofactors)之间相互作用的先验知识相结合来模拟细胞间通讯的概率。
在推断细胞间通信网络后,CellChat提供进一步数据探索、分析和可视化的功能。
library(CellChat)library(patchwork)options(stringsAsFactors=FALSE)#reticulate::use_python("/Users/suoqinjin/anaconda3/bin/python",required=T)1PartI:CellChat对象的数据输入、处理、和初始化CellChat需要两个用户输入:一个是细胞的基因表达数据,另一个是用户分配的细胞标签。
对于基因表达数据矩阵,行名为genes,列名为cells。Normalizeddata(e.g.,library-sizenormalizationandthenlog-transformedwithapseudocountof1)作为CellChat分析的输入。如果用户提供countdata,我们提供一个normalizeData函数来计算librarysize,然后进行log-transformed。
对于cellgroup信息,需要一个带有行名的dataframe作为CellChat的输入。
这里使用选项A导入数据
按照以下选项之一创建一个CellChat对象:选项A:使用数字基因表达矩阵和细胞标签信息作为输入,选项B:使用Seurat对象作为输入,选项C:使用SingleCellExperiment对象作为输入,选项D:使用AnnData对象作为输入
注意:如果用户加载了先前计算的CellChat对象(版本<0.5.0),请通过updateCellChat更新对象。
这里使用选项A创建一个CellChat对象
cellchat<-createCellChat(object=data.input,meta=meta,group.by="labels")如果在创建CellChat对象时未添加cellmeta数据信息,用户也可以在之后使用addMeta添加它,并使用setIdent设置默认细胞身份。
cellchat<-addMeta(cellchat,meta=meta)cellchat<-setIdent(cellchat,ident.use="labels")#set"labels"asdefaultcellidentitylevels(cellchat@idents)#showfactorlevelsofthecelllabelsgroupSize<-as.numeric(table(cellchat@idents))#numberofcellsineachcellgroup1.3设置配体-受体相互作用数据库在用户可以使用CellChat推断细胞间通讯之前,他们需要设置配体-受体(ligand-receptor)相互作用数据库,并确定过表达的配体(ligands)或受体(receptors)。
在分析人类样本时,请使用数据库CellChatDB.human;在分析小鼠样本时,请使用数据库CellChatDB.mouse。CellChatDB将配体-受体对分类为不同类型,包括“分泌信号(SecretedSignaling)”、“细胞外基质-受体(ECM-Receptor)”、“细胞-细胞接触(Cell-CellContact)”和“非蛋白质信号(Non-proteinSignaling)”。默认情况下,“非蛋白质信号(Non-proteinSignaling)”不被使用。
CellChatDB<-CellChatDB.human#useCellChatDB.mouseifrunningonmousedatashowDatabaseCategory(CellChatDB)
CellChat还提供了将基因表达数据投射到protein-proteininteraction(PPI)网络上的功能。具体来说,adiffusionprocess用于根据在高置信度实验验证的PPI网络中定义的邻居基因的表达值来平滑基因的表达值。此功能在分析具有浅测序深度的单细胞数据时非常有用,因为投射减少了信号基因的丢失效应,特别是对于配体/受体亚基的零表达。人们可能会担心这种diffusionprocess可能引入artifact,但是,它只会引入非常微弱的通信。默认情况下CellChat使用rawdata(i.e.,object@data.signaling)而不是projecteddata。为了使用projecteddata,用户需要在运行computeCommunProb()之前运行projectData,并在运行computeCommunProb时设置raw.use=TRUE。
#对信号基因的表达数据取子集以节省计算成本cellchat<-subsetData(cellchat)#Thisstepisnecessaryevenifusingthewholedatabasefuture::plan("multisession",workers=4)#doparallelcellchat<-identifyOverExpressedGenes(cellchat)cellchat<-identifyOverExpressedInteractions(cellchat)execution.time=Sys.time()-ptmprint(as.numeric(execution.time,units="secs"))#将基因表达数据投射到PPI上(可选:运行时,用户应在函数computeCommunProb()中设置raw.use=FALSE以使用投射数据)#cellchat<-projectData(cellchat,PPI.human)2PartII:推断细胞间通信网络CellChat通过为每个交互分配一个概率值并执行排列测试来推断具有生物学意义的细胞间通讯。CellChat通过使用质量作用定律将基因表达与信号配体、受体及其辅因子之间相互作用的先前已知知识相结合来模拟细胞间通讯的可能性。
注意:推断的配体-受体对数目明显取决于计算每个cellgroup的平均基因表达量的方法。默认情况下,CellChat使用一种称为“trimean”的统计上稳健的均值方法,与其他方法相比,它产生的交互更少。然而,我们发现CellChat在预测更强的相互作用方面表现良好,这对于缩小相互作用范围以进行进一步的实验验证非常有帮助。在computeCommunProb中,我们提供了使用其他方法,例如5%和10%truncatedmean来计算平均基因表达的选项。值得注意的是,“trimean”近似于25%truncatedmean,这意味着如果一组中表达细胞的百分比小于25%,则平均基因表达为零。要使用10%truncatedmean,用户可以设置type="truncatedMean"和trim=0.1。为了确定trim的合适值,CellChat提供了一个函数computeAveExpr可以帮助检查感兴趣的信号基因的平均表达,例如computeAveExpr(cellchat,features=c("CXCL12","CXCR4"),type="truncatedMean",trim=0.1)。因此,如果在研究的生物过程中未预测到已知的信号通路,用户可以尝试使用较低的trim值来改变计算每个cellgroup平均基因表达量的方法。
在分析未分类的单细胞转录组时,假设丰富的细胞群往往比稀有细胞群发送更强的信号,CellChat还可以在概率计算中考虑每个细胞群中细胞比例的影响。用户可以设置population.size=TRUE。
ptm=Sys.time()cellchat<-computeCommunProb(cellchat,type="triMean")这项分析的关键参数是type,用于计算每个细胞组平均基因表达量的方法。默认情况下,type="triMean",会产生较少但较强的相互作用。当设置type="truncatedMean"时,应为trim分配一个值,这会产生更多的相互作用。请在上面详细查看计算每个细胞组平均基因表达量的方法。
用户可以在某些细胞组中仅有少数细胞时过滤掉细胞间通讯。默认情况下,每个细胞组中用于细胞间通讯的最小细胞数为10。
cellchat<-filterCommunication(cellchat,min.cells=10)2.2将推断的细胞通信网络提取为一个数据框CellChat提供了一个函数subsetCommunication来轻松访问感兴趣的推断细胞间通信。例如,
注意:每个配体-受体对和每个信号通路的推断细胞间通讯网络分别存储在slot\`net\`和\`netP\`中。
cellchat<-computeCommunProbPathway(cellchat)2.4计算聚合的细胞间通信网络我们可以通过计算连接数或汇总通信概率来计算聚合的细胞间通信网络。用户还可以通过设置sources.use和targets.use来计算cellgroups子集之间的聚合网络。
cellchat<-aggregateNet(cellchat)execution.time=Sys.time()-ptmprint(as.numeric(execution.time,units="secs"))CellChat还可以可视化聚合的细胞间通信网络。例如,使用circleplot显示任意两个cellgroups之间的交互次数或总交互强度(权重)。
ptm=Sys.time()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="Numberofinteractions")netVisual_circle(cellchat@net$weight,vertex.weight=groupSize,weight.scale=T,label.edge=F,title.name="Interactionweights/strength")
由于复杂的细胞-细胞通信网络,我们可以检查从每个cellgroup发送的信号。在这里我们还控制参数edge.weight.max以便我们可以比较不同网络之间的边权重(edgeweights)。
mat<-cellchat@net$weightpar(mfrow=c(3,4),xpd=TRUE)for(iin1: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),label.edge=F,title.name=rownames(mat)[i])}
在推断细胞间通信网络后,CellChat为进一步的数据探索、分析和可视化提供了各种功能。具体有:
Hierarchyplot:用户应该定义vertex.receiver,它是一个numericvector,给出了cellgroups的索引作为层次结构图左侧的目标。该层次图由两个部分组成:左侧部分显示了对某些感兴趣的细胞群(即定义的vertex.receiver)的自分泌和旁分泌信号,右侧部分显示了对数据集中其余细胞群的自分泌和旁分泌信号。因此,层次结构图提供了一种信息丰富且直观的方式来可视化感兴趣的细胞群之间的自分泌和旁分泌信号通信。例如,在研究成纤维细胞和免疫细胞之间的细胞间通讯时,用户可以将vertex.receiver定义为所有成纤维细胞群。
Chorddiagram:CellChat提供了两个函数netVisual_chord_cell和netVisual_chord_gene用于可视化不同目的和不同层次的细胞间通信。netVisual_chord_cell用于可视化不同细胞群之间的细胞间通讯(弦图中每个扇区为一个细胞群),netVisual_chord_gene用于可视化由多个配体-受体或信号通路介导的细胞间通讯(其中弦图中的每个部分都是配体、受体或信号通路)。
Explnationsofedgecolor/weight,nodecolor/size/shape:在所有可视化图中,edgecolors与作为发送者的sources一致,edgeweights与交互强度成正比。edge线越粗表示信号越强。在Hierarchyplot和Circleplot中,圆形大小与每个cellgroup中的细胞数成正比。在hierarchyplot中,实心圆和空心圆分别代表source和target。在Chorddiagram中,内部较细的条形颜色表示从相应外部条形接收信号的目标。内部条形大小与目标接收到的信号强度成正比。这样的内部小节有助于解释复杂的和弦图。请注意,对于某些cellgroup,存在一些没有任何和弦的内部小节,请忽略它,因为这是circlize包尚未解决的问题。
这里我们以一个信号通路的输入为例。所有显示重要通信的信号通路都可以通过cellchat@netP$pathways访问。
pathways.show<-c("CXCL")#Hierarchyplot#Herewedefine`vertex.receive`sothattheleftportionofthehierarchyplotshowssignalingtofibroblastandtherightportionshowssignalingtoimmunecellsvertex.receiver=seq(1,4)#anumericvector.netVisual_aggregate(cellchat,signaling=pathways.show,vertex.receiver=vertex.receiver)#Circleplotpar(mfrow=c(1,1))netVisual_aggregate(cellchat,signaling=pathways.show,layout="circle")
#Chorddiagrampar(mfrow=c(1,1))netVisual_aggregate(cellchat,signaling=pathways.show,layout="chord")
#Heatmappar(mfrow=c(1,1))netVisual_heatmap(cellchat,signaling=pathways.show,color.heatmap="Reds")#>Doheatmapbasedonasingleobject
#Chorddiagramgroup.cellType<-c(rep("FIB",4),rep("DC",4),rep("TC",4))#groupingcellclustersintofibroblast,DCandTCcellsnames(group.cellType)<-levels(cellchat@idents)netVisual_chord_cell(cellchat,signaling=pathways.show,group=group.cellType,title.name=paste0(pathways.show,"signalingnetwork"))#>Plottheaggregatedcell-cellcommunicationnetworkatthesignalingpathwaylevel
netAnalysis_contribution(cellchat,signaling=pathways.show)
pairLR.CXCL<-extractEnrichedLR(cellchat,signaling=pathways.show,geneLR.return=FALSE)LR.show<-pairLR.CXCL[1,]#showoneligand-receptorpair#Hierarchyplotvertex.receiver=seq(1,4)#anumericvectornetVisual_individual(cellchat,signaling=pathways.show,pairLR.use=LR.show,vertex.receiver=vertex.receiver)#>[[1]]#CircleplotnetVisual_individual(cellchat,signaling=pathways.show,pairLR.use=LR.show,layout="circle")
#ChorddiagramnetVisual_individual(cellchat,signaling=pathways.show,pairLR.use=LR.show,layout="chord")
在实际使用中,用户可以使用'for...loop'来自动保存所有推断的网络,以便使用netVisual进行快速探索。netVisual支持svg、png和pdf格式的输出。
#Accessallthesignalingpathwaysshowingsignificantcommunicationspathways.show.all<-cellchat@netP$pathways#checktheorderofcellidentitytosetsuitablevertex.receiverlevels(cellchat@idents)vertex.receiver=seq(1,4)for(iin1:length(pathways.show.all)){#VisualizecommunicationnetworkassociatedwithbothsignalingpathwayandindividualL-RpairsnetVisual(cellchat,signaling=pathways.show.all[i],vertex.receiver=vertex.receiver,layout="hierarchy")#Computeandvisualizethecontributionofeachligand-receptorpairtotheoverallsignalingpathwaygg<-netAnalysis_contribution(cellchat,signaling=pathways.show.all[i])ggsave(filename=paste0(pathways.show.all[i],"_L-R_contribution.pdf"),plot=gg,width=3,height=2,units='in',dpi=300)}3.2可视化由多个配体-受体或信号通路介导的细胞间通讯CellChat还可以使用函数netVisual_bubble(选项A)和netVisual_chord_gene(选项B)显示由配体-受体对和信号通路介导的所有显著相互作用,以及用户从某些细胞组到其他细胞组提供的相互作用。
我们还可以使用netVisual_bubble显示从一些cellgroups到其他cellgroups的所有重要相互作用(L-Rpairs)。
#(1)showallthesignificantinteractions(L-Rpairs)fromsomecellgroups(definedby'sources.use')toothercellgroups(definedby'targets.use')netVisual_bubble(cellchat,sources.use=4,targets.use=c(5:11),remove.isolate=FALSE)#>Comparingcommunicationsonasingleobject
#(2)showallthesignificantinteractions(L-Rpairs)associatedwithcertainsignalingpathwaysnetVisual_bubble(cellchat,sources.use=4,targets.use=c(5:11),signaling=c("CCL","CXCL"),remove.isolate=FALSE)#>Comparingcommunicationsonasingleobject
#(3)showallthesignificantinteractions(L-Rpairs)basedonuser'sinput(definedby`pairLR.use`)pairLR.use<-extractEnrichedLR(cellchat,signaling=c("CCL","CXCL","FGF"))netVisual_bubble(cellchat,sources.use=c(3,4),targets.use=c(5:8),pairLR.use=pairLR.use,remove.isolate=TRUE)#>Comparingcommunicationsonasingleobject
#showallthesignificantinteractions(L-Rpairs)fromsomecellgroups(definedby'sources.use')toothercellgroups(definedby'targets.use')#showalltheinteractionssendingfromInflam.FIBnetVisual_chord_gene(cellchat,sources.use=4,targets.use=c(5:11),lab.cex=0.5,legend.pos.y=30)
#showalltheinteractionsreceivedbyInflam.DCnetVisual_chord_gene(cellchat,sources.use=c(1,2,3,4),targets.use=8,legend.pos.x=15)
#showallthesignificantinteractions(L-Rpairs)associatedwithcertainsignalingpathwaysnetVisual_chord_gene(cellchat,sources.use=c(1,2,3,4),targets.use=c(5:11),signaling=c("CCL","CXCL"),legend.pos.x=8)
#showallthesignificantsignalingpathwaysfromsomecellgroups(definedby'sources.use')toothercellgroups(definedby'targets.use')netVisual_chord_gene(cellchat,sources.use=c(1,2,3,4),targets.use=c(5:11),slot.name="netP",legend.pos.x=10)
注意:在生成绘图时请忽略注释,例如“注意:第一个链接端是从扇区‘MIF’绘制出来的。”。如果基因名称重叠,您可以通过减小值来调整参数small.gap。
plotGeneExpression(cellchat,signaling="CXCL",enriched.only=TRUE,type="violin")
plotGeneExpression(cellchat,signaling="CXCL",enriched.only=FALSE)execution.time=Sys.time()-ptm
为了便于解释复杂的细胞间通信网络,CellChat通过从图论、模式识别和流形学习中抽象出来的方法对网络进行定量测量。
用户可以在热图(选项A)和二维图(选项B)上可视化中心性分数。CellChat还可以回答某些细胞组的outgoing或incoming信号中哪些信号贡献最大的问题(选项C)。
ptm=Sys.time()#Computethenetworkcentralityscorescellchat<-netAnalysis_computeCentrality(cellchat,slot.name="netP")#theslot'netP'meanstheinferredintercellularcommunicationnetworkofsignalingpathways#Visualizethecomputedcentralityscoresusingheatmap,allowingreadyidentificationofmajorsignalingrolesofcellgroupsnetAnalysis_signalingRole_network(cellchat,signaling=pathways.show,width=8,height=2.5,font.size=10)
#Signalingroleanalysisontheaggregatedcell-cellcommunicationnetworkfromallsignalingpathwaysgg1<-netAnalysis_signalingRole_scatter(cellchat)#>Signalingroleanalysisontheaggregatedcell-cellcommunicationnetworkfromallsignalingpathways#Signalingroleanalysisonthecell-cellcommunicationnetworksofinterestgg2<-netAnalysis_signalingRole_scatter(cellchat,signaling=c("CXCL","CCL"))#>Signalingroleanalysisonthecell-cellcommunicationnetworkfromuser'sinputgg1+gg2
我们还可以回答哪些信号对某些cellgroups的outgoing或incoming信号贡献最大的问题。在这个热图中,颜色条代表一个信号通路在cellgroups间的相对信号强度(注:值经过row-scaled)。顶部有色条形图显示了通过总结热图中显示的所有信号通路来表示cellgroup的总信号强度。右侧的灰色条形图显示了通过总结热图中显示的所有cellgroup来表示信号通路的总信号强度。
#Signalingroleanalysisontheaggregatedcell-cellcommunicationnetworkfromallsignalingpathwaysht1<-netAnalysis_signalingRole_heatmap(cellchat,pattern="outgoing")ht2<-netAnalysis_signalingRole_heatmap(cellchat,pattern="incoming")ht1+ht2
#Signalingroleanalysisonthecell-cellcommunicationnetworksofinterestht<-netAnalysis_signalingRole_heatmap(cellchat,signaling=c("CXCL","CCL"))4.2识别全局通信模式以探索多种细胞类型和信号通路如何协调在一起除了探索各个通路的详细通信之外,一个重要的问题是多个cellgroups和信号通路如何协调发挥作用。CellChat采用模式识别方法来识别全局通信模式。
随着模式数量的增加,可能会出现冗余模式,从而难以解释通信模式。CellChat默认选择了五种模式。一般来说,模式的数量大于2才具有生物学意义。此外,CellChat还提供了一个函数selectK来推断模式的数量,它基于NMFR包中已经实现的两个指标,包括Cophenetic和Silhouette。这两个指标都基于共识矩阵的层次聚类来衡量特定数量模式的稳定性。对于一定数量的模式,合适数量的模式是Cophenetic和Silhouette值开始突然下降的模式。
这项分析可以针对outgoing(选项A)和incoming(选项B)信号模式进行。Outgoing模式显示了发送细胞(即cellsassignalsource)如何相互协调以及它们如何与特定信号通路协调以推动通讯。Incoming模式显示了目标细胞(即cellsassignalreceivers)如何相互协调以及它们如何与特定信号通路协调以响应入向信号。
Outgoingpatterns揭示了sendercells(i.e.cellsassignalsource)如何相互协调以及它们如何与某些信号通路协调以驱动通信。
加载通信模式分析所需的包
library(NMF)library(ggalluvial)在这里,我们运行selectK来推断模式的数量。
selectK(cellchat,pattern="outgoing")
当传出模式的数量为6时,Cophenetic和Silhouette值都开始突然下降。
nPatterns=6cellchat<-identifyCommunicationPatterns(cellchat,pattern="outgoing",k=nPatterns)
#riverplotnetAnalysis_river(cellchat,pattern="outgoing")#>Pleasemakesureyouhaveload`library(ggalluvial)`whenrunningthisfunction
#dotplotnetAnalysis_dot(cellchat,pattern="outgoing")
传入模式显示了目标细胞(i.e.cellsassignalreceivers)如何相互协调以及它们如何与某些信号通路协调以响应传入信号。
selectK(cellchat,pattern="incoming")
当传入模式的数量为3时,Cophenetic值开始下降。
nPatterns=3cellchat<-identifyCommunicationPatterns(cellchat,pattern="incoming",k=nPatterns)
#riverplotnetAnalysis_river(cellchat,pattern="incoming")#>Pleasemakesureyouhaveload`library(ggalluvial)`whenrunningthisfunction
#dotplotnetAnalysis_dot(cellchat,pattern="incoming")
此外,CellChat能够量化所有重要信号通路之间的相似性,然后根据它们的细胞通信网络相似性对它们进行分组。可以根据功能(Functional)或结构(Structural)相似性进行分组。
Functionalsimilarity:高度的功能相似性表明主要的发送者(senders)和接收者(receivers)相似,可以解释为两个信号通路或两个配体-受体对表现出相似和/或冗余的作用。功能相似性分析需要两个数据集之间具有相同的细胞群组成。
Structuralsimilarity:结构相似性用来比较它们的信号网络结构,不考虑发送者和接收者的相似性。
cellchat<-computeNetSimilarity(cellchat,type="functional")cellchat<-netEmbedding(cellchat,type="functional")#>Manifoldlearningofthesignalingnetworksforasingledatasetcellchat<-netClustering(cellchat,type="functional")#>Classificationlearningofthesignalingnetworksforasingledataset#Visualizationin2D-spacenetVisual_embedding(cellchat,type="functional",label.size=3.5)
netVisual_embeddingZoomIn(cellchat,type="functional",nCol=2)
cellchat<-computeNetSimilarity(cellchat,type="structural")cellchat<-netEmbedding(cellchat,type="structural")#>Manifoldlearningofthesignalingnetworksforasingledatasetcellchat<-netClustering(cellchat,type="structural")#>Classificationlearningofthesignalingnetworksforasingledataset#Visualizationin2D-spacenetVisual_embedding(cellchat,type="structural",label.size=3.5)
netVisual_embeddingZoomIn(cellchat,type="structural",nCol=2)
execution.time=Sys.time()-ptmprint(as.numeric(execution.time,units="secs"))5PartV:保存CellChat对象saveRDS(cellchat,file="cellchat_humanSkin_LS.rds")6PartVI:通过交互式CellChatExplorer探索细胞间通信对于单细胞转录组学的CellChat分析,请确保object@dr包含数据的低维空间,例如“umap”和“tsne”,以便生成信号基因的特征图。可以通过函数addReduction添加新的低维空间。