pdNeville {pdSpecEst}

Polynomial interpolation of curves (1D) or surfaces (2D) of HPD matrices

Description

pdNeville performs intrinsic polynomial interpolation of curves or surfaces of HPD matrices in the metric space of HPD matrices equipped with the affine-invariant Riemannian metric (see (Bhatia 2009)[Chapter 6] or (Pennec et al. 2006)) via Neville's algorithm based on iterative geodesic interpolation detailed in (Chau and von Sachs 2019) and in Chapter 3 and 5 of (Chau 2018).

Usage

pdNeville(P, X, x, metric = "Riemannian")

Arguments

P

for polynomial curve interpolation, a (d, d, N)-dimensional array corresponding to a length N sequence of (d, d)-dimensional HPD matrices (control points) through which the interpolating polynomial curve passes. For polynomial surface interpolation, a (d, d, N_1, N_2)-dimensional array corresponding to a tensor product grid of (d, d)-dimensional HPD matrices (control points) through which the interpolating polynomial surface passes.

X

for polynomial curve interpolation, a numeric vector of length N specifying the time points at which the interpolating polynomial passes through the control points P. For polynomial surface interpolation, a list with as elements two numeric vectors x and y of length N_1 and N_2 respectively. The numeric vectors specify the time points on the tensor product grid expand.grid(X$x, X$y) at which the interpolating polynomial passes trough the control points P.

x

for polynomial curve interpolation, a numeric vector specifying the time points (locations) at which the interpolating polynomial is evaluated. For polynomial surface interpolation, a list with as elements two numeric vectors x and y specifying the time points (locations) on the tensor product grid expand.grid(x$x, x$y) at which the interpolating polynomial surface is evaluated.

metric

the metric on the space of HPD matrices, by default metric = "Riemannian", but instead of the Riemannian metric this can also be set to metric = "Euclidean" to perform (standard) Euclidean polynomial interpolation of curves or surfaces in the space of HPD matrices.

Details

For polynomial curve interpolation, given N control points (i.e., HPD matrices), the degree of the interpolated polynomial is N - 1. For polynomial surface interpolation, given N_1 \times N_2 control points (i.e., HPD matrices) on a tensor product grid, the interpolated polynomial surface is of bi-degree (N_1 - 1, N_2 - 1). Depending on the input array P, the function decides whether polynomial curve or polynomial surface interpolation is performed.

Value

For polynomial curve interpolation, a (d, d, length(x))-dimensional array corresponding to the interpolating polynomial curve of (d,d)-dimensional matrices of degree N-1 evaluated at times x and passing through the control points P at times X. For polynomial surface interpolation, a (d, d, length(x$x), length(x$y))-dimensional array corresponding to the interpolating polynomial surface of (d,d)-dimensional matrices of bi-degree N_1 - 1, N_2 - 1 evaluated at times expand.grid(x$x, x$y) and passing through the control points P at times expand.grid(X$x, X$y).

Note

If metric = "Euclidean", the interpolating curve or surface may not be positive definite everywhere as the space of HPD matrices equipped with the Euclidean metric has its boundary at a finite distance.

The function does not check for positive definiteness of the input matrices, and may fail if metric = "Riemannian" and the input matrices are close to being singular.

References

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

Chau J (2018). Advances in Spectral Analysis for Multivariate, Nonstationary and Replicated Time Series. phdthesis, Universite catholique de Louvain.

Chau J, von Sachs R (2019). “Intrinsic wavelet regression for curves of Hermitian positive definite matrices.” Journal of the American Statistical Association. doi: 10.1080/01621459.2019.1700129.

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

See Also

pdPolynomial

Examples

### Polynomial curve interpolation
P <- rExamples1D(50, example = 'gaussian')$f[, , 10*(1:5)]
P.poly <- pdNeville(P, (1:5)/5, (1:50)/50)
## Examine matrix-component (1,1)
plot((1:50)/50, Re(P.poly[1, 1, ]), type = "l") ## interpolated polynomial
lines((1:5)/5, Re(P[1, 1, ]), col = 2) ## control points

### Polynomial surface interpolation
P.surf <- array(P[, , 1:4], dim = c(2,2,2,2)) ## control points
P.poly <- pdNeville(P.surf, list(x = c(0, 1), y = c(0, 1)), list(x = (0:10)/10, y = (0:10)/10))