hochberg:

Usage Arguments Examples

Usage

1
hochberg(x, x2 = NA, cil = NA, crit = NA, con = 0, tr = 0.2, alpha = 0.05, iter = 10000, SEED = TRUE)

Arguments

x
x2
cil
crit
con
tr
alpha
iter
SEED

Examples

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
##---- Should be DIRECTLY executable !! ----
##-- ==>  Define data, use random,
##--	or do  help(data=index)  for the standard data sets.

## The function is currently defined as
function (x, x2 = NA, cil = NA, crit = NA, con = 0, tr = 0.2, 
    alpha = 0.05, iter = 10000, SEED = TRUE) 
{
    x3 <- x2
    if (is.matrix(x)) 
        x <- listm(x)
    if (!is.list(x)) 
        stop("Data must be stored in list mode or in matrix mode.")
    J <- length(x)
    tempn <- 0
    svec <- NA
    for (j in 1:J) {
        temp <- x[[j]]
        temp <- temp[!is.na(temp)]
        tempn[j] <- length(temp)
        x[[j]] <- temp
        svec[j] <- winvar(temp, tr = tr)/(1 - 2 * tr)^2
    }
    tempt <- floor((1 - 2 * tr) * tempn)
    A <- sum(1/(tempt - 1))
    df <- J/A
    print(paste("If using the tables of Studentized range distribution,"))
    print(paste("the degrees of freedom are:", df))
    if (!is.list(x2) && !is.matrix(x2)) {
        x2 <- list()
        for (j in 1:J) x2[[j]] <- NA
    }
    if (is.na(cil)) 
        stop("To proceed, you must specify the maximum length of the confidence intervals.")
    if (is.na(crit)) {
        print("Approximating critical value")
        crit <- trange(tempn - 1, alpha = alpha, iter = iter, 
            SEED = SEED)
        print(paste("The critical value is ", crit))
    }
    if (con[1] == 0) {
        Jm <- J - 1
        ncon <- (J^2 - J)/2
        con <- matrix(0, J, ncon)
        id <- 0
        for (j in 1:Jm) {
            jp <- j + 1
            for (k in jp:J) {
                id <- id + 1
                con[j, id] <- 1
                con[k, id] <- 0 - 1
            }
        }
    }
    ncon <- ncol(con)
    avec <- NA
    for (i in 1:ncon) {
        temp <- con[, i]
        avec[i] <- sum(temp[temp > 0])
    }
    dvec <- (cil/(2 * crit * avec))^2
    d <- max(dvec)
    n.vec <- NA
    for (j in 1:J) {
        n.vec[j] <- max(tempn[j], floor(svec[j]/d) + 1)
        print(paste("Need an additional ", n.vec[j] - tempn[j], 
            " observations for group", j))
    }
    if (is.matrix(x2)) 
        x2 <- listm(x2)
    temp2 <- n.vec - tempn
    if (!is.list(x3) && !is.matrix(x3) && sum(temp2) > 0) 
        stop("No second stage data supplied; this function is terminating")
    if (length(x) != length(x2)) 
        warning("Number of groups in first stage data does not match the number in the second stage.")
    ci.mat <- NA
    if (!is.na(x2[1]) || sum(temp2) == 0) {
        xtil <- NA
        nvec2 <- NA
        for (j in 1:J) {
            nvec2[j] <- 0
            temp <- x2[[j]]
            if (!is.na(temp[1])) 
                nvec2[j] <- length(x2[[j]])
            if (nvec2[j] < n.vec[j] - tempn[j]) 
                warning(paste("The required number of observations for group", 
                  j, " in the second stage is ", n.vec[j] - tempn[j], 
                  " but only ", nvec2[j], " are available"))
            xtil[j] <- mean(c(x[[j]], x2[[j]]), tr = tr, na.rm = TRUE)
        }
        ci.mat <- matrix(0, ncol = 3, nrow = ncon)
        dimnames(ci.mat) <- list(NULL, c("con.num", "ci.low", 
            "ci.high"))
        for (ic in 1:ncon) {
            ci.mat[ic, 1] <- ic
            bvec <- con[, ic] * sqrt(svec/(nvec2 + tempn))
            A <- sum(bvec[bvec > 0])
            C <- 0 - sum(bvec[bvec < 0])
            D <- max(A, C)
            ci.mat[ic, 2] <- sum(con[, ic] * xtil) - crit * D
            ci.mat[ic, 3] <- sum(con[, ic] * xtil) + crit * D
        }
    }
    list(ci.mat = ci.mat, con = con)
  }

musto101/wilcox_R documentation built on May 23, 2019, 10:52 a.m.