[R] over/under flow

Martin Maechler maechler at stat.math.ethz.ch
Mon Jul 4 09:33:31 CEST 2005


>>>>> "William" == William H Asquith <wasquith at austin.rr.com>
>>>>>     on Sun, 3 Jul 2005 16:00:40 -0500 writes:

    William> Great, but to followup, how do I select the bounds (B,E) on the root 
    William> for an arbitrary machine?

    William> OVER = uniroot(function(x) lgamma(x)-log(.Machine$double.xmax), 
    William> c(B,E))$root

    William> If uniroot() is fast enough, is it appropriate for me to set B at say 1 
    William> and E at log(.Machine$double.max)?  Suggestions on do this the proper 
    William> "R way"?  Perhaps this . . .

    William> OVER = uniroot(function(x) lgamma(x)-log(.Machine$double.xmax), 
    William> c(1,log(.Machine$double.xmax)))$root

    William> I am working on a package so different machines will be involved thus 
    William> simple 171,172 might not be the best idea for the root?

Well,

1) We nowadays *require* compliance with
the usual ``IEEE arithmetic standard''  {put a bit too simplistically}
and I'd bet that those compiler / runtime library combinations
that can build R properly, would also give cutoffs in this same interval

2) lgamma(x) is very linear ``out there'' -- try
      curve(lgamma(x) - log(.Machine$double.xmax), 170, 180)
   --    so you can easily enlarge the interval a bit, e.g. to
   (170,175) or even to (100,200)

Martin


    William> On Jul 3, 2005, at 3:43 PM, Peter Dalgaard wrote:

    >> "William H. Asquith" <wasquith at austin.rr.com> writes:
    >> 
    >>> I am porting some FORTRAN to R in which an Inf triggers an if().  The
    >>> trigger is infinite on exp(lgamma(OVER)).  What is the canonical R
    >>> style of determining OVER when exp(OVER)== Inf?  The code structure
    >>> that I am
    >>> porting is best left intact--so I need to query R somehow to the value
    >>> of OVER that causes exp(lgamma(OVER)) to equal Inf.
    >>> 
    >>> On my system,
    >>> exp(lgamma(171)) is about first to equal Inf.
    >>> 
    >>> I asked similar question a few weeks ago on exp(OVER) and got the
    >>> answer back as log(.Machine$double.xmax).  I now have the lgamma
    >>> involved.  I think that answer is what is OVER such the
    >>> 
    >>> .Machine$double.xmax = lgamma(OVER),
    >> 
    >> Not quite... (see below)
    >> 
    >>> but I am not sure how to invert or solve for OVER
    >> 
    >> 
    >>> uniroot(function(x) lgamma(x)-log(.Machine$double.xmax), c(171,172))
    >> $root
    >> [1] 171.6244
    >> 
    >> $f.root
    >> [1] -1.462051e-07
    >> 
    >> $iter
    >> [1] 3
    >> 
    >> $estim.prec
    >> [1] 6.103516e-05
    >> 
    >> 
    >> -- 
    >> O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
    >> c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
    >> (*) \(*) -- University of Copenhagen   Denmark          Ph: (+45) 
    >> 35327918
    >> ~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk)                  FAX: (+45) 
    >> 35327907
    >> 

    William> ______________________________________________
    William> R-help at stat.math.ethz.ch mailing list
    William> https://stat.ethz.ch/mailman/listinfo/r-help
    William> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html


    William> !DSPAM:42c852c3295961260279199!




More information about the R-help mailing list