Debiased Sinkhorn barycenters

1. Introduction

• The paper focuses on fixed-support barycenters, which optimize masses on known supports and are faster in practice than free-support methods.Free-support methods jointly optimize support and mass when the barycenter support is unknown.Source: §1 · p. 1
• Entropy regularization makes OT scalable through fast, parallelizable Sinkhorn iterations, but it can blur barycenters and small regularization may cause numerical instability.Log-domain stabilization mitigates instability but considerably slows Sinkhorn iterations.Source: §1 · pp. 1–2
• Earlier debiasing approaches either leave the IBP framework, require gradient descent with repeated Sinkhorn loops, depend on small ε for deconvolution, or return to slower unregularized OT.These alternatives trade sharper solutions against computational efficiency or numerical practicality.Source: §1 · p. 2
• The choice of reference measures in the entropy term determines the induced barycenter bias, including blurring for Lebesgue and shrinking for product measures.The paper identifies the Sinkhorn divergence Sε as debiased relative to these reference-measure choices.Source: §1 · p. 2
• The paper proposes Sinkhorn divergence Sε, defined without specifying reference measures, to obtain debiased barycenters while retaining fast iterative computation.Its computational approach is a Sinkhorn-like fixed-point algorithm similar to IBP and supports differentiation with respect to barycentric weights.Source: §1 · p. 2
• The paper quantifies reference-measure bias for univariate Gaussian barycenters and establishes differentiability and convexity properties needed for its theoretical results.It then develops a fast debiased-barycenter algorithm and evaluates it numerically.Source: §1 · p. 2

Lebesgue reference and smoothing bias

• The weighted barycenter of univariate Gaussians under the Lebesgue reference is Gaussian with increased variance, demonstrating a smoothing bias.The result is illustrated with discrete histograms using a grid of 500 bins.Source: §2 · p. 3

Product reference and shrinking bias

• Entropy-regularized barycenters under the product reference can become degenerate: for sufficiently large ε, equal-variance Gaussian inputs yield a Dirac at their weighted mean.The product-reference barycenter may therefore exhibit a shrinking bias and collapse as regularization increases.Source: §2 · pp. 3–4

Reference measures and bias

• The reference measure determines the entropy bias: Lebesgue references induce blurring, product references induce shrinking, and Sε is debiased.These differences motivate using Sinkhorn divergences, which are defined independently of the reference measures.Source: §2 · p. 4

Sinkhorn divergence

• Sε is defined for arbitrary probability measures and retains convexity and differentiability under extended conditions for sub-Gaussian measures with unbounded supports.The extension assumes a negative-semidefinite cost, including costs of the form C(x,y)=∥x−y∥^d for 0<d≤2.Source: §2 · p. 4

Debiased barycenters

• For Gaussian inputs, the Sε barycenter is Gaussian and, when input variances are equal, has the same variance for every ε.Figure 3 compares the smoothing, shrinking, and debiased barycenters; the debiased result is independent of ε.Source: §2 · p. 4

Computational formulation

• The identity underlying Sε enables a Sinkhorn-like fixed-point algorithm, recovering an IBP-style computational structure for debiased barycenters.The method combines one IBP iteration with an update step and avoids the technical difficulties of directly using the product measure.Source: §2 · pp. 2–4

• Sε is differentiable on sub-Gaussian measures with unbounded supports, with gradients characterized through dual potentials.The proof uses quadratic bounds on potentials and dominated convergence to handle unbounded domains.Source: §3 · pp. 5–6
• Feasible-direction differentiability is sufficient to characterize convexity and first-order optimality over probability measures.The derivative notion avoids requiring ordinary Fréchet differentiability because the probability simplex has empty interior in the signed-measure space.Source: §3 · p. 6
• Sε is convex with respect to either argument on sub-Gaussian measures under the stated kernel assumption.The convexity proof is based on the dual formulation and the positivity of the kernel K.Source: §3 · pp. 5–6
• These properties yield global barycenter optimality conditions for both OT⊗ε and Sε.For Sε, the condition incorporates the autocorrelation potential hα⋆ into the gradient.Source: §3 · p. 6

• The debiased barycenter problem is formulated on a finite support using transport plans, constraint sets, and a convex objective.Minimizing over transport plans recovers an IBP problem with a modified kernel, while minimizing over d yields the fixed-point equation.Source: §4 · p. 7
• The method leverages the identity that makes Sε independent of the reference measures, allowing computation through the uniform-measure IBP machinery.Kernel operations can be used within the finite-grid formulation to compute the debiased barycenter.Source: §4 · p. 7
• For Gaussian barycenters, Algorithm 1 converges to the true barycenters and is reported as fast as IBP with a linear convergence rate.Figure 4 compares convergence for univariate Gaussian inputs N(−0.5, 0.1) and N(0.5, 0.1).Source: Figure 4 · p. 7, §4 · p. 7
• Algorithm 1 combines one IBP iteration with a debiasing update and has O(Kn^2) complexity per iteration, matching IBP's theoretical per-iteration complexity.The algorithm alternates transport scaling with an update for the debiasing variable d.Source: §4 · p. 7
• A convergence proof for Algorithm 1 is left for future work, although IBP convergence follows from alternating Bregman projections.The authors report that similar techniques were not successful for Algorithm 1.Source: §4 · p. 7

• The experiments compare the proposed method with sharp barycenters, free-support debiased barycenters, and unregularized Wasserstein barycenters.These baselines use accelerated gradient descent, Frank–Wolfe updates, or an accelerated interior-point method, respectively.Source: §5 · p. 8
• On nested ellipses, the proposed barycenters reduce blurring and remain computationally efficient relative to uniform-reference entropic barycenters.Figure 6 reports running times beneath the images; all six barycenters were computed on an Intel Core i5 3.1 GHz laptop.Source: Figure 6 · p. 8, §5 · p. 8
• Debiased Sε produces sharper 3D shape interpolations than OTUε while taking about the same time, approximately 5 seconds on a GPU.The shapes are represented on a 200 × 200 × 200 grid, and Gaussian convolutions accelerate kernel operations.Source: Figure 8 · p. 8, §5 · p. 8, Figure 7 · p. 8
• In MNIST barycentric embeddings, debiasing improves cross-validation accuracy and makes performance less sensitive to ε.The evaluation uses 500 MNIST images, a 50-image dictionary, 450 embedded test images, and 20 randomized dictionary selections.Source: §5 · p. 9, Figure 9 · p. 9

Conclusion

• Entropy bias can be a blur or a shrink depending on the reference measure, whereas Sinkhorn divergences yield debiased Gaussian barycenters.The paper also extends convexity and differentiability results to measures with non-compact supports.Source: Conclusion · p. 9

A. Convexity and Optimality condition

• Feasible-direction differentiability characterizes convexity through the first-order inequality for every pair of probability measures.The appendix proves both directions: the inequality implies convexity, and convexity implies the inequality.Source: A. Convexity and Optimality condition · p. 11
• For differentiable convex functions on probability measures, α⋆ minimizes F exactly when ⟨∇F(α⋆), α − α⋆⟩≥0.This condition is established using derivatives along feasible directions in P(Rd).Source: A. Convexity and Optimality condition · p. 11

Proof of Lemma 1

• The proof of Lemma 1 bounds perturbation ratios using suboptimal potentials and passes to the limit through pointwise convergence and dominated convergence.Source: Proof of Lemma 1 · p. 12

B.1. Differentiability and convexity of OTL

• Entropic OT is convex in each argument because its dual formulation is a supremum of linear functions, while the KL term is jointly convex.Source: §B.1 · p. 13
• Under positive-density, absolutely continuous sub-Gaussian assumptions, entropic OT is differentiable and its gradient is expressed through Sinkhorn dual potentials.Source: §B.1 · p. 13

C.1. Fixed point equations Lemmas

• For Gaussian barycenter variance equations, Lemma C.1 guarantees a positive solution, while Lemma C.2 gives existence exactly when ε′2 < ¯σ.Source: §C.1 · p. 14
• Under Theorem 3 assumptions, the relevant variance equation has a positive solution S⋆ lying between σ(0) and σ(K).Source: §C.1 · p. 14

C.2. Proofs of theorems 2 and 3

• The debiased barycenter is constructed by solving convex optimality systems with quadratic potentials and Gaussian measures parameterized by mean m and variance S.The coefficients are determined from the Sinkhorn equations and sufficient moment-based optimality conditions.Source: §C.2 · pp. 15, 19

C.2. Proofs of the theorems 2 and 3

• The variance solution lies between the smallest and largest input standard deviations, S⋆∈(σ(0), σ(K)).Source: §C.2 · p. 15

C.2. Proofs of theorems 2 and 3

• For the entropic Gaussian barycenter, the variance S is the unique non-zero solution of the fixed-point equation, and equal input variances preserve that common variance.Source: §C.2 · pp. 15, 19
• A unique positive variance solution exists for the debiased fixed-point equation when ε′2 < ¯σ, as guaranteed by Lemma C.2.Source: §C.2 · pp. 19–20
• When ε′2 ≥ ¯σ, the debiased barycenter construction instead verifies optimality of the Dirac measure at the weighted mean ¯µ.Source: §C.2 · pp. 20–21

C.3. Proof of theorem 1

• The entropic OT barycenter is characterized through first-order optimality conditions because entropic OT is convex and differentiable on positive-density sub-Gaussian measures.Source: §C.3 · p. 21
• The resulting theorem identifies the entropic barycenter as Gaussian with variance determined by that unique fixed-point solution.Source: §C.3 · p. 21
• Quadratic dual potentials and a Gaussian candidate reduce the proof to determining the candidate’s mean and variance from the Sinkhorn and optimality equations.Source: §C.3 · pp. 21–22
• Lemma C.1 guarantees that the fixed-point equation governing the barycenter variance has a unique positive solution S.Source: §C.3 · p. 24

D. The IBP algorithm

• The IBP algorithm computes the OT barycenter through alternating minimization over constraint sets, yielding Sinkhorn-like scaling iterations for transport-plan factors.Source: D. The IBP algorithm · p. 24

• The proposition proof introduces transport-plan constraint sets fixing the input marginals and shared barycenter marginal across plans.Source: Proof of proposition 5 · p. 25
• The debiased autocorrelation term is rewritten using its dual problem, enabling comparison with the IBP formulation.Source: Proof of proposition 5 · p. 25

F.1. Barycenters of nested ellipses

• Nested-ellipse barycenters use ten 60×60 simulated images, squared Euclidean ground cost on [0, 1]2, and the lowest numerically stable entropy parameter.The study compares uniform-reference OT, the proposed debiased divergence, iterative OT⊗ ε, Sharp barycenters, free-support barycenters, and unregularized Wasserstein barycenters.Source: §F.1 · p. 26, Figure 10 · p. 27
• The nested-ellipse comparison uses IBP for several methods, accelerated gradient descent for Sharp barycenters, Frank-Wolfe updates for free-support barycenters, and accelerated interior-point methods for W.Sharp barycenters required manual learning-rate tuning and were more prone to numerical instabilities; the W implementation used the authors’ online MATLAB code.Source: §F.1 · pp. 26–27

F.2. Barycenters of 3D shapes

• 3D-shape barycenters use smoothed, centered, rescaled tore and rabbit meshes sampled as normalized 200^3 histograms with 10^-10 added weights and ε = 0.01.Weighted barycenters are computed for w ∈ [0, 0.25, 0.5, 0.75, 1] using IBP and the proposed debiased Sinkhorn algorithm.Source: §F.2 · p. 27, Figure 11 · p. 27

F.3. OT barycentric embeddings

• The MNIST embedding experiment filters labels 0–4, uses 500 images, selects K = 50 dictionary samples, and embeds the remaining 450 by barycentric weight fitting.A 100-estimator Random Forest with scikit-learn defaults is trained on these embeddings.Source: §F.3 · p. 27