R/median.is.R

In cwendorf/dgb: DGB

# DGB
## Median Difference from Two (Independent) Samples

ci.median.is <- function(alpha, y1, y2) {
 # Computes confidence interval for a difference of
 # population medians in a 2-group design
 # Arguments:
 #   alpha: alpha level for 1-alpha confidence
 #   y1:    vector of scores for group 1
 #   y2:    vector of scores for group 2  
 # Values:
 #   confidence interval
 z <- qnorm(1 - alpha/2)
 n1 <- length(y1)
 y1 <- sort(y1)
 n2 <- length(y2)
 y2 <- sort(y2)
 median1 <- median(y1)
 median2 <- median(y2)
 a1 <- round((n1 + 1)/2 - sqrt(n1))
 if (a1 < 1) {a1 = 1}
 L1 <- y1[a1]
 U1 <- y1[n1 - a1 + 1]
 p <- pbinom(a1 - 1, size = n1, prob = .5)
 z0 <- qnorm(1 - p)
 se1 <- (U1 - L1)/(2*z0)
 a2 <- round((n2 + 1)/2 - sqrt(n2))
 if (a2 < 1) {a2 = 1}
 L2 <- y2[a2]
 U2 <- y2[n2 - a2 + 1]
 p <- pbinom(a2 - 1, size = n2, prob = .5)
 z0 <- qnorm(1 - p)
 se2 <- (U2 - L2)/(2*z0)
 diff <- median1 - median2
 se <- sqrt(se1^2 + se2^2)
 L <- diff - z*se
 U <- diff + z*se
 out <- t(c(median1, median2, diff, se, L, U))
 colnames(out) <- c("Median1", "Median2", "Median1-Median2", "SE", "LL", "UL")
 return(out)
}