[R] p-value for nonparamentric multiple comparison testing
Antonio Olinto
aolinto_r at bignet.com.br
Tue Nov 2 15:05:18 CET 2004
Dear R users,
I wrote a function to perform a nonparametric multiple comparison test. The
function below solves the example 11.10 (pg. 228) from ZarÂs Biostatistical
Analysis (3rd ed.). Nevertheless I couldnÂt find a function to get the p-value
associated to Q, so I still have to consult Table B.15, App106 (Critical values
of Q for Nonparamentric Multiple Comparison Testing). Is there any function in R
that returns this value?
Thanks in advance.
Antonio Olinto
==========
multcomp <- function(VarCat,VarNum) {
dat.multcomp <- data.frame(VarCat,VarNum)
names(dat.multcomp) <- c("VarCat","VarNum")
attach(dat.multcomp)
dat.multcomp$Rank <- rank(VarNum)
attach(dat.multcomp)
RankList <- aggregate(Rank,list(Rank=Rank),FUN=length)
t <- length(RankList$Rank)
st <- sum(ifelse(RankList[,2]>1,RankList[,2]^3-RankList[,2], 0))
LevCat <- levels(dat.multcomp$VarCat)
NLevCat <- aggregate(VarCat,list(LevCat=VarCat),FUN=length)
RLevCat <- aggregate(Rank,list(LevCat=VarCat),FUN=sum)
MLevCat <- aggregate(Rank,list(LevCat=VarCat),FUN=mean)
SampleSummary <- data.frame(LevCat,RLevCat[,2],NLevCat[,2],MLevCat[,2])
names(SampleSummary)<-c("Samples","RSum","N","RMean")
SampleSummary <- SampleSummary[order(SampleSummary$RMean,decreasing=T),]
NCat <- length(LevCat)
NComb <- choose(NCat,2)
N <- length(dat.multcomp$VarCat)
Results <- data.frame(rep(NA,NComb),rep(NA,NComb),rep(NA,NComb),rep(NA,NComb))
names(Results) <- c("Comparison","Difference","SE","Q")
l <- 1
for (i in 1:(NCat-1)) {
for (j in NCat:(i+1)) {
SE <- sqrt(((N*(N+1)/12)-(st/(12*(N-1))))*((1/SampleSummary[i,3])+
(1/SampleSummary[j,3])))
Dif <- SampleSummary[i,4]-SampleSummary[j,4]
Q=Dif/SE
Results[l,1] <- paste(SampleSummary[i,1],"vs",SampleSummary[j,1])
Results[l,2] <- round(Dif,4)
Results[l,3] <- round(SE,4)
Results[l,4] <- round(Q,4)
l <-l+1
}
}
print("Sample summary ranked by mean ranks")
print(SampleSummary)
print("")
print("Table of multiple comparisons")
print(Results)
}
-------------------------------------------------
WebMail Bignet - O seu provedor do litoral
www.bignet.com.br
More information about the R-help
mailing list