R/twwmcp.R

twwmcp <-
function(J,K,x,grp=c(1:p),p=J*K,tr=.2,alpha=.05,dif=FALSE){
#
#  For a J by K anova using quantiles with
#  repeated measures on both factors,
#  Perform all multiple comparisons for main effects
#  and interactions.
#
#  tr=.2. default trimming
#  bop=F means bootstrap option not used;
#  with bop=T, function uses usual medians rather
#  rather than a single order statistic to estimate median
#  in conjunction with bootstrap estimate of covariances
#  among the sample medians.
#
#  The R variable data is assumed to contain the raw
#  data stored in a matrix or in list mode.
#  When in list mode data[[1]] contains the data
#  for the first level of both factors: level 1,1.
#  data[[2]] is assumed to contain the data for level 1 of the
#  first factor and level 2 of the second: level 1,2
#  data[[K]] is the data for level 1,K
#  data[[K+1]] is the data for level 2,1, data[2K] is level 2,K, etc.
#
#  It is assumed that data has length JK, the total number of
#  groups being tested, but a subset of the data can be analyzed
#  using grp
#
Qa<-NA
Qab<-NA
if(is.list(x))x<-elimna(matl(x))
if(is.data.frame(x))x=as.matrix(x)
if(is.matrix(x))x<-elimna(x)
data<-x
if(is.matrix(data))data<-listm(data)
if(!is.list(data))stop("Data are not stored in list mode or a matrix")
if(p!=length(data)){
print("The total number of groups, based on the specified levels, is")
print(p)
print("The number of groups stored in x is")
print(length(data))
print("Warning: These two values are not equal")
}
if(p!=length(grp))stop("Apparently a subset of the groups was specified that does not match the total number of groups indicated by the values for J and K.")
tmeans<-0
temp<-con2way(J,K) # contrasts matrices stored in temp
Qa<-rmmcp(x,con=temp$conA,alpha=alpha,dif=dif,tr=tr)
# Do test for factor B
Qb<-rmmcp(x,con=temp$conB,alpha=alpha,dif=dif,tr=tr)
# Do test for factor A by B interaction
Qab<-rmmcp(x,con=temp$conAB,alpha=alpha,dif=dif,tr=tr)
list(Qa=Qa,Qb=Qb,Qab=Qab)
}
musto101/wilcox_R documentation built on May 23, 2019, 10:52 a.m.