Calculate the Riemannian Hessian of Karcher mean problem on positive definite matrices

Question

Consider a collection of positive definite matrices $\{A_1,...,A_n\}\in\mathbb{S}_{++}^d$, the Karcher mean of these matrices is given by (see (5.4) in [1]): $$ \min_{X\in\mathbb{S}_{++}^d} f(X):=\frac{1}{n}\sum_{i=1}^{n}\|\ln(X^{-1/2}A_i X^{-1/2})\|_F^2 $$ where $f$ is simply the sum of squared geodesic distances of $X$ to all the $A_i$s.

The Riemannian gradient of $f$ is given as (5.5) in [1]: $$ \mathrm{grad} f(X) = \frac{1}{n}\sum_{i=1}^{n} X^{1/2}\ln(X^{1/2}A_i^{-1}X^{1/2})X^{1/2} $$

My question is how to calculate the Riemannian Hessian (or Hessian-vector product) of this function $f$, if it enjoys closed form?

(A more generic question is that can the Riemannian Hessian of squared geodesic distance be expressed in closed form by operators on Riemannian manifolds, such as exponential and logarithm mappings?)

Any comment or reference are appreciated!

References: [1] O. P. Ferreira, M. S. Louzeiro, and L.F. Prudente. "Gradient method for optimization on Riemannian manifolds with lower bounded curvature." SIAM Journal on Optimization 29.4 (2019): 2517-2541.

Answer 1

I came up with an ad-hoc solution to calculate the Hessian-vector product below. I'm still looking forward to see if there are any comments on the more general question, i.e. calculating the Riemannian Hessian of the squared geodesic distance.

First denote $\mathrm{dist}(A, B):=\|\log(A^{-1/2}B A^{-1/2})\|_{F}$, which is the geodesic distance of two positive definite matrices (see (1), Chpater 6), then we just need to calculate the Euclidean Hessian vector product $\nabla^2 h(S)[V]$ of $h(S):=\mathrm{dist}(S, A)^2$ for any symmetric matrix $V$, and the Riemannian Hessian $\mathrm{H} h(S)[V]$ can be calculated by (see (2)) $$ \mathrm{H} h(S)[V] = S \nabla_S^2 h(S)[V] S +\mathrm{sym}(S \nabla_S h(S) V) $$ with $\mathrm{sym}(A)=(A+A^\top) / 2$.

Notice that (see (1), Chapter 6) $$ \nabla_S h(S) = S^{-1/2}\log(S^{1/2}A^{-1} S^{1/2}) S^{-1/2} = S^{-1} A^{1/2} \log(A^{-1/2} S A^{-1/2})A^{-1/2} $$ For any symmetric matrix $V$, we have $$ \langle\nabla_S h(S), V\rangle = \mathrm{tr} (V S^{-1} A^{1/2} \log(A^{-1/2} S A^{-1/2})A^{-1/2}) $$ To take the derivative of this, notice that $S$ appear twice in $\langle\nabla_S h(S), V\rangle$. Denote $\tilde{h}(S_1, S_2) = \mathrm{tr} (V S_1^{-1} A^{1/2} \log(A^{-1/2} S_2 A^{-1/2})A^{-1/2})$, we know that $$ \frac{\partial \tilde{h}}{\partial S_1} = -S_1^{-1} V A^{-1/2} \log(A^{-1/2} S_2 A^{-1/2})A^{1/2} S_1^{-1} $$

It remains to calculate ${\partial \tilde{h}} / {\partial S_2}$, which takes a form of $l(S):=\mathrm{tr}(C \log(P S Q))$. To the best of my knowledge, this function doesn't enjoy closed-form derivatives...

However I was inspired by this post to proceed my calculation. Denote $Y=PSQ$ and $L=\log(Y)$, we have $$ \begin{bmatrix} L & d L\\ 0 & L \end{bmatrix} = \log\left(\begin{bmatrix} Y & d Y\\ 0 & Y \end{bmatrix}\right) = \log\left(\begin{bmatrix} P & 0\\ 0 & P \end{bmatrix}\begin{bmatrix} S & d S\\ 0 & S \end{bmatrix}\begin{bmatrix} Q & 0\\ 0 & Q \end{bmatrix}\right) $$ therefore we get $$ d l = \left\langle \begin{bmatrix} 0 & C\\ 0 & 0 \end{bmatrix}, \log\left(\begin{bmatrix} P & 0\\ 0 & P \end{bmatrix}\begin{bmatrix} S & d S\\ 0 & S \end{bmatrix}\begin{bmatrix} Q & 0\\ 0 & Q \end{bmatrix}\right) \right\rangle $$ where the inner product is simply the Euclidean inner product. We can plug $d S$ as standard Euclidean basis to obtain an representation of ${d l} / {d S}$, which will take $\mathcal{O}(d^2)$ number of times to cover all the entries. Nevertheless this provides an implementable way to calculate the Euclidean and Riemannian Hessian.

To summarize, the Euclidean and Riemannian Hessian of $h$ are thus calculated as follows \begin{equation*} \begin{aligned} & \nabla^2 h(S)[V] = -S^{-1} V A^{-1/2} \log(A^{-1/2} S A^{-1/2})A^{1/2} S^{-1} + L \\ & \mathrm{H} (h(S))[V] = S \nabla^2 h(S)[V] S +\mathrm{sym}(S \nabla_S h(S) V) \end{aligned} \end{equation*} where each entry of matrix $L$ is calculated as follows: $$ L_{i, j} = \left\langle \begin{bmatrix} 0 & C\\ 0 & 0 \end{bmatrix}, \log\left(\begin{bmatrix} P & 0\\ 0 & P \end{bmatrix}\begin{bmatrix} S & E_{i,j}\\ 0 & S \end{bmatrix}\begin{bmatrix} Q & 0\\ 0 & Q \end{bmatrix}\right) \right\rangle $$ where $E_{i,j}\in\mathbb{R}^{d\times d}$ is the matrix which is all zero and only one at $i,j$-th entry, and \begin{equation*} \begin{aligned} & P = A^{-1/2}, \\ & Q = A^{-1/2}, \\ & C = A^{-1/2} V S^{-1} A^{1/2}. \end{aligned} \end{equation*}

References:

(1) Bhatia, Rajendra. Positive definite matrices. Princeton university press, 2009.

(2) Boumal, Nicolas. An introduction to optimization on smooth manifolds. Cambridge University Press, 2023.