社区首页 >专栏 >R语言ggplot2复现一下CELL论文中的基因共线性图

R语言ggplot2复现一下CELL论文中的基因共线性图

用户7010445

发布于 2024-06-18 09:35:51

30700

代码可运行

文章被收录于专栏：小明的数据分析笔记本小明的数据分析笔记本

运行总次数：0

代码可运行

论文

Reciprocal conversion between annual and polycarpic perennial flowering behavior in the Brassicaceae

https://www.sciencedirect.com/science/article/pii/S0092867424004732

今天的推文复现一下论文中的 FigureS3b

前提是已经通过分析确定了哪些基因存在共线性，现在想用ggplot2来作图

构造数据

（基因长度和基因间区的长度并没有按照真实比例来，论文中是否是按照真实比例来画的暂时没有看明白）

用来画基因的数据

从前向后

基因的起始位置

纵坐标

基因名

物种名

基因方向

基因类别

用来画共线性的数据

这个是手动整理的，如何用代码根据第一个数据得到这个数据暂时还没有想明白

作图代码

library(tidyverse)
library(readxl)
library(gggenes)
library(ggforce)

dat01<-read_excel("2024.data/20240618/dat01.xlsx")
dat01 %>% 
  mutate(X1=case_when(
    X5 == "species01" ~ X1,
    X5 == "species02" ~ X1 -5,
    X5 == "species03" ~ X1 + 10,
    X5 == "species04" ~ X1 + 5,
    X5 == "species05" ~ X1 + 13,
  ),
  X2=case_when(
    X5 == "species01" ~ X2,
    X5 == "species02" ~ X2 -5,
    X5 == "species03" ~ X2 + 10,
    X5 == "species04" ~ X2 + 5,
    X5 == "species05" ~ X2 + 13,
  )) -> dat01

dat01 %>% 
  group_by(X5) %>% 
  summarise(min=min(X1),
            max=min(X2))


dat01 %>% 
  group_by(X5) %>% 
  summarise(min=max(X1),
            max=max(X2))

dat02<-read_excel("2024.data/20240618/dat02.xlsx")

dat01 %>% 
  group_by(X5) %>% 
  summarise(min=min(X1),
            max=max(X2),
            X7=min(X3)) %>% 
  ungroup() -> segment.df

ggplot()+
  geom_segment(data=segment.df,
               aes(x=min-3,xend=max+3,y=X7,yend=X7),
               color="grey")+
  geom_diagonal_wide(data = dat02,
                     aes(x=x,y=y,group=group),
                     fill="grey",
                     strength = 0.5,
                     radius=unit(0,'mm'),
                     orientation = "y")+
  geom_gene_arrow(data=dat01,
                  aes(xmin = X1, 
                      xmax = X2, 
                      y = X3, 
                      fill = X7,
                      forward=X6),
                  arrowhead_height = unit(3, "mm"), 
                  arrowhead_width = unit(1.5, "mm"),
                  alpha=1)+
  theme_bw()+
  theme(panel.grid = element_blank(),
        panel.border = element_blank(),
        axis.ticks = element_blank(),
        axis.text.x = element_blank(),
        axis.title = element_blank(),
        legend.position = "none",
        axis.text.y = element_text(face="italic"))+
  scale_y_continuous(breaks = 1:5,
                     label=c("A","B","C","D","E"))+
  scale_fill_manual(values = c("#e59f01","#56b4e8",
                               "#009f73","#0072b1"))+
  geom_text(data=dat01,aes(x=X1,y=X3,label=X4),
            hjust=0,vjust=2,
            fontface="italic")