day 10 GSVA和CellChat

2024-07-01 16:03:59 浏览数 (1)

GSVA

单样本和多样本都适用

  • 输入数据
  • GSVA
  • 可视化
代码语言:R复制
rm(list = ls())
library(Seurat)
library(GSVA)
library(clusterProfiler)
load("../day7/scRNA.Rdata")
seu.obj = scRNA

#文献和教程一般使用平均值来做GSVA
exp  =  AverageExpression(seu.obj)[[1]] #exp =  AggregateExpression(seu.obj)[[1]]
exp  =  as.matrix(exp)
exp  =  exp[rowSums(exp)>0,] ;exp[1:4,1:4]

## GSVA 官网下载gmt文件 并读取
h_df = read.gmt("h.all.v2023.2.Hs.symbols.gmt")[,c(2,1)]
h_list = unstack(h_df)
ES  =  gsva(exp, h_list)
#ES = gsva(gsvaParam(exp,h_list,maxDiff = T)) #针对R4.4.0
ES[1:4,1:4]

#可视化
library(pheatmap)
pheatmap(ES, scale = "row",angle_col = "45",
         color = colorRampPalette(c("navy", "white", "firebrick3"))(50))

如果是多样本,我们可以多做一个组间比较

代码语言:R复制
#添加分组信息
seu.obj$group = ifelse(seu.obj$orig.ident %in% paste0("sample",1:3),"treat","control")

exp  =  AverageExpression(seu.obj,group.by = c("ident","group"))[[1]]
exp  =  as.matrix(exp)
exp  =  exp[rowSums(exp)>0,] 

h_df = read.gmt("h.all.v2023.2.Hs.symbols.gmt")[,c(2,1)]
h_list = unstack(h_df)
ES  =  gsva(exp, h_list);ES[1:4,1:4]

#可视化
library(pheatmap)
pheatmap(ES, scale = "row",angle_col = "45",
         color = colorRampPalette(c("navy", "white", "firebrick3"))(50),
         cluster_cols = F,
         gaps_col = seq(2,ncol(ES)-1,2)) 
         #是第二行开始空格,按两个两个空n到整个ES

CellChat

Cellchat的安装

一手教程,1,2

gitub上面下载zip,删掉src文件夹下的这三个文件

输入数据

代码语言:R复制
rm(list = ls())
if(!require(CellChat))devtools::install_local("CellChat-master/")
if(!require(presto))devtools::install_github('immunogenomics/presto')
library(CellChat)
library(ggplot2)
library(Seurat)
library(ggalluvial)

load("../day5-6/sce.Rdata")
scRNA = sce
DimPlot(scRNA)
scRNA = subset(scRNA,downsample = 100) #抽样节省计算

创建对象

代码语言:R复制
str(CellChatDB.human,max.level = 1) #mouse就换
table(CellChatDB.human$interaction$annotation)
scRNA$samples = scRNA$orig.ident
cellchat <- createCellChat(scRNA,
                           group.by = "ident",
                           assay = "RNA") 

cellchat@DB <- subsetDB(CellChatDB.human, 
                        search = "Secreted Signaling") 

cellchat <- subsetData(cellchat) 
dim(cellchat@data.signaling) 

开始计算

代码语言:R复制
cellchat <- identifyOverExpressedGenes(cellchat) # 识别过表达基因
cellchat <- identifyOverExpressedInteractions(cellchat) # 识别配体-受体对
cellchat <- projectData(cellchat, PPI.human)#慢 将配体、受体投射到PPI网络
cellchat <- smoothData(cellchat, adj = PPI.human) #v2版本的
cellchat <- computeCommunProb(cellchat) #慢 推测细胞通讯网络
cellchat <- computeCommunProbPathway(cellchat)
cellchat <- aggregateNet(cellchat)

选择通路并画图

 hierarchy plot hierarchy plot
代码语言:R复制
cellchat@netP$pathways
pathways.show <- "GALECTIN"
#1 hierarchy plot
groupSize <- as.numeric(table(cellchat@idents)) 
vertex.receiver = seq(1,nlevels(scRNA)/2);vertex.receiver
netVisual_aggregate(cellchat, signaling = pathways.show, layout = "hierarchy", vertex.receiver = vertex.receiver, vertex.weight  = groupSize)  

#2 circle plot
par(mfrow = c(1,1), xpd=TRUE,mar = c(2, 2, 2, 2))
netVisual_aggregate(cellchat, signaling = pathways.show, 
                    layout = "circle", 
                    vertex.receiver = vertex.receiver,
                    vertex.weight  = groupSize)

#3 chord plot
netVisual_aggregate(cellchat, signaling = pathways.show, layout = "chord", vertex.receiver = vertex.receiver, vertex.weight  = groupSize)

#4 heatmap
#纵坐标是发射端,横坐标是接收端,有颜色代表横纵坐标所指的两类细胞之间有通讯,颜色深浅代表通讯概率。右侧和上方的条形图是该行/列通讯概率之和
netVisual_heatmap(cellchat, signaling = pathways.show, color.heatmap = "Reds") 
#5 计算配体-受体对信号网络的贡献度
netAnalysis_contribution(cellchat, signaling = pathways.show) 



#6 热图展示 分析细胞在信号网络中角色:发送者、接收者、调解者和影响者。
cellchat <- netAnalysis_computeCentrality(cellchat, slot.name = "netP") # the slot 'netP' means the inferred intercellular communication network of signaling pathways
netAnalysis_signalingRole_network(cellchat, signaling = pathways.show, width = 12, height = 5, font.size = 10)


#7 气泡图
#可以分开,也可以合到一起 p<0.05的才会被画出来,颜色仍然是通讯概率,圈的大小是按照p值,p值越小圈越大
netVisual_bubble(cellchat, sources.use = 1, 
                 targets.use = 1:nlevels(scRNA), 
                 remove.isolate = FALSE) 

#从第一类细胞到全部细胞
netVisual_bubble(cellchat, sources.use = 1:nlevels(scRNA), 
                 targets.use = 1:nlevels(scRNA), 
                 remove.isolate = FALSE) #从全部细胞到全部细胞

细胞通讯模式和信号网络

  • 传出模式(outgoing),揭示了发射端细胞如何相互协调,以及它们如何与某些信号通路协调以驱动通信。
  • 传出模式(outgoing),显示接收端细胞如何相互协调,以及它们如何与某些信号通路协调以响应输入信号。
代码语言:R复制
library(NMF)
selectK(cellchat, pattern = "outgoing")

#选择合适的细胞通讯模式数量
#这两个指标都是评估聚类稳定性的,二者都突然下降的值对应的横坐标就是合适的聚类数,2到3突然下降就选2,3到4出现了下降就选3
nPatterns = 9 #根据上图选择的#identifyCommunicationPatterns函数识别通讯模式,并画出热图
cellchat <- identifyCommunicationPatterns(cellchat, pattern = "outgoing", k = nPatterns) #传出
nPatterns = 6
cellchat <- identifyCommunicationPatterns(cellchat, pattern = "incoming", k = nPatterns) #传入

可视化通路来源

代码语言:R复制
# 桑基图 用桑基图和气泡图展示每种细胞传入/传出的信号都是属于哪些通路的
p1 <- netAnalysis_river(cellchat, pattern = "outgoing")
p2 <- netAnalysis_river(cellchat, pattern = "incoming")
p1 p2
# 气泡图 气泡图的颜色是按照细胞类型来分配,大小按照每个通路对每个细胞类型的贡献程度分配
p3 <- netAnalysis_dot(cellchat, pattern = "outgoing",dot.size = 4)
p4 <- netAnalysis_dot(cellchat, pattern = "incoming",dot.size = 4)
p3 p4

0 人点赞