Riemannian Zeroth-Order Gradient Estimation with Structure-Preserving Metrics for Geodesically Incomplete Manifolds

AI summary

Overview

Research area: Riemannian optimization and zeroth-order (derivative-free) optimization, bridging non-smooth/black-box machine learning with differential geometry.

Technical level: Advanced. The paper assumes familiarity with Riemannian manifolds, exponential maps, retractions, sectional curvature, and stochastic gradient convergence proofs.

Scope: The paper builds a geodesically complete "structure-preserving" metric for manifolds whose canonical metric is geodesically incomplete, then supplies an intrinsic zeroth-order gradient estimator, an unbiased sampling routine, and SGD convergence guarantees that it validates on synthetic problems and a mesh optimization task.

What This Paper Is About

Riemannian zeroth-order methods approximate gradients using only function evaluations at perturbed points, but the perturbation is applied through the exponential map, which is only locally defined. When a manifold inherits a Euclidean metric from an ambient space, that metric may be geodesically incomplete, so a randomly sampled tangent vector can fall outside the exponential map's domain and the perturbed point becomes undefined. The paper's goal is to make Riemannian zeroth-order optimization work on such incomplete manifolds while still returning points that are stationary under the original metric.

Key Contributions

1. Structure-preserving metric construction. The authors define a metric g′ that is structure-preserving with respect to a given g (Definition 2.5) if it is geodesically complete, conformally equivalent to g, and preserves ε-stationarity. Theorem 2.6 shows such a g′ always exists for any smooth manifold (possibly non-compact) equipped with any Riemannian metric g, by modifying the classical Nomizu-Ozeki construction with an adjusted conformal coefficient.

2. Intrinsic zeroth-order gradient estimation. Rather than embedding the new metric in a Euclidean ambient space, the authors analyze the classical symmetric two-point estimator (Equation 3) purely intrinsically. Theorem 2.7 bounds its mean-squared error in terms of the dimension d and a uniform sectional curvature bound κ, under bounded third- and fourth-order derivatives (Assumption 2.3) and bounded sectional curvature (Assumption 2.4).

3. Unbiased sampling under general metrics. The paper shows that the common rescaling approach (draw a Gaussian vector, normalize to g-unit length) is biased under non-Euclidean metrics, and instead applies rejection sampling (Algorithm 1). Proposition 2.8 proves the output is exactly uniform over the set {v ∈ ℝ^d : vᵀAv = 1} determined by the metric.

4. Convergence guarantees and empirical validation. Theorem 2.9 establishes SGD convergence under a general geodesically complete metric, and Corollary 2.10 transfers this to the original Euclidean metric under stated conditions. Synthetic experiments and a mesh optimization experiment support the theory.

Main Findings

• Curvature enters the estimation error directly. For v sampled uniformly from the g-unit sphere, Theorem 2.7 gives an MSE bound of the form (1 + μ²κ²)/d · ‖∇f(p)‖²_p + O(μ²), where κ is a uniform upper bound on absolute sectional curvature and μ is the perturbation stepsize. The bound is required to hold for μ² ≤ min{1/(d−1), 1/2 + 6/d + 8/d²}.

• The flat case recovers the classical result. When κ = 0, the curvature term vanishes and the bound reduces to the standard approximation error for zeroth-order gradient estimation in Euclidean spaces.

• SGD rate under a general metric. Theorem 2.9 gives min over t ≤ T of ‖∇f(p_t)‖²_{p_t} ≤ C₁·d/(ηT) + C₂·η + C₃·d²μ², with η ≲ √(d/T) and μ² ≲ √(d/T). Choosing μ ≲ (1/d²)√(d/T) yields a rate of order √(d/T).

• Matching the best-known complexity under extra conditions. Corollary 2.10 states that if the original Euclidean metric g_E is geodesically complete, or if the set of ε-stationary points under g_E is compact, then at most T ≤ O(d/ε⁴) iterations are needed to reach min over t ≤ T of E[‖∇f(p_t)‖²_{p_t,g_E}] ≤ ε². The paper notes that Theorem 2.9 covers more general scenarios but with potentially worse complexity, and that the stationary-point preservation is asymmetric: an ε-stationary point under g′ need not be one under g.

• Sampling bias matters for stability. In the synthetic experiment on ℝ^d with metric g_A(u,v) = uᵀAv, over 16 independent runs with identical SGD hyperparameters, rejection sampling (Algorithm 1) consistently outperformed rescaling sampling, and rescaling even diverged on the logistic loss objective.

• Lower curvature gives more accurate estimation. In a second synthetic experiment comparing four conformally equivalent metrics on the probability simplex, MSE of the zeroth-order estimator (measured over 50,000 independent trials) decreased as sectional curvature decreased for both the KL objective against the simplex centroid and the squared Euclidean objective, matching the predicted trend.

• Mesh optimization converges more stably. On the Helmholtz equation with wave number k = 10, using a 200 × 200 ground-truth mesh and a 20 × 20 coarse mesh over 20,000 optimization steps, the unconstrained baseline showed unstable fluctuations (most notably around step 16,000), reversion stalled after 8,000 steps, and soft projection showed little improvement beyond 14,000 steps, while the structure-preserving approach reduced error throughout and reached the lowest final MSE.

• The reported text is truncated. The provided paper content cuts off mid-sentence in Section 3.3, so any later concluding remarks, additional experiments, or final quantitative tables are not reported here.

Methodology in Plain English

The authors start from a practical obstacle: to estimate a gradient without derivatives on a curved space, you nudge a point a small distance in a random direction and compare function values, but the "nudge" is only defined while you stay inside the manifold. Their fix has three parts.

First, they replace the problematic metric with a new one built by scaling the old metric pointwise by a positive smooth function. The scaling is chosen so that geodesics never run off the manifold, so the perturbation always lands at a valid point. Crucially, the construction is designed so that the stationary points of the new metric line up with stationary points of the old metric.

Second, because this new metric is generally not the Euclidean metric inherited from an ambient space, existing analyses do not apply. The authors therefore analyze the two-point estimator using only intrinsic manifold quantities — distances, gradients, and curvature measured in the metric itself. This yields an error bound whose leading correction term is driven by how curved the manifold is.

Third, they address how to draw random tangent directions. Uniform sampling on a sphere defined by a general metric is equivalent to uniformly sampling the ellipsoid {v : vᵀAv = 1}, and simply normalizing a Gaussian vector gives a biased distribution. They use rejection sampling instead, rejecting proposed points with a probability derived from the eigenvalues of A, and prove the accepted samples are exactly uniform.

Finally, they plug the estimator into standard Riemannian SGD with a retraction, require the learning rate and perturbation stepsize to shrink at specified rates, and derive a convergence rate. Transferring the guarantee back to the original metric requires an extra assumption, either that the original metric was complete all along or that the set of near-stationary points is compact, which bounds the conformal scaling factor uniformly.

Why This Matters

Research impact. Existing Riemannian zeroth-order analyses typically assume the manifold is viewed as a submanifold of Euclidean space and inherits the Euclidean metric, which silently requires that metric to be geodesically complete. This paper supplies a general-purpose repair: it extends zeroth-order guarantees from the Euclidean-metric special case to a much broader class of Riemannian metrics. It also makes the curvature-dependence of zeroth-order error explicit, which is a structural insight rather than a technical detail. The Nomizu-Ozeki-based construction gives the remedy a classical geometric foundation while the sampling analysis connects it to a practical numerical requirement.

Real-world applications named in the paper:

• Mesh optimization. Physical simulations such as CFD-GCN adjust vertex positions using an external PDE solver that lacks auto-differentiation. The valid configuration space excludes boundary configurations and is therefore geodesically incomplete under the Euclidean metric.
• Irrigation system layout design. Optimizing sprinkler coordinates to maximize coverage relies on non-differentiable simulation of spray overlap, pressure, and wind. Sprinklers cannot sit on the field boundary, so the valid configuration space is an incomplete open set.
• Covariance matrix estimation. Estimating covariance for PCA or Gaussian modeling means optimizing over the set of d × d symmetric positive definite matrices, an open convex cone that is geodesically incomplete.
• Black-box and non-differentiable pipeline optimization more broadly. Any objective involving a non-differentiable external solver or a black-box module is a candidate, since the framework needs only function evaluations.

Industry relevance. Derivative-free optimization over constrained, curved parameter spaces appears in physical simulation, engineering design, and statistics. The paper's specific value for practice is its treatment of validity constraints: rather than heuristically rejecting, projecting, or reverting invalid steps, the structure-preserving metric is designed so invalid steps are never proposed. The sampling result is also immediately practical, since it warns that a common shortcut — normalize a Gaussian to non-Euclidean unit length — can destabilize or even diverge training.

Future Directions

• Relaxing the smoothness assumptions. The proof of Theorem 2.9 notes that the L-smoothness assumption (Assumption 2.1) can be relaxed to an expected-smoothness condition, which is a natural next step to confirm and analyze.
• Closing the complexity gap for genuinely incomplete metrics. Theorem 2.9 covers general incomplete metrics but can give worse bounds than the geodesically complete case. Characterizing exactly when the O(d/ε⁴) bound transfers, beyond the two conditions in Corollary 2.10, remains open.
• Making the conformal construction numerically concrete. The construction relies on a positive smooth conformal function h, and the transfer argument requires h to admit a uniform upper bound. How to choose h effectively in practice for large or irregular meshes is not resolved in the reported content.
• Extending and scaling the empirical evaluation. The provided content reports a single mesh optimization problem, a synthetic sampling-bias study, and a curvature study. Broader benchmarking against other zeroth-order Riemannian methods and higher-dimensional applications (such as the covariance estimation setting) is a logical follow-up.

Target Audience

Researchers and graduate students working on Riemannian optimization, derivative-free or black-box optimization, and the theory of stochastic gradient methods on manifolds. It is also relevant to practitioners in scientific computing and physical simulation who need to optimize over constrained geometric domains without access to gradients, and who can tolerate an advanced, proof-heavy presentation. Readers without a background in Riemannian geometry will need the appendix preliminaries to follow the main results.

Authors’ abstract

In this paper, we study Riemannian zeroth-order optimization in settings where the underlying Riemannian metric $g$ is geodesically incomplete, and the goal is to approximate stationary points with respect to this incomplete metric. To address this challenge, we construct structure-preserving metrics that are geodesically complete while ensuring that every stationary point under the new metric remains stationary under the original one. Building on this foundation, we revisit the classical symmetric two-point zeroth-order estimator and analyze its mean-squared error from a purely intrinsic perspective, depending only on the manifold's geometry rather than any ambient embedding. Leveraging this intrinsic analysis, we establish convergence guarantees for stochastic gradient descent with this intrinsic estimator. Under additional suitable conditions, an $ε$-stationary point under the constructed metric $g'$ also corresponds to an $ε$-stationary point under the original metric $g$, thereby matching the best-known complexity in the geodesically complete setting. Empirical studies on synthetic problems confirm our theoretical findings, and experiments on a practical mesh optimization task demonstrate that our framework maintains stable convergence even in the absence of geodesic completeness.