今天推文重復(fù)的圖來自于 論文
Whole-genome resequencing of 445 Lactuca accessions reveals the domestication history of cultivated lettuce
這篇論文的數(shù)據(jù)是公開的皱碘，代碼也公開了一部分，那我們就可以按照他的提供的數(shù)據(jù)來試著復(fù)原一些論文中的圖了。
本來已經(jīng)重復(fù)到了論文中Fig1的2c，但是試著做局部放大的時候遇到了一些問題几颜，暫時還搞不定，看了其他的圖實現(xiàn)起來還有一定的難度殴玛。所以先挑一個相對比較簡單的箱線圖來模仿一下吧毕荐。 對應(yīng)的是論文中的Figure 2d
數(shù)據(jù)是論文附件的 Source Data Fig.2
首先是讀入數(shù)據(jù)
df<-readxl::read_excel("NG/41588_2021_831_MOESM5_ESM.xlsx",
sheet="Fig2d")
論文中提供的是寬格式數(shù)據(jù)，如果使用ggplot2作圖需要轉(zhuǎn)換成長格式冕碟，這里本來想嘗試一下tidyr包中的pivot_longer()函數(shù)了拦惋，幫助文檔沒有看明白。沒有搞定安寺，還是直接使用reshape2中的melt()函數(shù)吧
library(dplyr)
df %>%
mutate(new_col=paste(Group1,Group2,sep="_")) %>%
select(-c("Group1","Group2","Group3","Outgroup")) %>%
#reshape2::melt(var.ids=c("Group1")) %>%
#arrange(Group1,Group2) %>%
reshape2::melt(var.ids="new_col") -> df1
head(df1)
ggplot2 作圖
library(ggplot2)
library(stringr)
library(ggprism)
x_level<-paste(df$Group1,df$Group2,sep="_")
x_level
df1$group<-str_sub(df1$new_col,5,7)
df1$new_col<-factor(df1$new_col,
levels = x_level)
ggplot(df1,aes(x=new_col,y=value))+
stat_boxplot(geom = "errorbar",width=0.2)+
geom_boxplot(outlier.shape = 1,
aes(fill=group),
show.legend = F)+
scale_fill_manual(values = c("#e64b35",
"#4daf4a",
"#4dbbd5",
"#cab2d6",
"#b2df8a"))+
scale_x_discrete(labels=str_sub(x_level,1,3),
guide = "prism_offset")+
scale_y_continuous(limits = c(-0.021,0.085),
breaks = seq(-0.02,0.08,by=0.02))+
theme_prism(axis_text_angle = 90,
base_line_size = 0.1,
base_fontface = "plain",
base_family = "serif")+
labs(x=NULL,
y=expression(paste(italic("D")," statistic")))+
theme(plot.margin = unit(c(0.2,0.2,2,0.2),'cm'))+
geom_segment(x=1,xend=5,y=-0.04,yend=-0.04)+
annotate("text",x=3,y=-0.02,label="CAU",vjust=10)+
geom_segment(x=6,xend=9,y=-0.04,yend=-0.04)+
annotate("text",x=7.5,y=-0.02,label="SEU",vjust=10)+
geom_segment(x=10,xend=12,y=-0.04,yend=-0.04)+
annotate("text",x=11.5,y=-0.02,label="WEU",vjust=10)+
geom_segment(x=13,xend=14,y=-0.04,yend=-0.04)+
annotate("text",x=13.5,y=-0.02,label="EEU",vjust=10)+
geom_segment(x=14.5,xend=15.5,y=-0.04,yend=-0.04)+
annotate("text",x=15,y=-0.02,label="WAS",vjust=10)+
annotate("text",x=1,y=-0.02,label="P1",
hjust=2,vjust=5)+
annotate("text",x=1,y=-0.02,label="P2",
hjust=2,vjust=10)+
coord_cartesian(clip = "off")
箭頭指的地方如何用代碼縮短暫時不知道了厕妖，出圖收手動調(diào)整吧
最終結(jié)果
這里遇到的問題是
- 如何將箱線圖的垂直線改成虛線呢？
歡迎大家關(guān)注我的公眾號
小明的數(shù)據(jù)分析筆記本
小明的數(shù)據(jù)分析筆記本 公眾號 主要分享：1挑庶、R語言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡單小例子言秸；2软能、園藝植物相關(guān)轉(zhuǎn)錄組學(xué)、基因組學(xué)举畸、群體遺傳學(xué)文獻閱讀筆記查排；3、生物信息學(xué)入門學(xué)習(xí)資料及自己的學(xué)習(xí)筆記抄沮！
