你只是想畫個相關性熱圖!G镀鳌8卣妗!
生信分析中經(jīng)常要做基因與代謝物之間的相關性熱圖爽航,R中號稱繪制相關性熱圖的包有許多蚓让,但ggplot2始終是無敵的存在乾忱,運行以下代輕松搞定相關性熱圖。
library("ggplot2")
library("dplyr")
library("reshape")
library("psych")
library("ggtree")
library("aplot")
table1 <- mtcars[,1:5] %>% as.matrix()
table2 <- mtcars[,6:11] %>% as.matrix()
#讀入自己數(shù)據(jù)請執(zhí)行
table <- read.table("1.txt",header = T,sep="\t",check.names = F,row.names = 1)
##注:2張表的行名名稱與順序必須一致
#psych包的corr.test()計算相關性系數(shù)與P值,adjust = "fdr"進行多重矯正
pp <- corr.test(table1,table2,method="pearson",adjust = "fdr")
#抽提出相關性系數(shù)與P值
cor <- pp$r
pvalue <- pp$p
#根據(jù)P值用循環(huán)將P值做成*
display <- pvalue
l1 <- nrow(display);l2 <- ncol(display)
for(i in 1:l1){
for(k in 1:l2){
a <- as.numeric(display[i,k])
if(a <= 0.001){
a <- "***"
}
if( 0.001 < a && a <= 0.01){
a <- "**"
}
if(0.01 < a && a < 0.05){
a <- "*"
}
if(a >= 0.05){
a <- ""
}
display[i,k] <- a
}
}
#對數(shù)據(jù)進行格式轉換適用于ggplot2繪圖
heatmap <- melt(cor)
heatmap[,4] <- melt(pvalue)[,3]
heatmap[,5] <- melt(display)[,3]
names(heatmap) <- c("sample","gene","cor","pvalue","display")
#導出數(shù)據(jù)
write.table (heatmap,file ="heatmap.xls", sep ="\t", row.names = F)
#根據(jù)相關性系數(shù)矩陣历极,用ggtree繪制行聚類與列聚類樹
phr <- hclust(dist(cor)) %>%
ggtree(layout="rectangular", branch.length="none")
phc <- hclust(dist(t(cor))) %>% ggtree() + layout_dendrogram()
#ggplot2繪制熱圖
pp <- ggplot(heatmap,aes(gene,sample,fill=cor)) +
geom_tile()+
theme_minimal()+scale_fill_viridis_c()+
geom_text(aes(label=display),size=5,color="white")+
scale_y_discrete(position="right")+xlab(NULL) + ylab(NULL)+
theme(axis.text.x = element_text(angle = 90,hjust=0.8,vjust=0.6))+
theme(axis.text.x=element_text(family = "Times",
face = "plain",colour = "black",size=10))+
theme(axis.text.y=element_text(family = "Times",
face = "plain",colour = "black",size=10))+
theme(legend.text=element_text(face="plain",family = "Times",
colour = "black",size = 10))+
labs(fill = "expr")
#通過aplot包將聚類樹與熱圖拼接
pp %>% insert_left(phr, width=.1) %>%
insert_top(phc, height=.1)
heatmap.jpeg