这图怎么画| 多组箱线图+组间/内差异分析

rm(list = ls())
library(ggplot2) 
library(ggpubr)
library(reshape2) 
library(dplyr) 
# data example
mat <- matrix(runif(40*12, min=0, max=10), nrow=40, ncol=12)
rownames(mat) <- paste0("sample", 1:40)
colnames(mat) <- paste0("gene", 1:12)
mat
gene_group = data.frame(gene = paste0("gene", 1:12),
                        group = rep(paste0("group", 1:6), each = 2))
sample_group = data.frame(sample = paste0("sample", 1:40),
                          group = rep(c('c','w'),20))
# data preprocess
library(reshape2)
mat_long <- melt(mat, varnames = c("sample", "gene"), value.name = "expression")
df1 <- merge(mat_long, gene_group, by = "gene")
df2 <- merge(mat_long, sample_group, by = "sample")
## calculate mean value
df1_mean = df1 %>% 
  group_by(gene) %>% 
  summarise(mean_value=mean(expression))
df1_mean$x = 1:12
df2_mean <- df2 %>%
  group_by(gene, group) %>%
  summarise(mean = mean(expression), sd = sd(expression))
df2_mean$gene <- factor(df2_mean$gene)
# plot
## left boxplot
df1$facet = 'facet1'# for facet title
p1 <- ggplot(df1,aes(gene,expression))+
  geom_boxplot(aes(fill=group,color=group),outlier.shape = 18,size=0.6)+
  stat_boxplot(aes(color=group),geom = "errorbar", width=0.3,size=0.6)+
  geom_segment(df1_mean,
               mapping=aes(x=x-0.5,xend=x+0.5,y=mean_value,yend = mean_value),
               color="white",size=0.5)+
  geom_hline(yintercept = 5, linetype = 2, color = "grey60",linewidth=0.8)+
  ## compare with gene6
  stat_compare_means(label = "p.signif", method = "t.test",
                     ref.group = "gene1",size = 5)+
  coord_flip()+
  scale_fill_manual(values = c("#e3ac6d","#9d7bb8","#6caf83","#d9586e","#3c74bb","#f85b2b"))+
  scale_color_manual(values = c("#e3ac6d","#9d7bb8","#6caf83","#d9586e","#3c74bb","#f85b2b"))+
  theme_bw()+
  theme(legend.position = "none",
        panel.grid = element_blank(),
        axis.text = element_text(color = "black",size=10))+
  labs(y="exp",x='')+
  annotate("rect", xmin = 0, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.3,fill="#EAA9CF") +
  annotate("rect", xmin = 4.5, xmax = 8.5, ymin = -Inf, ymax = Inf, alpha = 0.3,fill="#D6BDFF") +
  annotate("rect", xmin = 8.5, xmax = 13, ymin = -Inf, ymax = Inf, alpha = 0.3,fill="#99DAFF")+
  facet_grid(~ facet)
p1
## right plot
df2_mean$facet = 'facet2' # for facet title
p2 <- ggplot(df2_mean,aes(x =gene, y = mean,
                          group = group,color = group)) +
  geom_point(position = position_dodge(0.8),shape=18,size=3)+
  geom_errorbar(aes(x =gene, y = mean,ymin = mean- sd, ymax = mean + sd,
                                       group = group,color = group), 
                width = 0,position = position_dodge(0.8),linewidth=0.5)+
  scale_color_manual(values = c("#7fc190","#efb684"))+
  stat_compare_means(df2, mapping = aes(gene,expression,group = group),
                     label = "p.signif", method = "t.test",size = 5)+
  coord_flip()+
  theme_bw()+
  theme(legend.position = c(0.9,0.95),
        legend.background = element_blank(),
        legend.key = element_blank(),
        panel.grid = element_blank(),
        axis.text.x = element_text(color = "black",size=10),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank())+
  labs(y='exp2',x='',color='')+
  annotate("rect", xmin = 0, xmax = 4.5, ymin = -Inf, ymax = Inf, alpha = 0.3,fill="#EAA9CF") +
  annotate("rect", xmin = 4.5, xmax = 8.5, ymin = -Inf, ymax = Inf, alpha = 0.3,fill="#D6BDFF") +
  annotate("rect", xmin = 8.5, xmax = 13, ymin = -Inf, ymax = Inf, alpha = 0.3,fill="#99DAFF")+
  facet_grid(~ facet)
p2
#merge
library(patchwork)
library(cowplot)
p1+p2+plot_layout(nrow= 1,
                  width = c(3, 2))
ggsave('boxplot.pdf',width = 10,height = 8)