GEORCE-FM: Joint Fréchet Mean & Geodesics

• GEORCE-FM is an efficient algorithm that jointly computes the Fréchet mean and associated geodesics by framing the problem as a discrete-time optimal control task in a local coordinate chart.
• The method circumvents repeated geodesic boundary value problems by employing closed-form joint updates, leading to global convergence and local quadratic convergence guarantees.
• Empirical results demonstrate that GEORCE-FM outperforms traditional methods in runtime and solution quality, with extensions to adaptive mini-batching and Finsler manifold settings.

The GEORCE-FM algorithm is an efficient and theoretically grounded method for simultaneously computing the Fréchet mean and associated geodesics on both Riemannian and Finsler manifolds. By framing the joint computation as a single discrete-time optimal control problem in a local coordinate chart, GEORCE-FM circumvents the traditional need to solve nested geodesic boundary value problems in each iteration of Fréchet mean estimation. The method accommodates adaptive extensions for large-scale data settings and provides rigorous guarantees of global convergence and local quadratic convergence. Empirical evidence demonstrates that GEORCE-FM dominates earlier manifold mean estimation techniques in both runtime and solution quality across classical and learned manifold models.

1. Fréchet Mean Estimation on Manifolds

The Fréchet mean $\mu$ for a set of points $a_1, \dots, a_N$ on a Riemannian manifold $(\mathcal{M}, g)$ is the minimizer of the sum of squared distances: $$\mu = \arg\min_{y \in \mathcal{M}} \sum_{i=1}^N \operatorname{dist}^2(y, a_i)$$ with $\operatorname{dist}$ the geodesic distance. Computing $\operatorname{dist}(y, a_i)$ generally requires solving a geodesic boundary value problem, for example by integrating the ODE

$$\ddot\gamma^k + \Gamma^k_{ij}(\gamma)\dot\gamma^i\dot\gamma^j = 0, \qquad \gamma(0)=y,\, \gamma(1)=a_i$$

or minimizing the energy functional

$\E_\gamma(y, a_i) = \frac12 \int_0^1 \dot\gamma(t)^\top G(\gamma(t)) \dot\gamma(t) dt$

where $G$ is the metric tensor in a local chart. Most existing algorithms, therefore, embed a costly geodesic solver within each Fréchet mean update, leading to compounding computational expense.

2. Joint Energy Formulation and Optimization

GEORCE-FM addresses this computational bottleneck by jointly optimizing the mean location $y$ and all associated geodesics’ discretizations. In local coordinates, each geodesic $\gamma_i$ from $a_i$ to $y$ is discretized by a set of $T$ points $x_{0,i}, \ldots, x_{T,i}$ with step vectors $u_{t,i} = x_{t+1,i} - x_{t,i}$, subject to boundary conditions $x_{0,i} = a_i$, $x_{T,i} = y$, and state updates $x_{t+1,i} = x_{t,i} + u_{t,i}$. The joint discrete energy minimized by GEORCE-FM is: $\min_{\substack{y,\x_{t,i}, u_{t,i}}} \sum_{i=1}^N\sum_{t=0}^{T-1} u_{t,i}^\top G(x_{t,i}) u_{t,i}$ with Lemma 1 showing that as $T \to \infty$, minimizing this energy over $y$ matches the original Fréchet mean problem. The algorithm thus eliminates the repeated inner geodesic solves typically required.

3. Algorithmic Workflow and Pseudocode

The GEORCE-FM algorithm operates in the following principal steps, all conducted within a single local chart (i.e., no need for chart transitions or global log/exponential map evaluations):

1. Initialization: Set $y^{(0)}$ (e.g., $a_1$), $x_{t,i}^{(0)} = a_i + (t/T)(y^{(0)} - a_i)$, $u_{t,i}^{(0)} = x_{t+1,i}^{(0)} - x_{t,i}^{(0)}$.
2. Metric Evaluation: For all $t, i$, set $G_{t,i} = G(x_{t,i})$.
3. Gradient Computation: For all $t>0, i$, compute $$\nu_{t,i} = \partial_x(u_{t,i}^\top G(x) u_{t,i})|_{x = x_{t,i}}$$.
4. Closed-Form Joint Update:
   • Compute

   $$W = \sum_i\left(\sum_t G_{t,i}^{-1}\right)^{-1},\quad V = \sum_i\left(\sum_t G_{t,i}^{-1}\right)^{-1}a_i - \frac12\sum_i\left(\sum_t G_{t,i}^{-1}\right)^{-1}\sum_t G_{t,i}^{-1}\sum_{j>t}\nu_{j,i}$$

• Update the mean: $y = W^{-1}V$.
• Compute dual variables and new controls for all $i, t$:

  $$\mu_{T-1,i} = \left(\sum_t G_{t,i}^{-1}\right)^{-1}(2(a_i-y) - \sum_t G_{t,i}^{-1}\sum_{j>t}\nu_{j,i})$$

  $$u_{t,i} = -\frac12 G_{t,i}^{-1}(\mu_{T-1,i} + \sum_{j>t} \nu_{j,i}),\quad x_{t+1,i} = x_{t,i} + u_{t,i}$$

1. Projected Line Search: Armijo backtracking line search is performed over the control updates to ensure sufficient energy decrease.
2. Mixing: Weighted averaging (with step size $\alpha$) of previous and proposed $(x_{t,i}, u_{t,i})$ ensures convergence and regularization.
3. Repeat until the decrease in sum energy $\sum_{i,t} u_{t,i}^\top G(x_{t,i})u_{t,i}$ meets the tolerance.

This algorithmic procedure matches Algorithm 2 and its Riemannian update rules as given in (Rygaard et al., 6 Nov 2025). The Finslerian extension involves similar machinery but adapts $G(x, u)$ to depend on both position and direction (see Section 4).

4. Extensions: Finsler Manifolds and Adaptive Mini-Batching

Finsler Geometry: In Finsler manifolds, the metric is defined via a fundamental tensor $$G_{ij}(x, v) = \frac12 \partial_{v^i} \partial_{v^j}[F^2(x, v)]$$ that is homogeneous and strongly convex in the tangent vector $v$. GEORCE-FM extends naturally to this setting by replacing the Riemannian $G(x)$ with $G(x, u)$, and updating the optimization objective accordingly: $$\min \sum_{i=1}^N\sum_{t=0}^{T-1} u_{t,i}^\top G(x_{t,i}, u_{t,i}) u_{t,i}$$ The update rules, line-search procedure, and convergence proofs carry over nearly unchanged.

Adaptive (Stochastic) GEORCE-FM: For large $N$, computational complexity is reduced by employing mini-batch stochastic updates: - At each outer iteration, a random subset $\mathcal{I}$ of size $n$ is sampled. - GEORCE-FM is run for several steps within $\mathcal{I}$ yielding estimates $(\tilde W, \tilde V)$. - The global mean iterate is formed using an exponential moving average in $W, V$, with decaying step size $\alpha$. - Under mild unbiasedness and correlation assumptions, the resulting sequence converges in expectation to a local minimizer (Proposition 7).

5. Theoretical Guarantees

Global Convergence: If the initial iterates satisfy the boundary constraints, and with the Armijo line search, GEORCE-FM generates a non-increasing sequence of objective values; energy is bounded below, so accumulation points are local minima (Prop. 5).

Local Quadratic Convergence: When the discrete energy $E(z)$ is strongly convex in a neighborhood of a unique minimizer $z^*=(x^*,u^*,y^*)$, and the line search accepts full steps ($\alpha=1$), quadratic convergence holds: $\|z^{(k+1)} - z^*\| \leq C \|z^{(k)} - z^*\|^2$ for some $C>0$, by a Taylor expansion argument under smoothness and convexity assumptions.

Stochastic Convergence (Adaptive Mini-Batching): When mini-batch updates are unbiased and satisfy positivity (Eq. (26)), and each mini-batch uses a line search, convergence in expectation is assured (Prop. 7).

6. Computational Complexity and Empirical Performance

Theoretical Complexity

• Per Iteration: Each GEORCE-FM iteration is $\mathcal{O}(N T d^3)$ due to $N\times T$ matrix inversions in $d$-dimensional charts.
• Comparison: This is comparable to standard Fréchet mean algorithms (which require $N$ geodesic solves per iteration at a similar cost) but benefits from amortized reuse of chart structures and avoids repeated boundary value problem solves.

Empirical Results

Benchmark	Baseline Iterations (ADAM/SGD)	GEORCE-FM Iterations	Acceleration (Runtime)	Accuracy
$\mathbb{S}^n$, $E(n)$	$100$–$500$	$3$–$5$	$10$–$30\times$	Comparable
Fisher–Rao statistics	$20+$	$5$–$6$	$5\times$	Comparable
Learned VAE manifolds	Minutes	Tens of seconds	$2\times$–$6\times$	Comparable
Adaptive ($N=1000$)	Full batch	$10\%$ geodesics	$10\times$	$\sim$Few %

• GEORCE-FM converges rapidly on synthetic spheres, ellipsoids, and information geometry manifolds.
• On real-data latent spaces (e.g., MNIST, CelebA decoder manifolds), GEORCE-FM computes both geodesics and Fréchet mean on the order of seconds.
• In Finsler and Randers-type geometries, GEORCE-FM attains lower energy solutions and $100$–$10,000\times$ improved runtimes versus BFGS, ADAM, and Newton-type baselines.

7. Significance, Limitations, and Directions

GEORCE-FM’s principal innovation is the joint optimization over mean location and geodesic fields, facilitated by closed-form local updates within a single chart. This avoids compounding the computational and numerical overhead of embedded geodesic solvers typical of prior approaches. Its design allows extensibility to asymmetric and velocity-dependent metrics (Finsler), as well as natural stochastic mini-batching for large-scale data.

Limitations include the reliance on local charts (potentially problematic on highly nontrivial manifold topologies) and cubic scaling in intrinsic dimension $d$ due to matrix inversions, although $d$ is typically modest.

• A plausible implication is that further acceleration may be possible using approximate inversion strategies or manifold preconditioning for very high-dimensional problems.
• The uniform chart assumption could be relaxed for more general manifold geometries via atlas-based or manifold learning techniques.

Recent empirical benchmarking on classical and learned manifolds supports the algorithmic dominance of GEORCE-FM in both runtime and accuracy for Fréchet mean estimation, with strong theoretical backing for its global and local convergence properties (Rygaard et al., 6 Nov 2025).