n_10;p_0.5;jjx_0:n;qbinom(pbinom(jjx,n,p),n,p) # This one works as expected n_100;p_0.5;jjx_0:n;qbinom(pbinom(jjx,n,p),n,p) # This one causes severe problems I cannot interrupt using ESC and I finally have to resort to the Windows Task manager to kill the R session. A friend of mine told me that he faced similar problems under Unix. --please do not edit the information below-- Version: platform = Windows arch = x86 os = Win32 system = x86, Win32 status = major = 1 minor = 0.0 year = 2000 month = February day = 29 language = R Windows NT 4.0 (build 1381) Service Pack 3 Search Path: .GlobalEnv, Autoloads, package:base --- D.Trenkler --- ************************************************************************ ********* Dietrich Trenkler (trenkler@oec.uni-osnabrueck.de) Statistik / Empirische Wirtschaftsforschung Universitaet Osnabrueck Rolandstrasse 8 Phone: +49(0) 541-969-2753 D-49069 Osnabrueck Fax : +49(0) 541-969-2745 GERMANY ************************************************************************ ********* -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-devel-request@stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
On Fri, 7 Apr 2000 dtrenkler@nts6.oec.uni-osnabrueck.de wrote:> n_10;p_0.5;jjx_0:n;qbinom(pbinom(jjx,n,p),n,p) # This one works as > expected > > n_100;p_0.5;jjx_0:n;qbinom(pbinom(jjx,n,p),n,p) # This one causes > severe problemsThe problem is in qbinom(1-epsilon,100,0.5) for epsilon between 10^-15 and 10^-17. It appears to go for ever.> pbinom(jjx,n,p)[89]-1[1] -1.110223e-16 An even simpler version is qbinom(1-1e-16,2,0.5) Presumably a check for zero needs some fuzz. -thomas Thomas Lumley Assistant Professor, Biostatistics University of Washington, Seattle -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-devel-request@stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
On Fri, 7 Apr 2000 dtrenkler@nts6.oec.uni-osnabrueck.de wrote:> n_10;p_0.5;jjx_0:n;qbinom(pbinom(jjx,n,p),n,p) # This one works as > expected > > n_100;p_0.5;jjx_0:n;qbinom(pbinom(jjx,n,p),n,p) # This one causes > severe problemsIn qbinom we calculate an estimate y based on a Cornish-Fisher expansion and then search up or down to the right number. However, for p=1-1e-16 or thereabouts we get y=83931130944911952 which is bigger than 1/machine epsilon. Thus y=y-1 and searching up or down doesnt work. A crude solution is to just treat this as p=1, but it would be better to search up or down in bigger steps. -thomas -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-devel-request@stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
>>>>> "TL" == Thomas Lumley <thomas@biostat.washington.edu> writes:TL> On Fri, 7 Apr 2000 dtrenkler@nts6.oec.uni-osnabrueck.de wrote: >> n_10;p_0.5;jjx_0:n;qbinom(pbinom(jjx,n,p),n,p) # This one works as >> expected >> >> n_100;p_0.5;jjx_0:n;qbinom(pbinom(jjx,n,p),n,p) # This one causes >> severe problems TL> The problem is in qbinom(1-epsilon,100,0.5) for epsilon between TL> 10^-15 and 10^-17. It appears to go for ever. >> pbinom(jjx,n,p)[89]-1 TL> [1] -1.110223e-16 TL> An even simpler version is TL> qbinom(1-1e-16,2,0.5) TL> Presumably a check for zero needs some fuzz. Note that we (I) also have a considerable improvement for dbinom() in the pipeline, from Clive Loader's post of a few weeks ago. But this may be a completely orthogonal issue... I will have a look over the weekend. Martin -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-devel-request@stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
I've committed a fix to the development branch. We now start searching at max(0,min(n,normal-approximation)) We still don't get qbinom(pbinom(x,100,0.5))==x for x>86, but that's to be expected. Thomas Lumley Assistant Professor, Biostatistics University of Washington, Seattle -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-devel-request@stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._