[onc]复杂热图之瀑布图

2024-01-31  本文已影响0人  花生学生信
##输入文件:
#01矩阵,每行一个基因,每列一个sample。01矩阵代表基因在这个sample里是否变异。
#
#pathway.txt,基因所在的pathway;
# 
# input_cli.csv,信息,例如人种、病史、治疗方式等。
# 
# 不同文件之间,基因名、sampleID必须一致。


library(ComplexHeatmap)


#读取基因变异数据
mygene <- read.csv("1matrix01.csv", header = T,row.names = 1, as.is = T)

#mygene<-read.table("clipboard",header = T,row.names = 1,as.is = T)

#查看部分数据的格式
mygene[1:3,1:4]


#用"mut"替换1,用""替换0
mygene[mygene == 1] <- "mut"
mygene[mygene == 0] <- ""

#读取基因所在的pathway
mypathway <- read.csv("2pathway.csv", header = T, row.names = 1, as.is = T)
head(mypathway)

#把pathway添加到mygene的最后一列
mygene$pathway <- mypathway[as.character(rownames(mygene)),]

#读取临床信息
mytype <- read.csv("3clib.csv")
head(mytype)

mytype$V2


#查看有多少种临床信息类型
unique(mytype$type_A)

unique(mytype$type_B)

unique(mytype$type_C)


##设置变异的颜色、形状
#用fill = 设置mutation以及背景用什么颜色
alter_fun = list(
  background = function(x, y, w, h) {
    grid.rect(x, y, w-unit(0.5, "mm"), h-unit(0.5, "mm"), gp = gpar(fill = NA, col = NA)) #不要背景色
  },
  mut = function(x, y, w, h) {
    grid.rect(x, y, w-unit(0.5, "mm"), h-unit(0.5, "mm"), gp = gpar(fill = "#008000", col = NA)) #mut是绿色
  })

#mutation bar plot的颜色,跟瀑布图一致
col = c("mut" = "#008000")

#定义足够多的颜色
mycol <- c("#223D6C","#D20A13","#FFD121","#088247","#11AA4D","#58CDD9","#7A142C","#5D90BA","#431A3D","#91612D","#6E568C","#E0367A","#D8D155","#64495D","#7CC767")


###从临床数据中删除没有变异的sample
#参数remove_empty_columns = TRUE会删掉没有变异的sample,我们需要把临床数据中相应的sample也删掉。

#为了提取瀑布图中的sample,修改了oncoprint函数,保存为oncoPrint_plus.R,搜索#Ya#查看修改的位置。

#注意:在后面所有画图命令中,column_order和remove_empty_columns这两个参数,都要跟p1的参数一致。

#提取瀑布图中画出的sample
svg("1.onco.svg",width = 30,height = 15)
p1 <- oncoPrint(mygene[1:(ncol(mygene)-1)], get_type = function(x) x,
                alter_fun = alter_fun, col = col,
                remove_empty_columns = TRUE, #删除没有突变的sample
                #column_order = NULL, #不按突变频率给sample排序
                #row_order = NULL, #不按突变频率给基因排序
                row_split = mygene$path)
p1
dev.off()


matrix <- p1@matrix
sampleOrder <-  data.frame(p1@column_order)
rownames(sampleOrder) <- p1@column_names_param$labels
sampleOrder$oriOrder <- row.names(sampleOrder)
sampleOrder <- sampleOrder[order(as.numeric(sampleOrder[,1]),decreasing=F),]
rownames(mytype) <- mytype$sample

#在临床数据中,只保留瀑布图中画出的sample
mytype <- mytype[sampleOrder$oriOrder,]

mytype


##现在,这个mytype就跟瀑布图中的sample一致了。
mytype[3]

#画临床数据的heatmap
my_annotation = HeatmapAnnotation(df = data.frame(mytype[3]),
                                  col = list(Chr = c("chr01" = mycol[5],"chr06" = mycol[2], "chr02" = mycol[4], "chr11" = mycol[1],"chr07" = mycol[3],"chr08" = mycol[6])))


library(ggplot2)

mytype[2]

#画瀑布图
png("2.lubby.png")
p2 <- oncoPrint(mygene[1:(ncol(mygene)-1)], get_type = function(x) x,
                alter_fun = alter_fun, col = col,
                remove_empty_columns = TRUE,#删除没有突变的sample
                #column_order = NULL, #不按突变频率给sample排序
                #row_order = NULL, #不按突变频率给基因排序
                show_pct = FALSE, #左侧不显示百分比
                bottom_annotation = my_annotation,#把临床信息画在下面
                row_split = mygene$path, #按照pathway分开画
                show_heatmap_legend = FALSE) #不显示突变的图例

p2
dev.off()

我放一张结果图 只能放一半
上一篇下一篇

猜你喜欢

热点阅读