permute.residuals:

In brentp/bootstrip: Resampling residuals of a linear model to determine significance.

Usage

permute.residuals(mat, mod, mod0, iterations = 100, p_samples = 1, mc.cores = 12)

Arguments

mat
mod
mod0
iterations
p_samples
mc.cores

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 (mat, mod, mod0, iterations = 100, p_samples = 1, mc.cores = 12) 
{
    stopifnot(nrow(mod) == ncol(mat))
    reduced_lm = lmFit(mat, mod0)
    reduced_residuals = residuals(reduced_lm, mat)
    reduced_fitted = fitted(reduced_lm)
    fit = lmFit(mat, mod)
    size = p_samples * nrow(mod)
    coef.name = setdiff(colnames(mod), colnames(mod0))
    beta.orig = coefficients(fit)[, coef.name]
    rm(reduced_lm, fit)
    gc()
    nc = ncol(reduced_residuals)
    beta.list = mclapply(1:iterations, function(ix) {
        if (p_samples < 1) {
            sub_ids = sample.int(nc, size = size)
        }
        else {
            sub_ids = 1:nc
        }
        mat_sim = reduced_fitted[, sub_ids] + reduced_residuals[, 
            sample(sub_ids)]
        coefficients(lmFit(mat_sim, mod[sub_ids, ]))[, coef.name]
    }, mc.cores = mc.cores)
    beta = matrix(NA, nrow(mat), ncol = iterations)
    for (i in 1:iterations) {
        beta[, i] = beta.list[[i]]
    }
    df = data.frame(pvalue = unlist(lapply(1:nrow(beta), function(i) {
        row = beta[i, ]
        val = ecdf(row)(beta.orig[i])
        max(min(val, 1 - val) * 2, 1/length(row))
    })), beta.orig = beta.orig)
    rownames(df) = rownames(mat)
    df
  }

brentp/bootstrip documentation built on May 13, 2019, 5:11 a.m.