ixxmu / mp_duty

抓取网络文章到github issues保存
https://archives.duty-machine.now.sh/
110 stars 30 forks source link

绘图练习 | 分组散点图 / 蜜蜂图 #5560

Closed ixxmu closed 3 weeks ago

ixxmu commented 3 weeks ago

https://mp.weixin.qq.com/s/z044AxGC12Ax9qq2fr4Djg

ixxmu commented 3 weeks ago

绘图练习 | 分组散点图 / 蜜蜂图 by 被炸熟的虾

昨天我们记录了ggplot2包的position参数使用方法(对下图感兴趣但没看昨天笔记的建议了解一下),今天借助一个图片复现来练习一下:

Gut microbiota modulates weight gain in mice after discontinued smoke exposure. Nature. 2021 Dec; PMID: 34880502.

绘图思路:

1、为了避免遮盖,第一步先添加灰色背景;

2分组散点图:使用geom_point函数绘制,既要分组也要设置随机抖动,所以需要设置position_jitterdodge()函数;

3均值\中位数横线:使用geom_crossbar()函数绘制,注意转换函数设置。只需要分组,设置position_dodge()就可以。

全文画图数据都可以在附表中找到(Fig1B

Step1、数据整理与添加背景

library(readxl)
exprSet <- as.data.frame(read_excel("./Data2.xlsx",sheet = 1))
exprSet2 <- reshape2::melt(exprSet)
head(exprSet2)
  # Group variable     value
# 1    NS  Smoking 0.4837348
# 2    NS  Smoking 0.2165572
# 3    NS  Smoking 0.5785105
# 4    NS  Smoking 0.7888895
# 5    NS  Smoking 0.5367586
# 6    NS  Smoking 0.4482413
exprSet2$Group <- factor(exprSet2$Group,
                         levels = c("NS","SMK","NS+abx","SMK+abx"),
                         labels = c("Non-\nSMK","SMK","Non-\nSMK+\nabx","SMK+\nabx"))       
colnames(exprSet2) <- c("Group2","Group","Value")
library(ggplot2)
mytheme <- theme_classic() +
           theme(legend.title = element_blank(),
                 legend.text = element_text(size = 20),
      legend.position = "top",
      plot.margin = margin(t = 5,  # 顶部边缘距离
                                      r = 10,  # 右边边缘距离
                                      b = 5,  # 底部边缘距离
                                      l = 5), # 左边边缘距离
                 axis.text = element_text(size = 20),
      axis.title = element_text(size = 20))

生成一个空图,添加分类背景图层

library(ggplot2)
p1 <- ggplot(exprSet2,aes(x = Group,y = Value)) +
      scale_y_continuous(expand = c(0,0),limits = c(-1,4),breaks = seq(-1,4,1)) + 
      labs(x = NULL,y = "Weight Change(%)",title = NULL) + 
      mytheme 
p2 <- p1 + 
      annotate("rect", xmin = 1.5, xmax = Inf, fill = "grey",
               ymin = -1, ymax = Inf, alpha = .9

此时的背景位置有点偏,不用担心

Step2、添加图形主体

添加散点:

p3 <- p2 + 
      # 添加散点
      geom_point(aes(fill = Group2),
                 show.legend = F,color = "black",
                 position = position_jitterdodge(seed = 12345,
                                      dodge.width = 0.8,
                                                 jitter.width = 0.15),
                 shape = 21,size = 4)

均值\中位数横线:

p4 <- p3 + 
      # 均值\中位数横线
    geom_crossbar(aes(fill = Group2),
                    stat = "summary",
                  fun = median,
                  position = position_dodge(width = 0.8),
                    size = 0.5,width = 0.5,
                    color = 'black')
      # stat_summary(fun = median,
                 # position = position_dodge(width = 0.8),
                   # geom = "crossbar",
                   # size = 0.5,width = 0.5,
                   # color = 'black')

修改颜色并添加额外水平线:

p5 <- p4 + 
      # 修改颜色
      scale_color_manual(values = c("#4d97cd","#db6968","#008343","#e8c559")) +
      scale_fill_manual(values = c("#4d97cd","#db6968","#008343","#e8c559")) + 
    # 添加额外水平线
      geom_hline(yintercept = 0,
                 size = 1,color = "black",lty = "dashed")

Step3、添加显著性检验

library(rstatix)
library(tidyverse)
stat.test <- exprSet2 %>% 
             group_by(Group) %>% 
             t_test(Value ~ Group2) %>%
             adjust_pvalue() %>% 
             add_significance("p.adj") %>% 
             add_xy_position(x = "Group")
data.frame(stat.test)
stat.test$y.position <- stat.test$y.position - 0.5
p5 + stat_pvalue_manual(stat.test,label = "p.adj.signif",hide.ns = T,#label = "p" or label = "p.adj"
                        tip.length = 0.01,label.size = 8)

手动添加(因为作者只标注了部分):

tmp <- stat.test[c(1,6,7,11),]
data.frame(tmp)
#      Group   .y.          group1    group2 n1 n2 statistic       df        p
#1   Smoking Value       Non-\nSMK       SMK 37 38  8.657194 52.58912 1.08e-11
#2   Smoking Value Non-\nSMK+\nabx SMK+\nabx 40 39  9.813350 55.89030 9.20e-14
#3 Cessation Value       Non-\nSMK       SMK 33 34 -9.476267 62.06300 1.13e-13
#4 Cessation Value             SMK SMK+\nabx 34 34  5.275416 59.57884 1.94e-06
#      p.adj p.adj.signif y.position                     groups x xmin xmax
#1 8.640e-11         ****   2.062760             Non-\nSMK, SMK 1  0.7  0.9
#2 1.104e-12         ****   3.627320 Non-\nSMK+\nabx, SMK+\nabx 1  1.1  1.3
#3 1.130e-12         ****   2.567760             Non-\nSMK, SMK 2  1.7  1.9
#4 1.164e-05         ****   3.819408             SMK, SMK+\nabx 2  1.9  2.3
library(ggsignif)
p5 + geom_signif(annotations = tmp$p.adj.signif,
                 y_position = c(2.5,2.8,3.1,3.6),
                 xmin = tmp$xmin,
                 xmax = tmp$xmax,
                 size = 1.2,textsize = 8,
                 tip_length = 0)


被炸熟的虾
自己的摸索,发现问题麻烦告诉作者,光速回来改正