論文

Graph pangenome captures missing heritability and empowers tomato breeding

https://www.nature.com/articles/s41586-022-04808-9#MOESM8

s41586-022-04808-9.pdf

沒有找到論文里的作圖的代碼，但是找到了部分做圖數(shù)據(jù)，我們可以用論文中提供的原始數(shù)據(jù)模仿出論文中的圖

今天的推文重復(fù)一下論文中的 Figure3a Figure3b Figure3c 頻率分布直方圖,堆積柱形圖,散點(diǎn)圖

image.png

頻率分布直方圖代碼

library(readxl)
fig3a<-read_excel("data/20220711/41586_2022_4808_MOESM7_ESM.xlsx",
                  sheet = "Fig3a",
                  skip = 1)
head(fig3a)
dim(fig3a)

library(ggplot2)
library(latex2exp)

ggplot(data=fig3a,aes(x=h2))+
  geom_histogram(aes(fill=type),
                 bins = 100)+
  scale_fill_manual(values = c("#d1edcd","#94c2db","#fdbeb8"),
                    label=c(TeX(r"(\textit{h}${^2}$ (all variants)        )"),
                            TeX(r"(\textit{h}${^2}$ (leading variants))"),
                            TeX(r"(\textit{h}${^2}$ (local variants)    )")),
                    name="")+
  theme_bw()+
  theme(panel.grid = element_blank(),
        legend.position = c(0.8,0.8))+
  scale_x_continuous(expand = expansion(mult = c(0,0)),
                     breaks = seq(0,1,0.2))+
  scale_y_continuous(expand = expansion(mult = c(0,0)))+
  labs(y="Counts",
       x=TeX(r"(\textit{h}${^2}$)"))+
  geom_vline(xintercept = 0.27,lty="dashed",color="#94c2db")+
  geom_vline(xintercept = 0.37,lty="dashed",color="#fdbeb8")+
  geom_vline(xintercept = 0.62,lty="dashed",color="#d1edcd") -> p1

x<-c(0.27,0.37,0.62)

for (i in 1:3){
  p1<-p1+
    annotate(geom = "text",x=x[i],y=80,label=x[i],hjust=0)
}
p1

image.png

堆積柱形圖

fig3b<-read_excel("data/20220711/41586_2022_4808_MOESM7_ESM.xlsx",
                  sheet = "Fig3b")
head(fig3b)
dim(fig3b)

fig3b$var2<-factor(fig3b$var2,
                   levels = c("MLM","LASSO","Overlapping"))

library(tidyverse)
fig3b %>% 
  group_by(var1) %>% 
  summarise(y=stack.bar.label.position(value),
            y_label=value) %>% 
  ungroup() -> df.label


stack.bar.label.position<-function(x){
  x<-rev(x)
  new.x<-vector()
  
  for (i in 1:length(x)){
    if (i == 1){
      new.x<-append(new.x,x[i]/2)
    }
    
    else{
      new.x<-append(new.x,sum(x[1:i-1])+x[i]/2)
    }
  }
  return(new.x)
}


ggplot(data=fig3b,aes(x=var1,y=value))+
  geom_bar(stat="identity",
           position = "stack",
           aes(fill=var2))+
  scale_fill_manual(values = c("#5ba555","#2baae1","#c6dcf0"),
                    name="",
                    label=c("MLM unique (11)",
                            "LASSO unique (1,249)",
                            "Overlapping (538)"))+
  theme_classic()+
  theme(legend.position = c(0.8,0.8))+
  geom_text(data=df.label,
            aes(x=var1,y=y,label=y_label)) -> p2
p2

image.png

最后的散點(diǎn)圖

fig3c<-read_excel("data/20220711/41586_2022_4808_MOESM7_ESM.xlsx",
                  sheet = "Fig3c",
                  skip = 1)
head(fig3c)
dim(fig3c)



ggplot(data=fig3c %>% filter(Type=="MLM"),
       aes(x=pos,y=-log10(pvalue)))+
  geom_point(aes(shape=Variant,color=Variant,size=Variant))+
  scale_color_manual(values = c("#868686","#b8275a"))+
  theme_classic()+
  scale_x_continuous(labels = function(x)
    {sprintf("%0.2f",x/1000000)})+
  labs(x="Chr3 (Mb)",
       y=TeX(r"(-log${_1}{_0}$$\left[$\textit{P}$\right]$)"))+
  geom_hline(yintercept = 6,lty="dashed")+
  ggtitle("MLM")+
  theme(legend.position = "none")+
  scale_y_continuous(limits = c(0,10),
                     breaks = c(0,5,10)) -> p3.1

ggplot(data=fig3c %>% filter(Type=="LASSO"),
       aes(x=pos,y=-log10(pvalue)))+
  geom_point(aes(shape=Variant,color=Variant),
             size=3)+
  scale_color_manual(values = c("#b8275a"))+
  scale_shape_manual(values = 17)+
  theme_classic()+
  scale_x_continuous(labels = function(x)
  {sprintf("%0.2f",x/1000000)},
  limits = c(42.90*1000000,43*1000000))+
  labs(x="Chr3 (Mb)",
       y=TeX(r"(-log${_1}{_0}$$\left[$\textit{P}$\right]$)"))+
  geom_hline(yintercept = 6,lty="dashed")+
  ggtitle("LASSO")+
  theme(legend.position = "none")+
  scale_y_continuous(breaks = c(0,5,10))+
  geom_text(aes(label=ID),hjust=1.2) -> p3.2
library(patchwork)
p3.1/p3.2

image.png

最終的拼圖

p1+p2 + (p3.1/p3.2)

image.png

示例數(shù)據(jù)和代碼可以自己到論文中獲取，或者給本篇推文點(diǎn)贊，點(diǎn)擊在看，然后留言獲取

歡迎大家關(guān)注我的公眾號(hào)

小明的數(shù)據(jù)分析筆記本

小明的數(shù)據(jù)分析筆記本 公眾號(hào) 主要分享：1、R語言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡(jiǎn)單小例子；2、園藝植物相關(guān)轉(zhuǎn)錄組學(xué)、基因組學(xué)、群體遺傳學(xué)文獻(xiàn)閱讀筆記；3、生物信息學(xué)入門學(xué)習(xí)資料及自己的學(xué)習(xí)筆記！