[R] permutation test on paired samples
Holger Taschenberger
Holger.Taschenberger at mpi-bpc.mpg.de
Mon Jul 16 16:11:18 CEST 2012
Henric,
thanks a lot for the clarification. I guess your response also
implies that there is no 'simple' equivalent or replacement for the
(exactRankTests) function perm.test(y,x,paired = TRUE,...) in the
package coin. Perhaps it will be added in the future. Meanwhile I can
use your 'recipe'.
--Holger
On Sun, 15 Jul 2012 19:36:01 +0200
"Henric (Nilsson) Winell" <nilsson.henric at gmail.com> wrote:
> Holger,
>
> Thanks for providing a reproducible example. However, since your
> space key only works sporadically, the below is a little hard to read... ;)
>
> On 2012-07-12 20:26, Holger Taschenberger wrote:
> > Hi,
> >
> > I'm trying to run a permutation test on paired samples.
> >
> > First I tried the package "exactRankTests":
> >
> > require("exactRankTests")
> > x <- c(1.83,0.50,1.62,2.48,1.68,1.88,1.55,3.06,1.30)
> > y <- c(0.878,0.647,0.598,2.05,1.06,1.29,1.06,3.14,1.29)
>
> The relevant output missing here is
>
> > wilcox.test(x,y,paired = TRUE,alternative = "greater")
>
> Wilcoxon signed rank test
>
> data: x and y
> V = 40, p-value = 0.01953
> alternative hypothesis: true location shift is greater than 0
>
> > perm.test(y,x,paired = TRUE,exact = TRUE,alternative = "greater")
>
> 1-sample Permutation Test (scores mapped into 1:m using rounded
> scores)
>
> data: y and x
> T = 41, p-value = 0.003906
> alternative hypothesis: true mu is greater than 0
>
>
> Firstly, you've interchanged the 'x' and 'y' in the second call.
> Secondly, and more important, the output says that "(scores mapped into
> 1:m using rounded scores)". In this case this can easily be avoided,
> and note the interchange of 'x' and 'y' to match your 'wilcox.test'
> call, using:
>
> > yy <- 1000 * y
> > xx <- 1000 * x
> > perm.test(xx, yy, paired = TRUE, exact = TRUE,
> + alternative = "greater")
>
> 1-sample Permutation Test
>
> data: xx and yy
> T = 4114, p-value = 0.01367
> alternative hypothesis: true mu is greater than 0
>
>
> So, now that we've computed the correct p-value, let's see how to obtain
> this using the 'coin' package:
>
> >
> > Then I wanted to use the package 'coin':
> >
> > require("coin")
> > x <- c(1.83,0.50,1.62,2.48,1.68,1.88,1.55,3.06,1.30)
> > y <- c(0.878,0.647,0.598,2.05,1.06,1.29,1.06,3.14,1.29)
> > xydat <- data.frame(y = c(y,x),x = gl(2,length(x)),block =
> factor(rep(1:length(x),2)))
>
> The relevant output missing here is
>
> > wilcoxsign_test(y ~ x | block,data = xydat,alternative =
> "greater",distribution = exact())
>
> Exact Wilcoxon-Signed-Rank Test
>
> data: y by x (neg, pos)
> stratified by block
> Z = 2.0732, p-value = 0.01953
> alternative hypothesis: true mu is greater than 0
>
> > oneway_test(y ~ x | block,data = xydat,alternative =
> "greater",distribution = exact())
>
> Exact 2-Sample Permutation Test
>
> data: y by x (1, 2)
> stratified by block
> Z = -2.1948, p-value = 0.6982
> alternative hypothesis: true mu is greater than 0
>
>
> Using 'oneway_test' in this way does *not* correspond to a paired test.
> The "raw" scores version of the Wilcoxon signed-rank test can be
> constructed using
>
> > diff <- x - y
> > y <- as.vector(t(cbind(abs(diff) * (diff < 0),
> + abs(diff) * (diff >= 0))))
> > x <- factor(rep(c("neg", "pos"), length(diff)),
> + levels = c("pos", "neg"))
> > b <- gl(length(diff), 2)
> >
> > oneway_test(y ~ x | b, alternative = "greater", distr = "exact")
>
> Exact 2-Sample Permutation Test
>
> data: y by x (pos, neg)
> stratified by b
> Z = 2.1948, p-value = 0.01367
> alternative hypothesis: true mu is greater than 0
>
>
> And, as you can see, this is equal to the 'perm.test' result.
>
>
> HTH,
> Henric
>
>
>
> >
> > While the results of the Wilcoxon test are the same for both packages
> > are the same, those of the permutation test are very different. So,
> > obviously I'm doing something wrong here. Can somebody please help?
> >
> > Thanks a lot,
> > Holger
> >
> > ______________________________________________
> > R-help at r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
> >
--
Holger Taschenberger, PhD
Dept. of Membrane Biophysics
Max Planck Institute for Biophysical Chemistry
Am Fassberg 11
D-37077 Goettingen, Germany
Tel.: ++49-551-201-1668
Fax: ++49-551-201-1688
E-mail: <Holger.Taschenberger at mpi-bpc.mpg.de>
More information about the R-help
mailing list