regtestMC:

In musto101/wilcox_R: What the package does (short line)

Usage

regtestMC(x, y, regfun = tsreg, nboot = 600, alpha = 0.05, plotit = TRUE, grp = c(1:ncol(x)), nullvec = c(rep(0, length(grp))), xout = FALSE, outfun = outpro, SEED = TRUE, pr = TRUE, ...)

Arguments

x
y
regfun
nboot
alpha
plotit
grp
nullvec
xout
outfun
SEED
pr
...

Examples

##---- 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, y, regfun = tsreg, nboot = 600, alpha = 0.05, plotit = TRUE, 
    grp = c(1:ncol(x)), nullvec = c(rep(0, length(grp))), xout = FALSE, 
    outfun = outpro, SEED = TRUE, pr = TRUE, ...) 
{
    library(parallel)
    x <- as.matrix(x)
    p1 <- ncol(x) + 1
    p <- ncol(x)
    xy <- cbind(x, y)
    xy <- elimna(xy)
    x <- xy[, 1:p]
    y <- xy[, p1]
    if (xout) {
        if (pr) 
            print("Default for outfun is now outpro")
        m <- cbind(x, y)
        flag <- outfun(x, plotit = FALSE, ...)$keep
        m <- m[flag, ]
        x <- m[, 1:p]
        y <- m[, p1]
    }
    x <- as.matrix(x)
    if (length(grp) != length(nullvec)) 
        stop("The arguments grp and nullvec must have the same length.")
    if (SEED) 
        set.seed(2)
    data <- matrix(sample(length(y), size = length(y) * nboot, 
        replace = TRUE), nrow = nboot)
    data = listm(t(data))
    bvec = mclapply(data, regbootMC, x, y, regfun, mc.preschedule = TRUE)
    bvec = matl(bvec)
    grp <- grp + 1
    est <- regfun(x, y)$coef
    estsub <- est[grp]
    bsub <- t(bvec[grp, ])
    if (length(grp) == 1) {
        m1 <- sum((bvec[grp, ] - est)^2)/(length(y) - 1)
        dis <- (bsub - estsub)^2/m1
    }
    if (length(grp) > 1) {
        mvec <- apply(bsub, 2, FUN = mean)
        m1 <- var(t(t(bsub) - mvec + estsub))
        dis <- mahalanobis(bsub, estsub, m1)
    }
    dis2 <- order(dis)
    dis <- sort(dis)
    critn <- floor((1 - alpha) * nboot)
    crit <- dis[critn]
    test <- mahalanobis(t(estsub), nullvec, m1)
    sig.level <- 1 - sum(test > dis)/nboot
    if (length(grp) == 2 && plotit) {
        plot(bsub, xlab = "Parameter 1", ylab = "Parameter 2")
        points(nullvec[1], nullvec[2], pch = 0)
        xx <- bsub[dis2[1:critn], ]
        xord <- order(xx[, 1])
        xx <- xx[xord, ]
        temp <- chull(xx)
        lines(xx[temp, ])
        lines(xx[c(temp[1], temp[length(temp)]), ])
    }
    list(test = test, crit = crit, sig.level = sig.level, nullvec = nullvec, 
        est = estsub)
  }

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