摘要:本笔记展示了如何应用 CellChat 来通过定量对比和联合流形学习识别不同生物条件下的主要信号变化。


我们通过将 CellChat 应用于来自两种生物条件的 scRNA-seq 数据:非病变(NL,正常)和病变(LS,病态)的人类皮肤上的细胞,从患有特应性皮炎的患者中,来展示 CellChat 的多样化功能。这两个数据集(条件)在联合聚类后具有相同的细胞群体组成。如果不同数据集之间的细胞群体组成略有或极大不同,请查看另一个相关教程。

CellChat 采用自上而下的方法,即从大局着手,然后在信号机制上进行更详细的细化,以识别不同层次的信号变化,包括细胞间通信的一般原则和功能失调的细胞群体 / 信号通路 / 配体 - 受体。

# 准备工作

# 加载所需的库

r
ptm = Sys.time()
library(CellChat)
载入需要的程辑包:dplyr


载入程辑包:'dplyr'

The following objects are masked from 'package:stats':

    filter, lag

The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union

载入需要的程辑包:igraph


载入程辑包:'igraph'

The following objects are masked from 'package:dplyr':

    as_data_frame, groups, union

The following objects are masked from 'package:stats':

    decompose, spectrum

The following object is masked from 'package:base':

    union

载入需要的程辑包:ggplot2
r
library(patchwork)

# 创建目录以保存图形

r
data.dir <- './comparison'
dir.create(data.dir)
setwd(data.dir)

# 加载每个数据集的 CellChat 对象,然后合并在一起

用户需要分别对每个数据集运行 CellChat,然后将不同的 CellChat 对象合并在一起。如果您有使用较早版本(<1.6.0)获得的 CellChat 对象,请执行 updateCellChat

r
cellchat.NL <- readRDS("cellchat_humanSkin_NL.rds")
cellchat.LS <- readRDS("cellchat_humanSkin_LS.rds")
object.list <- list(NL = cellchat.NL, LS = cellchat.LS)
cellchat <- mergeCellChat(object.list, add.names = names(object.list))
Merge the following slots: 'data.signaling','images','net', 'netP','meta', 'idents', 'var.features' , 'DB', and 'LR'.
r
cellchat
An object of class CellChat created from a merged object with multiple datasets 
 593 signaling genes.
 7563 cells. 
CellChat analysis of single cell RNA-seq data! 
r
execution.time = Sys.time() - ptm
print(as.numeric(execution.time, units = "secs"))
[1] 4.650729
r
# 用户现在可以导出合并的 CellChat 对象和两个单独对象的列表以供以后使用
save(object.list, file = "cellchat_object.list_humanSkin_NL_LS.RData")
save(cellchat, file = "cellchat_merged_humanSkin_NL_LS.RData")

# 第一部分:预测细胞间通信的一般原则

CellChat 从大局着手,以预测细胞间通信的一般原则。在比较多个生物条件下的细胞间通信时,它可以回答以下生物学问题:

  • 细胞间通信是否增强

  • 哪些细胞类型之间的相互作用显著改变

  • 从一个条件到另一个条件,主要源和目标如何变化

# 比较交互作用总数和交互强度

为了回答细胞间通信是否增强的问题,CellChat 比较了来自不同生物条件的推断细胞间通信网络的总交互作用数和交互强度。

r
ptm = Sys.time()
gg1 <- compareInteractions(cellchat, show.legend = F, group = c(1,2))
gg2 <- compareInteractions(cellchat, show.legend = F, group = c(1,2), measure = "weight")
gg1 + gg2

# 比较不同细胞群体间的交互作用数量和强度

为了识别哪些细胞群体之间的相互作用显示出显著变化,CellChat 比较了不同细胞群体间的交互作用数量和强度。

# 不同细胞群体间的交互作用数量或强度差异

可以使用圆形图可视化两个数据集之间的细胞间通信网络的交互作用数量或强度差异,其中红色(或蓝色)边缘代表第二个数据集与第一个数据集相比增加(或减少)的信号。

r
par(mfrow = c(1,2), xpd=TRUE)
netVisual_diffInteraction(cellchat, weight.scale = T)
netVisual_diffInteraction(cellchat, weight.scale = T, measure = "weight")

我们还可以使用热图更详细地显示交互作用数量或强度的差异。顶部的彩色条形图代表热图中显示的值的列之和(进入的信号)。右侧的彩色条形图代表行之和的值(发出的信号)。在颜色条中,红色(或蓝色)代表第二个数据集与第一个数据集相比增加(或减少)的信号。

r
gg1 <- netVisual_heatmap(cellchat)
gg2 <- netVisual_heatmap(cellchat, measure = "weight")
gg1 + gg2

差异网络分析只适用于成对的数据集。如果有更多的数据集用于比较,我们可以直接展示每个数据集中任何两个细胞群体之间的交互作用数量或强度。

为了更好地控制不同数据集中推断网络的节点大小和边缘权重,我们计算了每个细胞群体的最大细胞数和所有数据

集中的最大交互作用数量(或交互作用权重)。

r
weight.max <- getMaxWeight(object.list, attribute = c("idents","count"))
par(mfrow = c(1,2), xpd=TRUE)
for (i in 1:length(object.list)) {
    netVisual_circle(object.list[[i]]@net$count, weight.scale = T, label.edge= F, edge.weight.max = weight.max[2], edge.width.max = 12, title.name = paste0("Number of interactions - ", names(object.list)[i]))
}

# 不同细胞类型间的交互作用数量或强度差异

为了简化复杂的网络并洞察细胞类型级别的细胞间通信,我们可以根据定义的细胞群体聚合细胞间通信。这里我们将细胞群体分类为三种细胞类型,然后重新合并 CellChat 对象列表。

r
group.cellType <- c(rep("FIB", 4), rep("DC", 4), rep("TC", 4))
group.cellType <- factor(group.cellType, levels = c("FIB", "DC", "TC"))
object.list <- lapply(object.list, function(x) {mergeInteractions(x, group.cellType)})
cellchat <- mergeCellChat(object.list, add.names = names(object.list))
Merge the following slots: 'data.signaling','images','net', 'netP','meta', 'idents', 'var.features' , 'DB', and 'LR'.

然后我们可以展示每个数据集中任何两种细胞类型之间的交互作用数量或强度。

r
weight.max <- getMaxWeight(object.list, slot.name = c("idents", "net", "net"), attribute = c("idents","count", "count.merged"))
par(mfrow = c(1,2), xpd=TRUE)
for (i in 1:length(object.list)) {
    netVisual_circle(object.list[[i]]@net$count.merged, weight.scale = T, label.edge= T, edge.weight.max = weight.max[3], edge.width.max = 12, title.name = paste0("Number of interactions - ", names(object.list)[i]))
}

同样,我们也可以使用圆形图显示任何两种细胞类型之间的交互作用数量或强度差异。红色(或蓝色)边缘代表第二个数据集与第一个数据集相比增加(或减少)的信号。

r
par(mfrow = c(1,2), xpd=TRUE)
netVisual_diffInteraction(cellchat, weight.scale = T, measure = "count.merged", label.edge = T)
netVisual_diffInteraction(cellchat, weight.scale = T, measure = "weight.merged", label.edge = T)

# 在二维空间中比较主要的源和目标

通过在二维空间中比较传出和传入的交互强度,可以轻松识别在不同数据集之间发送或接收信号发生显著变化的细胞群体。

r
num.link <- sapply(object.list, function(x) {rowSums(x@net$count) + colSums(x@net$count)-diag(x@net$count)})
weight.MinMax <- c(min(num.link), max(num.link)) # 控制不同数据集中的点大小
gg <- list()
for (i in 1:length(object.list)) {
    gg[[i]] <- netAnalysis_signalingRole_scatter(object.list[[i]], title = names(object.list)[i], weight.MinMax = weight.MinMax)
}
Signaling role analysis on the aggregated cell-cell communication network from all signaling pathways
Signaling role analysis on the aggregated cell-cell communication network from all signaling pathways
r
patchwork::wrap_plots(plots = gg)

从散点图可以看出,Inflam.DC 和 cDC1 在 LS 中相比于 NL 成为了主要的源和目标。LS 中纤维细胞群体也成为主要的源。

此外,我们可以识别 Inflam.DC 和 cDC1 在 NL 和 LS 之间的具体信号变化。

# 识别与一个细胞群体相关的信号变化

r
gg1 <- netAnalysis_signalingChanges_scatter(cellchat, idents.use = "Inflam. DC", signaling.exclude = "MIF")
Visualizing differential outgoing and incoming signaling changes from NL to LS

The following `from` values were not present in `x`: 0

The following `from` values were not present in `x`: 0, -1
r
gg2 <- netAnalysis_signalingChanges_scatter(cellchat, idents.use = "cDC1", signaling.exclude = c("MIF"))
Visualizing differential outgoing and incoming signaling changes from NL to LS

The following `from` values were not present in `x`: 0, 2

The following `from` values were not present in `x`: 0, -1
r
patchwork::wrap_plots(plots = list(gg1,gg2))

r
execution.time = Sys.time() - ptm
print(as.numeric(execution.time, units = "secs"))
[1] 2.031146

# 第二部分:识别具有独特网络架构和相互作用强度的改变信号

CellChat
随后可以识别具有较大(或较小)差异的信号网络、信号群组,以及基于其在多个生物条件下的细胞间通信网络的保守和特定环境的信号路径。

# 基于它们的功能 / 结构相似性识别具有较大(或较小)差异以及信号群组的信号网络

CellChat
执行联合流形学习和分类,以推断通信网络基于它们的功能和拓扑相似性。注意:此分析适用于多于两个数据集。

功能相似性:高度的功能相似性表明主要的发送者和接收者是相似的,可以解释为两个信号路径或两个配体 - 受体对表现出相似和 / 或冗余的角色。注意:功能相似性分析不适用于具有不同细胞类型组成的多个数据集。

结构相似性:用于比较它们的信号网络结构的结构相似性,并不考虑发送者和接收者的相似性。注意:结构相似性分析适用于具有相同细胞类型组成或极大不同细胞类型组成的多个数据集。

我们可以基于功能相似性运行流形和分类学习分析,因为两个数据集具有相同的细胞类型组成。

# 基于功能相似性识别信号群组

r
ptm = Sys.time()
cellchat <- computeNetSimilarityPairwise(cellchat, type = "functional")
Compute signaling network similarity for datasets 1 2 
r
cellchat <- netEmbedding(cellchat, type = "functional")
Manifold learning of the signaling networks for datasets 1 2 
r
cellchat <- netClustering(cellchat, type = "functional")
Classification learning of the signaling networks for datasets 1 2 
r
# 在二维空间中可视化
netVisual_embeddingPairwise(cellchat, type = "functional", label.size = 3.5)
2D visualization of signaling networks from datasets 1 2 

![](https://image.marswh.top/unnamed-chunk-13-1.png)

# 基于结构相似性识别信号群组

r
cellchat <- computeNetSimilarityPairwise(cellchat, type = "structural")
cellchat <- netEmbedding(cellchat, type = "structural")
cellchat <- netClustering(cellchat, type = "structural")
# 在二维空间中可视化
netVisual_embeddingPairwise(cellchat, type = "structural", label.size = 3.5)

# 计算并可视化在学习到的联合流形中的路径距离

我们可以基于它们在共享的二维空间中的欧几里得距离识别具有较大(或较小)差异的信号网络。较大的距离意味着两个数据集之间通信网络在功能或结构相似性方面的较大差异。注意:我们只计算两个数据集之间重叠的信号路径的距离。这里不考虑仅在一个数据集中识别的信号路径。如果有超过三个数据集,可以通过在函数 rankSimilarity 中定义 comparison 进行成对比较。

r
rankSimilarity(cellchat, type = "functional")
Compute the distance of signaling networks between datasets 1 2 

![](https://image.marswh.top/unnamed-chunk-15-1.png)

# 识别和可视化保守的和特定环境的信号路径

通过比较每个信号路径的信息流 / 互动强度,我们可以识别信号路径,(i) 关闭,(ii) 减少,(iii) 打开或 (iv) 增加,通过改变它们在一种条件下的信息流与另一种条件相比。

# 比较每个信号通路的整体信息流

通过比较每个信号通路的信息流,我们可以识别保守和特定于上下文的信号通路。信息流定义为推断网络中所有细胞群之间的通信概率之和(即网络中的总权重)。

这个条形图可以以堆叠模式或非堆叠模式绘制。基于 NL 和 LS 皮肤之间推断网络内整体信息流差异,将重要的信号通路进行排名。在 NL 皮肤中富集的顶级信号通路呈红色,而在 LS 皮肤中富集的呈绿色。

# 比较与每个细胞群体相关的外向(或内向)信号

上述分析总结了外向和内向信号的信息。我们还可以比较两个数据集之间的外向(或内向)信号模式,从而识别表现出不同信号模式的信号通路 / 配体 - 受体。

我们可以结合来自不同数据集的所有已识别信号通路,从而进行侧面比较,包括外向信号、内向信号和通过汇总外向和内向信号在一起的整体信号。注意: rankNet 也显示了整体信号的比较,但它没有显示特定细胞群体中的信号强度。

# 第三部分:识别上调和下调的信号配体 - 受体对

# 通过比较通信概率识别功能失调的信号

我们可以通过设置 comparison 函数在 netVisual_bubble 中,比较一些细胞群体到其他细胞群体介导的配体 - 受体对的通信概率。

此外,我们可以通过指定 max.datasetmin.dataset 在函数 netVisual_bubble 中,识别相对于另一个数据集在一个数据集中上调(增加)和下调(减少)的信号配体 - 受体对。增加的信号意味着这些信号在一个数据集中的通信概率(强度)比另一个数据集高。

# 通过使用差异表达分析识别功能失调的信号

通过比较两个数据集中每对 L-R 和每对细胞群体之间的通信概率来识别上调和下调的信号,这是上述方法的一种。另一种方法,我们可以基于差异表达分析(DEA)识别上调和下调的信号配体 - 受体对。具体来说,我们对每个细胞群体在两种生物条件(即 NL 和 LS)之间进行差异表达分析,然后根据发送细胞中配体和接收细胞中受体的折叠变化获得上调和下调的信号。

请注意,用户可能会观察到相同的 LR 对同时出现在上调和下调的结果中,因为条件之间的 DEA 是针对每个细胞群体进行的。为了忽略细胞群体信息进行条件间的 DEA,用户可以在 CellChat v2.1.1 中将 group.DE.combined 设置为 TRUE。

由于 net.upnet.down 中的信号基因可能是复杂的多亚基体,我们可以进一步解构以获得单个信号基因。

我们然后使用气泡图或和弦图可视化上调和下调的信号配体 - 受体对。

# 第四部分:使用层次图、圆形图或和弦图可视化比较细胞 - 细胞通信

与单个数据集的 CellChat
分析类似,我们可以使用层次图、圆形图或和弦图可视化细胞 - 细胞通信网络。

边缘颜色 / 权重,节点颜色 / 大小 / 形状:在所有可视化图中,边缘颜色与作为发送者的源一致,边缘权重与交互强度成正比。更粗的边线表示更强的信号。在层次图和圆形图中,圆的大小与每个细胞群体中的细胞数量成正比。在层次图中,实心和空心圆分别代表源和目标。在和弦图中,内部较薄的条颜色代表接收来自相应外部条的目标的信号。内部条的大小与目标接收的信号强度成正比。这样的内部条有助于解读复杂的和弦图。请注意,一些细胞群体的内部条没有任何和弦,请忽略它,因为这是
circlize 包尚未解决的问题。

# 第五部分:比较不同数据集之间的信号基因表达分布

我们可以使用 Seurat 包装器函数 plotGeneExpression 绘制与 L-R 对或信号通路相关的信号基因的基因表达分布。

# 保存合并的 CellChat 对象

r
save(object.list, file = "cellchat_object.list_humanSkin_NL_LS.RData")
save(cellchat, file = "cellchat_merged_humanSkin_NL_LS.RData")
r
sessionInfo()
R version 4.1.2 (2021-11-01)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Linux Mint 21.2

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0

locale:
 [1] LC_CTYPE=zh_CN.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=zh_CN.UTF-8        LC_COLLATE=zh_CN.UTF-8    
 [5] LC_MONETARY=zh_CN.UTF-8    LC_MESSAGES=zh_CN.UTF-8   
 [7] LC_PAPER=zh_CN.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=zh_CN.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] patchwork_1.1.3     CellChat_2.1.1      Biobase_2.54.0     
[4] BiocGenerics_0.40.0 ggplot2_3.4.4       igraph_1.5.1       
[7] dplyr_1.1.3        

loaded via a namespace (and not attached):
  [1] backports_1.4.1       circlize_0.4.15       systemfonts_1.0.4    
  [4] NMF_0.26              plyr_1.8.9            lazyeval_0.2.2       
  [7] BiocParallel_1.28.3   listenv_0.9.0         ggnetwork_0.5.12     
 [10] gridBase_0.4-7        digest_0.6.33         foreach_1.5.2        
 [13] htmltools_0.5.6       ggalluvial_0.12.5     fansi_1.0.4          
 [16] magrittr_2.0.3        cluster_2.1.2         doParallel_1.0.17    
 [19] sna_2.7-2             ComplexHeatmap_2.10.0 globals_0.16.2       
 [22] matrixStats_1.0.0     svglite_2.1.3         colorspace_2.1-0     
 [25] ggrepel_0.9.3         textshaping_0.3.6     xfun_0.40            
 [28] crayon_1.5.2          jsonlite_1.8.7        iterators_1.0.14     
 [31] glue_1.6.2            registry_0.5-1        gtable_0.3.4         
 [34] GetoptLong_1.0.5      car_3.1-2             future.apply_1.11.0  
 [37] shape_1.4.6           abind_1.4-5           scales_1.2.1         
 [40] rngtools_1.5.2        rstatix_0.7.2         Rcpp_1.0.11          
 [43] viridisLite_0.4.2     xtable_1.8-4          clue_0.3-65          
 [46] reticulate_1.34.0     stats4_4.1.2          htmlwidgets_1.6.2    
 [49] httr_1.4.7            FNN_1.1.3.2           RColorBrewer_1.1-3   
 [52] ellipsis_0.3.2        pkgconfig_2.0.3       farver_2.1.1         
 [55] sass_0.4.7            utf8_1.2.3            tidyselect_1.2.0     
 [58] labeling_0.4.3        rlang_1.1.1           reshape2_1.4.4       
 [61] later_1.3.1           munsell_0.5.0         tools_4.1.2          
 [64] cachem_1.0.8          cli_3.6.1             generics_0.1.3       
 [67] statnet.common_4.9.0  broom_1.0.5           evaluate_0.22        
 [70] stringr_1.5.0         fastmap_1.1.1         yaml_2.3.7           
 [73] ragg_1.2.5            knitr_1.45            purrr_1.0.2          
 [76] pbapply_1.7-2         future_1.33.0         mime_0.12            
 [79] compiler_4.1.2        rstudioapi_0.15.0     plotly_4.10.2        
 [82] png_0.1-8             ggsignif_0.6.4        tibble_3.2.1         
 [85] bslib_0.5.1           stringi_1.7.12        RSpectra_0.16-1      
 [88] lattice_0.20-45       Matrix_1.6-1.1        vctrs_0.6.3          
 [91] pillar_1.9.0          lifecycle_1.0.3       BiocManager_1.30.22  
 [94] jquerylib_0.1.4       GlobalOptions_0.1.2   BiocNeighbors_1.12.0 
 [97] data.table_1.14.8     cowplot_1.1.1         irlba_2.3.5.1        
[100] httpuv_1.6.11         R6_2.5.1              promises_1.2.1       
[103] network_1.18.2        IRanges_2.28.0        parallelly_1.36.0    
[106] codetools_0.2-18      rjson_0.2.21          withr_2.5.1          
[109] presto_1.0.0          S4Vectors_0.32.4      parallel_4.1.2       
[112] grid_4.1.2            tidyverse_2.0.0       tidyr_1.3.0          
[115] coda_0.19-4           rmarkdown_2.25        carData_3.0-5        
[118] ggpubr_0.6.0          shiny_1.7.5          
更新于