R/sampleSizeWilcoxTwoSample.R

In felix-hof/biostatUZH: Misc Tools of the Department of Biostatistics, EBPI, University of Zurich

#' Sample size for two sample Wilcoxon (Mann-Whitney) test
#' 
#' Compute the sample size required to test the hypothesis that two sets of independent samples,
#' \eqn{X_1, ..., X_m} and \eqn{Y_1, ..., Y_n} come from the
#' same population against that \eqn{Y}'s is to be larger than \eqn{X}'s.
#'  
#' @param delta Probability \eqn{P(Y_i > X_j)}.
#' @param sig.level Significance level with default 0.05.
#' @param power Power with default 0.9.
#' @param ratio Ratio of the sample size of group 1 to the total samples size:
#' \code{ratio} \eqn{= m / (m + n)}, where \eqn{m} and \eqn{m} are the sample
#' sizes of group 1 and 2, respectively.
#' A \code{ratio} of 0.5 implies equal group sizes.
#' @param two.sided If \code{TRUE}, the test is two-sided.
#' Otherwise, it is one-sided.
#' @return data.frame with columns: 
#' \item{delta}{\code{delta} from input.}
#' \item{sig.level}{\code{sig.level} from input.}
#' \item{power}{\code{power} from input.}
#' \item{ratio}{\code{ratio} from input.}
#' \item{two.sided}{\code{two.sided} from input.}
#' \item{m}{Sample size of the group 1.}
#' \item{n}{Sample size of the group 2.}
#' @author Kaspar Rufibach
#' @references Nother, G.E. (1987). Sample Size Determination for Some Common
#' Nonparametric Tests \emph{JASA}, \bold{82}, 644--647.
#' @keywords htest
#' @examples
#' 
#' sampleSizeWilcoxTwoSample(delta = 0.75)
#' sampleSizeWilcoxTwoSample(delta = 0.85)
#' sampleSizeWilcoxTwoSample(delta = 0.85, power = 0.8)
#' sampleSizeWilcoxTwoSample(delta = 0.85, sig.level = 0.01)
#' sampleSizeWilcoxTwoSample(delta = 0.85, two.sided = FALSE)
#'
#' tmp <- lapply(X = seq(0.55, 0.95, .1),
#'               FUN = function(delta) sampleSizeWilcoxTwoSample(delta = delta, ratio = .25))
#' do.call(what = rbind, args = tmp)
#' 
#' @export
sampleSizeWilcoxTwoSample <- function(delta = 0.75, sig.level = 0.05, power = 0.9,
                                      ratio = 0.5, two.sided = TRUE){
    stopifnot(is.numeric(delta), length(delta) == 1,
              is.finite(delta),
              0 < delta, delta < 1,
              is.numeric(sig.level), length(sig.level) == 1,
              is.finite(sig.level), 
              0 < sig.level, sig.level < 1,
              is.numeric(power), length(power) == 1,
              is.finite(power),
              0 < power, power < 1,
              is.numeric(ratio), length(ratio) == 1,
              is.finite(ratio),
              0 < ratio, ratio < 1,
              is.logical(two.sided), length(two.sided) == 1,
              is.finite(two.sided))

    za <- if(two.sided) qnorm(1 - sig.level / 2) else qnorm(1 - sig.level)
    zb <- qnorm(power)
    N <- (za + zb)^2 / (12 * ratio * (1 - ratio) * (delta - 1/2)^2)    
    data.frame("delta" = delta,
               "sig.level" = sig.level,
               "power" = power,
               "ratio" = ratio,
               "two.sided" = two.sided, 
               "m" = ceiling(ratio * N),
               "n" = ceiling(N * (1 - ratio)))
}