######--------ssGSEA---------##############
#加载18条神经通路合集
load("./score/ssGSEA/geneset_list_18.RData")
geneset_list
select_pathway_name <- c("G_protein_coupled_ACh_receptor_signalling","Response_to_ACh","ACh_receptor_binding")#个性化修改# ssGSEA计算
expression_matrix <- as.matrix(MERGE@assays$RNA@data)# 提取表达矩阵ssgsea_scores <- gsva(expr = expression_matrix, gset.idx.list = geneset_list[select_pathway_name], method = "ssgsea", kcdf = "Gaussian", abs.ranking = TRUE
)ssgsea_scores[1:3,1:4]# 将ssGSEA分数添加到 Seurat 对象
ssgsea_scores <- t(ssgsea_scores) # 转置,使每行对应细胞,每列对应基因集
MERGE <- AddMetaData(MERGE, metadata = as.data.frame(ssgsea_scores))
head(MERGE@meta.data)# 画图
# 1.山脊图
score_columns <- c("G_protein_coupled_ACh_receptor_signalling","Response_to_ACh","ACh_receptor_binding")for (score_name in score_columns) {p <- ggplot(MERGE@meta.data, aes_string(x = paste0(score_name), y = "cluster", fill = "cluster")) +geom_density_ridges(scale = 1.5, alpha = 0.7, color = "white") +scale_fill_manual(values = color_epi) +theme_bw() +theme(axis.text.x = element_text(size = 12, face = "plain"),axis.text.y = element_text(size = 12, face = "plain"),axis.title.x = element_blank(), # 移除x轴标题axis.title.y = element_text(size = 14, face = "plain"),plot.title = element_text(size = 16, face = "bold", hjust = 0.5), # 设置标题样式panel.grid.major = element_blank(),panel.grid.minor = element_blank(),legend.position = "none") +labs(title = paste0(score_name, " Score"), # 主标题y = "", fill = NULL)# 保存图片ggsave(paste0("ssGSEA_", score_name, "_ridge.png"), p, height = 4, width = 7)
}