R语言作图——热图聚类及其聚类结果输出

代码

不多说了，做个记录，代码如下。

library(pheatmap)
library(RColorBrewer)
# args = commandArgs(TRUE)
betafile = "twist_common_panel_434.csv"
infofile = "twist_common_panel_434.txt"
title = "twist_common_panel"
# readbeta <- function(data){
#   dat <- read.table(data, header = T,stringsAsFactors = F)
#   dat$pos <- paste(dat$chr, dat$start, sep = "_")
#   dat$chr = dat$start = NULL
#   rownames(dat) <- dat$pos
#   dat$pos <- NULL
#   dat
# }
# d <- readbeta(betafile)
# d<-read.table(betafile,header=T,stringsAsFactors=FALSE,row.names = 1)
# d <- as.matrix(t(d))
d <- read.csv(betafile, row.names = 1)
info <- read.table(infofile,header = F,sep = "\t")
colnames(info) <- c("id","group")
rownames(info) <- info$id
# table(rownames(d) == info$id)
#排序
loc = match(rownames(d), rownames(info))
dd <- d[loc,]
col_anno <- as.data.frame(info[,2])
rownames(col_anno) <- rownames(info)
colnames(col_anno) <- 'group'
tiff(paste0(title,".tif"),width = 2000,height = 1500)
pheatmap(dd,# cellwidth = 45,# cellheight = 12,fontsize = 18,# border = 'white',cluster_rows = F,cluster_cols = F,# annotation_col = col_anno,annotation_row = col_anno,angle_col = "0",main = paste0(title, " pheatmap"))
dev.off()result <- pheatmap(dd,# cellwidth = 45,# cellheight = 12,fontsize = 18,# border = 'white',cluster_rows = T,cluster_cols = T,# annotation_col = col_anno,annotation_row = col_anno,angle_col = "0",main = paste0(title, " pheatmap"))
col_oder=result$tree_col$order # 保存热图列顺序（序号）
row_oder=result$tree_row$order # 保存热图行顺序（序号）
cn_new <- colnames(dd)[col_oder] # 保存热图的列名
rn_new <- rownames(dd)[row_oder] # 保存热图的行名## 生成两个新的与原来的总表长宽相同的数据框（暂时都填写0）
new_dd <- matrix(rep(0,ncol(dd)*nrow(dd)),nrow = nrow(dd),ncol = ncol(dd))
out <- matrix(rep(0,ncol(dd)*nrow(dd)),nrow = nrow(dd),ncol = ncol(dd))## 将数据读入数据框new_dd,行顺序已经按照热图的顺序重排
for (i in 1:ncol(dd)){new_dd[,i]=dd[,i][row_oder]
}
## 将数据读入数据框out,列顺序已经进一步按照热图的顺序重排
for (i in 1:ncol(dd)){out[,i]= new_dd[,col_oder[ i ]]
}
# 将热图的行名和列名导入到排序后的表达量总表中
rownames(out)=rn_new
colnames(out)=cn_newwrite.table(out,"twist大panel聚类.txt",sep="\t",quote = F) #输出重排后的表达量表