-------- Original Message --------
Subject: Re: p-value
Date: Sat, 5 Sep 2009 06:13:48 -0700 (PDT)
From: Michael Collyer <[email protected]>
To: [email protected]
References: <[email protected]>
Samor,
Your resampling procedure is not correct for the type of test you want
to do. You are resampling from each group with replacement (i.e.,
bootstrap). This might be fine for empirically calculating e.g.,
standard errors, for your group means but it does not create random
outcomes under a null hypothesis of no difference in means from two
populations. To do that you need a resampling procedure that creates
random combinations of two groups sampled from the same population. You
can do that by randomly shuffling all of your values into either of the
two groups.
Also when you do this procedure, you need not calculate t-values but
just the difference in means.
Here is an R-script which follows your initial set-up, which I found
gives me the same answer as the t-test. I use the absolute difference
in means, which is like doing a 2-tailed t-test.
x1 <- rgamma(10,3,0.5)
x2 <- rgamma(10,2,1)
t.test(x1,x2)
# permutation test
ptest<-function(x1,x2,p){
x1.mean<-mean(x1)
x2.mean<-mean(x2)
n1<-length(x1)
n2<-length(x2)
dif<-abs(x1.mean-x2.mean)
X<-c(x1,x2)
cnt<-1
for(i in 1:p){
X.random<-sample(X,replace=F)
x1.random<-X.random[1:n1]
x2.random<-X.random[-(1:n1)]
x1.r.mean<-mean(x1.random)
x2.r.mean<-mean(x2.random)
dif.random<-abs(x1.r.mean-x2.r.mean)
if(dif.random>=dif) cnt<-cnt+1
}
pval<-cnt/(p+1)
result<-c(dif,pval) # function will produce two values: absolute
difference in means and P-value
}
ptest1<-ptest(x1,x2,99)
ptest1
Michael Collyer
Assistant Professor
Department of Biology
Stephen F. Austin State University
PO Box 13003 SFA Station
Nacogdoches, TX 75962-3003
Phone (936) 468-2322 Fax (936) 468-2056
Email [email protected]
--
Replies will be sent to the list.
For more information visit http://www.morphometrics.org
