From thomas@biostat.washington.edu  Wed Jul 15 01:33:58 1998
From: thomas@biostat.washington.edu (Thomas Lumley)
Date: Tue, 14 Jul 1998 17:33:58 -0700 (PDT)
Subject: Are post-hoc tests being developed for R?
In-Reply-To: <35ABAA17.74C1B359@spud.me.wustl.edu>
Message-ID: <Pine.SOL.3.95.980714172634.13636A-100000@sol>


Here is code for three p-value based post hoc tests: bonferroni, holm, and
hochberg. Holm is better than Bonferroni under all circumtances. Hochberg
is more powerful than Holm, but under some slightly pathological
situations it can exceed the nominal Type I error rate. 

I haven't looked at this code for a long time, but I think it works
correctly. Refeerences for the methods are in "Adjusted p-values for
simultaneous inference" by S. Paul Wright in Biometrics some time in
the early 90s (I'm in the wrong country at the moment so I don't have the
full reference).

	-thomas
-------------------------

   p.adjust.holm <-function(p,n=length(p)) {  
        ##n<-length(p)
        r<-rank(p)
        index<-order(p)
        qi<-p*(n+1-r)
        for (i in 2:length(p)) {
                qi[index[i]]<-max(qi[index[i]], qi[index[i-1]])
        }
        list(adjp=pmin(qi,1),p=p,method="Holm")
        }
  p.adjust.hochberg <-function(p) {
        n<-length(p)
        r<-rank(p)
        index<-order(p)
        qi<-p*(n+1-r)
        for (i in (n-1):1) {
                qi[index[i]]<-min(qi[index[i]], qi[index[i+1]])
        }
        list(adjp=qi,p=p,method="Hochberg")
        }
  p.adjust.bonferroni<-function(p,n=length(p)){ 
        list(adjp=pmin(p*n,1),p=p,method="Bonferroni")
        }

p.adjust<-function(p,method=c("hochberg","holm","bonferroni"),...){
  how<-pmatch(method[1],c("hochberg","holm","bonferroni"))
  if (is.na(how)) stop(paste("Don't know method:",method))
  m<-match.call()
  m[[1]]<-as.name(paste("p.adjust",c("hochberg","holm","bonferroni")[how],sep="."))
  m$method<-NULL
  eval(m,sys.parent())
}

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._


Are post-hoc tests being developed for R?
Thomas Lumley