R devel - Feb 2005 - Re: [R-SIG-Mac] Bug running pbinom() in R-GUI?

On Feb 10, 2005, at 7:38 PM, George W. Gilchrist wrote:
> Today I was running a graduate level stats lab using R and we 
> encountered a
> major problem while using the current build of the Cocoa GUI:
>
>> From the GUI:
>> system.time(pbinom(80, 1e5, 806/1e6))
> [1] 14.37  4.94 30.29  0.00  0.00
>>
>
>> From the command line on the same machine:
>> system.time(pbinom(80, 1e5, 806/1e6))
> [1] 0.02 0.00 0.02 0.00 0.00
>>
>
> I've tried the CRAN version and the latest build of the R-GUI and both
> deliver the same terrible performance. It seems as if this only occurs 
> for
> certain arguments, but we saw this on ~15 different machines today. 
> And it
> seems to only be when running the Cocoa GUI. No problems at all with 
> this
> under Windoze. Any ideas?
The cause is pbeta_raw calling (indirectly via R_CheckUserInterrupt) 
R_ProcessEvents for every iteration - and for small p the number of 
iterations goes really high (e.g. for the case above n=99919). 
R_ProcessEvents is not a cheap operation, because it enters system 
event loop and processes any pending events. A quick fix could be the 
following patch, which checks for break only after several iterations 
(including the first one, just in case this is part of a sequence that 
may need user interaction).

Index: src/nmath/pbeta.c
==================================================================---
src/nmath/pbeta.c   (revision 33148)
+++ src/nmath/pbeta.c   (working copy)
@@ -139,7 +139,8 @@
             for(i= 1; i <= n; i++) {
  #ifndef MATHLIB_STANDALONE
                 /* for now, at least allow this:*/
-               R_CheckUserInterrupt();
+               if ((i&1023)==1)
+                       R_CheckUserInterrupt();
  #endif
                 if (p1 <= 1 && term / eps <= finsum)
                     break;

after this patch has been applied I get in the GUI:

 > system.time(pbinom(80, 1e5, 806/1e6))
[1] 0.02 0.00 0.08 0.00 0.00

Cheers,
Simon

The problem rather is that if R_CheckUserInterrupt is so expensive, we 
need to redesign it, for it should not be (and is not on other platforms 
as the comments below show).  See the comment in src/main/errors.c.

One idea might be for the GUI to set a flag that R_CheckUserInterrupt 
consults (which is what happens on Windows via another thread).

I agree that pbeta seems far too enthusiastic on checking, but the 
`weight' of R_CheckUserInterrupt needs to be comparable across platforms 
or this will recur.  Looks to me like it at least 100x more expensive on 
OS X than anywhere else.  (I timed 1us on Windows, and the timings below 
suggest 140us.)

(I think we need to understand why there is a check in that loop, rather 
than just choose some arbitrary frequency: looks to me that it is simple 
enough to check every million since pbeta(0.5, 1e6, 1e6) takes less than 
0.1s on my machine.  But integer overflow of n cuts in before it gets 
really slow, and I suspect that by n=1e6 it is better to use an 
approximation than a sum.)


On Fri, 11 Feb 2005, Simon Urbanek wrote:
> On Feb 10, 2005, at 7:38 PM, George W. Gilchrist wrote:
>
>> Today I was running a graduate level stats lab using R and we
encountered a
>> major problem while using the current build of the Cocoa GUI:
>> 
>>> From the GUI:
>>> system.time(pbinom(80, 1e5, 806/1e6))
>> [1] 14.37  4.94 30.29  0.00  0.00
>>> 
>> 
>>> From the command line on the same machine:
>>> system.time(pbinom(80, 1e5, 806/1e6))
>> [1] 0.02 0.00 0.02 0.00 0.00
>>> 
>> 
>> I've tried the CRAN version and the latest build of the R-GUI and
both
>> deliver the same terrible performance. It seems as if this only occurs
for
>> certain arguments, but we saw this on ~15 different machines today. And
it
>> seems to only be when running the Cocoa GUI. No problems at all with
this
>> under Windoze. Any ideas?
>
> The cause is pbeta_raw calling (indirectly via R_CheckUserInterrupt) 
> R_ProcessEvents for every iteration - and for small p the number of 
> iterations goes really high (e.g. for the case above n=99919). 
> R_ProcessEvents is not a cheap operation, because it enters system event
loop
> and processes any pending events. A quick fix could be the following patch,
> which checks for break only after several iterations (including the first 
> one, just in case this is part of a sequence that may need user
interaction).
>
> Index: src/nmath/pbeta.c
> ==================================================================> ---
src/nmath/pbeta.c   (revision 33148)
> +++ src/nmath/pbeta.c   (working copy)
> @@ -139,7 +139,8 @@
>            for(i= 1; i <= n; i++) {
> #ifndef MATHLIB_STANDALONE
>                /* for now, at least allow this:*/
> -               R_CheckUserInterrupt();
> +               if ((i&1023)==1)
> +                       R_CheckUserInterrupt();
> #endif
>                if (p1 <= 1 && term / eps <= finsum)
>                    break;
>
> after this patch has been applied I get in the GUI:
>
>> system.time(pbinom(80, 1e5, 806/1e6))
> [1] 0.02 0.00 0.08 0.00 0.00
>
> Cheers,
> Simon
>
> ______________________________________________
> R-devel@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
>
-- 
Brian D. Ripley,                  ripley@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