R 2.3.0 (alpha) on FreeBSD 6.1 fails make check-all

I tried it.  That is, in qhyper.c I replaced

p *= 1 - 64*DBL_EPSILON;

with

p *= 1 - 1000*DBL_EPSILON;

and recompiled.  2.3.0 alpha now passses make check-all and returns
sensible results on both of my FreeBSD 6.1 machines.

Rhyper <- 1:20
Phyper         <- phyper   (Rhyper, m = 40, n = 30, k = 20)
qhyper(Phyper, m = 40, n = 30, k = 20)

[1]  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20

Thanks very much for your guidance, this has been educational.  Is
there anything else that I should check?




On Tue, Mar 28, 2006 at 01:53:09PM +0100, Prof Brian Ripley wrote:

qgamma.c comtains the line

    p *= 1 - 64*DBL_EPSILON;

which seems to be too small nowadays.  I find for example

p <-  phyper(7, m = 40, n = 30, k = 20)

qhyper(p*(1+5e-16), m = 40, n = 30, k = 20)

[1] 7

qhyper(p*(1+6e-16), m = 40, n = 30, k = 20)

[1] 8

so it seems far more sensitive than the tolerance indicates is intended.
Now, phyper has been changed since that tolerance was added, and that
may be part of the explanation.

Can you try a larger tolerance, e.g 1000?

A better fix is probably not to use a tolerance but rather to do a
downwards search using phyper.


On Tue, 28 Mar 2006, Andrew Robinson wrote:

This is from 2.2.1:

Phyper <- phyper(7, m = 40, n = 30, k = 20)
print(Phyper, digits=16)

[1] 0.01796062766370490

This is from 2.3.0:

Phyper <- phyper(7, m = 40, n = 30, k = 20)
print(Phyper, digits=16)

[1] 0.01796062766370491

qhyper(Phyper, m = 40, n = 30, k = 20)

[1] 8

Then if I save Phyper from 2.2.1, and load it in 2.3.0,

qhyper(Phyper, m = 40, n = 30, k = 20)

[1] 7

so it appears that the difference between 0.01796062766370490 and
0.01796062766370491 is important here.

(2.3.0)

qhyper(0.01796062766370490, m = 40, n = 30, k = 20)

[1] 7

qhyper(0.01796062766370491, m = 40, n = 30, k = 20)

[1] 8



Comparing the full sample:

This is from 2.2.1:

Rhyper <- scan()

1: 16 11 11 15 11 13 13 12 13 10 10  7 11 14 13  9 14 13 13 11
21:
Read 20 items

Phyper         <- phyper   (Rhyper, m = 40, n = 30, k = 20)
print(Phyper, digits=16)

[1] 0.99744478362220412 0.51295659016196715 0.51295659016196715
[4] 0.98680472393309626 0.51295659016196715 0.86627412777488222
[7] 0.86627412777488222 0.71494883441868140 0.86627412777488222
[10] 0.30864259597126753 0.30864259597126753 0.01796062766370490
[13] 0.51295659016196715 0.95139460528774533 0.86627412777488222
[16] 0.15132082044442866 0.95139460528774533 0.86627412777488222
[19] 0.86627412777488222 0.51295659016196715



This is from 2.3.0:

print(Phyper, digits=16)

[1] 0.99744478362220401 0.51295659016196726 0.51295659016196726
[4] 0.98680472393309615 0.51295659016196726 0.86627412777488211
[7] 0.86627412777488211 0.71494883441868129 0.86627412777488211
[10] 0.30864259597126747 0.30864259597126747 0.01796062766370491
[13] 0.51295659016196726 0.95139460528774533 0.86627412777488211
[16] 0.15132082044442868 0.95139460528774533 0.86627412777488211
[19] 0.86627412777488211 0.51295659016196726


The 14's (elements 14 and 17) are identical, everything else is
slightly different.


I also got this interesting result from 2.2.1:

Rhyper <- 1:20
Phyper         <- phyper   (Rhyper, m = 40, n = 30, k = 20)
qhyper(Phyper, m = 40, n = 30, k = 20)

[1]  1  2  3  4  6  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20


Repeating the operation for 2.3.0:

Rhyper <- 1:20
Phyper         <- phyper   (Rhyper, m = 40, n = 30, k = 20)
qhyper(Phyper, m = 40, n = 30, k = 20)

[1]  1  2  3  4  6  6  8  8  9 10 11 12 13 14 15 16 17 18 19 20


(I also tried on WinXP R 2.2.1, and got the expected 1:20 back.)


On Tue, Mar 28, 2006 at 12:02:11PM +0100, Prof Brian Ripley wrote:

Then we need you to dig in to find out why.  I'd start by seeing if Phyper
had changed (either by printing it to 16dp or by saving it from 2.2.1 and
reloading into 2.3.0).

On Tue, 28 Mar 2006, Andrew Robinson wrote:

I get:

Rhyper <- scan()

1: 16 11 11 15 11 13 13 12 13 10 10  7 11 14 13  9 14 13 13 11
21:
Read 20 items

Phyper         <- phyper   (Rhyper, m = 40, n = 30, k = 20)
qhyper(Phyper, m = 40, n = 30, k = 20)

[1] 16 11 11 15 11 13 13 12 13 10 10  8 11 14 13  9 14 13 13 11


The 12th element (8) differs from the input (7).

Phyper <- phyper (7, m = 40, n = 30, k = 20)
qhyper(Phyper, m = 40, n = 30, k = 20)

[1] 8


If I do this using 2.2.1 then the input and the output are identical.




On Tue, Mar 28, 2006 at 10:27:08AM +0100, Prof Brian Ripley wrote:

On Tue, 28 Mar 2006, Andrew Robinson wrote:

You're welcome.  You are correct.  d-p-q-r-tests.Rout.fail
shows:

All.eq(Rhyper,          qhyper   (Phyper, m = 40, n = 30, k = 20))

[1] "Mean scaled  difference: 0.08333333"

Yes, please run the lines below, e.g.

Rhyper <- scan()
16 11 11 15 11 13 13 12 13 10 10  7 11 14 13  9 14 13 13 11

Phyper         <- phyper   (Rhyper, m = 40, n = 30, k = 20)
qhyper(Phyper, m = 40, n = 30, k = 20)

and tell us what answer you get.

Let me know if/how I can further assist.

Andrew



On Tue, Mar 28, 2006 at 09:03:48AM +0100, Prof Brian Ripley wrote:

Thanks for checking.

Please look in d-p-q-r-tests.Rout.fail and see what immediately
preceeds
the line

[1] "Mean scaled  difference: 0.08333333"

Some experimentation suggests it is

All.eq(Rhyper,	  qhyper   (Phyper, m = 40, n = 30, k = 20))

If so, we have

Rhyper <- scan()
16 11 11 15 11 13 13 12 13 10 10  7 11 14 13  9 14 13 13 11

Phyper	  <- phyper   (Rhyper, m = 40, n = 30, k = 20)

and those have been checked.  So the error would appear to be in

qhyper(Phyper, m = 40, n = 30, k = 20)

and indeed a mean scaled difference of 1/12 is plausible, since the
mean
of Rhyper is 12. So I deduce that your platform has a problem in
qhyper,
but please cross-check.

If so, this is strange as the only recent change to qhyper.c (or
things
I
can see it uses such as lfastchoose) is cosmetic.

Can you confirm the diagnosis is correct so far?


On Tue, 28 Mar 2006, Andrew Robinson wrote:

Hi Developers,

The alpha, compiles successfully, but it is failing make check-all
(on
two seperate machines, both FreeBSD 6.1).

Here is the version string:

platform       i386-unknown-freebsd6.1
arch           i386
os             freebsd6.1
system         i386, freebsd6.1
status         alpha
major          2
minor          3.0
year           2006
month          03
day            27
svn rev        37584
language       R
version.string Version 2.3.0 alpha (2006-03-27 r37584)



Here is the error message from make check-all

comparing 'd-p-q-r-tests.Rout' to './d-p-q-r-tests.Rout.save'
...706c706
< [1] "Mean scaled  difference: 0.08333333"
---

[1] TRUE

gmake[3]: *** [d-p-q-r-tests.Rout] Error 1
gmake[3]: Leaving directory `/usr/local/beta/R-alpha/tests'
gmake[2]: *** [test-Specific] Error 2
gmake[2]: Leaving directory `/usr/local/beta/R-alpha/tests'
gmake[1]: *** [test-all-basics] Error 1
gmake[1]: Leaving directory `/usr/local/beta/R-alpha/tests'
gmake: *** [check-all] Error 2




I have checked d-p-q-r-tests.Rout.fail for any obvious problems - I
found some warnings, viz.

pgamma(1,Inf,scale=Inf) == 0

[1] TRUE

## Also pgamma(Inf,Inf) == 1 for which NaN was slightly more

appropriate

all(is.nan(c(pgamma(Inf,  1,scale=Inf),

+              pgamma(Inf,Inf,scale=Inf))))
[1] TRUE
Warning messages:
1: NaNs produced in: pgamma(q, shape, scale, lower.tail, log.p)
2: NaNs produced in: pgamma(q, shape, scale, lower.tail, log.p)

x0 <- -2 * 10^-c(22,10,7,5)
stopifnot(pbinom(x0, size = 3, prob = 0.1) == 0,

+           dbinom(x0, 3, 0.1) == 0) # d*() warns about non-integer
Warning messages:
1: non-integer x = -0.000000
2: non-integer x = -0.000020

## very small negatives were rounded to 0 in R 2.2.1 and earlier

I hope that this is helpful.  Thanks are due to Peter Dalgaard for
guidance.  So, thanks Peter :).

Cheers

Andrew

--
Brian D. Ripley,                  ripley at 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

--
Brian D. Ripley,                  ripley at 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

--
Brian D. Ripley,                  ripley at 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

--
Brian D. Ripley,                  ripley at 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