I am using R to add some p-values to boxplots. I would like to specifically compare pairs of boxplots using a wilcoxon rank sum test. I tried using the function stat_compare_means() from ggpubr but unfortunately there is no option to add adjusted p-values. I searched online for a bit and noticed a function from the same package called stat_pvalue_manual() which would allow me to put the p-values in. However, I have not been able to figure out the package.

I am using the tooth growth data from R as an example:

ToothGrowth$dose=as.factor(ToothGrowth$dose)

I have the following boxplot:

ToothGrowth%>%ggplot(aes(x=dose, y=len, fill=supp))+
          geom_boxplot()

I would like to do a wilcoxon sum rank test on the data and label each comparison between OJ and VC.

Below is the stats test, which when I check it corresponds to the three comparisons I want between OJ and VC:

 stat.test <- ToothGrowth%>%group_by(dose)%>%wilcox_test(len~supp, p.adjust.method = 'BH')%>%
          mutate(y.position = c(29, 35, 39))

However, when I try to add this to my boxplot I get an error:

ToothGrowth%>%ggplot(aes(x=dose, y=len, fill=supp))+
          geom_boxplot()+
          stat_pvalue_manual(stat.test, label = 'p')

Error in FUN(X[[i]], ...) : object 'supp' not found

What I think is happening is that it is using the aes() from ggplot, but when I set inherit.aes=false within the stat_pvalue_manual function OJ and VC to the x-axis and does a weird comparison between the two groups instead of the other boxplots. I have been trying for hours to figure out what I am doing wrong but I have had no luck so far.