#9 Precision loss in NonCentralChiSquare

[Roby Joehanes]

This precision loss might be caused by the long double problem.

Output of test_noncentralchisq():

-- non central Chi^2 :

df =0.1Error: xx=38.1000, df=0.100000, ncp=1, pchisq = 0.999999952483711000, qchisq = 38.1000000126858540 != xx
[1] FALSE
df =0.5Error: xx=38.5000, df=0.500000, ncp=1, pchisq = 0.999999940533911500, qchisq = 38.5000000021106100 != xx
[1] FALSE
df =1.5Error: xx=39.5000, df=1.50000, ncp=1, pchisq = 0.999999898134906100, qchisq = 39.5000000027148560 != xx
[1] FALSE
df =4.7Error: xx=39.7000, df=4.70000, ncp=1, pchisq = 0.999998397115614200, qchisq = 39.7000000003507700 != xx
Error: xx=42.7000, df=4.70000, ncp=1, pchisq = 0.999999535394586300, qchisq = 42.7000000008072100 != xx
[1] FALSE
df =10.0[1] TRUE
df =20.0[1] TRUE
df =50.0[1] TRUE
df =100.0Error: xx=1.00000e-05, df=100.000, ncp=1, pchisq = 0.00000000000000000, qchisq = 0.00000000000000000 != xx
[1] FALSE
[1] TRUE
[1] TRUE

[Roby Joehanes]

The last one is particularly devastating. The results of qchisq have to be equal to xx.

[Roby Joehanes]

Vis-a-vis comparison:

options(digits=22)
pchisq(38.1, 0.1, 1) = 0.999999952483710874418
JDistlib = 0.999999952483711000
qchisq(pchisq(38.1, 0.1, 1), 0.1, 1) = 38.10000000167171663179
JDistlib = 38.1000000126858540

pchisq(38.5, 0.5, 1) = 0.9999999405339115021718
JDistlib = 0.999999940533911500
qchisq(pchisq(38.5, 0.5, 1), 0.5, 1) = 38.50000000297760038848
JDistlib = 38.5000000021106100

pchisq(39.5, 1.5, 1) = 0.9999998981349061955015
JDistlib = 0.999999898134906100
qchisq(pchisq(39.5, 1.5, 1), 1.5, 1) = 39.50000000224288498885
JDistlib = 39.5000000027148560

pchisq(39.7, 4.7, 1) = 0.9999983971156142237291
JDistlib = 0.999998397115614200
qchisq(pchisq(39.7, 4.7, 1), 4.7, 1) = 39.70000000002929141374
JDistlib = 39.7000000003507700

pchisq(42.7, 4.7, 1) = 0.9999995353945863030276
JDistlib = 0.999999535394586300
qchisq(pchisq(42.7, 4.7, 1), 4.7, 1) = 42.70000000051550159696
JDistlib = 42.7000000008072100

pchisq(1e-5, 100, 1) = 0
JDistlib = 0
qchisq(pchisq(1e-5, 100, 1), 100, 1) = 0
JDistlib = 0

[Roby Joehanes]

The last one seems to be an underflow bug, at least to me. Let's see what R developers think.
https://bugs.r-project.org/bugzilla3/show_bug.cgi?id=15635

[Roby Joehanes]

Martin said that the fix is committed to SVN r65384, but I can't find it at r65385. Lag maybe? The fix is on the verge of getting released.

[Roby Joehanes]

Not found in 3.1.0 release either.

[Roby Joehanes]

Fixes found in R_devel_2014-04-10. Synced and incorporated. However, the long double issue remains standing.