You crossed the value ncp = 80: the help page warned you of cancellation there, and at p = 1 - 1e-12, that is what you got (badly).

On Sun, 11 Oct 2009, witkov...@savba.sk wrote:

Full_Name: Viktor Witkovsky
Version: 2.9.2
OS: Windows XP
Submission from: (NULL) (78.98.89.227)


Hello,

I have found strange behavior of the function qchisq (the non-central qchisq is
based on inversion of pchisq, which is further based on pgamma). The function
gives wrong results without any warning. For example:

qchisq(1e-12,1,8.94^2,lower.tail=FALSE) gives 255.1840972465858 (notice that
here the correct value should be 255.1841334848075),
but
qchisq(1e-12,1,8.95^2,lower.tail=FALSE) gives 1249.681320136836
Here, the correct value should be 255.5037231613135.

So, it seems that qchisq is inaccurate for small probability values and larger
non/centrality parameter.

Actually, for probability values near 1, and it does say so in the documentation. But see

qchisq(1e-12, 1, 81, lower.tail=FALSE)
[1] 1258.412
qchisq(1-1e-12, 1, 81)
[1] 257.1488

which I suggest gives you a workaround -- we'll look into giving a warning from the code.


I am using the precompiled binary version of R, under Windows XP.

              _
platform       i386-pc-mingw32
arch           i386
os             mingw32
system         i386, mingw32
status
major          2
minor          9.2
year           2009
month          08
day            24
svn rev        49384
language       R
version.string R version 2.9.2 (2009-08-24)

______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


--
Brian D. Ripley,                  rip...@stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595

______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel

Reply via email to