trq.fit:

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

Usage

trq.fit(x, y, a1 = 0.1, a2, z, int = T, method = "primal", tol = 1e-04)

Arguments

x
y
a1
a2
z
int
method
tol

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, a1 = 0.1, a2, z, int = T, method = "primal", 
    tol = 1e-04) 
{
    if (missing(a2)) 
        a2 <- a1
    if (a1 < 0 | a2 < 0) 
        stop("trimming proportion negative.")
    if (a1 + a2 - 1 > tol) 
        stop("trimming proportion greater than 1.")
    if (method != "primal" & method != "dual") 
        stop("invalid method: should be 'primal' or 'dual'.")
    x <- as.matrix(x)
    if (missing(z)) {
        ny = length(y)
        xx = cbind(rep(1, ny), x)
        z <- rq.fit.br(xx, y, tau = -1)
        print(z)
    }
    p <- z$sol[1, ]
    q <- matrix(z$sol[-c(1:3), ], nrow(z$sol) - 3, ncol(z$sol))
    n <- nrow(z$dsol)
    s <- NULL
    if (length(dimnames(x)[[2]]) == 0) 
        dimnames(x) <- list(NULL, paste("X", 1:(nrow(q) - 1), 
            sep = ""))
    if (int) {
        x <- cbind(1, x)
        dimnames(x)[[2]][1] <- "Intercept"
    }
    xbar <- apply(x, 2, "mean")
    xxinv <- solve(t(x) %*% x)
    if (abs(a1 + a2 - 1) <= tol) {
        i <- sum(p < a1)
        s$coef <- q[, i]
        names(s$coef) <- dimnames(x)[[2]]
        s$resid <- y - x %*% s$coef
        PI <- 3.14159
        x0 <- qnorm(a1)
        d0 <- (1/sqrt(2 * PI)) * exp(-(x0^2/2))
        d0 <- ((4.5 * d0^4)/(2 * x0^2 + 1)^2)^0.2
        d <- d0 * (length(s$resid) - length(s$coef))^(-0.2)
        if (d > min(a1, 1 - a1)) 
            d <- min(a1, 1 - a1)
        s$d <- d
        i <- sum(p < a1 + d)
        j <- sum(p < a1 - d)
        shat <- as.numeric(xbar %*% t(q[, i] - q[, j]))/(2 * 
            d)
        s$int <- int
        s$v <- a1 * (1 - a1) * shat^2
        s$cov <- s$v * xxinv
    }
    else {
        p1 <- p[-1]
        f <- 1/(1 - a1 - a2)
        d <- pmax((pmin(p1, 1 - a2) - c(a1, pmax(p1[-length(p1)], 
            a1))), 0)
        if (method == "primal") {
            s$coef <- q[, 1:length(p1)] %*% d * f
            s$resid <- y - x %*% s$coef
            s$int <- int
        }
        else {
            i <- max(1, sum(p < a1))
            g <- (z$dsol[, i + 1] - z$dsol[, i])/(p[i + 1] - 
                p[i])
            wa <- z$dsol[, i] + (a1 - p[i]) * g
            j <- sum(p < 1 - a2)
            g <- (z$dsol[, j + 1] - z$dsol[, j])/(p[j + 1] - 
                p[j])
            wb <- z$dsol[, j] + (1 - a2 - p[j]) * g
            wt <- wa - wb
            if (min(wt) < 0) 
                warning("some weights negative!")
            s <- lsfit(x, y, abs(wt), int = F)
        }
        mu <- xbar %*% s$coef
        v <- d %*% (z$sol[2, 1:length(d)] - mu)^2
        k <- qrq(z, c(a1, a2)) - mu
        v <- v + a1 * k[1]^2 + a2 * k[2]^2 + (a1 * k[1] + a2 * 
            k[2])^2
        names(s$coef) <- dimnames(x)[[2]]
        s$v <- as.vector(f^2 * v)
        s$cov <- s$v * xxinv
    }
    class(s) <- "trq"
    s
  }

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