前往小程序,Get更优阅读体验!
立即前往
首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >上手即用,分组统计检验直方图的绘图脚本分享

上手即用,分组统计检验直方图的绘图脚本分享

作者头像
用户1075469
发布2022-05-25 14:55:41
5540
发布2022-05-25 14:55:41
举报
文章被收录于专栏:科技记者

发文章,写论文,分组统计检验直方图是最常见和最实用的,你是否还在烦恼如果把图画好,帮你解决困难啦!这里分享下同事新鲜写就的绘图脚本,自带了示例数据,可以一键出图,助力你的科研和学习。

脚本和文件在这里:

链接: https://pan.baidu.com/s/1A5Opyo-rB4iXpxcdqahNQg?pwd=wmbd 提取码: wmbd 复制这段内容后打开百度网盘手机App,操作更方便哦

先来看下他提供的示例数据长啥样:

准备起来还是比较容易的,就6列。

再来看下他的脚本:

三个函数组成,代码写的还是比较规范的,管道用的66的,虽然赋值符合和我不是一个风格,这个萝卜青菜啦!

代码语言:javascript
复制
# 读入表达矩阵
readFlie=function(input,type,row=T,header=T){
  # input 为读入文件的路径,type为读入文件的类型,格式为‘.txt’或‘.csv’,row=T,将文件的第一列设置为列名
  library(data.table,quietly = TRUE)
  if(type=='txt'){
    dat = fread(input,header = header,sep='\t',stringsAsFactors = F,check.names = F)
    if(row){
      dat = as.data.frame(dat,stringsAsFactors = F)
      rownames(dat) = dat[,1]
      dat = dat[,-1]
    }else{
      dat = as.data.frame(dat,stringsAsFactors = F)
    }
  }else{
    dat = fread(input,header = header,sep=',',stringsAsFactors = F,check.names = F)
    if(row){
      dat = as.data.frame(dat,stringsAsFactors = F)
      rownames(dat) = dat[,1]
      dat = dat[,-1]
    }else{
      dat = as.data.frame(dat,stringsAsFactors = F)
    }
  }
  return(dat)
}

# 保存图片,只支持ggplot对象
savePlots=function(path,plot,type=c('pdf','png','tiff','svg')[1],width=10,height=8,dpi=300){
  # path表示保存图片路径,需要加上相应的文件扩展名称,plot表示ggplot对象,type表示要保存的图片格式
  if(type=='pdf'){
    ggplot2::ggsave(filename = path,plot = plot,width = width,height = height,device = 'pdf')
  }else if(type=='svg'){
    svglite::ggsave(filename = path,plot = plot,width = width,height = height,device = 'svg')
  }else if(type=='png'){
    ggplot2::ggsave(filename = path,plot = plot,width = width,height = height,device = 'png',dpi = dpi)
  }else{
    ggplot2::ggsave(filename = path,plot = plot,width = width,height = height,device = 'tiff',dpi = dpi)
  }
}
# 分组条形图
PlotGroup_barplot=function(data,xlab = '', ylab = 'Value',legend = 'Group',
                             color =NULL,show_compare= T,cmp_show_compare = T,
                             test_method = c('wilcox.test','t.test','kruskal.test','anova')[1],
                             cmp_test_method = c('wilcox.test','t.test')[1],
                             orientation = c('vertical','horizontal','reverse')[1],
                             p_show = c('p.signif','p.format')[1],
                             add = c( 'mean_se', 'mean_sd', 'mean_ci', 'mean_range', 'median','mean', 'median_iqr')[1]){
# data是一个数据框,包含三列,第一列是样本名称(X轴坐标),第二列是样本分组信息,第三列为值(Y轴坐标)
  # xlab,ylab和labs 分别自定义X轴名称,Y轴名称和图例标题名称
  # colors 自定义颜色,默认为NULL,函数会根据RColorBrewer包进行配色
  # show_compare 是否展示同组样本之间所有组别之间差异比较
  # test_method 同组样本之间所有组别之间差异比较方法
  # cmp_show_compare T or F 是否展示同组样本之间两两进行比较,样本分组必须大于2
  # cmp_test_method 同组样本之间两两进行比较方法
  # p_show 图片上差异显著性展示方式,*: p <= 0.05 **: p <= 0.01 ***: p <= 0.001 ****: p <= 0.0001
  # add 在条形图上添加统计值如样本的标准差,标准误,详细参数见https://www.rdocumentation.org/packages/ggpubr/versions/0.4.0/topics/add_summary 
  library(ggplot2)
  library(ggpubr)
  library(rstatix)
  library(dplyr)
  if(is.null(colors)){
    library(RColorBrewer)
    color = c(brewer.pal(9,'Set1'),brewer.pal(8,'Set2'),brewer.pal(12,'Set3'),
              brewer.pal(9,'Pastel1'),brewer.pal(8,'Pastel2'),brewer.pal(8,'Accent'),
              brewer.pal(8,'Dark2'),brewer.pal(12,'Paired')) %>% unique()
  }
  names(data) = c('Sample','Group',"Value")
  # 绘图
  p = ggbarplot(data, x = "Sample", y = "Value",fill = "Group", xlab = xlab,ylab = ylab,add = add,width = 0.5,legend.title = legend,palette = color,position = position_dodge(0.6),ggtheme = theme_bw()) + scale_y_continuous(expand = expansion(mult = c(0, 0.1)))
  if(show_compare){
    p = p + stat_compare_means(aes(group = Group), label = p_show,method = test_method )
  }
  # 组内样本两两之间差异比较
  if(cmp_show_compare & length(unique(data$Group))>2){
    if(cmp_test_method=='wilcox.test'){
      stat.test = data %>%
        group_by(Sample) %>%
        wilcox_test(Value ~ Group,p.adjust.method = "none") 
      stat.test = stat.test %>%
        add_xy_position(x = "Sample", fun = "mean_sd", dodge = 0.8) %>% dplyr::filter(p<0.05)
    }else{
      stat.test = data %>%
        group_by(Sample) %>%
        t_test(Value ~ Group,p.adjust.method = "none") 
      stat.test = stat.test %>%
        add_xy_position(x = "Sample", fun = "mean_sd", dodge = 0.8) %>% dplyr::filter(p<0.05)
    }
    if(p_show == 'p.signif'){
      p = p + stat_pvalue_manual(stat.test, label = "p.adj.signif", tip.length = 0.005,bracket.nudge.y = -2)
    }else{
      p = p + stat_pvalue_manual(stat.test, label = "p.adj", tip.length = 0.005,bracket.nudge.y = -2)
    }
  }
  return(p)
}
df = readFlie(input = 'C:\\Users\\zd200\\Desktop\\不常用\\直方图实例文件.txt',type = 'txt',row = F)
library(tidyr)
library(dplyr)
dt = df %>% pivot_longer(cols = -Group,names_to  = 'gene')
data1 = dt[,c(2,1,3)]
data2 = dt[,c(2,1,3)] %>% filter(Group!='Med')
p1 = PlotGroup_barplot(data = data1,test_method = 'kruskal.test',add = 'mean_sd')
p1
p2 = PlotGroup_barplot(data = data2,test_method = 'wilcox.test',add = 'mean_sd')
p2
savePlots(path = './PlotGroup_barplot1.pdf',plot = p1,type = 'pdf',width = 6,height = 4)
savePlots(path = './PlotGroup_barplot2.pdf',plot = p2,type = 'pdf',width = 6,height = 4)

最后,看下图的效果

配色还不错,颜值在线,你可以根据自己的需求微调。

怎样,快试试吧!

本文参与 腾讯云自媒体同步曝光计划,分享自微信公众号。
原始发表:2022-05-07,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 微因 微信公众号,前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体同步曝光计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
目录
  • 先来看下他提供的示例数据长啥样:
  • 再来看下他的脚本:
  • 最后,看下图的效果
领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档