R devel - Nov 2013 - R interface to C API Rf_logspace_{add,sub}?

Is there an R-language interface to the R API C-language functions
Rf_logspace_add()
and Rf_logspace_sub()?   I don't see one but I may not looking under the
right name.

Various packages have functions which do that same sort
of thing (log(exp(x)+exp(y)) and log(exp(x)-exp(y)) without unnecessary
floating point errors).  They have names like
    matrixStats::logSumExp(lx, na.rm=FALSE, ...) (the ... is ignored by the
function)
    sna::logSub(x,y), logSum(x), logMean(x)
    BTYD::addLogs(loga, logb) and subLogs(loga, logb)
Googling for logSumExp (the Python name) indicates that many know this as
"LSE"
(the "Log-Sum-Exp trick").

I've seen several instances in R-help recently where user code could be
made more accurate if these were available.  They can be written in pure
R code, but would be slower because they would require something like
pmin() and pmax() to work.  I think they belong in the same package as log1p()
and expm1().  I think the add and subtraction versions are the most useful
and the sum version would be very nice to have.

Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com

>>>>> William Dunlap <wdunlap at tibco.com>
>>>>>     on Thu, 7 Nov 2013 21:32:46 +0000 writes:
    > Is there an R-language interface to the R API C-language
    > functions Rf_logspace_add() and Rf_logspace_sub()?  I
    > don't see one but I may not looking under the right name.

    > Various packages have functions which do that same sort of
    > thing (log(exp(x)+exp(y)) and log(exp(x)-exp(y)) without
    > unnecessary floating point errors).  They have names like
    > matrixStats::logSumExp(lx, na.rm=FALSE, ...) (the ... is
    > ignored by the function) sna::logSub(x,y), logSum(x),
    > logMean(x) BTYD::addLogs(loga, logb) and subLogs(loga,
    > logb) Googling for logSumExp (the Python name) indicates
    > that many know this as "LSE" (the "Log-Sum-Exp
trick").

    > I've seen several instances in R-help recently where user
    > code could be made more accurate if these were available.
    > They can be written in pure R code, but would be slower
    > because they would require something like pmin() and
    > pmax() to work.  I think they belong in the same package
    > as log1p() and expm1().  I think the add and subtraction
    > versions are the most useful and the sum version would be
    > very nice to have.

I agree.... As a matter of fact,
Marius Hofert and I have been using pure R "vectorized", i.e.,
working on a *matrix* argument, returning a vector,
versions of these in our package 'copula' 
 (--> copula/R/special-func.R ), 
not exported though,
however somewhat documented and explained on page 49 of

 M Hofert, M M?chler, AJ McNeil (2013)
 Archimedean Copulas in High Dimensions:
 Estimators and Numerical Challenges Motivated by Financial Applications
 Journal de la Soci?t? Fran?aise de Statistique 154 (1), 25-63

p.49 is also attached.

The clue is that the  lssum() function is really an n-term
generalization of logspace_add() and logspace_sub()
by being able to specifying the *signs* and lsum() is a special
case with all signs being +1.

Before jumping ahead and providing a faster R interface in
"base" (or rather "utils" ?), I think we should think twice
before embarking on a 2-argument version only rather than a
vectorized version as our lsum() and lssum() provide.

Martin Maechler, ETH Zurich


-------------- next part --------------
A non-text attachment was scrubbed...
Name: HoMaMcN2013-p49.pdf
Type: application/pdf
Size: 117551 bytes
Desc: not available
URL:
<https://stat.ethz.ch/pipermail/r-devel/attachments/20131111/bc142e21/attachment.pdf>
-------------- next part --------------



    > Bill Dunlap Spotfire, TIBCO Software wdunlap tibco.com