R：在箱线图/GGPLOT 上添加 p 值（emmeans 的 p 值与“emmeans”函数形成对比）答案

【问题标题】：R：在箱线图/GGPLOT 上添加 p 值（emmeans 的 p 值与“emmeans”函数形成对比）
【发布时间】：2021-11-20 14:31:42
【问题描述】：

我拟合了一个线性混合模型（带有 lmer 函数）。 “比例”作为我的响应变量，“阶段”作为固定因子（有 5 个水平），主题作为随机因子。

LMM.fit<-lmer(Proportion~Stage+(1|Subject),data=tab)

       Chisq Df Pr(>Chisq)    
Stage 290.07  4  < 2.2e-16 ***

使用 emmeans 函数进行事后成对比较，以确定不同阶段发生的差异：

emmeans(LMM.fit, pairwise ~ Stage)$contrasts

 contrast        estimate     SE     df t.ratio p.value
 Stage1 - Stage2    0.178 0.0505 105087   3.528  0.0038
 Stage1 - Stage3    0.601 0.0363 565787  16.529  <.0001
 Stage1 - Stage4    0.438 0.0496 224258   8.833  <.0001
 Stage1 - Stage5    0.272 0.0497 295762   5.472  <.0001
 Stage2 - Stage3    0.422 0.0514 103631   8.216  <.0001
 Stage2 - Stage4    0.260 0.0615  95267   4.224  0.0002
 Stage2 - Stage5    0.094 0.0620 107696   1.517  0.5516
 Stage3 - Stage4   -0.163 0.0496 229001  -3.279  0.0092
 Stage3 - Stage5   -0.328 0.0491 314661  -6.683  <.0001
 Stage4 - Stage5   -0.166 0.0599 190576  -2.769  0.0446

Degrees-of-freedom method: kenward-roger 
P value adjustment: tukey method for comparing a family of 5 estimates

我基本上想在上面显示的箱线图上添加 emmeans 结果中显示的 p 值（在同一图中两个两个之间的所有组之间）。我知道有函数stat_pvalue_manual()，但我很想知道如何将它与 emmeans contrasts output

一起使用

【问题讨论】：

标签： r plot boxplot p-value emmeans

【解决方案1】：

这就是我设法做到的：

source("http://news.mrdwab.com/install_github.R")
install_github("mrdwab/overflow-mrdwab") 
install_github("mrdwab/SOfun")
library(SOfun)

LMM.fit<-lmer(Proportion~Stage+(1|Subject),data=tab)

p.val.test<-pwpm(emmeans(LMM.fit,  ~ Stage),means = FALSE, flip = TRUE,reverse = TRUE)
p.val.test<-sub("[<>]", "", p.val.test)
p.matx<-matrix(as.numeric((p.val.test)),nrow = length(p.val.test[,1]),ncol = length(p.val.test[,1])) #if your factor has 5 levels ncol and nrow=5
rownames(p.matx) <- colnames(p.matx) <-colnames(p.val.test)
p.matx[upper.tri(p.matx, diag=FALSE)] <- NA
stat.test<-subset(melt(p.matx),!is.na(value))
names(stat.test)<-c("group1","group2","p.adj")
stat.test[stat.test$p.adj<=0.001,"p.adj.signif"]<-"***"
stat.test[stat.test$p.adj>0.001 & stat.test$p.adj<=0.01,"p.adj.signif"]<-"**"
stat.test[stat.test$p.adj>0.01 & stat.test$p.adj<=0.05,"p.adj.signif"]<-"*"
stat.test[ stat.test$p.adj>0.05,"p.adj.signif"]<-"ns"
stat.test<-mc_tribble(stat.test) # copy & paste the result of this line before the line of ggboxplot!!

stat.test <- tribble(
        ~group1, ~group2, ~p.adj, ~p.adj.signif,
        "Stage2","Stage1",0.0038,"**",
        "Stage3","Stage1",1e-04,"***",
        "Stage4","Stage1",1e-04,"***",
        "Stage5","Stage1",1e-04,"***",
        "Stage3","Stage2",1e-04,"***",
        "Stage4","Stage2",2e-04,"***",
        "Stage5","Stage2",0.5516,"ns",
        "Stage4","Stage3",0.0092,"**",
        "Stage5","Stage3",1e-04,"***",
        "Stage5","Stage4",0.0446,"*")

bxp<-ggboxplot(tab, x = "Stage", y = "Proportion",color = "Stage", palette = "jco",add = "jitter",bxp.errorbar=T) 
bxp + stat_pvalue_manual(stat.test,
                         y.position = max(tab$Proportion)+sd(tab$Proportion),
                         color = "midnightblue")

PS：this website 帮助（“GGPUBR：如何将在其他地方生成的 P 值添加到 GGPLOT”）。

【讨论】：