On Riemannian Optimization over Positive Definite Matrices with the Bures-Wasserstein Geometry
Authors: Andi Han
Bamdev Mishra
Pratik Jawanpuria
Junbin Gao
Venue: NeurIPS 2021 Poster (2021)
Keywords: Riemannian optimization
Riemannian manifold
Symmetric Positive Definite
Bures-Wasserstein
Affine-Invariant
geodesic convexity
non-negative curvature

Abstract: In this paper, we comparatively analyze the Bures-Wasserstein (BW) geometry with the popular Affine-Invariant (AI) geometry for Riemannian optimization on the symmetric positive definite (SPD) matrix manifold. Our study begins with an observation that the BW metric has a linear dependence on SPD matrices in contrast to the quadratic dependence of the AI metric. We build on this to show that the BW metric is a more suitable and robust choice for several Riemannian optimization problems over ill-conditioned SPD matrices. We show that the BW geometry has a non-negative curvature, which further improves convergence rates of algorithms over the non-positively curved AI geometry. Finally, we verify that several popular cost functions, which are known to be geodesic convex under the AI geometry, are also geodesic convex under the BW geometry. Extensive experiments on various applications support our findings. 


On Riemannian Optimization over Positive Definite Matrices with the Bures- Wasserstein Geometry
Andi Han', Bamdev Mishra?, Pratik Jawanpuria', Junbin Gao!
'The University of Sydney, Australia "Microsoft, India
han, junbin.gas @sydney.edu.au (bamdevm, jawanurra)) Saterriaa@mmmmooooffi @microsoft.com.
Abstract
this paper, we comparatively analyze the Weres-WWsssrsttinn (BW) geometry with the popular Affine-Invarian (Al) geometry for Riemannian optimization on the symmetric positive definite (SPD) matrix manifold. Our study begins with observation that the BW metric has linear dependence contrast to the quadratic dependence the AI metric. We build on this to show that the BW metric is more suitable and robust choice for several Riemannian optimization problems over ill-conditioned SPD matrices. We show that the BW geometry has non-negative curvature, which further improves convergence rates of algorithms over the non non-positively curved Al geometry. Finally, we verify that several popular cost functions, which known be geodesic convex under the Al geometry, are also geodesic convex under the BW geometry. Extensive experiments on various applications support our findings.
SPD matrices in
Introduction
Learning symmetric positive definite (SPD) matrices is fundamental problem in various machine
learning applications, including metric and kernel learning [69
24
38
13
66), medical imaging [57
56). natural language processing [37 40), computer vision 28
35
34
multi-task learning (39] [4
domain adaptation [48, 41), modeling time-varying data
18 object detection (70), and quantum
mechanics [53 47 among others.
The set of SPD matrices of size n x n, defined as $ := [X X p2Xn X = X, and X > 0}.
has smooth manifold structure with richer geometry than the Euclidean space. When endowed with metric (inner product structure), the set of SPD matrices becomes a Riemannian manifold [T]. Hence, numerous existing works [57
36 34 46 56 have studied and employed the Riemannian
optimization framework for learning over the space of SPD matrices 2 [6).
Several Riemannian metrics on $ have been proposed such as the Affine-Invarian [57 the Log-Fuclidean 7 |X], the Log-Det 64 20] and the Log-Cholesky [46), name few. One can additionally obtain different families of Riemannian metrics on S by appropriate parameterizations
based on the principles invariance and symmetry [23 20
67 [6). However, to the best of
knowledge, systematic study comparing the different metrics for optimizing generic cost functions defined on $ is missing. Practically, the Affine-Invarian (AI) metric seems to be the most widely used metric in Riemannian first- and second-orde algorithms (e.g.. steepest descent, conjugate gradients, and trust regions) as it is the only Riemannian SPD metric available in manifold optimization toolboxes, such Manopt [17], Manopt.jl [10), Pymanopt [68]. ROPTLIB [32), and McTorch [50]. Moreover, many interesting problems in machine learning are found to be geodesic convex (generalization of Euclidean convexity) under the AI metrie, which allows fast convergence
of optimization algorithms 76 31
35th Conference on Neural Information Processing Systems (NcurIPS 2021).

---

Recent works have studied the Bures- Wasserstein (BW) distance on SPD matrices 49 12 72 It well-known result that the Wasserstein distance between two multivariate Gaussian densities is function of the BW distance between their covariance matrices. Indeed, the BW metric is Riemannian metric. Under this metric, the necessary tools for Riemannian optimization, including the Riemannian gradient and Hessian expressions, can be efficiently computed [49]. Hence, it is promising candidate for Riemannian optimization on S
In this work, we theoretically and empirically analyze the quality of optimization with the BW geome- try and show that it is a viable alternative to the default choice of AI geometry. Our analysis discusses the classes of cost functions (e.g,, polynomial) for which the BW metric has better convergence rates than the AI metric. We also discuss cases (e.g., log-det) where the reverse is true. In particular, our contributions are as follows.
We observe that the BW metric has linear dependence on SPD matrices while the Al metric has a quadratic dependence. We show this impacts the condition number of the Riemannian Hessian and makes the BW metric more suited to learning ill-conditioned SPD matrices than the Al metric.
In contrast the non-positively curved AI geometry, the BW geometry is shown to be non-negatively curved, which leads to tighter trigonometry distance bound and faster convergence rates for optimization algorithms.
For both metrics, we analyze the convergence rates of Riemannian steepest descent and trust region methods and highlight the issues arising from the differences in the curvature and condition number of the Riemannian Hessian.
We verify that common optimization problems that are geodesic convex under the Al metric are also geodesic convex under the BW metric.
We support our analysis with extensive experiments on applications such as weighted least squares, trace regression, metric learning, and Gaussian mixture model.
2 Preliminaries
The fundamental ingredients for Riemannian optimization are Riemannian metric, exponential map, Riemannian gradient, and Riemannian Hessian. We refer readers 2 for a general treatment on Riemannian optimization.
Riemannian metric is a smooth, bilinear, and symmetric positive definite function on the tangent space T.M for any x E M. That is, T.M x T.M + R, which is often written as an inner product (z.) The induced norm of a tangent vector u E T.M is given by |ul|= /(u, u),. A geodesic on a manifold 7 [0. 1] + .M is defined as a locally shortest curve with zero acceleration. For any M, u E M, the exponential map, Exp,, T. + M is defined such that there exists a geodesic curve with 7(0) = x, 7(1) = Exp, (u) and (0))) = u.
First-order geometry and Riemannian steepest descent. Riemannian gradient a differentiable function f M + IR at x, denoted as grad (z). is tangent vector that satisfies for any u E T,M, (gradf(x) (z), u)xx f.,((r), where fxx) is the directional derivative of fxx along u. The Riemannian steepest descent method [71 generalizes the standard gradient descent in the Euclidean
space to Riemannian manifolds by ensuring that the updates are along the geodesic and stay on the manifolds. That is, I++1=
Expz, (-1, g (x:)) for some step size T/t-
Second-order geometry and Riemannian trust region. Second-order methods such as trust re- gion and cubic regularized Newton methods are generalized to Riemannian manifolds 33 They make use of the Riemannian Hessian, Hess f(x) T + T,M, which a linear operator that is defined as the covariant derivative of the Riemannian gradient. Both the trust region and cubic regu- larized Newton methods are Hessian-free in the sense that only evaluation of the Hessian acting on a [u] is required. Similar the Euclidean counterpart, the Riemannian trust region method approximates the Newton step by solving a subproblem, i.e..
tangent vector, i.e., Hess (x)
1 (H,, [u], u) mzx (u) f(xe) + (zr) u)z. + 2
min "c MET ullz, ull,, SA

---

Table 1: Riemannian optimization ingredients for Al and BW geometries.
Invarintt
Bures- WarsssWassrrttiinn
R.Metric
(u, V) tr(X 1exx--U0 UV)
9bw (U, V) 2tt((x
E.Exp
Expai, (U) X/2 exp UUX
ExPow: x (U) = X U Lx [
R. G.Grdieett
grad.. XV)(XX
grad bw (X) = 4(V[(x)x);
R.Hessian
Hess, / [X] = XX
Hess / [X]
4(V² [XX]s +
(UV (UV](X)X}s
2
1 Lx grad bw f(X) is
where Ha, .T. M is symmetric and linear operator that approximates the Riemannian Hessian. A is the radius of trust region, which may be increased or decreased depending how model value mg, (u) changes. The subproblem is solved iteratively using a truncated conjugate gradient algorithm. The next iterate is given by 21+1 Exp,, (u) with the optimized t..
Next, the eigenvalues and the condition number of the Riemannian Hessian are defined as follows. which we use for analysis in Section 3
The minimum and maximum eigenvalues of Hess are defined as Amin = min (Hess f (x)] u) and Amax max//ull =1 (Hesss ([)ssf( u).. The condition number of H (x) is defined as ((Hess(rr) == Amax/ Xmin.
Definition 1.
Function classes on Riemannian manifolds. For analyzing algorithm convergence, we require the definitions for several important function classes on Riemannian manifolds, including geodesic convexity and smoothness. Similarly, we require the definition for geodesic convex sets that generalize (Euclidean) convex sets to manifolds [65 (4].
Definition (Geodesic convex set 665 74)..
A set X C M is geodesic convex for any x E X, the distance minimizing geodesic 7 joining the two points lies entirely in X.
Indeed, this notion is well-defined for any manifold because a sufficiently small geodesic ball is always geodesic convex.
Definition 3 (Geodesic convexity 65 74
Consider geodesic convex set X C M. A function f: X I is called geodesic convex for any x, y E X. the distance minimizing geodesic y joining I and satisfies Vt E [0, 1], f()())) < (1 tzf(2) + tf(y). Function f is strictly geodesic convex the equality holds only when 0,
Definition 4 (Geodesic strong convexity and smoothness [65] 33
Under the same settings in
Definition 3 A twice-continuously diffcrentiable function f X + is called geodesic p-strongly
convex if for any distance minimizing geodesic "y in X' with ||0( 1, satisfies d'(()(()) 2 p.
dba
for some μ 0. Function is called geodesic L-smooth if f(((()) L, for some L > 0.
de2
3 Comparing BW with AI for Riemannian optimization
This section starts with an observation of a veraar-verru--uudaattte dependency between the two metrics. From this observation, we analyze the condition number of the Riemannian Hessian. Then, we further compare the sectional curvature of the two geometries. Together with the differences in the condition number, this allows us compare the convergence rates of optimization algorithms on the two geometries. We conclude this section by showing geodesic convexity of several generic cost functions under the BW geometry. The proofs for this section are our extended report (7].
AI and BW geometries on SPD matrices. When endowed with a Riemannian metric the set of SPD matrices of size becomes Riemannian manifold M (8)++ g). The tangent space at X is M = {U U E Raxx and UT = U). Under the and BW metrics, the Riemannian exponential map. Riemannian gradient, and Hessian are compared in Table I where we denote {A}s == (A + A T)/2 and exp(A) ) as the matrix exponential of A. [x[U] the solution to the
matrix linear system Lx U + XC U and is known the Lyapunov operator. We use

---

Vf(X) and fX) to represent the first-order and second-order derivatives, i.e., the Euclidean gradient and Hessian, respectively. The derivations in Table can be found in T2. the rest of the paper, we use Mai and Mbw to denote the SPD manifolds under the two metrics. From Table I the computational costs for evaluating the Al and BW ingredients are dominated by the matrix exponential/inversion operations and the Lyapunov operator L computation, respectively. Both at most cost (n²)) which implies a comparable per-iteration cost of optimization algorithms between the two metric choices. This claim is validated in Section 4
key observation. From Table I the Affine-Invariant metric on the SPD manifold can be rewritten as for any U, V E MxM,
(U, V)ai trai-tXX--- -xx--X-X----- - -ve(U)) (x @ X) -ve((V)
(1)
where vec(U) and vec(V) are the vectorizations of U and V. respectively. Note that we omit the subscript X for inner product (...) to simplify the notation. The specific tangent space where the inner product is computed should be clear from contexts.
The Bures Wasserstein metric is rewritten as, for any U. Tx.M
1 (U, h/ [tr(Cx - veo (X) + ~ (V),
where X e & I + I & X is the Kronecker sum.
Remark 1. Comparing Eq. [ and (2) reveals that the BW metric has a linear dependence on X while the AI metric has quadratic dependence. This suggests that optimization algorithms under the Al metric should be more sensitive to the condition number of X compared to the BW metric.
The above observation serves as key motivation for further analysis.
3.1 Condition number of Riemannian Hessian at optimality
Throughout the rest of the paper, we make the following assumptions.
Assumption 1. (a). f is least twice continuously differentiable with a non-degenerate local minimizer X" (b). The subset X' C.M (usually as a neighbourhood of a center point) we consider throughout this paper is totally normal, i.e., the exponential map is diffeomorphimm.
Assumption T is easy to satisfy. Particularly, Assumption bb) is guaranteed for the SPD manifold under the Al metric because its geodesic is unique. Under the BW metric, for center point X, we can choose the neighbourhood such that X' = [Expx (U): I + [x[U] E 34. as in 49 In other words, X' assumed to be unique-geodesic under both the metrics.
We now formalize the impact of the linear-veruus-uaaraee dependency, highlighted in Remark At a local minimizer X* where the Riemannian gradient vanishes, we first simplify the expression for the Riemannian Hessian in Table
On Mais (Hess (X" [U], U ( RRn xx² the matrix representation of the Euclidean Hessian v² and U Tx-Mai: The max- imum eigenvalue of He (X*) is then given by Xinax ma vec(U) H(X*)vee(U)) where ||], vec(U) (X* @ X*) vee U). This is a generalized eigenvalue problem with the solution to be the maximum eigenvalue of (X* 2 (X*). Similarly, Xin. corresponds the minimum eigenvalue of (X* @ (X*))
vec(U) H( (ee(U), with H(X)
(X*)
On Mbws (Hessow f (X* [U], bw = vec(U))H(("")))))) "((e))vec(UU and the norm is | vec(U) (X* + X*) - (U). Hence, the minimum/maximum eigenvalue of Hesstw ((x*) equals the minimum/maximum eigenvalue of (X* @ XX*)
Let E := x: ((X))s <((" (X")) and Kip == x(Hess) f(X*)) x((X" e X X())( The following lemma bounds these two condition numbers.
Lemma 1. For a local minimizer X of fXX), the condition number of Hess (X*) satisfies
<(X* ((X()))(H(((*)) SKL. S <(X*))
x(X*)/<<H(XX Skfw I <((") <<(X'<(((XX")

---

It is clear that Now I Sa when R x(((X*)) - VR(X).. This is true for linear, quadratic, higher-ordet
polynomial functions and in general holds for several machine learning optimization problems on the SPD matrices (discussed in Section 4
Case 1 (Condition number for linear and quadratic optimization). For linear function f(X)
Ir(XA), A, its Euclidean Hessian matrix is H(X)
O. For quadratic function f(X) tr(XAXB)
with A, B E + HX) = & B + B A. Therefore, x(H((*)) constant and for @
ill-conditioned X", we have x(H(X))) S Vx(X"), which leads to HAi 2 Kbww
Case 2 (Condition number for higher-order polynomial optimization). For an integer a 2 3. consider function [(X = tr(X*) with derived H(X) - a E (Xl Xa-l-2). We get Mai a' Σi= ((*)) & (X*) and Kiw = axx + x') (X' X x)) is apparent that X = (x(( )9) while K (({((*)a-)) Hence, for ill-conditioned X K 2 kiw
One counter-example where Kiw 2 N is the log-det function.
Case 3 (Condition number for det optimization), For the log-det function f(X) = - log det(X). its Euclidean Hessian is [U = X --- -UX- and H(X) = X x At local minimizer X", H [U] U with K = 1. While on Mbw: we have kiw = x((X* @ X*)((*)-
(X*) = x((x===1(+0)) (X*)=)) = x(X*). Therefore, Kii S kibw'
3.2 Sectional curvature and trigonometry distance bound
To study the curvature of Mbb: we first show in Lemma 2 the existence of a matching geodesic between the Wasserstein geometry of zero-centered non-degenerate Gaussian measures and the BW geometry of SPD matrices. Denote the manifold of such Gaussian measures under the 22-Wasserstein distance as (2), W2) with Σ E S
Lemma 2. For any X, Y E $4 geodesic between (X) and No(Y) on (No(E), W2) is given by (o)((t)), where 7(t is the geodesic between X and Y on Mbw-
The following lemma builds on a result from the Wasserstein geometry 6 and uses Lemma 2 to analyze the sectional curvature of Mbw-
Lemma 3. Mow is an Alexandrov space with non non-negative sectional curvature.
It is well-known that Mai is a non-positively curved space 221 56 while, in Lemma we show that Mbw is non-negatively curved. The difference affects the curvature constant in the trigonometry distance bound of Alexandrov space (77.. This bound is crucial in analyzing convergence for optimization algorithms on Riemannian manifolds [77] 76). In Section [333 only local convergence to minimizer X is analyzed. Therefore, it suffices to consider neighbourhood f2 around X In such a compact set, the sectional curvature is known to be bounded and we denote the lower bound as
The following lemma compares the trigonometry distance bounds under the AI and BW geometries. This bound was originally introduced for Alexandrov space with lower bounded sectional curva- ture (77] The result for non-negatively curved spaces has been applied in many work [76] 60 25 though without a formal proof. We show the proof in (27), where it follows from the Toponogov comparison theorem [52) on the unit hypersphere and Assumption
[
Lemma 4. Let X. Y, Z E 22 which forms a geodesic triangle on M. Denote x d(Y, Z),
d(X, Z), d(X, Y) as the geodesic side lengths and let 0 be the angle between sides y and such
that cos ()) = (Expx) -1 (Y), Expx -1 (Z)) (yz). Then, have
a s <3² + 22 2yz cos(O),
VIKTID tambl / (/KIDD) diameter bound of i.e. maxx, d(X, Y) S D. Under the BW meiric, = Sow = 1.
where Ç is a curvature constant. Under the Al metric, Ç = Cai
with D as the
It is clear that Cai fow 1, which leads a tighter bound under the BW metric.
3.3 Convergence analysis
We now analyze the local convergence properties of the Riemannian steepest descent and trust region methods under the two Riemannian geometries. Convergence is established in terms of the
5

---

Riemannian distance induced from the geodesics. We begin by presenting a lemma that shows in neighbourhood of X*, the second-order derivatives of o Expx are both lower and upper bounded. Lemma 5. In a totally normal neighbourhood Qt around non-degenerate local minimizer X for ărz
any X E 22 it satisfies that X S d² [Fxxx
a.Xmaxe for some 2 and | 1. Xmxx Xmin > O are the largest and smallest eigenvalues of Hesss(X*)
For simplicity of the analysis, we assume such an is universal under both the Riemannian geometries. We, therefore, can work with a neighbourhood f with diameter uniformly bounded by D, where we can choose D == min Dai Dow) } such that is universal.
One can readily check that under Lemma 5 the function f is both y-geodesic strongly convex and L-geodesic smooth in S where = / a and L = a We now present the local convergence analysis of the two algorithms, which are based on results [.
Theorem (Local convergence of Riemannian steepest descent). Under Assumption 7 and consider a non degeeereee local minimizer X". For a neighbourhood 52 3 X* with diameter bounded by D on two Riemannian geometries Mai. Mow, running Riemannian steepest descent from Xo E St with a fixed step size 1j adiixx yields for
min C a² - ))²2
d (X, X S Sa
Theorem 2 (Local convergence of Riemannian trust region). Under the same settings as in Theorem
7 assume further in Ωt it holds that (1) | - ( S f(lgrad((x))) and (2) V[V
Expx, (U) ( a Expx, (0)] S pl U|] some É, universal on Mais Mow. Then running
Riemannian trust region from Xo E 2 yield.
d(),,")) & (2V0+ ()x (X.-.X").
Theorems and and Sshw show that Mbw has a clear advantage compared to Mai for learning ill-conditioned SPD matrices where "bw, I Ki For first-order algorithms, Mow has an additional benefit due to its non non-negative sectional curvature. As Ças Çbw the convergence rate degrades on Mai- Although the convergence is presented in Riemannian distance, it can be readily converted to function value gap by noticing (X,, X*) S f(X,) f(x*) S 1 +d X*). Additionally, we note that these local convergence results hold regardless of whether the function is geodesic convex or not, and similar comparisons also exist for other Riemannian optimization methods.
3.4 Geodesic convexity under BW metric for cost functions of interest
Finally we show geodesic convexity of common optimization problems on Miw- Particularly, we verify that linear, quadratic, log-det optimization, and also certain geometric optimization problems, that are geodesic convex under the Al metric, are also geodesic convex under the BW metric.
Proposition 1. For any A SiL where St {Z functions (X( = Lr XA), [x(X) (r(XAX) and (x(X) = - log det (X) are geodesic convex on Mbw-
Raxn Z
Z, and Z > 0), the
Based on the result in Proposition T we also prove geodesic convexity of a reparameterized version of the Gaussian density estimation and mixture model Mbw (discussed in Section Similar claims on Mai can be found in [31].
We further show that monotonic functions on sorted eigenvalues are geodesic convex on Mbw This an analogue of (65 Theorem 2.3] on Mai-
Proposition 2. Let X+ S4 + R", be the decreasingly sorted eigenvalue map and h R, + IR be an increasing and convex function. Then s(X) = Σ (X)) for 1 S n is geodesic conver on Mbw. Examples of such functions include fi (X) = tr(exp(X)) and fz(X) = tr(X*), 2
Experiments
In this section, we compare the empirical performance of optimization algorithms under different Riemannian geometries for various problems. In addition to Al and BW, we also include the Log-Euee Log-Euclidean (LE) geometry (7) in our experiments.

---

10
Db 204
1000
1003
10:10
300
10000
2000
4002
Dense, LowCN (RTR:iter)
(b) Dense, HighCN (RTR: iter)
(c) Sparse, HighCN (d) Sparse, HighCN (RTR:iter)
(c) Sparse, HighCN (RSD: time)
(RsD:
Figure I: Weighted least squares problem.
The LE geometry explores the the linear space symmetrie matrices where the matrix exponential acts global diffeomorphism from the space to $ The LE metric is defined as
(U, Vi - tr(Du log Dx log(X))
(3)
for any U, V E M, where Du log(X) is the directional derivative of matrix logarithm at X along U. Following (69 (5 for deriving various Riemannian optimization ingredients under the LE metric 3. we consider the parameterization X = cxp(S), where S E S", i.e., the space of n x n symmetric matrices. Equivalently, optimization the SPD manifold with the LE metric is identified with optimization on S" and the function of interest becomes f(exp(S) for S E S While the Riemannian gradient be computed efficiently by exploiting the directional derivative of the matrix exponential [4), deriving the Riemannian Hessian is tricky and we rely on finiteddifference Hessian approximations 15]
We present convergence mainly in terms of the distance to the solution X whenever applicable. The distance is measured in the Frobenius norm, i.e., X, -X* When X is not known, convergence is shown in the modified Euclidean gradient norm ((x)(ll This is comparable across different metrics the optimality condition X xXX" initialize the algorithms with the identity matrix for the Al and BW metrics and zero matrix for the
0 arises from problem structure itself (43. We
LE metric (i.e., the matrix logarithm of the identity).
We mainly present the results on the Riemannian trust region (RTR) method, which is the method of choice for Riemannian optimization. Note for RTR, the results are shown against the cumulative sum of inner iterations (which are required to solve the trust region subproblem at every iteration). We also include the Riemannian steepest descent (RSD) and Riemannian stochastic gradient (RSGD) [14] methods for some examples, The experiments are conducted in Matlab using the Manopt toolbox [1 on 10500500 3.1 GHz CPU processor.
In our extended report 277, we include additional experiments comparing convergence objective function values for the three geometries. We also present results for the Riemannian conjugate gradient method, and results with different initializations (other than the identity and zero matrices) to further support our claims.
The code can be found httts::/ com/ andaaadddjjnnamae http://tthu...ocmaaaaya/////////
Weighted least squares. We first consider the weighted least squares problem with the symmetric positive definite constraint. The optimization problem is minxes;, f(X) = | - 13. which is encountered in for example, SPD matrix completion (63] where A is a sparse matrix. The Euclidean gradient and Hessian are (x(X) (A O X B) A and v² Hence, at optimal X the Euclidean Hessian in matrix representation is H
A o U c A, respectively.
diag(vec( A)).
We experiment with two choices of A, i.e. A = 1 1, (Dense) and A as a random sparse matrix (Sparse). The former choice for A leads to well-conditioned H(X*) while the latter choice leads to an ill-conditioned H(X*). Also note that when A = 1,, 1, K x(x*)² and Ki =(X*).
We generate X as SPD matrix with size n 50 and exponentially decaying eigenvalues. We consider two cases with condition numbers x(X*) = 10 (LowCN) and 10 (HighCH). The matrix B generated as = A & X*. F compares both RSD and RTR for different metrics. When A is either dense or sparse, convergence is significantly faster on Mhw than both Mai and Mia- The advantage of using Mhw becomes more prominent the setting when condition number of X* high. Figure [[) shows that Mbw is also superior in terms of runtime.
7

---

10
-01
Oitince 50308
Dieence 20
2030
101
5,920
10,000
25,500
42,002
65,800
5,500
(a) LYA: Ex1Full
(b) LYA: Ex1Low
(c) LYA: Ex2Full
(d) LYA: Ex2Low
600
Eutldes gnd gndrem rm
Eut ppd rm
600
10
200
203
400
20
(e) TR: SynPull
(f) TR: SynLow
(g) DML: glass
(h) DML: phoneme
Figure 2: Lyapunov equation (a, b, c, d), trace regression (e, f), and metric learning (g, h) problems.
Lyapunov equations. Continuous Lyapunov matrix equation, AX XA = C with X commonly employed in analyzing optimal control systems and differential equations 59
E ST
45)
When A is stable, i... (A) > O and C E 84 the solution X > 0 and is unique 44]. When Ce $4 and is low rank, X S is also low rank. We optimize the following problem for solving
the Lyapunov equation (73.. i.e., minxes" f(X)= = tr(XAx) tr(XC). The Euclidean gradient and Hessian are respectively V ((()
AX XA C and [x] AU UA with H(X)= = A e A. At optimal X", the condition number x xHH(")) = k(A).
We experiment with two settings for the matrix A, i... A as the Laplace operator on the unit square where we generate 7 interior points so that n = 49 (Ex1), and is a particular Toeplitz matrix with n 50 (Ex2). The generated A matrices are ill- ill-conditioned. The above settings correspond to Examples 7.1 and 7.3 in [45]. Under each setting, X* is set to be either full or low rank. The matrix C' is generated as C AX X'A. The full rank X is generated from the full-rank Wishart distribution while the low rank X" is diagonal matrix with r = 10 ones and n r zeros in the diagonal. We label the four cases Ex1Full, ExiLow, Ex2Full, and Ex2Low. The results shown in Figures 2 (d), where we observe that in all four cases, the BW geometry outperforms both AI and LE geometries.
Trace regression. We consider the regularization-free trace regression model 62 for estimat- ing covariance and kernel matrices [6T] 19]. The optimization problem is milxees, f(X)
1, ..... matrices. Thus, we have ff(x) E (a. Xa; 1o4-yya y,)A, and V [X) (a. Ua A.,
1 Ee (y: tr x))? where A; a,a;
rank-one measurement
2m
We create X as a rank-r Wishart matrix and {A,} as rank-one Wishart matrices and generate y; = tr( tr( aC; with 6 ~ NO, 1), 0.1. We consider two choices, (m, r) (1000, 50, 50) and (1000, 50, 10), which are respectively labelled as SynFull and SynLow. From Figures &()).
also observe that convergence to the optimal solution is faster for the BW geometry.
Metric learning. Distance metric learning (DML) aims to learn a distance function from sam- ples and popular family of such distances is the Mahalanobis distance, i.e. dM(x,y) V(x- y) M(x y) for any x,y E Rd.. The distance is parameterized by symmetric posi- tive semi-definite matrix M. We refer readers to this survey [66] for more discussions on this topic. We particularly consider a logistic discriminant learning formulation (24]. Given a train- ing sample {x;, Vi ):) denote the link Lij = if ys Yj and by = 0 otherwise. The ob- jective is given by min f(M) = - Ei, (big log + (1 - b4g) log(1 - P())), with Po
+ exp(dm x,))) -1, We derive the matrix Hessian H(M) Er, Pi Pty)
(x, ) @ (x, X,) (x; - -x,) Notice & (H(M*)) depends on M only through the constants Piy- Thus, the condition number will not be much affected by s(M*).
We consider two real datasets, glass and phoneme, from the Keel database [5]. The number of
classes is denoted c. The statistics of these two datasets are (N,d, d, c) (241, 9, 7) for glass
8

---

101
Disnne 3 ad4i0o
muDi ELL
Dee
3.5
3.48
1/48
20
2022
1002
2000
3000
(a) LD: LowCN (RTR) (b) LD: HighCN (RTR)
(c) GMM: (RTR)
(d) GMK: (RSGD)
(c) GMM: (RSGD)
Figure 3: Log-det maximization (a, b) and Gaussian mixture model (c, d. e) problems.
(5404, 5, 2) for phoneme. In Figures &h), we similarly see the advantage of using the BW metric compared the other two metries that behave similarly.
Log-det maximization. As discussed in Section 3.1 log-det optimization is one instance where Nbw 2 Kaj We first consider minimizing negative doeeeetermnnant along with a linear function stud- ied in 75 That is, for some C E @ the objective is minxes", f(X) tr(XC)
log det(X). C X and V f² X --- This problem is geodesic convex under both AI and BW metrics. We generate X the same way as in the example of weighted least square with n 50 and C (X*) We consider two cases with condition number en 10 (LowCN) and 103 (HighCN). As expeted, we observe faster convergence of Al and LE metrics over the BW metric in Figures (a) This even more evident when the condition number increases.
The Euclidean gradient and Hessian are given by f(X)
Gaussian mixture model. Another notable example of log-det optimization is the Gaussian density estimation and mixture model problem. Following [3T] we consider a reformulated problem on aug- mented samples y: = : = I...., N where x; E Rd are the original samples. The density is parameterized by the augmented covariance matrix E Rd+1, Notice that the log log-likelihood of Gaus- sian is geodesic convex on Mais but not on Mbw. We, therefore, define S = Σ and the reparame- terized log- likelihood is Pv PN(Y; :Y Σ Ei=1 log ((2m))-4/2 exp( 1/2) det(S) 1/2 exp(- by; Sy;)). which is now geodesic convex on Maw due to Proposition 1 Hence, we can solve the problem of Gaussian mixture model similar in (31]
Here, we test on a dataset included in the MixEst package [30] The dataset has 1580 samples R² with 3 Gaussian components. In Figure 3c). we observe a similar pattern with RTR as in the log-det example. We also include performance of RSGD, which is often preferred for large scale problems. We set the batch size to be 50 and consider decaying step size, with the best initialized step size shown in Figures 3 (c)) Following (8), the algorithms are initialized with kmeans++, We find that the AI geometry still maintains its advantage under the stochastic setting.
5 Conclusion and discussion
In this paper, we show that the less explored Bures- Wasserstein (BW) geometry for SPD matrices should often be the preferable choice than the Affine-Invarian geometry for optimization, par- ticularly for learning ill-conditioned matrices. A systematic analysis shows that the BW metric preserves geodesic convexity of some popular cost functions and leads better convergence rates in optimization.
We also theoretically discuss 'counter-example' of log-det optimization where the (AI) geometry enjoys better second-order conditioning and validate our findings empirically. This issue addressed in our recent work (26), where propose generalized Bures- W (GBW) geometry that is built on generalization of the Lyapunov operator in the metric:
(U, g) vee (U) (X M M M - vvc((V)
(4)
2
where U and are symmetric matrices and M is a given SPD matrix. When M I, the metric reduces to the special BW metric (2). The use of the parameter M in (4) allows great flexibility for optimization algorithms. For one, choosing a particular M allows preconditioning the Hessian by locally approximating the AI geometry. This leads to improved convergence for log-det optimization.

---

Our comparisons between Al and BW geometries are based on optimization over generic cost functions. For specific problems, however, there may exist other alternative metrics that potentially work better. This is an interesting research direction to pursue. We also remark that optimization is not the only area where the Al and BW geometries can compared. would be useful qualitatively compare the two metrics for other learning problems on SPD matrices, such as barycenter learning. In addition, kernel methods have been studied on the SPD manifold [42 29 78 that embed SPD matrices to high dimensional feature space, known as the Reproducing Kernel Hilbert Space (RKHS). Such representations are used for subsequent learning tasks such as clustering or classification. But only positive definite kernel provides a valid RKHS. We show in appendix Section A that the induced Gaussian kernel based on the BW distance is a positive definite kernel unlike the case for the AI metric. This difference highlights potential advantage of the BW metric for representation learning on SPD matrices.
A Positive definite kernel from BW geometry
In this section, show the existence of a positive definite Gaussian kernel on Miw- Existing works (55] [2 have studied the Wasserstein distance kernel. First, we present the definition of positive (resp. negative) definite function in [9]
A function f xx - IR is called positive definite and only if k symmetric and for all integers m 2 2, {011, Tm S X and {Ci...., Cm C R, satisfies Σx,==1 c.cy (I;, 2z) 2 0. A function f is called negative definite if and only if under the same conditions, satisfies Σnj=1 ( x3)) S 0 with E 0.
Definition 5, Consider X' to be nonempty
The following theorem shows the Gaussian kernel induced from the BW distance is positive definite on the SPD manifold.
Theorem 3. The induced Gaussian kernel k(.s-) == exp( -app--ppom d.w (1 (2(2)) is positive definite.
Proof of Theorem 3 From Theorem 4.3 in 42 definite. Indeed for any m 2 2, [X1,.., } C MMww {ei,..., Cm. S R with E C, 0, we have
it suffices to prove the BW distance is negative th
m
GC Gejaf
Am (X, X,)
m
m
m
172
m C (xx((1)) 2 Σ
ctr(X,) +
C2
cxttr(X) 1/2 xx,x //2//31112 i==1
j=1
i=1
1=1
1=1
Σ
1/2
x,X!
12xxx1121111 1/231/2 0.
=
This shows &yyn is negative definite, and therefore, the exponentiated Gaussian kernel is positive definite.
Broader Impact
This paper is of a theoretical nature and does not foresee any immediate societal impacts.
References
[1] P-A Absil, Christopher Baker, and Kyle A Gallivan. Trust-region methods on Riemannian manifolds. Foundations of Computational Mathematics, 7(3):303- -3380, 2007.
[2] P-A Absil, Robert Mahony, and Rodolphe Sepulchre. Optimization algorithms on matrix manifolds. Princeton University Press, 2009.
[3] Naman Agarwal, Nicolas Boumal, Brian Bullins, and Coralia Cartis. Adaptive regularization with cubics manifolds. Mathematical Programming, pages 1-50, 2020.
10

---

[4] Awad H Al-Mohy and Nicholas J Higham. Computing the Fréchet derivative of the matrix exponential, with an application condition number estimation. SIAM Journal on Matrix Analysis and Applications, 30(4): 1639-16577 2009.
[5] Jesús Alcalá-Fdez Luciano Sánchez, Salvador Garcia, Maria Jose del Jesus, Sebastian Ventura, Josep Maria Garrell, José Otero, Cristóbal Romero, Jaume Bacardit, Victor M Rivas, et al. Keel: software tool to assess evolutionary algorithms for data mining problems. Soft Computing, 13(3):307-318, 2009.
[6] Luigi Ambrosio, Nicola Gigli, and Giuseppe Savaré. Gradient flows: metric spaces and in the space of probability measures. Springer Science & Business Media, 2008.
[7] Vincent Arsigny, Pierre Fillard, Xavier Pennec, and Nicholas Ayache. Geometric means in novel vector space structure on symmetric positive-defintte matrices. SIAM journal on matrix analysis and applications, 29(1):22-3477, 2007.
[8] David Arthur and Sergei Vassilvitskii. k-means++: The advantages of careful seeding. Technical report, Stanford, 2006.
[9] Christian Berg, Jens Peter Reus Christensen, and Paul Ressel. Harmonic analysis on semigroups: theory of positive definite and related functions, volume 100. Springer, 1984.
[10] Ronny Bergmann. Optimisation on Manifolds in Julia, 2019. https: /ttpsb//giihub..... thub com/ kellertuer/Manopt..]
[11] Rajendra Bhatia. Positive definite matrices. Princeton university press, 2009.
[12] Rajendra Bhatia, Tanvi Jain, and Yongdo Lim. On the Bures- Wasserstein distance between positive definite matrices. Expositiones Mathematicae, 37(2):165-199 19,
[13] Mukul Bhutani, Pratik Jawanpuria, Hiroyuki Kasai, and Bamdev Mishra. Low-rank geometric mean metric learning. ICML workshop on Geometry in Machine Learning (GiMLi), 2018
[14] Silvère Bonnabel. Stochastic gradient descent on Riemannian manifolds. IEEE Transactions on Automatic Control, 58(9):2117-2299 2013.
[15] Nicolas Boumal. Riemannian trust regions with finite-difference Hessian approximations are globally convergent. In International Conference on Geometric Science of Information, pages 467-475. Springer, 2015.
[16] Nicolas Boumal. An introduction to optimization on smooth manifolds. Available online, May,
2020.
[17] Nicolas Boumal, Bamdev Mishra, P-A Absil, and Rodolphe Sepulchre. Manopt, a Matlab toolbox for optimization manifolds. The Journal of Machine Learning Research, 15(1):1455- 1459, 2014.
[18] Daniel A Brooks, Olivier Schwander, Frédéric Barbaresco, Jean-Y Yves Schneider, and Matthieu Cord, Exploring complex time-series representations for Riemannian machine learning of radar data. IEEE International Conference on Acoustics, Speech and Signal Processing, pages 3672-3676. IEEF, 2019.
[19] T Tony Cai and Anru Zhang. ROP: Matrix recovery via rank-one projections. Annals of Statistics, 43(1): (3(1):102-133) 2015.
[20] Zeineb Chebbi and Maher Moakher. Means of Hermitian positive-definite matrices based on the log-determinant a-divergence function. Linear Algebra and its Applications, 436(7):1872-188) 2012.
[21] Calin Cruceru, Gary Becigneul, and Octavian-Fine -ugen Ganea. Computationally tractable Rie- mannian manifolds for graph embeddings. In AAAI Conference Artificial Intelligence, pages 7133-7141, 2021.
[22] Henri De Plaen, Michaël Fanuel, and Johan AK Suykens. Wasserstein exponential kernels, In International Joint Conference Neural Networks, pages 1-6. IEEE, 2020.
11

---

[23] Ian L. Dryden, Alexey Koloydenko, and Diwei Zhou. Non-Euclidean statistics for covariance matrices, with applications to diffusion tensor imaging. The Annals of Applied Statisties, 3(3):1102-1123, 2009.
[24] Matthieu Guillaumin, Jakob Verbeek, and Cordelia Schmid. Is that you? Metric learning approaches for face identification. In International Conference on Computer Vision, pages 498-505. IEEE, 2009.
[25] Andi Han and Junbin Gao. Improved variance reduction methods for Riemannian non-convex optimization. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2021.
[26] Andi Han, Bamdev Mishra, Pratik Jawanpuria, and Junbin Gao. Generalized Bures- Wasserstein geometry for positive definite matrices. arXiv:2 10. 10464, 2021.
[27] Andi Han, Bamdev Mishra, Pratik Jawanpuria, and Junbin Gao. On Riemannian optimization over positive definite matrices with the Bures- Wasserstein geometry. urXiv: (tt::2006.0222 (02600 2021.
[28] Mehrtash T Harandi, Mathieu Salzmann, and Richard Hartley. From manifold to manifold: Geometry-aware dimensionality reduction for SPD matrices. In European conference on computer vision, pages 17-32. Springer, 2014.
[29] Mehrtash T Harandi, Conrad Sanderson, Richard Hartley, and Brian Lovell Sparse coding and dictionary learning for symmetric positive definite matrices: A kernel approach. In European Conference on Computer Vision, pages 216-229. Springer, 2012.
[30] Reshad Hosseini and Mohamadreza Mash'al al. MixEst: An estimation toolbox for mixture models. orXiv:1507.06065. 2015.
[31] Reshad Hosseini and Suvrit Sra. An alternative to EM for Gaussian mixture models: batch and stochastic Riemannian optimization. Mathematical Programming, 181(1):187-223, 2020.
[32] Wen Huang. Pierre-Antoine Absil, Kyle A. Gallivan, and Paul Hand. Roptlib: an object-oriented C++ library for optimization on Riemannian manifolds. Technical Report ESU16-14.2, Florida State University, 2016.
[33] Wen Huang, Kyle A Gallivan, and P-A Absil. A Broyden class of quasi-Newton methods for Riemannian optimization. SIAM Journal Optimization, 25(3):1660-166) 1685, 2015.
[34] Zhiwu Huang and Lue Van Gool. A Riemannian network for SPD matrix learning. In AAAI Conference on Artificial Intelligence, 2017.
[35] Zhiwu Huang, Ruiping Wang, Xianqiu Li, Wenxian Liu, Shiguang Shan, Lue Van Gool, and Xilin Chen. Geometry-aware similarity learning SPD manifolds for visual recognition. IEEE Trunsuctions on Circuits and Systems for Video Technology, 28(10):2533-252, 2017.
[36] Zhiwu Huang, Ruiping Wang, Shiguang Shan, Xianqiu Li, and Xilin Chen. -uellden met- ric learning on symmetric positive definite manifold with application to image set classification. In International Conference on Machine Learning, pages 720-729. PMLR, 2015.
[37] Pratik Jawanpuria, Arjun Balgovind, Anoop Kunchukutran, and Bamdev Mishra. Learning multilingual word embeddings in latent metric space: A geometric approach. Transactions of the Association for Computational Linguistics, 7:107-120, 2019.
[38] Pratik Jawanpuria, Maksim Lapin, Matthias Hein, and Bernt Schiele. Efficient output kernel learning for multiple tasks. In Conference on Neural Information Processing Systems, 2015.
[39] Pratik Jawanpuria and Bamdev Mishra. A unified framework for structured low-rank matrix learning. In ICML, 2018
[40] Pratik Jawanpuria, Satya Dev V, Anoop Kunchukuttan, and Bamdev Mishra. Learning geometric word meta-embeddings. In ACL Workshop on Representation Learning for NLP, 2020.
[41] Pratik Jawanpuria, Satya Dev N V, and Bamdev Mishra. Efficient robust optimal transport with application to multi-label classification. In IEEE CDC, 2021.
12

---

[42] Sadeep Jayasumana, Richard Hartley, Mathieu Salzmann, Hongdong Li, and Mehrtash Harandi. Kernel methods on the Riemannian manifold of symmetric positive definite matrices. In Conference Computer Vision and Pattern Recognition, pages 73-80, 2013.
[43] Michel Journée, Francis Bach, P-A Absil, and Rodolphe Sepulchre. Low-rank optimization the cone of positive semidefinite matrices. SIAM Journal on Optimization, 20(5):2327-2355 2010.
[44] Peter Lancaster. Explicit solutions of linear matrix equations. SIAM Review, 12(4):544-566, 1970.
[45] Yiding Lin and Valeria Simoncini. A new subspace iteration method for the algebraic riccati equation. Numerical Linear Algebra with Applications, 22(1):26 2015.
[46] Zhenhua Lin. Riemannian geometry of symmetric positive definite matrices via cholesky decomposition. SIAM Journal on Matrix Analysis and Applications, 40(4):1353-133 2019.
[47] Ilia A. Luchnikov, Alexander Ryzhov, Sergey N. Filippov, and Henni Ouerdane. QGOpt: Riemannian optimization for quantum technologies. SciPost Phys., 10:79, 2021.
[48] Sridhar Mahadevan, Bamdev Mishra, and Shalini Ghosh. A unified framework for domain adaptation using metric learning on manifolds. In ECML-PKDD, 2019.
[49] Luigi Malagò, Luigi Montrucchio, and Giovanni Pistone. Wasserstein Riemannian geometry of Gaussian densities. Information Geometry, 1()):137-179, 2018.
[50] Mayank Meghwanshi, Pratik Jawanpuria, Anoop Kunchukuttan, Hiroyuki Kasai, and Bamdev Mishra. Mctorch, a manifold optimization library for deep learning. arXiv:1810.01811, 2018.
[51] Gilles Meyer, Silvère Bonnabel, and Rodolphe Sepulchre. Regression on fixed-rank positive semidefinite matrices: a Riemannian approach. The Journal of Machine Learning Research, 12:593-625, 625, 2011.
[52] Wolfgang Meyer. Toponogov's theorem and applications. Lecture Notes, Trieste, 1989.
[53] Bamdev Mishra, Hiroyuki Kasai, and Pratik Jawanpuria. Riemannian optimization the simplex of positive definite matrices. NeurlPS workshop on Optimization for Machine Learning, 2020.
[54] Madhav Nimishakavi, Pratik Jawanpuria, and Bamdev Mishra. A dual framework for low-rank tensor completion. In NeurIPS, 2018.
[55] Jung Hun Oh, Maryam Pouryahya, Aditi Iyer, Aditya P Apte, Allen Tannenbaum, and Joseph O Deasy. Kernel Wasserstein distance. arXiv:1905.09314 2019.
[56] Xavier Pennec. Manifold-valued image processing with SPD matrices. In Riemannian Geomet- ric Statistics Medical Image Analysis, pages 75-134. Elsevier, 2020.
[57] Xavier Pennec, Pierre Fillard, and Nicholas Ayache. A Riemannian framework for tensor computing. International Journal of computer vision, 66(1):41-66 2006.
[58] Minh Ha Quang, Marco San Biagio, and Vittorio Murino. Log-hilbert-schmiit metric between positive definite operators on hilbert spaces. In Conference on Neural Information Processing Systems, pages 388-396, 2014.
[59] Rothschild and A Jameson. Comparison of four numerical algorithms for solving the Liapunov matrix equation. International Journal of Control, 11(2):181-198, 1970.
[60] Hiroyuki Sato, Hiroyuki Kasai, and Bamdev Mishra. Riemannian stochastic variance re- duced gradient algorithm with retraction and vector transport. SIAM Journal on Optimization, 29(2):1444-1472 2019.
[61] Bernhard Schölkopf and Alexander Smola, Learning with kernels: support vector machines, regularization, optimization, and beyond. MIT press, 2002.
13

---

[62] Martin Slawski, Ping Li, and Matthias Hein. Regularization-Fie free estimation in trace regression with symmetric positive semidefinite matrices. In Conference on Neural Information Processing Systems, 2015.
[63] Ronald L. Smith. The positive definite completion problem revisited. Linear Algebra and its Applications, 429(7): (497)::44-1155, 2008
[64] Suvrit Sra. A new metric on the manifold of kernel matrices with application to matrix geometric means. In Conference Neural Information Processing Systems, volume 25, pages 144 152, 2012.
[65] Suvrit Sra and Reshad Hosseini. Conic geometric optimization on the manifold of positive definite matrices. SIAM Journal on Optimization, 25(1):713-739, 2015.
[66] Juan Luis Suárez, Salvador Gareía, and Francisco Herrera. A tutorial on distance metric learn- ing: Mathematical foundations, algorithms, experimental analysis, prospects and challenges. Neurocomputing. 425:300-322, 2021.
[67] Yann Thanwerdas and Xavier Pennec. Is affineiinvariance well defined on SPD matrices? principled continuum of metrics. In International Conference Geometric Science of Information, pages 502-510. Springer, 2019.
[68] James Townsend, Niklas Koep, and Sebastian Weichwald. Pymanopt: A Python toolbox for optimization manifolds using automatic differentiation. The Journal of Machine Learning Research, 17()):75554452, 2016.
[69] Koji Tsuda, Gunnar Rätsch, and Manfred K. Warmuth. Matrix exponentiated gradient updates for on-line learning and Bregman projection. Journal of Machine Learning Research, 6(34):995- 1018, 2005.
[70] Oncel Tuzel, Fatih Porikli, and Peter Meer. Pedestrian detection via classification on Riemannian manifolds. IEEE Transactions on Puttern Analysis and Machine Intelligence, 30(10): 2008.
[71] Constantin Udriste. Convex functions and optimization methods Riemannian manifolds, volume 297. Springer Science Business Media, 2013.
[72] Jesse van Oostrum. Bures- Wassersesnn geometry. arXiv:2001.08055 2020.
[73] Bart Vandereycken and Stefan Vandewalle. A Riemannian optimization approach for computing low-rank solutions of Lyapunov equations. SIAM Journal on Matrix Analysis and Applications, (15))2553-2559 2010.
[74] Nisheeth K Vishnoi. Geodesic convex optimization: Differentiation on manifolds, geodesics, and convexity. 1rXiv:1106.06373. 2018
[75] Chengjing Wang, Defeng Sun, and Kim-Chuar Toh. Solving log-determinant optimization problems by Newton-CG primal proximal point algorithm. SIAM Journal Optimization, 2994)3399,3 2010.
[76] Hongyi Zhang, Sashank J Reddi, and Suvrit Sra. Riemannian SVRG: Fast stochastic opti- mization on Riemannian manifolds. In Conference on Neural Information Processing Systems, volume 29, pages 4502-4600, 2016.
[77] Hongyi Zhang and Suvrit Sra. First-order methods for geodesically convex optimization. In Conference on Learning Theory, pages 1617-1638. PMLR, 2016.
[78] Jianjia Zhang, Lei Wang, Luping Zhou, and Wanqing Li. Learning discriminative Stein kernel for SPD matrices and its applications. IEEE Transactions on Neural Networks and Learning Systems, 27(5):1020-1033 2015.
14


Reviews:
[{"review_id": "yFj39RNV8E_", "reviewer": "Reviewer_crTV", "summary": "This paper investigates the advantages of the Bures-Wasserstein (BW) metric over the widely used Affine Invariant (AI) one, for the manifold of symmetric positive definite (SPD) matrices. Rewriting the two metrics in a compatible manner reveals that for a wide variety of objectives, the BW metric is a linear function of the input matrix whereas the AI has a quadratic dependence. One exception is the log-det function which is included in the analysis. The paper experimentally validates that for ill-conditioned problems BW is a better choice compared to AI.", "questions": "", "limitations": "", "rating": 8, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "PROS\n\n- SPD matrices are used ubiquitously in computer vision. Hence, the fundamental contribution of this paper is expected to have a significant impact in the field -- at least for the portion which consider optimization on SPDs. \n- I like the principled approach of the paper. It starts with an analysis of the problem and later fills in some theoretical gaps in the literature. \n- The experiments include various relevant problems and reflect what the theory stands for: BW does indeed converge better.\n- The content, albeit involved, is clearly presented with adequate references. I would be happy to see this paper at NeurIPS.\n\nCONS\n\nWhile I have not gone through all the derivations, I have not spotted major weaknesses of the paper. Here are some minor points of concern:\n- The fact that BW provides a valid RKHS is very important. In fact the non-positive definiteness of the kernel has already caused nuisances for computer vision people. See for example:\n* Birdal, T., Arbel, M., Simsekli, U., & Guibas, L. J. (2020). Synchronizing probability measures on rotations via optimal transport. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition (pp. 1569-1579). *\nTherefore, I would be happy to see further discussions on the connections with kernel methods here.\n- Can we extend this study and compare with other metrics? Are there other metrics which also show linear dependence but somewhat inferior to BW? \n- Log-det seems to be just one counter-example. I feel like there are many other functions of this type. Could we have a general characterization of when affinte-invariant metric is a good choice and when not? Maybe a general recipe towards choosing a metric?\n- While the analysis is conducted for an analytical Hessian, in practice the Hessian is approximated using a finite differencing scheme. Does that create a significant gap with the theory? Can we actually study/quantify this gap? Why would the paper not use auto-differentiation to be exact?\n- The paper assumes the availability of an exponential map. For SPDs this is okay. But for conducting a similar analysis on other manifolds, a (maybe first-order) retraction should be considered. Is it possible to extend the analysis to the case of retractions? \n- I suggest that the findings here are incorporated into Manopt. ", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "This paper has a theoretical nature and does not have direct negative impacts on society. The limitations are explained in the discussions section.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4"}, "extra_text": {}}, {"review_id": "jokr7x_nLo", "reviewer": "Reviewer_w1Wx", "summary": "This paper provides several results about the Bures-Wasserstein  (BW) geometry on the positive definite matrices. A thorough comparison with the Affine Invariant (AI) metric is made. The main focus of the paper is on optimization: is it faster to optimize a function under the AI or under the BW metric? The authors argue that when the solution of the problem $X^*$ is ill-conditioned, then algorithms derived under BW metric have faster convergence theoretical convergence. This is due to the conditioning of the Riemannian Hessian at the optimum, which is better for the BW metric when $X^*$ has poor conditioning. The authors also demonstrate that a variety of functions that are convex for the AI geometry are also convex for the BW geometry.\nExperiments on an array of synthetic problems validate that both Riemannian gradient descent and Riemannian trust region methods are faster in the BW geometry.", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "**Originality** :  This paper derives some novel results about the Bures-Wasserstein geometry, and provides some insights comparing it to the affine-invariant metric. Such comparisons are important because both metrics are widely used in practice. To the best of my knowledge, insights regarding the conditioning (lemma 1), distance bounds (lemma 4) or convexity (prop 1) are novel.\n- Lemma 2 is already known, see Example 1.7 in [1]\n- I think that lemma 4 and Lemma 5 in [2] are similar.\n\n**Quality** : I think that comparing in a thorough manner the convergence properties of optimization algorithms under different metrics is an important and interesting problem. However, I have some concerns regarding the results of the authors. In particular, I think that the authors paint a picture that is too biased in favor the BW metric, and that there are some flaws in the analysis. \n\n- One of the most important result of the paper is that the conditioning of the Hessian under the BW metric, $\\kappa_{BW}$, is smaller than the conditioning under the AI metric, $\\kappa_{AI}$, when $\\kappa(H) \\leq \\sqrt{\\kappa(X^*)}$, where $H$ is the Euclidean Hessian of the function, and $X^*$ is the solution. This means that algorithms using the BW geometry will be faster when $X^*$ is poorly conditioned, not when $X^*$ is well conditioned. Some parts of the text seem to skip this important condition (e.g. L163). Also, when $\\kappa(H) \\leq \\sqrt{\\kappa(X^*)}$, then lemma 1 gives $\\kappa_{BW} \\geq \\kappa(X^*)/ \\kappa(H) \\geq \\kappa(H)$. Therefore, BW has itself a worse conditioning that using simply the Euclidean metric ! This fact should be mentioned, and this also raises the question of using plain Euclidean algorithms (for instance, projected gradient descent) in the experiments.\n\n- The constant $\\alpha$ in Lemma 5 plays a key role in the convegence speed of the algorithms. However, I think it is important to insist that $\\alpha$ depends on the metric used, since it is defined with the Exp. As a consequence, it is impossible to conclude that the rates in Thm.1 and 2 are better or worse for the BW/AI metric, even when $\\kappa_{BW} < \\kappa_{AI}$. In order to draw a conclusion, a study of the behavior of $\\alpha$ under the two metrics should be carried.\n\n- I also have concerns regarding the experiments. I replicated the least squares experiment with success. However, I tried another  related cost function, where here there is a matrix product instead of Hadamard product $f(X) = \\frac12 \\|AX - B\\|$, where $B = AX^*$ with $X^*$ ill-conditioned, and I found that the two algorithms performed quite comparably, even when $X^*$ is poorly conditioned. I suspect that this is due to the constant $\\alpha$ in the bounds. Also, I think that an experiment based on symmetric matrices taken from a real dataset would strengthen the point of the paper.\n\n\n\n\n**Clarity**: the paper is well written and pleasant to read. It is also well organized.\n\n**Significance** : as I already stated, I think that a clear and accurate comparison of the two metrics is a significant result that is of broad interest to any researcher working with PSD matrices.\n\n**Misc** :\n- Another field where PSD matrices with the AI/BW matrices are extensively used is EEG/MEG processing, with works such as [3]\n- It should be mentioned that the paper assumes that the Euclidean gradient of $f$ is always symmetric.\n- The actual BW and AI distances between two PSD matrices is never explicitly written in the paper, while these are the most used quantities. Including them in Table 1 would be worthwhile.\n- Following the definition of conditioning of the authors, the conditioning of a linear problem is ill-defined.\n- In lemma 1, are the inequalities tight?\n- In Lemma 5, what is $t$?\n- In fig.1 with gradient descent, how was chosen the step size? It is of crucial importance.\n- In fig1. c, d, e, using a log scale for the y axis would make things more readable. \n\n[1]: McCann, Robert J. \"A convexity principle for interacting gases.\" Advances in mathematics 128.1 (1997): 153-179.\n\n[2]: Hongyi Zhang and Suvrit Sra. First-order methods for geodesically convex optimization. In\nConference on Learning Theory, pages 1617–1638. PMLR, 2016.\n\n[3]: Barachant, Alexandre, et al. \"Multiclass brain–computer interface classification by Riemannian geometry.\" IEEE Transactions on Biomedical Engineering 59.4 (2011): 920-928.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The authors seem over optimistic about the performance of BW metric compared to AI metric.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "10"}, "extra_text": {}}, {"review_id": "MiVLzij7_QM", "reviewer": "Reviewer_VX14", "summary": "This paper studies the Bures-Wasserstein (BW) geometry of symmetric positive-definite (SPD) matrices and compares it to the commonly-used affine-invariant (AI) geometry in the context of Riemannian optimization. It is shown theoretically that the BW geometry has several advantages. Particularly, it is shown that (i) the BW metric is linearly dependent on the SPD matrix, (ii) the condition number of the Riemannian Hessian of the BW metric is smaller than that of the AI metric, and (iii) the curvature constant of the BW metric is smaller than that of the AI metric. Based on these properties, it is shown that the convergence of Riemannian steepest descent and Riemannian trust-region based on the BW metric is faster than based on the AI metric. The theoretical results are demonstrated on simulations of six problems: weighted least squares, Lyapunov equations, trace regression, metric learning, log-det maximization, and Gaussian mixture model.", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The Riemannian geometry of SPD matrices has become popular recently, and the choice of geometry/metric is of paramount importance. This paper provides a systematic way to choose the metric from an optimization standpoint. As I outline below, I believe the paper has few drawbacks. However, despite its weaknesses, I believe this paper presents a solid contribution to the increasing literature and practice of Riemannian geometry of SPD matrices.\n\nDetailed comments:\n- The argument that a smaller condition number and/or smaller constant curvature is necessarily better is problematic. In terms of optimization, it is indeed shown for example that a smaller condition number translates to better convergence. However, if this is the only criterion, then one could claim that the Euclidean metric outperforms both the AI and the BW metrics. In other words, the paper claims that “simpler” metrics are better for optimization, but it does not consider the effect of the metric on the induced geometry. \n- To make my comment above more concrete, I would love to see, in addition to the simulations, results on some real data, where the contribution of using Riemannian geometry is established, and where the choice of the metric also determines the “quality” of the analysis.\n- In the experiments (line 265), the Frobenius (Euclidean) norm is used as a measure of distance to the desired solution. This departs from the entire message of the paper advocating for Riemannian geometry and could be particularly problematic when measuring the convergence rate.\n- In many figures (e.g., Fig 1 (a)+(b)), it appears as if not all the graphs have reached convergence. Why do the graphs have different lengths?  \n- line 185: the geodesic gamma(t) is undefined (a missing reference to Prop 3 in the appendix).\n- line 200: it would be better to explicitly define “geodesic triangle”, although it is implied from the context.\n\n**After rebuttal:**\nI thank the authors for their response to my comments.\nAfter reading the other reviews and the authors' responses, my two concerns regarding the Euclidean metric and norm still remain. That said, I am increasing my score and recommend accepting the paper, since I believe it presents analysis and results that will be interesting to the community of people working on the Riemannian geometry of SPD matrices.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The limitations of the work could be better presented and discussed, especially since the paper considers only optimization properties. I do not foresee any potential negative societal impact.", "review_text": "", "extra_scores": {"ethical_concerns": "I do not have any ethical concerns.", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4"}, "extra_text": {}}, {"review_id": "MTrsy6fpAo5", "reviewer": "Reviewer_LtgQ", "summary": "The authors consider the Riemannian optimization with Bures-Wasserstein (BW) geometry for symmetric positive definite (SPD) matrix manifold. The authors compare the proposed approach with its counterpart with the popular Affine-Invariant (AI) geometry.\n\nThe authors illustrate that the proposed approach (with BW geometry) is more suitable than its counterpart with AI geometry when SPD matrices are ill-conditioned.\n", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Riemannian optimization with SPD matrices is important for several approaches in machine learning. The authors consider the new raising Bures-Wasserstein geometry for this Riemannian optimization.\n\nOverall, the paper is easy to follow and well-organized. The authors also evaluate the proposed approaches on several problems in machine learning to illustrate the theoretical findings.\n\nBy analyzing the condition number (Lemma 1), the authors show that BW geometry is more suitable for Riemannian optimization with SPD matrices than its counterpart with the popular AI geometry when SPD matrices are ill-conditioned. \n\nI have some concerns as follow:\n+ In Remark 1, from Equ. (1), why BM metric has a linear dependence on X? (but not its inverse?)\n+ For Lemma 5, is there any assumption (e.g., the strongly convex and smooth) for the $f$ function?\n+ It seems that results about convergence in Theorem 1 and 2 are straightforward results from the literature? It is better in case the authors give more discussions about what new results are here, comparing with the literature?\n+ For Figure 2, row 1, (b) and (d), why the AI and LE has the same behaviors for some first iteration, but LE converges better with more iterations?\n+ Although the authors perform several experiments, most of them are with synthesized data. It is better in case the authors consider more illustrations with real-world problem/data (such as with the metric learning).\n\n~~~~~~~~\nAfter the rebuttal:\n\nI thank the authors for the response. I agree with the other reviewers that it is better to give general properties for the loss function when AI geometry is better (currently, only 1 example is listed in the submission). Overall, I keep leaning on the positive side.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "Yes", "review_text": "", "extra_scores": {"ethical_concerns": "N/A", "needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "8 hours"}, "extra_text": {}}]

Author rebuttal:

Meta-review: The authors consider the Riemannian optimization with Bures-Wasserstein geometry for symmetric positive definite matrix manifold. Since the reviewers agree that the content of the paper will interest many researchers working on this subject, we conclude that the paper is worthy of acceptance. I think that the comment on the Frobenius norm by a reviewer is important. I expect that the authors add some comments on that issue in the revision. 

Decision: Accept (Poster)