Weiszfeld-Type Algorithms Overview

• Weiszfeld-type algorithms are iterative fixed-point methods that minimize weighted distances to anchor points in both Euclidean and non-Euclidean spaces.
• They extend the classical Fermat–Weber formulation by incorporating convex constraints, d.c. decompositions, and subgradient strategies to handle singularities.
• Their robust convergence properties and low per-iteration cost make them effective for applications in location science, robust PCA, and rotation averaging.

Weiszfeld-type algorithms are iterative fixed-point methods designed to solve the geometric median or Fermat–Weber type location problems, in which the objective is to find a point minimizing a sum (or more general combination) of norms or distances—commonly weighted Euclidean distances—to a collection of anchor points, convex sets, or other geometric objects. Originally introduced by E. Weiszfeld in 1937 for the unconstrained, finite-dimensional Fermat–Weber problem, these algorithms and their variants have since been rigorously analyzed and extended to constrained, non-Euclidean, and nonconvex settings, as well as to applications such as robust PCA and rotation averaging. Their central appeal lies in their descent property, modest per-iteration complexity, and strong theoretical guarantees under convexity and nondegeneracy assumptions.

1. Classical Formulation and Theoretical Foundations

The prototypical setting is the weighted Fermat–Weber problem: given anchor points $a_1, \dots, a_m$ in $\mathbb{R}^n$ and positive weights $w_1, \dots, w_m$, find

$$\min_{x\in\mathbb{R}^{n}} F(x), \qquad F(x) = \sum_{i=1}^{m} w_i \|x-a_i\|.$$

When the $a_i$ are not collinear, $F$ is strictly convex and possesses a unique minimizer. The necessary and sufficient optimality condition for $x^*$ outside of the anchor set is

$$\sum_{i=1}^m w_i \frac{x^* - a_i}{\|x^* - a_i\|} = 0.$$

This yields the fixed-point equation

$$x^* = \frac{\sum_{i=1}^{m} w_i \frac{a_i}{\|x^*-a_i\|}}{\sum_{i=1}^{m} w_i \frac{1}{\|x^*-a_i\|}},$$

which defines the standard Weiszfeld mapping. The resulting iterative scheme is globally convergent as established by Kuhn (1973), so long as all iterates avoid the anchor points; in the convex-analytic framework, the objective strictly decreases at each iteration until convergence to the unique minimizer (Mordukhovich et al., 2013).

2. Extensions: Constraints, Generalized Sets, and Modified Iterations

Weiszfeld-type methods have been generalized to handle convex constraints (i.e., location restricted to a closed convex $C \subset \mathbb{R}^n$ or in a real Hilbert space), and to models with distances to convex sets and possibly negative weights.

Constrained Weiszfeld Problem: For closed convex $C$, the constrained objective

$$\min_{x \in C} f(x) := \sum_{i=1}^{m} w_i \|x-a_i\|$$

is approached by projecting each Weiszfeld step onto $C$: $$x_{k+1} = P_C \left( \frac{\sum_{i=1}^m w_i a_i / \|x_k - a_i\|}{\sum_{i=1}^m w_i / \|x_k - a_i\|} \right ),$$ where $P_C$ denotes the nearest-point projection. Under standard assumptions, the sequence remains in $P_C(\mathrm{co}\{a_i\})$ and converges strongly to the unique solution (Nguyen, 2018, Torres, 2012).

Sum-of-Distances to Convex Sets, D.C. Decomposition: When objective terms include weighted distances to convex sets with possible negative weights, the function becomes a difference-of-convex (d.c.) function. In this context, the DCA (difference-of-convex algorithm) framework is used; the inner step at each iteration involves a generalized Weiszfeld mapping acting on sets, followed by projection. This algorithm is shown to produce DCA-critical points under either compactness or coercivity assumptions (An et al., 2014).

Extended Power Weiszfeld for $q$-norms and Singularity Handling: The cost

$$C_q(y) = \sum_{i=1}^m \xi_i \|y - x_i\|^q, \quad 1 \leq q < 2$$

leads to a generalized Weiszfeld iteration with weights $\|y - x_i\|^{q-2}$, but at $y = x_i$ for $q < 2$ the expression becomes singular. Recent methods circumvent this using subgradient descent at the singular points to guarantee global convergence and superlinear local convergence in special cases (Lai et al., 2024).

3. Algorithmic Structure, Convergence, and Robustness

The canonical Weiszfeld-type algorithm has the structure:

1. Initialize $x_0$ (avoiding anchors if required).
2. For $k=0,1,2,\ldots$:
   1. Compute the Weiszfeld map (possibly modified).
   2. Project onto the feasible set (if constrained).
   3. Check for termination (based on iterate difference or, for singularities, subgradient step).

3. Output the limit point.

Table: Convergence Properties Across Variants

Variant	Descent Property	Global Conv.	Rate
Unconstrained (Euclidean)	Strict	Yes	Sublinear
Constrained (Hilbert space, convex $C$)	Strict	Yes	Strong
Convex-set D.C. version	Sufficient	Yes (DCA-critical)	Linear (empir.)
Power-$q$ ($1 \leq q<2$); QPWAWS	Strict	Yes	Superlinear (at sing. min)
Sphere (unit vector, robust PCA)	Strict (KL)	Yes (KL, whole seq)	Quadratic type

Strict descent and boundedness of iterates ensure global convergence or convergence to a critical point under the relevant convexity (or KL property) and continuity conditions (Mordukhovich et al., 2013, Nguyen, 2018, Neumayer et al., 2019, Lai et al., 2024).

4. Adaptations for Non-Euclidean, Nonconvex, and Structural Problems

Weiszfeld-type methods have been successfully adapted for robust subspace estimation, rotation averaging, and similar geometric objectives where both the feasible set and the cost landscape are nonlinear or manifold-valued:

• Robust PCA: Principal direction estimation can be approached by minimizing the sum of distances from projected data points onto lines or subspaces. The resulting objective is Lipschitz and subdifferential at anchor sets. Iterations employ Picard updates with a Weiszfeld-like flavor on the sphere and handle non-differentiability using subgradients and one-sided derivative conditions at nondifferentiable anchors (Neumayer et al., 2019).
• Rotation Averaging: For $\mathrm{SO}(3)$ data, the L$_1$-median of geodesic (or chordal) distances is found by tangent-space Weiszfeld updates and retractions, with initialization via entrywise medians and optional in-iteration outlier rejection for robustness. This algorithm matches or exceeds the performance of other L$_1$-rotation averaging methods and displays $O(N)$ complexity per iteration (Lee et al., 2020).

5. Practical Implementation and Numerical Performance

Implementation of Weiszfeld-type algorithms requires careful handling of edge cases:

• Division by zero: If an iterate lands on a data point or set, classical weights become undefined. Remedies include staying at the anchor, taking a convex combination step, or employing subgradient descent/perturbation (Torres, 2012, Lai et al., 2024).
• Projection step: Convex constraint sets (including general Hilbert space constraints) allow projection steps by nearest-point maps, preserving nonexpansiveness (Nguyen, 2018).
• Complexity: Each iteration typically requires $O(mn)$ operations for points (or more for convex sets), with per-iteration cost scaling linearly with the number of anchors and ambient dimension (Mordukhovich et al., 2013, An et al., 2014).
• Numerical evidence: For facility location, robust PCA, or rotation averaging, Weiszfeld-type methods consistently demonstrate strict descent, monotonic convergence, and favorable objective values relative to general-purpose nonlinear solvers. They exhibit rapid convergence in practice—often within tens of iterations for moderate-sized problems—and are robust to nondifferentiability and moderate ill-posedness (e.g., in robust statistics) (Torres, 2012, Neumayer et al., 2019, Lee et al., 2020, Lai et al., 2024).

6. Applications and Research Directions

Weiszfeld-type algorithms are foundational across several fields:

• Location Science: Optimal facility location, constrained service placement, and coverage design (Mordukhovich et al., 2013, Nguyen, 2018).
• Data Science & Robust Statistics: Geometric median computation, robust subspace estimation (robust PCA), $q$-norm location models for machine learning (Neumayer et al., 2019, Lai et al., 2024).
• Computer Vision & Robotics: L$_1$-rotation averaging for multi-view geometry, robust initialization and mean computation for rotations or transformations (Lee et al., 2020).
• Convex & Nonconvex Optimization: D.C. programs involving sums of distances to convex sets or objects, including cases where some weights are negative or some sets are unbounded (An et al., 2014).

Recent works have established full convergence even in singular regimes ($1Lai et al., 2024, Nguyen, 2018).

7. Key Properties and Theoretical Guarantees

• Existence and uniqueness: For noncollinear points (or analogous nondegeneracy for sets), the minimizer is unique under strict convexity.
• Stability: Solution mappings are Lipschitz-continuous with respect to perturbations in anchor locations, with precise subgradient descriptions (Nguyen, 2018).
• Descent and fixed-point characterization: Every nonstationary step produces a strictly lower objective; KKT points and fixed points of the Weiszfeld mapping coincide under suitable regularity (Torres, 2012).
• Robustness: Subgradient or desingularization strategies resolve singularities and ensure escape from non-minimizing bad anchors (Lai et al., 2024).
• Generality: The methodology extends to locations over convex sets, spheres, and manifolds, and to objectives involving distances to sets or subspaces, provided projection and evaluation can be efficiently performed (An et al., 2014, Neumayer et al., 2019).

Weiszfeld-type algorithms thus provide a versatile, theoretically sound approach for a spectrum of location, median, and robust estimation problems, with extensions addressing practical and theoretical challenges in modern optimization and data analysis.