Dear Gerrit,
cause for the bug within subroutine qwilcox are the following two
lines:
cumprob <- cumsum(dwilcox(x, m, n))
cumprob[length(cumprob)] <- 1
As a result of cumulated round-off-errors, the last elements of
cumprob could be >1 (circa 1+4.4e-016).
As a quick fix, please insert the line #2
cumprob <- cumsum(dwilcox(x, m, n)) #1
cumprob<- ifelse(cumprob>1,1,cumprob) #2
cumprob[length(cumprob)] <- 1 #3
Regards,
Herwig
--
Dr. Herwig Meschke
Wissenschaftliche Beratung
Hagsbucher Weg 27
D-89150 Laichingen
phone +49 7333 210 417 / fax +49 7333 210 418
email HerwigMeschke@t-online.de
> Dear list,
>
> I've encountered the following strange behaviour when computing a quantile
> of Wilcoxons rank sum distribution using qwilcox():
>
> > qwilcox( 0.95, 36, 36)
> Problem in cut.default(x, ..1): breaks must be given in ascending order
> and contain no NA's
> Use traceback() to see the call stack
>
>
> Similar sample sizes (even or odd, and at least the ones, I've tried)
> don't seem to present any problem, e. g.:
>
> > qwilcox( 0.95, 34, 34)
> [1] 1307
> > qwilcox( 0.95, 35, 35)
> [1] 1383
> > qwilcox( 0.95, 37, 37)
> [1] 1540
> > qwilcox( 0.95, 38, 38)
> [1] 1621
>
>
> I've checked the archive, but didn't find anything. Does anybody have an
> idea about what is happening above, and even better: how to work around
> this problem (bug?) ? Thank you.
>
> Best regards -- Gerrit
>
> PS: S-PLUS Version 6.0 Release 1 for Sun SPARC, SunOS 5.6 : 2000.
>
> -------------------------------------------------------------------------
> Dr. Gerrit Eichner Mathematical Institute of the
> gerrit.eichner@math.uni-giessen.de Justus-Liebig-Univ. Giessen
> Tel: +49-(0)641-99-32104 Arndtstr. 2, 35392 Giessen, Germany
> Fax: +49-(0)641-99-32029 http://www.math.uni-giessen.de/Stochastik
> -------------------------------------------------------------------------
>
>
|