pdMean {pdSpecEst}

Weighted Karcher mean of HPD matrices

Description

pdMean calculates an (approximate) weighted Karcher or Frechet mean of a sample of (d,d)-dimensional HPD matrices intrinsic to a user-specified metric. In the case of the affine-invariant Riemannian metric as detailed in e.g., (Bhatia 2009)[Chapter 6] or (Pennec et al. 2006), the weighted Karcher mean is either approximated via the fast recursive algorithm in (Ho et al. 2013) or computed via the slower, but more accurate, gradient descent algorithm in (Pennec 2006). By default, the unweighted Karcher mean is computed.

Usage

pdMean(M, w, metric = "Riemannian", grad_desc = FALSE, maxit = 1000,
  reltol)

Arguments

M

a (d,d,S)-dimensional array corresponding to a sample of (d,d)-dimensional HPD matrices of size S.

w

an S-dimensional nonnegative weight vector, such that sum(w) = 1.

metric

the distance measure, one of 'Riemannian', 'logEuclidean', 'Cholesky', 'Euclidean' or 'rootEuclidean'. Defaults to 'Riemannian'.

grad_desc

if metric = "Riemannian", a logical value indicating if the gradient descent algorithm in (Pennec 2006) should be used, defaults to FALSE.

maxit

maximum number of iterations in gradient descent algorithm, only used if grad_desc = TRUE and metric = "Riemannian". Defaults to 1000

reltol

optional tolerance parameter in gradient descent algorithm, only used if grad_desc = TRUE and metric = "Riemannian". Defaults to 1E-10.

Note

The function does not check for positive definiteness of the input matrices, and (depending on the specified metric) may fail if matrices are close to being singular.

References

Bhatia R (2009). Positive Definite Matrices. Princeton University Press, New Jersey.

Ho J, Cheng G, Salehian H, Vemuri B (2013). “Recursive Karcher expectation estimators and recursive law of large numbers.” Artificial Intelligence and Statistics, 325–332.

Pennec X (2006). “Intrinsic statistics on Riemannian manifolds: Basic tools for geometric measurements.” Journal of Mathematical Imaging and Vision, 25(1), 127–154.

Pennec X, Fillard P, Ayache N (2006). “A Riemannian framework for tensor computing.” International Journal of Computer Vision, 66(1), 41–66.

See Also

Mid, pdMedian

Examples

## Generate random sample of HPD matrices
m <- function(){
 X <- matrix(complex(real=rnorm(9), imaginary=rnorm(9)), nrow=3)
 t(Conj(X)) %*% X
}
M <- replicate(100, m())
z <- rnorm(100)
## Generate random weight vector
w <- abs(z)/sum(abs(z))
## Compute weighted (Riemannian) Karcher mean
pdMean(M, w)