R/wwtrimbt.R

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

wwtrimbt <-
function(J, K, x, tr = 0.2, JK = J*K, con = 0,
 alpha = 0.05, grp =c(1:JK), nboot = 599,SEED = TRUE, ...){
        #
        # A bootstrap-t for a within-by-within omnibus tests
        #  for all main effects and interactions
        #
        #  The R variable x is assumed to contain the raw
        #  data stored in list mode or in a matrix.
        #  If in list mode, x[[1]] contains the data
        #  for the first level of both factors: level 1,1.
        #  x[[2]] is assumed to contain the data for level 1 of the
        #  first factor and level 2 of the second: level 1,2
        #  x[[K]] is the data for level 1,K
        #  x[[K+1]] is the data for level 2,1, x[[2K]] is level 2,K, etc.
        #
        #  If the data are in a matrix, column 1 is assumed to
        #  correspond to x[[1]], column 2 to x[[2]], etc.
        #
        #  When in list mode x is assumed to have length JK, the total number
        #  groups being tested, but a subset of the data can be analyzed
        #  using grp
        #
if(is.data.frame(x))x=as.matrix(x)
        if(is.matrix(x)) {
                y <- list()
                for(j in 1:ncol(x))
                        y[[j]] <- x[, j]
x=y
}
ncon=ncol(con)
 p <- J*K
JK=p
if(p>length(x))stop("JK is less than the Number of groups")
JK=J*K
        data <- list()
xx=list()
        for(j in 1:length(x)) {
xx[[j]]=x[[grp[j]]] # save input data
data[[j]] = xx[[j]] - mean(xx[[j]], tr = tr)
#                # Now have the groups in proper order.
        }
        if(SEED)set.seed(2)
        # set seed of random number generator so that
        #             results can be duplicated.
        bsam = list()
        bdat = list()
aboot=NA
bboot=NA
cboot=NA
abboot=NA
test.stat=wwtrim(J,K,xx,tr=tr)
nv=length(x[[1]])
        for(ib in 1:nboot) {
bdat[[j]] = sample(nv, size = nv, replace =T)
for(k in 1:JK) bsam[[k]] = data[[k]][bdat[[j]]]
temp=wwtrim(J,K,bsam,tr=tr)
aboot[ib]=temp$Qa
bboot[ib]=temp$Qb
abboot[ib]=temp$Qab
}
pbA=NA
pbB=NA
pbAB=NA
pbA=mean(test.stat$Qa[1,1]<aboot)
pbB=mean(test.stat$Qb[1,1]<bboot)
pbAB=mean(test.stat$Qab[1,1]<abboot)
list(p.value.A=pbA,p.value.B=pbB,p.value.AB=pbAB)
}