Metric MDS & Wasserstein Geometry

• Metric MDS with Wasserstein Geometry is a framework that integrates optimal transport and Riemannian geometry to produce faithful low-dimensional embeddings of structured, non-Euclidean data.
• It employs explicit formulas and spectral decomposition techniques to compute Wasserstein distances for data types such as SPD(n) matrices and von Mises–Fisher distributions.
• The approach guarantees global geodesic convexity and computational stability, enabling scalable and interpretable manifold learning across various applications.

Metric Multidimensional Scaling (MDS) with Wasserstein Geometry refers to the extension of classical MDS methodologies to data endowed with the (optimal transport) Wasserstein metric, especially in the context of structured statistical objects such as covariance matrices (SPD(n)), distributions on spheres, empirical distributions, and parametric families like von Mises–Fisher (vMF) laws. This framework leverages explicit Riemannian geometry induced by the Wasserstein metric to enable faithful, interpretable, and globally convex low-dimensional embeddings of non-Euclidean data, with rigorous spectral and algorithmic underpinnings.

1. Theoretical Underpinnings: Wasserstein Geometry on Structured Spaces

The 2-Wasserstein metric, defined via optimal transport, induces a Riemannian structure on several statistical manifolds. For the manifold SPD(n) of $n \times n$ real symmetric positive-definite matrices, the 2-Wasserstein metric coincides with the geodesic distance between zero-mean Gaussian distributions with these covariances:

$$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$

where $(P^{1/2} Q P^{1/2})^{1/2}$ is evaluated via eigendecomposition (Luo et al., 2020, Malagò et al., 2018). The induced geometry is globally geodesically convex, with non-negative but non-constant sectional curvature, and the absence of cut-locus or conjugate points. The exponential and logarithm maps, as well as geodesic paths, admit closed-form expressions involving matrix square roots and solutions to Sylvester (or Lyapunov) equations, all computable in $O(n^3)$ time due to dominance by eigendecompositions.

For vMF distributions on the hypersphere $S^{d-1}$, a Wasserstein-like distance in the high-concentration regime reduces to a sum of squared angular (great-circle) distance and a variance-like concentration term:

$$d_W^2(\mathrm{vMF}(\mu_1, \kappa_1), \mathrm{vMF}(\mu_2, \kappa_2)) = \arccos^2(\langle\mu_1, \mu_2\rangle) + (d-1) \left(\frac{1}{\sqrt{\kappa_1}} - \frac{1}{\sqrt{\kappa_2}}\right)^2$$

This formula directly encodes the intrinsic geometry of the parameter space $S^{d-1} \times (0,\infty)$ with product Riemannian structure (You et al., 19 Apr 2025).

2. Metric MDS Workflow with Wasserstein Distances

The classical metric MDS workflow involves constructing the pairwise distance matrix $D_{ij}$ corresponding to $d_W(P_i, P_j)$ for SPD(n) matrices, or their analogues for other statistical objects, followed by double-centering and eigendecomposition to yield a Euclidean embedding:

1. Collect data points $\{P_i\}$ (e.g., SPD(n) matrices, distributions).
2. Build $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$0 using explicit formulas.
3. Form the double-centered Gram matrix $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$1, $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$2.
4. Perform spectral decomposition $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$3. Embedding $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$4 is given by $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$5 for the top $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$6 positive eigenvalues (Luo et al., 2020, Malagò et al., 2018).
5. For non-Euclidean or approximately metric distances, stress-majorization (“SMACOF”) approaches or nearest PSD projections of $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$7 may be used (Haviv et al., 2024, Lim et al., 2022).

Due to the global geodesic convexity and absence of cut-locus in Wasserstein–SPD(n), the pairwise distance matrix corresponds to unique minimal geodesics and admits stable Euclidean approximation via MDS.

3. Spectral Properties and Generalized cMDS Theory

Generalized classical MDS (cMDS) extends to arbitrary metric measure spaces $$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$8 using the kernel

$$d_W(P, Q)=\left\{\mathrm{tr}(P)+\mathrm{tr}(Q)-2\,\mathrm{tr}[(P^{1/2} Q P^{1/2})^{1/2}]\right\}^{1/2}$$9

The associated self-adjoint operator $(P^{1/2} Q P^{1/2})^{1/2}$0 possesses a discrete spectrum with positive and negative eigenvalues. The sum of negatives, $(P^{1/2} Q P^{1/2})^{1/2}$1, quantifies the "non-flatness" or $(P^{1/2} Q P^{1/2})^{1/2}$2 distortion relative to a true Euclidean embedding. For finite samples, this equals the sum of negative eigenvalues of the Gram matrix. Stability with respect to Gromov–Wasserstein perturbations is furnished by explicit Lipschitz bounds (Lim et al., 2022).

In the Wasserstein–SPD(n) case, the distance matrix is generally Euclidean in the Riemannian sense, so negative eigenvalues are small or vanish, ensuring high-quality low-dimensional embeddings with stable, convex stress landscapes.

4. Algorithmic and Computational Considerations

For SPD(n) matrices, the dominant cost arises from matrix square roots and Sylvester/Lyapunov solves: each distance requires $(P^{1/2} Q P^{1/2})^{1/2}$3 operations, yielding $(P^{1/2} Q P^{1/2})^{1/2}$4 for $(P^{1/2} Q P^{1/2})^{1/2}$5 samples. Final Gram decomposition is $(P^{1/2} Q P^{1/2})^{1/2}$6. For large-scale empirical distributions, scalable algorithmic solutions are needed.

Wasserstein Wormhole employs a transformer-based autoencoder to learn Euclidean embeddings that preserve pairwise OT (Wasserstein) or Sinkhorn-regularized distances in embedding space. The loss matches $(P^{1/2} Q P^{1/2})^{1/2}$7, where $(P^{1/2} Q P^{1/2})^{1/2}$8 are computed Sinkhorn divergences. The encoding-decoding paradigm, in conjunction with SGD and mini-batching, attains linear scaling in the inference phase ($(P^{1/2} Q P^{1/2})^{1/2}$9 for $O(n^3)$0 distributions, $O(n^3)$1 embedding dimension), drastically improving over $O(n^3)$2 complexity of all-pairs OT computations (Haviv et al., 2024).

For vMF mixtures, Wasserstein-like distances can be computed in $O(n^3)$3 time via direct evaluation of the angular and variance terms (You et al., 19 Apr 2025).

5. Practical Applications and Empirical Validation

Wasserstein metric MDS has broad applications in computational geometry, machine learning, single-cell biology, and structured statistical analysis:

• SPD(n): Embedding clouds of covariance matrices, with faithful low-dimensional representations that preserve statistical geometry—used in AI, computer vision, and hypothesis testing (Luo et al., 2020, Malagò et al., 2018).
• Empirical Distributions: Wasserstein Wormhole achieves high-fidelity embeddings for large-scale point cloud data, supporting efficient OT distance computation, barycenter estimation, and interpolation in latent space (Haviv et al., 2024).
• vMF Distributions: Embeddings derived from Wasserstein-like geometry cleanly separate clusters differing in mean direction and concentration, as demonstrated on synthetic and high-dimensional data (You et al., 19 Apr 2025).
• General Metric Spaces: The cMDS framework, extended to arbitrary metric measure spaces, is theoretically well-behaved under Gromov–Wasserstein perturbations and provides spectral invariants quantifying intrinsic non-flatness (Lim et al., 2022).

6. Geometry-Induced Convexity, Stability, and Interpretability

The global geodesic convexity and non-negative sectional curvature of SPD(n) with $O(n^3)$4 guarantee uniqueness and stability of geodesic paths and absence of ambiguities in the metric MDS embedding. The spectral invariants from cMDS (sum of negative eigenvalues) measure O(T)-induced non-flatness and embedding distortion. In the Wasserstein Wormhole approach, the theoretical lower and upper bounds on the minimal stress achievable by embedding a non-Euclidean distance matrix are rigorously related to the negative spectrum of the “centered distance” matrix (Haviv et al., 2024). This facilitates principled control and evaluation of embedding error in practical settings.

7. Extensions and Recent Advances

Recent work extends the Wasserstein metric MDS paradigm to data families with intrinsic geometry (e.g., the vMF family), introducing metrics that respect both angular and variance-like factors, and providing explicit algorithms for mixture reduction and high-dimensional inference (You et al., 19 Apr 2025). The combination of Riemannian geometry, spectral theory, and scalable neural architectures (e.g., Wasserstein Wormhole transformers) greatly broadens the practical domain of metric MDS with Wasserstein geometry.

A plausible implication is that as computational and statistical tools further mature, MDS with Wasserstein geometry will underpin a wide class of structure-preserving representations for non-Euclidean, high-dimensional, and distributional data. Its rigorous geometric foundations and scalability uniquely position it for interpretable, stable, and theoretically principled manifold learning.