[R绘图]oncoprint的一些细节优化

目的：参考临床里的变异0/1信息，绘制特定一个TE在水稻群体里的分布

如图：

我们在这里将pathway换成水稻里的亚群（group）、type换成具体的染色体

原始数据只要准备对应的TE分布和样本对应Group即可，并将染色体分布整和成新表

image.png

library(ComplexHeatmap)
#读取基因变异数据
mygene <- read.csv("karma.csv", header = T,row.names = 1)
#查看部分数据的格式
mygene[1:3,1:4]

数据格式，用0/1的矩阵表示TE的有无

Te_name <- "karma"  # 定义变量Te_name
mygene[mygene == 1] <- Te_name  # 将mygene中值为1的元素替换为Te_name的值
mygene[mygene == 0] <- ""

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

即对应的亚群

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

head(mygene)

#读取染色体信息
mytype <- read.csv("3clib.csv")
head(mytype)

即TE对应的Chr，后面会出教程整理

#查看有多少种TE位点及chr信息类型
unique(mytype$Pos)
unique(mytype$Chr)

##设置变异的颜色、形状
#用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的颜色，跟瀑布图一致
Te_name
col = c(mut = "#FF0000")
#col[Te_name] <- "#FF0000"
col
#定义足够多的颜色
mycol <- c("#223D6C","#D20A13","#FFD121","#088247","#11AA4D","#58CDD9","#7A142C","#5D90BA","#431A3D","#91612D","#6E568C","#E0367A","#D8D155","#64495D","#7CC767")

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

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

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


dev.new()
p1<- oncoPrint(mygene[1:(ncol(mygene)-1)], get_type = function(x) x,
                alter_fun = alter_fun, col = col,
               alter_fun_is_vectorized =TRUE,
                remove_empty_columns = TRUE, #删除没有突变的sample
                column_order = 1:ncol(mygene), #不按突变频率给sample排序
                row_order = NULL, #不按突变频率给基因排序
                row_split = mygene$path)
p1

p1

matrix <- p1@matrix
matrix
sampleOrder <-  data.frame(p1@column_order)
sampleOrder

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

##heatmap
my_annotation = HeatmapAnnotation(df = data.frame(mytype[3]),
                                  col = list(Chr = c("chr01" = mycol[5], "chr02" = mycol[4], "chr03" = mycol[3],"chr04" = mycol[2],"chr05" = mycol[1],"chr06" = mycol[6],"chr07" = mycol[7],"chr08" = mycol[8],"chr09" = mycol[9],"chr10" = mycol[10],"chr11" = mycol[11],"chr12" = mycol[12])))

p2