论文
Single-cell profiling of vascular endothelial cells reveals progressive organ-specific vulnerabilities during obesity
https://www.nature.com/articles/s42255-022-00674-x#Sec58
s42255-022-00674-x.pdf
https://github.com/Osynchronika/sc_EC_obesity_atlas
大部分 作图的数据都有，可以试着用论文中提供的数据复现一下论文中的图
今天的推文我们试着复现一下论文中的Figure2中的热图，figure2中有3个热图，按照复杂程度排序是 figure2m figure2f 和figure2o
今天的推文我们先复现figure2o
论文中提供的数据如下
数据中有很多缺失值，看论文中的配色 我猜是把缺失值替换成0了，我不太确定这种处理方式是否可以
加载要用到的R包
library(readxl)
library(tidyverse)
library(stringr)
library(ggtree)
library(ggh4x)读取数据
dat01<-read_excel("data/20230207/figure2o.xlsx",
sheet = "Sheet3",na="NA")
dat01这里我把第一列的列名除了第一列修改成A01:A21，因为最上层的两层分组需要人为构造数据，列名修改成A01:A21在excel里好输入些
colnames(dat01)<-c("gene_name",paste0("A",str_pad(1:21,2,side = "left",pad = "0")))首先是做最左侧的树图
因为有了树图后要根据树图的前后顺序调整热图的y轴的顺序
dat01 %>%
pivot_longer(!gene_name) %>%
mutate(value=replace_na(value,0)) %>%
pivot_wider() %>%
column_to_rownames("gene_name") %>%
dist() %>%
hclust() -> dat01.hclust
ggtree(dat01.hclust)+
geom_tiplab()
y_levels<-c("Fabp5","Fabp4","Cd36","Fabp1","Dbi","Lpl")
ggtree(dat01.hclust) -> p2
p2
然后是热图的代码
dat01 %>%
mutate(gene_name=factor(gene_name,
levels = y_levels )) %>%
pivot_longer(!gene_name) %>%
mutate(value=replace_na(value,0)) %>%
ggplot(aes(x=name,y=gene_name))+
geom_tile(aes(fill=value),color="black")+
scale_fill_gradient2(low="blue",mid="white",high = "red",
midpoint = 0,
breaks=c(-0.5,-0.25,0,0.25,0.5),
labels=c("<-0.5","",0,"",">0.5"))+
scale_y_discrete(position = "right")+
scale_x_discrete(expand = expansion(mult = c(0,0)))+
theme_bw()+
theme(legend.position = "bottom",
axis.text.x = element_blank(),
axis.text.y = element_text(face="italic"),
axis.ticks = element_blank(),
axis.title = element_blank(),
panel.border = element_blank(),
panel.grid = element_blank())+
guides(fill=guide_colorbar(title.position = "top",
title.hjust = 0.5,
barwidth = 10)) -> p1
p1
最后是分组标记的代码
我人为构造的数据如下
两层的分组标记：一层是热图，另外一层是分面的形式，分面的形式可以添加好多层
dat<-read_excel("data/20230207/figure2o.xlsx",
sheet = "Sheet2")
dat %>% pull(group02) %>% unique()
dat %>%
mutate(group=factor(group02,levels = dat %>% pull(group02) %>% unique())) -> dat
strip.fill<-c("#ee716b","#c99e0d","#16b5ea","#5ab033","#37b28f","#a08ec3","#d96aa6")
ggplot(data=dat,aes(x=x,y=1))+
geom_tile(aes(fill=group01),color="black")+
geom_text(aes(label=group01),color="white")+
facet_nested(.~group02,space = "free",scales = "free",
strip=strip_nested(background_x =
elem_list_rect(fill=strip.fill,
by_layer_x=FALSE)))+
scale_x_discrete(expand = expansion(mult = c(0,0)))+
scale_y_discrete(expand = expansion(mult = c(0,0)))+
scale_fill_manual(values = c("#dd0c1b","#48ad96","#4d6eb4"))+
theme(panel.spacing = unit(0,"lines"),
strip.background = element_rect(fill=c("white"),color="black"),
axis.text = element_blank(),
axis.title = element_blank(),
axis.ticks = element_blank(),
strip.text.x = element_text(margin = margin(0.5,0,0.5,0, "cm")),
legend.position = "none",) -> p3
p3
这里遇到一个问题是
facet_nested(.~group02,space = "free",scales = "free",
strip=strip_nested(background_x =
elem_list_rect(fill=strip.fill,
by_layer_x=FALSE)))代码里scales=free不能和coord_equal()组合使用，暂时不知道如何解决这个问题
最后是组合图
library(patchwork)
(plot_spacer() + p3 + plot_layout(widths = c(1,8)))/(p2+p1+plot_layout(widths = c(1,8)))+
plot_layout(heights = c(1,4))最终结果如下
示例数据和代码可以给推文点赞，然后点击在看，最后留言获取
欢迎大家关注我的公众号
小明的数据分析笔记本
小明的数据分析笔记本 公众号 主要分享：1、R语言和python做数据分析和数据可视化的简单小例子；2、园艺植物相关转录组学、基因组学、群体遗传学文献阅读笔记；3、生物信息学入门学习资料及自己的学习笔记！
微信公众号好像又有改动，如果没有将这个公众号设为星标的话，会经常错过公众号的推文，个人建议将 小明的数据分析笔记本 公众号添加星标，添加方法是
点开公众号的页面，右上角有三个点
点击三个点，会跳出界面
直接点击 设为星标 就可以了
