regpre:

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

Usage

regpre(x, y, regfun = lsfit, error = absfun, nboot = 100, adz = TRUE, mval = round(5 * log(length(y))), model = NULL, locfun = mean, pr = FALSE, xout = FALSE, outfun = out, STAND = TRUE, plotit = TRUE, xlab = "Model Number", ylab = "Prediction Error", SEED = TRUE, ...)

Arguments

x
y
regfun
error
nboot
adz
mval
model
locfun
pr
xout
outfun
STAND
plotit
xlab
ylab
SEED
...

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 = lsfit, error = absfun, nboot = 100, 
    adz = TRUE, mval = round(5 * log(length(y))), model = NULL, 
    locfun = mean, pr = FALSE, xout = FALSE, outfun = out, STAND = TRUE, 
    plotit = TRUE, xlab = "Model Number", ylab = "Prediction Error", 
    SEED = TRUE, ...) 
{
    if (pr) {
        print("By default, least squares regression is used, ")
        print("But from Wilcox, R. R. 2008, Journal of Applied Statistics, 35, 1-8")
        print("Setting regfun=tsreg appears to be a better choice for general use.")
        print("That is, replace least squares with the Theil-Sen estimator")
        print("Note: Default for the argument error is now absfun")
        print(" meaning absolute error is used")
        print("To use squared error, set error=sqfun")
    }
    x <- as.matrix(x)
    d <- ncol(x)
    p1 <- d + 1
    temp <- elimna(cbind(x, y))
    x <- temp[, 1:d]
    y <- temp[, d + 1]
    x <- as.matrix(x)
    if (xout) {
        x <- as.matrix(x)
        if (!STAND) 
            flag <- outfun(x, plotit = FALSE, ...)$keep
        if (STAND) 
            flag <- outpro(x, STAND = TRUE, plotit = FALSE)$keep
        x <- x[flag, ]
        y <- y[flag]
        x <- as.matrix(x)
    }
    if (is.null(model)) {
        if (d <= 5) 
            model <- modgen(d, adz = adz)
        if (d > 5) 
            model[[1]] <- c(1:ncol(x))
    }
    mout <- matrix(NA, length(model), 5, dimnames = list(NULL, 
        c("apparent.error", "boot.est", "err.632", "var.used", 
            "rank")))
    if (SEED) 
        set.seed(2)
    data <- matrix(sample(length(y), size = mval * nboot, replace = TRUE), 
        nrow = nboot)
    bid <- apply(data, 1, idb, length(y))
    for (imod in 1:length(model)) {
        nmod = length(model[[imod]]) - 1
        temp = c(nmod:0)
        mout[imod, 4] = sum(model[[imod]] * 10^temp)
        if (sum(model[[imod]] == 0) != 1) {
            xx <- x[, model[[imod]]]
            xx <- as.matrix(xx)
            if (sum(model[[imod]] == 0) != 1) 
                bvec <- apply(data, 1, regpres1, xx, y, regfun, 
                  mval, ...)
            if (sum(model[[imod]] == 0) != 1) 
                yhat <- cbind(1, xx) %*% bvec
            if (sum(model[[imod]] == 0) == 1) {
                bvec0 <- matrix(0, nrow = p1, ncol = nboot)
                for (it in 1:nboot) {
                  bvec0[1, it] <- locfun(y[data[it, ]])
                }
                yhat <- cbind(1, x) %*% bvec0
            }
            bi <- apply(bid, 1, sum)
            temp <- (bid * (yhat - y))
            diff <- apply(temp, 1, error)
            ep0 <- sum(diff/bi)/length(y)
            aperror <- error(regfun(xx, y, ...)$resid)/length(y)
            regpre <- 0.368 * aperror + 0.632 * ep0
            mout[imod, 1] <- aperror
            mout[imod, 3] <- regpre
            temp <- yhat - y
            diff <- apply(temp, 1, error)
            mout[imod, 2] <- sum(diff)/(nboot * length(y))
        }
        if (sum(model[[imod]] == 0) == 1) {
            mout[imod, 3] <- locpre(y, error = error, est = locfun, 
                SEED = SEED, mval = mval)
        }
    }
    mout[, 5] = rank(mout[, 3])
    if (plotit) 
        plot(c(1:nrow(mout)), mout[, 3], xlab = xlab, ylab = ylab)
    list(estimates = mout)
  }

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