gamplotv2:

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

Usage

gamplotv2(x, y, sop = FALSE, pyhat = FALSE, eout = FALSE, xout = FALSE, outfun = out, plotit = TRUE, varfun = pbvar, xlab = "X", ylab = "", zlab = "", theta = 50, phi = 25, expand = 0.5, SCALE = FALSE, cor.fun = pbcor, ADJ = FALSE, nboot = 20, pr = TRUE, SEED = TRUE, ticktype = "simple")

Arguments

x
y
sop
pyhat
eout
xout
outfun
plotit
varfun
xlab
ylab
zlab
theta
phi
expand
SCALE
cor.fun
ADJ
nboot
pr
SEED
ticktype

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, sop = FALSE, pyhat = FALSE, eout = FALSE, xout = FALSE, 
    outfun = out, plotit = TRUE, varfun = pbvar, xlab = "X", 
    ylab = "", zlab = "", theta = 50, phi = 25, expand = 0.5, 
    SCALE = FALSE, cor.fun = pbcor, ADJ = FALSE, nboot = 20, 
    pr = TRUE, SEED = TRUE, ticktype = "simple") 
{
    if (ADJ) {
        if (SEED) 
            set.seed(2)
    }
    if (pr) {
        if (!ADJ) {
            print("To get adjusted estimates of strength of association, use ADJ=T")
            print("The strength of association is estimated under independence")
            print(" and then rescaled")
        }
    }
    library(akima)
    library(mgcv)
    x <- as.matrix(x)
    np <- ncol(x)
    np1 <- np + 1
    if (ncol(x) > 4) 
        stop("x should have at most four columns of data")
    m <- elimna(cbind(x, y))
    if (xout && eout) 
        stop("Can't have xout=eout=T")
    if (eout) {
        flag <- outfun(m)$keep
        m <- m[flag, ]
    }
    if (xout) {
        flag <- outfun(x, plotit = FALSE)$keep
        m <- m[flag, ]
    }
    x <- m[, 1:np]
    x = as.matrix(x)
    y <- m[, np1]
    if (!sop) {
        if (ncol(x) == 1) 
            fitr <- fitted(gam(y ~ x[, 1]))
        if (ncol(x) == 2) 
            fitr <- fitted(gam(y ~ x[, 1] + x[, 2]))
        if (ncol(x) == 3) 
            fitr <- fitted(gam(y ~ x[, 1] + x[, 2] + x[, 3]))
        if (ncol(x) == 4) 
            fitr <- fitted(gam(y ~ x[, 1] + x[, 2] + x[, 3] + 
                x[, 4]))
    }
    if (sop) {
        if (ncol(x) == 1) 
            fitr <- fitted(gam(y ~ s(x[, 1])))
        if (ncol(x) == 2) 
            fitr <- fitted(gam(y ~ s(x[, 1]) + s(x[, 2])))
        if (ncol(x) == 3) 
            fitr <- fitted(gam(y ~ s(x[, 1]) + s(x[, 2]) + s(x[, 
                3])))
        if (ncol(x) == 4) 
            fitr <- fitted(gam(y ~ s(x[, 1]) + s(x[, 2]) + s(x[, 
                3]) + s(x[, 4])))
    }
    last <- fitr
    if (plotit) {
        if (ncol(x) == 1) {
            plot(x, fitr, xlab = xlab, ylab = ylab)
        }
        if (ncol(x) == 2) {
            iout <- c(1:length(fitr))
            nm1 <- length(fitr) - 1
            for (i in 1:nm1) {
                ip1 <- i + 1
                for (k in ip1:length(fitr)) if (sum(x[i, ] == 
                  x[k, ]) == 2) 
                  iout[k] <- 0
            }
            fitr <- fitr[iout >= 1]
            mkeep <- x[iout >= 1, ]
            fitr <- interp(mkeep[, 1], mkeep[, 2], fitr)
            persp(fitr, theta = theta, phi = phi, expand = expand, 
                xlab = "x1", ylab = "x2", zlab = "", scale = scale, 
                ticktype = ticktype)
        }
    }
    top = varfun(last)
    ep = top/varfun(y)
    if (ep >= 1) 
        ep = cor.fun(last, y)$cor^2
    eta = sqrt(ep)
    st.adj = NULL
    e.adj = NULL
    if (ADJ) {
        x = as.matrix(x)
        val = NA
        n = length(y)
        data1 <- matrix(sample(n, size = n * nboot, replace = TRUE), 
            nrow = nboot)
        data2 <- matrix(sample(n, size = n * nboot, replace = TRUE), 
            nrow = nboot)
        for (i in 1:nboot) {
            temp = gamplotv2.sub(x[data1[i, ], ], y[data2[i, 
                ]], plotit = FALSE)
            val[i] = temp$Explanatory.power
        }
        vindt = median(val)
        v2indt = median(sqrt(val))
        st.adj = (sqrt(ep) - max(c(0, v2indt)))/(1 - max(c(0, 
            v2indt)))
        e.adj = (ep - max(c(0, vindt)))/(1 - max(c(0, vindt)))
        st.adj = max(c(0, st.adj))
        e.adj = max(c(0, e.adj))
    }
    eta = as.matrix(eta)
    ep = as.matrix(ep)
    dimnames(eta) = NULL
    dimnames(ep) = NULL
    eta = eta[1]
    ep = ep[1]
    list(Strength.Assoc = eta, Explanatory.power = ep, Strength.Adj = st.adj, 
        Explanatory.Adj = e.adj)
  }

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