Become an expert in R — Interactive courses, Cheat Sheets, certificates and more!
Get Started for Free

logSumExp

Accurately computes the logarithm of the sum of exponentials


Description

Accurately computes the logarithm of the sum of exponentials, that is, log(sum(exp(lx))). If lx = log(x), then this is equivalently to calculating log(sum(x)).

Usage

logSumExp(lx, idxs = NULL, na.rm = FALSE, ...)

Arguments

lx

A numeric vector. Typically lx are log(x) values.

idxs

A vector indicating subset of elements to operate over. If NULL, no subsetting is done.

na.rm

If TRUE, any missing values are ignored, otherwise not.

...

Not used.

Details

This function, which avoid numerical underflow, is often used when computing the logarithm of the sum of small numbers (|x| << 1) such as probabilities.

This is function is more accurate than log(sum(exp(lx))) when the values of x = exp(lx) are |x| << 1. The implementation of this function is based on the observation that

log(a + b) = [ la = log(a), lb = log(b) ] = log( exp(la) + exp(lb) ) = la + log ( 1 + exp(lb - la) )

Assuming la > lb, then |lb - la| < |lb|, and it is less likely that the computation of 1 + exp(lb - la) will not underflow/overflow numerically. Because of this, the overall result from this function should be more accurate. Analogously to this, the implementation of this function finds the maximum value of lx and subtracts it from the remaining values in lx.

Value

Returns a numeric scalar.

Benchmarking

This method is optimized for correctness, that avoiding underflowing. It is implemented in native code that is optimized for speed and memory.

Author(s)

Henrik Bengtsson

References

[1] R Core Team, Writing R Extensions, v3.0.0, April 2013.
[2] Laurent El Ghaoui, Hyper-Textbook: Optimization Models and Applications, University of California at Berkeley, August 2012. (Chapter 'Log-Sum-Exp (LSE) Function and Properties')
[3] R-help thread logsumexp function in R, 2011-02-17. https://stat.ethz.ch/pipermail/r-help/2011-February/269205.html

See Also

To compute this function on rows or columns of a matrix, see rowLogSumExps().

For adding two double values in native code, R provides the C function logspace_add() [1]. For properties of the log-sum-exponential function, see [2].

Examples

## EXAMPLE #1
lx <- c(1000.01, 1000.02)
y0 <- log(sum(exp(lx)))
print(y0) ## Inf

y1 <- logSumExp(lx)
print(y1) ## 1000.708


## EXAMPLE #2
lx <- c(-1000.01, -1000.02)
y0 <- log(sum(exp(lx)))
print(y0) ## -Inf

y1 <- logSumExp(lx)
print(y1) ## -999.3218


## EXAMPLE #3
## R-help thread 'Beyond double-precision?' on May 9, 2009.

set.seed(1)
x <- runif(50)

## The logarithm of the harmonic mean
y0 <- log(1 / mean(1 / x))
print(y0)  ## -1.600885

lx <- log(x)
y1 <- log(length(x)) - logSumExp(-lx)
print(y1)  ## [1] -1.600885

# Sanity check
stopifnot(all.equal(y1, y0))

matrixStats

Functions that Apply to Rows and Columns of Matrices (and to Vectors)

v0.58.0
Artistic-2.0
Authors
Henrik Bengtsson [aut, cre, cph], Constantin Ahlmann-Eltze [ctb], Hector Corrada Bravo [ctb], Robert Gentleman [ctb], Jan Gleixner [ctb], Peter Hickey [ctb], Ola Hossjer [ctb], Harris Jaffee [ctb], Dongcan Jiang [ctb], Peter Langfelder [ctb], Brian Montgomery [ctb], Hugh Parsonage [ctb]
Initial release

We don't support your browser anymore

Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.