statnet.common (version 4.3.0)

logspace.utils: Utilities for performing calculations on logarithmic scale.

Description

A small suite of functions to compute sums, means, and weighted means on logarithmic scale, minimizing loss of precision.

Usage

log_sum_exp(logx, use_ldouble = FALSE)
log_mean_exp(logx, use_ldouble = FALSE)
lweighted.mean(x, logw)
lweighted.var(x, logw)

Arguments

logx

Numeric vector of \(\log(x)\), the natural logarithms of the values to be summed or averaged.

use_ldouble

Whether to use long double precision in the calculation. If TRUE, 's C built-in logspace_sum() is used. If FALSE, the package's own implementation based on it is used, using double precision, which is (on most systems) several times faster, at the cost of precision.

x

Numeric vector of \(x\), the (raw) values to be summed or averaged. For lweighted.mean, x may also be a matrix, in which case the weighted mean will be computed for each column of x.

logw

Numeric vector of \(\log(w)\), the natural logarithms of the weights.

Value

The functions return the equivalents of the following R expressions, but faster and with less loss of precision:

log_sum_exp(logx)

log(sum(exp(logx)))

log_mean_exp(logx)

log(mean(exp(logx)))

lweighted.mean(x,logw)

sum(x*exp(logw))/sum(exp(logw)) for x scalar and colSums(x*exp(logw))/sum(exp(logw)) for x matrix

lweighted.var(x,logw)

crossprod(x*exp(logw/2))/sum(exp(logw))

Examples

# NOT RUN {
logx <- rnorm(1000)
stopifnot(all.equal(log(sum(exp(logx))), log_sum_exp(logx)))
stopifnot(all.equal(log(mean(exp(logx))), log_mean_exp(logx)))

x <- rnorm(1000)
logw <- rnorm(1000)
stopifnot(all.equal(m <- sum(x*exp(logw))/sum(exp(logw)),lweighted.mean(x, logw)))
stopifnot(all.equal(sum((x-m)^2*exp(logw))/sum(exp(logw)),
                    lweighted.var(x, logw), check.attributes=FALSE))

x <- cbind(x, rnorm(1000))
stopifnot(all.equal(m <- colSums(x*exp(logw))/sum(exp(logw)),
                    lweighted.mean(x, logw), check.attributes=FALSE))
stopifnot(all.equal(crossprod(t(t(x)-m)*exp(logw/2))/sum(exp(logw)),
                    lweighted.var(x, logw), check.attributes=FALSE))
# }