❝本节来介绍如何灵活使用「rstatix」来进行数据统计分析,以前都是使用R内置数据来进行展示这次使用论文中的原始数据来进行展示;
加载R包
代码语言:javascript复制library(tidyverse)
library(rstatix)
library(ggpubr)
library(ggprism)
数据清洗
代码语言:javascript复制df <- read_tsv('F1-a.txt') %>% pivot_longer(-`MUFA-PI / total PI [%]`) %>%
group_by(name) %>%
summarise(value_mean=mean(value),sd=sd(value),se=sd(value)/sqrt(n())) %>%
mutate(group=case_when(name=="w/o" ~ "black",name=="TNFα" ~ "grey",
name=="I3M" ~ "grey",TRUE ~ "green"))
代码语言:javascript复制# A tibble: 10 × 5
name value_mean sd se group
<fct> <dbl> <dbl> <dbl> <chr>
1 CHX 19.6 2.45 1.22 green
2 ETO 32.7 2.46 1.23 green
3 I3M 50.1 3.51 1.76 grey
4 MC 39.8 4.85 2.43 green
5 Serum 22.4 3.36 1.68 green
6 STS 21.7 1.16 0.578 green
7 TNFα 52.7 5.97 2.98 grey
8 TPG 35.1 11.0 5.49 green
9 VAL 31.1 5.60 2.80 green
10 w/o 45.1 3.29 1.65 black
定义因子
代码语言:javascript复制df$name <- factor(df$name,levels=read_tsv('F1-a.txt') %>%
pivot_longer(-`MUFA-PI / total PI [%]`) %>%
select(name) %>% distinct() %>% pull())
统计分析
代码语言:javascript复制❝此处还是通过整合数据得到最高点位置信息,如果使用代码生成则会存在问题 ❞
stat.test <- read_tsv('F1-a.txt') %>% pivot_longer(-`MUFA-PI / total PI [%]`) %>%
t_test(data =., value ~ name, ref.group = "w/o") %>%
mutate(p.adj.signif = replace_na(p.adj.signif,""),across("p.adj.signif",str_replace,"ns","")) %>%
select(group1,group2,p.adj,p.adj.signif) %>%
left_join(.,df,by=c("group2"="name")) %>%
mutate(y.position=value_mean se 0.3)
定义主题
代码语言:javascript复制theme_niwot <- function(){
theme_test()
theme(axis.title = element_blank(),
axis.ticks.x=element_blank(),
panel.grid.major.y = element_line(color = "#DAE1E7"),
panel.grid.major.x = element_blank(),
plot.margin = unit(rep(0.2,4),"cm"),
axis.text = element_text(size = 10, color = "#22292F"),
axis.text.y = element_text(margin = margin(r = 5)),
axis.text.x = element_text(margin = margin(t = 5)),
legend.position = "non")
}
数据可视化
代码语言:javascript复制df %>% ggplot(.,aes(name,value_mean))
geom_errorbar(aes(ymax = value_mean se, ymin = value_mean - se),width = 0.2,color = "grey30")
geom_col(width=0.5,aes(fill=group),color="grey50")
add_pvalue(stat.test,label = "p.adj.signif",label.size=6,
coord.flip = TRUE, remove.bracket = TRUE)
scale_y_continuous(expand=c(0,0),limits = c(0,57))
theme_niwot()
scale_fill_brewer(palette="Blues")