從這個系列開始一睁,師兄就帶著大家從各大頂級期刊中的Figuer入手厅翔,從仿照別人的作圖風(fēng)格到最后實(shí)現(xiàn)自己游刃有余的套用在自己的分析數(shù)據(jù)上呢簸!這一系列絕對是高質(zhì)量!還不趕緊點(diǎn)贊+在看漓概,學(xué)起來漾月!
示例數(shù)據(jù)和代碼獲取
話不多說,直接上圖胃珍!
讀圖
這張圖理解起來沒什么復(fù)雜的梁肿,就是一個分組提琴圖蜓陌,然后將兩個組的小提琴分別顯示一半,這樣更方便讀者直觀比較吩蔑。本小節(jié)我們介紹兩種實(shí)現(xiàn)方法钮热,一種是基于
gghalves
包中的geom_half_violin
函數(shù),另一種是借助github
大佬編寫的geom_split_violin
函數(shù)烛芬。
效果展示
由于本次使用的數(shù)據(jù)分布并不是很好隧期,所以提琴的形狀并不是很美觀,但是圖形的外觀和細(xì)節(jié)都基本復(fù)現(xiàn)了原文赘娄。本次復(fù)現(xiàn)完全在R語言中進(jìn)行仆潮,請大家放心食用!
數(shù)據(jù)構(gòu)建
####################### 分半提琴圖 ####################
library(ggplot2)
library(gghalves)
library(tidyverse)
# 讀取測試數(shù)據(jù):此數(shù)據(jù)集來源于GSE142651遣臼,隨機(jī)挑選25個基因:
data <- read.csv("data.csv")
data <- data[sample(1:nrow(data), 10),]
# 寬數(shù)據(jù)轉(zhuǎn)長數(shù)據(jù):
data_new <- data %>%
pivot_longer(cols = !X,
names_to = "Samples",
values_to = "Values")
colnames(data_new)[1] <- "Genes"
# 添加分組信息:
data_new$group <- str_split(data_new$Samples, "_", simplify = T)[,4]
# 查看數(shù)據(jù)
head(data_new)
# # A tibble: 6 x 4
# Genes Samples Values group
# <chr> <chr> <dbl> <chr>
# 1 MCM5 Chip91481_r20_c71_Untreated 7.84 Untreated
# 2 MCM5 Chip91481_r47_c21_Untreated 5.12 Untreated
# 3 MCM5 Chip91484_r0_c62_Untreated 5.67 Untreated
# 4 MCM5 Chip91481_r16_c70_Untreated 5.12 Untreated
# 5 MCM5 Chip91484_r0_c35_Treated 6.67 Treated
# 6 MCM5 Chip91484_r37_c38_Untreated 5.12 Untreated
繪圖代碼
geom_half_violin函數(shù)
# 繪圖:
ggplot()+
geom_half_violin(
data = data_new %>% filter(group == "Treated"),
aes(x = Genes,y = Values),colour="white",fill="#1ba7b3",side = "l"
)+
geom_half_violin(
data = data_new %>% filter(group == "Untreated"),
aes(x = Genes,y = Values),colour="white",fill="#dfb424",side = "r"
)+
theme_bw()+
xlab("")+
ylab("log2(CPM)")+
geom_point(data = data_new, aes(x = Genes,y = Values, fill = group),
stat = 'summary', fun=mean,
position = position_dodge(width = 0.2))+
stat_summary(data = data_new, aes(x = Genes,y = Values, fill = group),
fun.min = function(x){quantile(x)[2]},
fun.max = function(x){quantile(x)[4]},
geom = 'errorbar', color='black',
width=0.01,size=0.5,
position = position_dodge(width = 0.2))+
stat_compare_means(data = data_new, aes(x = Genes,y = Values, fill = group),
# 修改顯著性標(biāo)注:
symnum.args=list(cutpoints = c(0, 0.001, 0.01, 0.05, 1),
symbols = c("***", "**", "*", "-")),
label = "p.signif",
label.y = max(data_new$Values),
hide.ns = F)+
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.justification = "right")
ggsave("violin_plot.pdf", height = 5, width = 10)
方法二
# 方法二:使用geom_split_violion函數(shù):
# 函數(shù)來源:https://github.com/tidyverse/ggplot2/blob/eecc450f7f13c5144069705ef22feefe0b8f53f7/R/geom-violin.r#L102
GeomSplitViolin <- ggproto("GeomSplitViolin", GeomViolin,
draw_group = function(self, data, ..., draw_quantiles = NULL) {
data <- transform(data, xminv = x - violinwidth * (x - xmin), xmaxv = x + violinwidth * (xmax - x))
grp <- data[1, "group"]
newdata <- plyr::arrange(transform(data, x = if (grp %% 2 == 1) xminv else xmaxv), if (grp %% 2 == 1) y else -y)
newdata <- rbind(newdata[1, ], newdata, newdata[nrow(newdata), ], newdata[1, ])
newdata[c(1, nrow(newdata) - 1, nrow(newdata)), "x"] <- round(newdata[1, "x"])
if (length(draw_quantiles) > 0 & !scales::zero_range(range(data$y))) {
stopifnot(all(draw_quantiles >= 0), all(draw_quantiles <=
1))
quantiles <- ggplot2:::create_quantile_segment_frame(data, draw_quantiles)
aesthetics <- data[rep(1, nrow(quantiles)), setdiff(names(data), c("x", "y")), drop = FALSE]
aesthetics$alpha <- rep(1, nrow(quantiles))
both <- cbind(quantiles, aesthetics)
quantile_grob <- GeomPath$draw_panel(both, ...)
ggplot2:::ggname("geom_split_violin", grid::grobTree(GeomPolygon$draw_panel(newdata, ...), quantile_grob))
}
else {
ggplot2:::ggname("geom_split_violin", GeomPolygon$draw_panel(newdata, ...))
}
})
geom_split_violin <- function(mapping = NULL, data = NULL, stat = "ydensity", position = "identity", ...,
draw_quantiles = NULL, trim = TRUE, scale = "area", na.rm = FALSE,
show.legend = NA, inherit.aes = TRUE) {
layer(data = data, mapping = mapping, stat = stat, geom = GeomSplitViolin,
position = position, show.legend = show.legend, inherit.aes = inherit.aes,
params = list(trim = trim, scale = scale, draw_quantiles = draw_quantiles, na.rm = na.rm, ...))
}
ggplot(data_new, aes(x = Genes,y = Values, fill = group))+
geom_split_violin(trim = T,colour="white")+
geom_point(stat = 'summary',fun=mean,
position = position_dodge(width = 0.2))+
scale_fill_manual(values = c("#1ba7b3","#dfb424"))+
stat_summary(fun.min = function(x){quantile(x)[2]},
fun.max = function(x){quantile(x)[4]},
geom = 'errorbar',color='black',
width=0.01,size=0.5,
position = position_dodge(width = 0.2))+
stat_compare_means(data = data_new, aes(x = Genes,y = Values),
# 修改顯著性標(biāo)注:
symnum.args=list(cutpoints = c(0, 0.001, 0.01, 0.05, 1),
symbols = c("***", "**", "*", "-")),
label = "p.signif",
label.y = max(data_new$Values),
hide.ns = F)+
theme_bw()+
xlab("")+
ylab("log2(CPM)")+
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
#legend.key = element_rect(fill = c("#1ba7b3","#dfb424")),
legend.justification = "right")
ggsave("violin_plot2.pdf", height = 5, width = 10)
結(jié)果展示
示例數(shù)據(jù)和代碼獲取
以上就是本期的全部內(nèi)容啦性置!歡迎點(diǎn)贊,點(diǎn)在看揍堰!師兄會盡快更新哦鹏浅!制作不易,你的打賞將成為師兄繼續(xù)更新的十足動力个榕!
往期文章
- 跟著Nature Medicine學(xué)作圖--箱線圖+散點(diǎn)圖
- 跟著Nature Communications學(xué)作圖--漸變火山圖
- 跟著Nature Communications學(xué)作圖--氣泡圖+相關(guān)性熱圖
- 跟著Nature Communications學(xué)作圖 -- 復(fù)雜提琴圖
- 跟著Nature Communications學(xué)作圖 -- 復(fù)雜熱圖
- 跟著Nature Communications學(xué)作圖--復(fù)雜散點(diǎn)圖
- 跟著Nature Communications學(xué)作圖 -- 復(fù)雜百分比柱狀圖