Thanks once again Peter, I understand your points -- I fiddled and googled and read some more and found an eas(ier) route:
library(coin) t <- wilcoxsign_test(a ~ b, alternative = "two.sided", distribution = exact()) # This is equivalent to paired wilcox.test pval <- pvalue(t) sweet.zscore <- statistic(t) # signed and everything # correct for multiple testing if you are doing a number of the above... I would never have got here without your guidance - so kudos to you. JP On 5 May 2011 14:52, peter dalgaard <pda...@gmail.com> wrote: > > On May 5, 2011, at 10:58 , JP wrote: > >> On 4 May 2011 15:32, peter dalgaard <pda...@gmail.com> wrote: >> >>> >>> On May 4, 2011, at 15:11 , JP wrote: >>> >>>> Peter thanks for the fantastically simple and understandable explanation... >>>> >>>> To sum it up... to find the z values of a number of pairwise wilcox >>>> tests do the following: >>>> >>>> # pairwise tests with bonferroni correction >>>> x <- pairwise.wilcox.test(a, b, alternative="two.sided", >>>> p.adj="bonferroni", exact=F, paired=T) >>> >>> >>> You probably don't want the bonferroni correction there. Rather >>> p.adj="none". You generally correct the p values for multiple testing, not >>> the test statistics. >>> >> >> Oh, I see thanks... of course since I have 5 groups (samples) and 10 >> comparisons I still have to correct when quoting p values... >> >> >>> (My sentiment would be to pick apart the stats:::wilcox.test.default >>> function and clone the computation of Z from it, but presumably >>> backtracking from the p value is a useful expedient.) >>> >> >> Should this be so onerous for the user [read non-statistician] ? > > > My main reservation is that if the p value gets very small or close to 1 (for > 1-sided tests), then it might not be possible to reconstruct the underlying Z > (qnorm(pnorm(Z)) breaks down around Z = -37 and Z=+8). If you get sensible > values, there's probably not a problem. > > >> >>>> # what is the data structure we got back >>>> is.matrix(x$p.value) >>>> # p vals >>>> x$p.value >>>> # z.scores for each >>>> z.score <- qnorm(x$p.value / 2) >>>> >>> >>> Hmm, you're not actually getting a signed z out of this, you might want to >>> try alternative="greater" and drop the division by 2 inside qnorm(). (If >>> the signs come out inverted, I meant "less" not "greater"...) >>> >> >> But I need a two sided test (changing the alternative would change the >> hypothesis!)... do I still do this? > > For the z's, yes. > > Hmm, the asymmetry of qnorm(pnorm(...)) indicates that you might want to be a > bit more careful. > >> All my z values are negative.... >> >> Is this correct? > > No, that's the problem. The sign of Z should depend on whether the signed > ranks are predominantly positive or negative. If you do a two-sided p-value > and divide by two, then by definition you get something less than .5 and > qnorm of that will be negative. So if the 2-sided p is 0.04, the one-sided p > will be either 0.02 or 0.98, depending on the alternative and on whether the > V statistic is above or below its expectation. Notice that this corresponds > to Z statistics of opposite sign: > >> qnorm(.02) > [1] -2.053749 >> qnorm(.98) > [1] 2.053749 > > -- > Peter Dalgaard > Center for Statistics, Copenhagen Business School > Solbjerg Plads 3, 2000 Frederiksberg, Denmark > Phone: (+45)38153501 > Email: pd....@cbs.dk Priv: pda...@gmail.com > > > ______________________________________________ R-help@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.
