Nothing
#' R function for visually displaying Mann-Whitney test's results
#'
#' The function allows to perform Mann-Whitney test, and to display the test's results in a plot
#' along with two boxplots. For information about the test, and on what it is actually testing, see
#' for instance the interesting article by R M Conroy, "What hypotheses do "nonparametric" two-group
#' tests actually test?", in The Stata Journal 12 (2012): 1-9.\cr
#'
#' The returned boxplots display the
#' distribution of the values of the two samples, and jittered points represent the individual
#' observations.
#'
#' At the bottom of the chart, a subtitle arranged on three lines reports relevant
#' statistics:\cr -test statistic (namely, U) and the associated z and p value;\cr -Probability of
#' Superiority value (which can be interpreted as an effect-size measure, as discussed in:
#' https://nickredfern.wordpress.com/2011/05/12/the-mann-whitney-u-test/);\cr -another measure of
#' effect size, namely r (see
#' https://stats.stackexchange.com/questions/124501/mann-whitney-u-test-confidence-interval-for-effect-size),
#' whose thresholds are indicated in the last line of the plot's subtitle.\cr
#'
#' The function may also
#' return a density plot (coupled with a rug plot at the bottom of the same chart) that displays the
#' distribution of the pairwise differences between the values of the two samples being compared.
#' The median of this distribution (which is represented by a blue reference line in the same chart)
#' corresponds to the Hodges-Lehmann estimator.
#'
#' @param x Object storing the values of the first group being compared.
#' @param y Object storing either the values of the second group being compared or a grouping
#' variable with 2 levels.
#' @param xlabl If y is not a grouping variable, user may want to specify here the name of the x
#' group that will show up in the returned boxplots (default is "x").
#' @param ylabl If y is not a grouping variable, user may want to specify here the name of the y
#' group that will show up in the returned boxplots (default is "y").
#' @param strip Logical value which takes FALSE (by default) or TRUE if the user wants jittered
#' points to represent individual values.
#' @param notch Logical value which takes FALSE (by default) or TRUE if user does not or do want to
#' have notched boxplots in the final display, respectively; it is worth noting that overlapping
#' of notches indicates a not significant difference at about 95 percent confidence.
#' @param omm It stands for overall mean and median; takes FALSE (by default) or TRUE if user
#' wants the mean and median of the overall sample plotted in the chart (as a dashed RED line and
#' dotted BLUE line respectively).
#' @param outl Logical value which takes FALSE or TRUE (by default) if users want the boxplots to
#' display outlying values.
#' @param HL Logical value that takes TRUE or FALSE (default) if the user wants to display the
#' distribution of the pairwise differences between the values of the two samples being compared;
#' the median of that distribution is the Hodges-Lehmann estimator.
#'
#' @keywords mwPlot
#'
#' @export
#'
#' @importFrom plyr count
#' @importFrom coin wilcox_test
#'
#' @examples
#' #create a toy dataset
#' mydata <- data.frame(values=c(rnorm(30, 100,10),rnorm(30, 80,10)),
#' group = as.factor(gl(2, 30, labels = c("A", "B"))))
#'
#' # performs the test, displays the test's result, including jittered points, notches,
#' #overall median and mean, and the Hodges-Lehmann estimator
#' mwPlot(x=mydata$values, y=mydata$group, strip=TRUE, omm=TRUE, notch=TRUE, HL=TRUE)
#'
mwPlot <- function (x,y,xlabl="x",ylabl="y", strip=FALSE,notch=FALSE,omm=FALSE, outl=TRUE, HL=FALSE){
options(scipen=999)
if (is.numeric(y)==FALSE) {
data <- data.frame(value=x, group=y)
} else {data <- data.frame(value=c(x,y), group=c(rep(xlabl, length(x)), rep(ylabl, length(y))))
}
res <- wilcox.test(data[,1] ~ data[,2], conf.int=TRUE)
U <- wilcox.test(data[,1] ~ data[,2])$statistic
p <- ifelse(res$p.value < 0.001, "< 0.001", ifelse(res$p.value < 0.01, "< 0.01", ifelse(res$p.value < 0.05, "< 0.05",round(res$p.value, 3))))
print(paste("p-value=",res$p.value))
samples.size <- plyr::count (data[,2]) #requires the plyr package
PS <- round(U/(samples.size[1,2] * samples.size[2,2]), 3)
z <- round(wilcox_test(data[,1] ~ data[,2])@statistic@teststatistic,3) #requires the coin package
r <- round(abs(z/sqrt(samples.size[1,2] + samples.size[2,2])), 3)
graphics::boxplot(data[,1] ~ data[,2], data = data, notch = notch, outline=outl)
chart.title="Box Plots"
if (strip==TRUE) {
stripchart(data[,1] ~ data[,2], vertical = TRUE, data = data, method = "jitter", add = TRUE, pch = 16, col="#00000088", cex = 0.5)
chart.title="Jittered Box Plots"
} else {
}
title(main=chart.title, sub=paste("Mann-Whitney U=", U, ", z=",z, ", p=", p, "; Probability of Superiority=", PS, "; r=", r, "\nP{value(group to the left) > value(group to the right)}=", PS,"\nEffect size thresholds [r]: small (0.10), medium (0.30), large (0.50)"), cex.sub=0.8)
if (omm==TRUE) {
abline(h=mean(data[,1]), lty=2, col="red")
abline(h=median(data[,1]), lty=3, col="blue")
} else {
}
if (HL==TRUE) {
unstacked.data <- utils::unstack(data)
diff <- outer(unstacked.data[[1]], unstacked.data[[2]],"-")
m <- round(median(diff), 3)
graphics::plot(stats::density(diff), main="Pairwise differences distribution", xlab="", sub=paste("difference in location", "\nmedian (Hodges-Lehmann estimator):", m))
polygon(stats::density(diff), col="grey")
rug(diff, col="red")
abline(v=m, lty=2, col="blue")
} else {
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.