Skip to contents

Calculate a generalized logarithmic mean / extended mean.

Usage

extended_mean(r, s)

generalized_logmean(r)

logmean(a, b, tol = .Machine$double.eps^0.5)

Arguments

r, s

A finite number giving the order of the generalized logarithmic mean / extended mean.

a, b

A strictly positive numeric vector.

tol

The tolerance used to determine if a == b.

Value

generalized_logmean() and extended_mean() return a function:

function(a, b, tol = .Machine$double.eps^0.5){...}

This computes the component-wise generalized logarithmic mean of order r, or the extended mean of orders r and s, of a and b.

logmean() returns a numeric vector, the same length as max(length(a), length(b)), giving the component-wise logarithmic mean of a and b.

Details

The function extended_mean() returns a function to compute the component-wise extended mean of a and b of orders r and s. See Bullen (2003, p. 393) for a definition. This is also called the difference mean, Stolarsky mean, or extended mean-value mean.

The function generalized_logmean() returns a function to compute the component-wise generalized logarithmic mean of a and b of order r. See Bullen (2003, p. 385) for a definition, or https://en.wikipedia.org/wiki/Stolarsky_mean. The generalized logarithmic mean is a special case of the extended mean, corresponding to extended_mean(r, 1)(), but is more commonly used for price indexes.

The function logmean() returns the ordinary component-wise logarithmic mean of a and b, and corresponds to generalized_logmean(1)().

Both a and b should be strictly positive. This is not enforced, but the results may not make sense when the generalized logarithmic mean / extended mean is not defined. The usual recycling rules apply when a and b are not the same length.

By definition, the generalized logarithmic mean / extended mean of a and b is a when a == b. The tol argument is used to test equality by checking if abs(a - b) <= tol. The default value is the same as all.equal(). Setting tol = 0 tests for exact equality, but can give misleading results when a and b are computed values. In some cases it's useful to multiply tol by a scale factor, such as max(abs(a), abs(b)). This often doesn't matter when making price indexes, however, as a and b are usually around 1.

Note

generalized_logmean() can be defined on the extended real line, so that r = -Inf / Inf returns pmin()/pmax(), to agree with the definition in, e.g., Bullen (2003). This is not implemented, and r must be finite as in the original formulation by Stolarsky (1975).

References

Balk, B. M. (2008). Price and Quantity Index Numbers. Cambridge University Press.

Bullen, P. S. (2003). Handbook of Means and Their Inequalities. Springer Science+Business Media.

Stolarsky, K. B. (1975). Generalizations of the Logarithmic Mean. Mathematics Magazine, 48(2): 87-92.

See also

transmute_weights() uses the extended mean to turn a generalized mean of order \(r\) into a generalized mean of order \(s\).

Other means: generalized_mean(), lehmer_mean(), nested_mean()

Examples

x <- 8:5
y <- 1:4

#---- Comparing logarithmic means and generalized means ----

# The arithmetic and geometric means are special cases of the
# generalized logarithmic mean

all.equal(generalized_logmean(2)(x, y), (x + y) / 2)
#> [1] TRUE
all.equal(generalized_logmean(-1)(x, y), sqrt(x * y))
#> [1] TRUE

# The logarithmic mean lies between the arithmetic and geometric means
# because the generalized logarithmic mean is increasing in r

all(logmean(x, y) < (x + y) / 2) &
  all(logmean(x, y) > sqrt(x * y))
#> [1] TRUE

# The harmonic mean cannot be expressed as a logarithmic mean, but can
# be expressed as an extended mean

all.equal(extended_mean(-2, -1)(x, y), 2 / (1 / x + 1 / y))
#> [1] TRUE

# The quadratic mean is also a type of extended mean

all.equal(extended_mean(2, 4)(x, y), sqrt(x^2 / 2 + y^2 / 2))
#> [1] TRUE

# As are heronian and centroidal means

all.equal(
  extended_mean(0.5, 1.5)(x, y),
  (x + sqrt(x * y) + y) / 3
)
#> [1] TRUE
all.equal(
  extended_mean(2, 3)(x, y),
  2 / 3 * (x^2 + x * y + y^2) / (x + y)
)
#> [1] TRUE

#---- Approximating the logarithmic mean ----

# The logarithmic mean can be approximated as a convex combination of
# the arithmetic and geometric means that gives more weight to the
# geometric mean

approx1 <- 1 / 3 * (x + y) / 2 + 2 / 3 * sqrt(x * y)
approx2 <- ((x + y) / 2)^(1 / 3) * (sqrt(x * y))^(2 / 3)

approx1 - logmean(x, y) # always a positive approximation error
#> [1] 1.932965e-02 3.260257e-03 3.420021e-04 3.852275e-06
approx2 - logmean(x, y) # a negative approximation error
#> [1] -6.436118e-02 -1.212079e-02 -1.336412e-03 -1.537117e-05

# A better approximation

correction <- (log(x / y) / pi)^4 / 32
approx1 / (1 + correction) - logmean(x, y)
#> [1] -8.576372e-04  1.064231e-04  2.148563e-05  2.877344e-07

#---- Some identities ----

# A useful identity for turning an additive change into a proportionate
# change

all.equal(logmean(x, y) * log(x / y), x - y)
#> [1] TRUE

# Works for other orders, too

r <- 2

all.equal(
  generalized_logmean(r)(x, y)^(r - 1) * (r * (x - y)),
  (x^r - y^r)
)
#> [1] TRUE

# Some other identities

all.equal(
  generalized_logmean(-2)(1, 2),
  (harmonic_mean(1:2) * geometric_mean(1:2)^2)^(1 / 3)
)
#> [1] TRUE

all.equal(
  generalized_logmean(0.5)(1, 2),
  (arithmetic_mean(1:2) + geometric_mean(1:2)) / 2
)
#> [1] TRUE

all.equal(
  logmean(1, 2),
  geometric_mean(1:2)^2 * logmean(1, 1 / 2)
)
#> [1] TRUE

#---- Integral representations of the logarithmic mean ----

logmean(2, 3)
#> [1] 2.466303

integrate(function(t) 2^(1 - t) * 3^t, 0, 1)$value
#> [1] 2.466303
1 / integrate(function(t) 1 / (2 * (1 - t) + 3 * t), 0, 1)$value
#> [1] 2.466303