R中优雅的可视化WGCNA模块热图

2023-10-24 14:01:51 浏览数 (1)

欢迎关注R语言数据分析指南

加载R包

代码语言：javascript复制
library(tidyverse)
library(ggthemes)
library(magrittr)
library(WGCNA)
library(linkET)

导入WGCNA结果

代码语言：javascript复制
load("df.Rdata")
代码语言：javascript复制
gene_exp <- read.table(file ='heatmap.txt',row.names=1) %>% dplyr::select(1:24)

# 转置
datExpr0 <- t(gene_exp)

# 缺失数据及无波动数据过滤
gsg <- goodSamplesGenes(datExpr0,minFraction = 1/2) #基因的缺失数据比例阈值
datExpr <- datExpr0[gsg$goodSamples, gsg$goodGenes]

WGCNA绘制模块热图

代码语言：javascript复制
MEs2 <- MEs %>% dplyr::select(1:20)

moduleTraitCor <- cor(MEs2,design %>% dplyr::select(1:6),
                      use = "p",method = 'pearson')

# 计算 Pvalue
moduleTraitPvalue <- corPvalueStudent(moduleTraitCor,nrow(datExpr))
sizeGrWindow(10,6)

# 连接相关性和 pvalue
textMatrix <- paste(signif(moduleTraitCor, 2), "n(",
                    signif(moduleTraitPvalue, 1), ")", sep = "");
dim(textMatrix) <- dim(moduleTraitCor)
# heatmap 画图
par(mar = c(5,7,3,1))

labeledHeatmap(Matrix = moduleTraitCor,
                     xLabels = names(design %>% dplyr::select(1:6)),
                     yLabels = names(MEs2),
                     ySymbols = names(MEs2),
                     colorLabels = FALSE,
                     colors = blueWhiteRed(50),
                     textMatrix = textMatrix,
                     setStdMargins = FALSE,
                     cex.text = 0.5,
                     zlim = c(-1,1),
                     main = paste("Module-trait relationships"))

绘制模块热图

代码语言：javascript复制
link_cor <- mantel_test(design %>% select(1,2,3,4,5,6),MEs2,
                        spec_select =list(CAR15=1,CAR5=2,
                                      Cd_CK=3,Cd_treatment=4,
                                      CK_elo=5,CK_mat=6)) %>% 
  mutate(rd=cut(r,breaks=c(-Inf,0.2,0.6,Inf),
                labels=c("<0.2","0.2-0.6",">=0.6")),
         pd=cut(p,breaks=c(-Inf,0.01,0.5,Inf),
                labels=c("<0.01","0.01-0.05",">=0.05")))

qcorrplot(correlate(MEs2),type = "upper", diag = FALSE)  
  geom_square()  
  geom_couple(aes(colour = pd, size = rd), data = link_cor, 
              curvature = nice_curvature())  
  scale_fill_gradientn(colours = rev(RColorBrewer::brewer.pal(11, "RdBu"))) 
  scale_size_manual(values = c(0.5, 1, 2))  
  scale_colour_manual(values = color_pal(3))  
  guides(size = guide_legend(title = "Mantel's r",override.aes = list(colour = "grey35"), order = 2),
         colour = guide_legend(title = "Mantel's p", override.aes = list(size = 3), order = 1),
         fill = guide_colorbar(title = "Pearson's r",order = 3)) 
  theme(plot.margin = unit(c(0,0,0,0),units="cm"),
        legend.background = element_blank(),
        legend.key = element_blank())
数据分析 select 可视化 连接 数据

0 人点赞