(本次操作有配套的视频教程,在果子学生信公众号回复“果子对角线”自行获取,和代码一起以project的形式分享)
果子唠嗑
这是一篇引起血战的帖子,一开始发表在Y叔的公众号上,然后就出现了
ggplot2版本,excel版本,complexheatmap版本,base plot版本。
这篇帖子实现的略显复杂,但是,这个帖子是我做事风格的见证:
用有限的工具,尽量去做无限的事情。就像钢琴的88键,键有数,律无穷。
比如,我只熟悉R语言,那么批量读取文件,数据调整,绘图都只用R语言完成,尽管有很多更好的工具。
比如,我只会ggplot2,所以当出现以下需求的时候,我都尽量用ggplot2去完成,因为没有选择。
这篇帖子最终不是教大家如何画一张图,画一张图没有意义,而是展示R语言里面的一些基本技能,包括
- 如何创建数据框?
- 如何融合多个数据框?
- 如何实现数据的宽长转换?
- 如何调整x轴label的顺序?
- 多张图如何拼接?
- ggplot2怎么用?
- 多边形怎么画?
- 如何把大需求拆解成可完成的小步骤?
- R语言中批量操作神器lapply怎么用?
既然有这么多技能,估计初学者也消化不了,所以我就变废为宝,录了一个视频教程,来讲解这些技能。
教程以R语言project的形式发布,里面包含了原始数据,代码,生成的图片,以及手把手教程。
有需要的在果子学生信公众号回复"果子对角线"自助获取。
以下是正文
最近有人在群里提问,下面的热图该怎么画。
这张热图,在每一个格子里面用对角线一分为二,呈现了两类信息。图片来自于这篇文章的图1D。
通常,遇到这种需求,我都会直接检索,看看有没有现成的R包可用。可以提取检索词检索一下,这里主要信息应该是"对角线 热图 R语言",
查阅字典后,关键词就变成了英文,“diagonal heatmap r”。在浏览器中检索后得到以下结果。
然后有人给出了使用base plot画图的方案
d=data.frame(p1=rep(LETTERS[1:2],times=2,each=2),
p2=rep(LETTERS[4:5],times=4),
value=c(.2,.4,.56,.32,.7,.16,.12,.71),
group=rep(c("A1","B1"),each=4))
x=as.numeric(d$p1)
y=as.numeric(d$p2)
plot(1,xlim=c(1,length(unique(x))+1),ylim=c(1,length(unique(y))+1),
type="n",bty="n",xaxt="n",yaxt="n",xlab="",ylab="")
for(i in 1:nrow(d)) {
if(d$group[i]=="A1") polygon(x[i]+c(0,1,1),y[i]+c(0,0,1),col=gray(d$value[i]))
if(d$group[i]=="B1") polygon(x[i]+c(0,1,0),y[i]+c(0,1,1),col=gray(d$value[i]))
}
axis(1,at=sort(unique(x))+.5,labels=levels(d$p1),lty=0)
axis(2,at=sort(unique(y))+.5,labels=levels(d$p2),lty=0)
复制到R语言后运行得到如下结果
这意味着,只要用基本功扎实,复现原图是没有问题的。
大鱼海棠应该就是用的这个方法。
base plot版本的对角线热图
但是,我不会base plot,我只会用ggplot2, 检索了一圈也没有现成的解决方案,所以我就准备自己画。
自己画我就开始拆分这张图,一个热图的每一个格子由上下两个三角形构成,我可以先画上三角形,再画下三角形,然后给他们分别配色就可以。
(果子:此时,我还不会写图层,尽管有这样的思想,但是没有趁手的工具)
三角形可以通过ggplot2中内置的geom_polygon
来画。
## 下三角数据
cone1 = data.frame(x = c(1,2,2),
y = c(1,1,2))
library(dplyr)
library(ggplot2)
## 下三角画图
cone1 %>%
ggplot(aes(x=x, y=y)) +
geom_polygon(fill="#213c18")
上三角也是一样的
## 上三角数据
cone2 = data.frame(x = c(1,1,2),
y = c(1,2,2))
## 上三角画图
cone2 %>%
ggplot(aes(x=x, y=y)) +
geom_polygon(fill = "#668c6f")
通过叠加图层,我们可以把两个三角形放在一起,这就是热图的一个小格子
cone1 %>%
ggplot(aes(x=x, y=y)) +
geom_polygon(fill="#213c18")+
geom_polygon(data=cone2, fill = "#668c6f")
解决了基本的元素,我就可以批量产生数据。
下面写了一个函数,批量算出给定行和列的三角形数据。
这里用上了,我很喜欢的lapply和do.call函数
trianle <- function(a,b,type="up"){
## 单个下三角数据的函数
trianle_down <- function(a,b){
data.frame(x=c(0,1,1)+a,
y= c(0,0,1)+b,
group=paste0(a,"_",b),
stringsAsFactors = F)
}
## 单个上三角数据的函数
trianle_up <- function(a,b){
data.frame(x=c(0,0,1)+a,
y= c(0,1,1)+b,
group=paste0(a,"_",b),
stringsAsFactors = F)
}
### 批量产生上三角的数据
if(type=="up"){
data <- do.call(rbind,lapply(1:b, function(i){
do.call(rbind,lapply(1:a,trianle_up,i))
}))
}
### 批量产生下三角的数据
if(type=="down"){
data <- do.call(rbind,lapply(1:b, function(i){
do.call(rbind,lapply(1:a,trianle_down,i))
}))
}
return(data)
}
因为我们的测试数据是33行,20列,所以用这个函数来产生对应个数的三角形
updata <- trianle(33,20,"up")
downdata <- trianle(33,20,"down")
updata的数据是这样的,每一个框线里面都是一个三角形的数据,group那一列是为了geom_polygon
画图需要。
可以作图来测试一下上三角形
library(ggplot2)
updata %>%
ggplot(aes(x=x, y=y)) +
geom_polygon(aes(group=group,fill=group))+
theme(legend.position = "none")
下三角也是一样的
downdata %>%
ggplot(aes(x=x, y=y)) +
geom_polygon(aes(group=group,fill=group))+
theme(legend.position = "none")
现在我比较安心了,只要调整原始数据的格式,然后把他跟三角形的数据匹配起来就行。匹配就依靠group那一列
读入第一个数据, 行是基因,列是癌症类型
dd <- data.table::fread("easy_input_amp.txt",data.table = F)
接下来宽数据变长
这里用到的技能我也很喜欢,可以参考这个帖子
我喜欢的gather快要被淘汰了,迎面走来了更好用的宽长转换工具
library(dplyr)
library(tidyr)
data1 <- dd %>%
## 数据变长
pivot_longer(-1,
names_to = "type",
values_to = "exp") %>%
rename(gene=V1)
现在的数据是这个样子的
如果想要改变未来热图基因的顺序,可以使用因子的水平来控制
这个技能在这里讲过
因子(factor)就像贤内助,让你始终分清主次,拨开云雾。
data1$gene <- factor(data1$gene,
levels = c("ALKBH5","FTO",
"HNRNPA2B1","HNRNPC", "IGF2BP1", "IGF2BP2","IGF2BP3","RBMX","YTHDC1","YTHDC2","YTHDF1","YTHDF2","YTHDF3",
"METTL14","METTL3","RBM15","RBM15B","VIRMA","WTAP","ZC3H13"),
ordered = F)
接下来就是给这个数据增加能够和三角形group匹配的列,
这一步的操作有trick,比如,因子转为数值,是需要先变成字符串再变成字符的,但是我们这里直接as.numeric
可以获取他的levels数值
data1 <- data1 %>%
mutate(a=as.numeric(as.factor(type)),
b=as.numeric(gene)) %>%
mutate(group = paste0(a,"_",b)) %>%
as.data.frame() %>%
inner_join(updata,by="group")
最终获取了完整的作图数据
有了数据,作图就很简单了
library(ggplot2)
data1 %>%
ggplot(aes(x=x, y=y)) +
geom_polygon(aes(group=group,fill=exp))+
theme(legend.position = "none")
第二个数据也是同样的处理
dd <- data.table::fread("easy_input_del.txt",data.table = F)
data2 <- dd %>%
pivot_longer(-1,
names_to = "type",
values_to = "del") %>%
rename(gene=V1)
data2$gene <- factor(data2$gene,
levels = c("ALKBH5","FTO",
"HNRNPA2B1","HNRNPC", "IGF2BP1", "IGF2BP2","IGF2BP3","RBMX","YTHDC1","YTHDC2","YTHDF1","YTHDF2","YTHDF3",
"METTL14","METTL3","RBM15","RBM15B","VIRMA","WTAP","ZC3H13"),
ordered = F)
data2 <- data2 %>%
mutate(a=as.numeric(as.factor(type)),
b=as.numeric(as.factor(gene))) %>%
mutate(group = paste0(a,"_",b)) %>%
as.data.frame() %>%
inner_join(downdata,by="group")
data2 %>%
ggplot(aes(x=x, y=y)) +
geom_polygon(aes(group=group,fill=del))+
theme(legend.position = "none")
现在只要把两个图叠加在一起就行了,但是这里面有个问题。这两个图是分别用两个不同的变量来控制颜色的。
但是如果画在一起,ggplot2是不允许的,但是我记得Y叔发过的大部分帖子,我记得他曾经解决过这种情况,就去翻了一下。
答案在这个帖子里面
连续型和离散型数据一起画热图,外加分开配色和图例
这是我们常见的使用gheatmap做图的代码,为了实现不同数据类型的热图放在一起,让我们祭出
ggnewscale
包,只需要用new_scall_fill()
,就可以把后面的fill映射和前面的分离开,然后我们就可以愉快地,再加一条gheatmap再画一个热图。
这样我的问题就完美解决了。
安装和加载 R包
options("repos"=c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/"))
if(!require("ggnewscale")) install.packages("ggnewscale",update = F,ask = F)
library(ggnewscale)
分别抽离出基因名称和亚型名称,用来加在图上
## x轴文字
xlabels <- data1 %>%
filter(y==1) %>%
distinct(type,.keep_all = T) %>%
arrange(x) %>%
pull(type)
## y轴文字
ylabels <- data1 %>%
filter(x==1) %>%
distinct(gene,.keep_all = T) %>%
arrange(y) %>%
pull(gene) %>%
as.character()
指定三个颜色,来自小丫画图的144号作品
red <- "#AB221F"
blue <- "#3878C1"
nake <- "#FFFADD"
画图开始
p1 <- ggplot(data1,aes(x=x, y=y)) +
## 画出上三角
geom_polygon(data = data1,aes(x=x, y=y, group=group,fill=exp),color="black") +
## 上三角用表达值来配色
scale_fill_gradientn(colours = c(nake,red)) +
## 神技能,清空aes
new_scale_fill() +
## 画出下三角
geom_polygon(data = data2,aes(x=x, y=y, group=group,fill=del),color="black") +
## 下三角用表达值来配色
scale_fill_gradientn(colours = c(nake,blue))+
## 调整x轴
scale_x_continuous(limits = c(1, NA), expand = c(0,0),breaks = c(1:33)+0.5,labels=xlabels) +
## 调整y轴
scale_y_continuous(limits = c(1, NA), expand = c(0,0),breaks = c(1:20)+0.5,labels=ylabels) +
## 定主题
theme(axis.title.x=element_blank(),
axis.title.y=element_blank(),
axis.text.x = element_text(angle = 90, hjust = 1))+
## legend置于底部
theme(legend.position="bottom")+
theme(plot.margin = margin(0.5,0.01,0.5,0.01, "cm"))
p1
现在来画侧边的标识,也是自己设计四边形的数据就行了
data <- data.frame(
x = rep(c(1,1,2,2),3),
y = c(1,3,3,1,3,14,14,3,14,21,21,14),
type= rep(c("Eraser","Reader","Writer"),each=4)
)
p2 <- ggplot(data,aes(x=x, y=y)) +
geom_polygon(aes(x=x, y=y, group=type,fill=type),color="black",alpha = 0.5)+
scale_fill_manual(values = c(red,nake,blue))+
## 打标签
geom_text(data=data.frame(x=1.5,y=c(2,8.5,17.5)),
aes(label=c("Eraser","Reader","Writer")),
angle = 90,
size=4)+
scale_x_continuous(limits = c(1, NA), expand = c(0,0))+
scale_y_continuous(limits = c(1, NA), expand = c(0,0)) +
theme(axis.title=element_blank(),
axis.text = element_blank(),
axis.ticks = element_blank())+
theme(legend.position="none")+
theme(plot.margin = margin(0.5,0.1,0.5,0.01, "cm"))
p2
现在我们用cowplot
把他们拼接起来
library(cowplot)
plot_grid(p1,p2,align = "h", axis = "tb",nrow = 1, rel_widths = c(33, 1))
如果用patchwork就更加简单
options("repos"=c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/"))
if(!require("patchwork")) install.packages("patchwork",update = F,ask = F)
library(patchwork)
p1 + p2 +
plot_layout(widths = c(33, 1))
虽然解决了问题,但是我觉得还是不够优雅,我还在学习图层中。Y叔有一篇帖子是让我夜不能寐的
扪心自问,meme几何?
这一年多来,很多次都想去回答这里面的问题,但是次次无功而返。好在近期洲更给了初步的回答,
答: 扪心自问,meme几何?
洲更说他对于图层顿悟了,就跟我问Y叔,他给出的答案是一样的,顿悟。但是我还没有,我还在找寻属于我自己的答案。学习图层,对我个人而言并不是为了画图,我就那么点数据,以后的知识以及足够去呈现了,我只是想从高手的视角去看待作图这件事,就跟Y叔给我启发最大的话那样:
图就是数据,数据就是图。
果子备注
这还没玩,下一篇帖子就讲一下我在画图上的突破。