logspace.utils {statnet.common} | R Documentation |
A small suite of functions to compute sums, means, and weighted means on logarithmic scale, minimizing loss of precision.
log_sum_exp(logx, use_ldouble = FALSE)
log_mean_exp(logx, use_ldouble = FALSE)
lweighted.mean(x, logw)
lweighted.var(x, logw)
logx |
Numeric vector of |
use_ldouble |
Whether to use |
x |
Numeric vector of |
logw |
Numeric vector of |
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))
Pavel N. Krivitsky
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))