Tight Nonasymptotic Local Convergence of Sinkhorn-Knopp

The paper Tight Nonasymptotic Local Convergence of Sinkhorn-Knopp revisits the Sinkhorn-Knopp (SK) algorithm for the matrix scaling problem. Despite extensive literature on global convergence, the local linear convergence behavior was less understood. The authors address this gap by providing the first nonasymptotic local analysis of SK that matches the rate obtained from existing asymptotic Jacobian-based arguments.

Abstract: "We revisit the Sinkhorn-Knopp (SK) algorithm for the matrix scaling problem. Despite extensive literature on the global convergence of SK and its variants, its local linear convergence behavior remains less understood. We address this gap by providing the first nonasymptotic local analysis of SK that matches the rate obtained from existing asymptotic Jacobian-based arguments. We show that under certain connectivity conditions, SK is a polynomial-time algorithm for doubly stochastic matrix scaling. With the developed tools, we showcase the local suboptimality of SK and provide accelerated variants. Finally, for dense matrices, we improve the complexity of existing first-order matrix scaling algorithms from O(n/ε(2/3)) to O(n(7/4)/√ε)."

Main Contributions:

• "We provide the first tight nonasymptotic local linear convergence analysis of SK, establishing an O(c + 1/σ2 log(1/ε)) iteration complexity. In particular, by explicitly computing the problem-dependent constant c, we show that SK is a polynomial-time algorithm for the doubly-stochastic matrix scaling problem when σ2 is treated as a fixed constant. Our analysis also provides a general template for proving nonasymptotic local convergence rates of alternating minimization algorithms."

• Building on these tools, we derive accelerated variants of SK that achieve improved nonasymptotic global and local complexity bounds for the matrix scaling problem.

Local Convergence Analysis:

The paper establishes that the asymptotic convergence rate is given by 1 − σ2:= λm−1(P −1/2 A⋆ Q−1 (A⋆)⊤ P −1/2), where P:= diag(p), Q:= diag(q). The authors note this Jacobian-based argument only yields an asymptotic guarantee, motivating their nonasymptotic analysis.

The key result (Theorem 3.2) states: Suppose A1 holds and ε ∈ (0, 16∥p∥1]. Then SK finds an ε-approximate scaling in Kε = O(min ∥p∥1D2/σ24, 1/(sσ22) + 1/σ2 log(1/ε)) iterations.

Polynomiality of SK:

Corollary 3.1 states: "Let p = q = 1/n 1n (case 1 in Lemma 3.1) and suppose σ2 > 0 is a fixed constant. Then SK finds an ε-approximate scaling in O(n3 + log(n/ε)) iterations. If A > 0 is positive (case 3 in Lemma 3.1), then the complexity is further improved to O(n + log(n/ε))."

Tightness:

Proposition 3.2 (informal) states: There exists (a family of) matrix scaling instances (A, p, q) such that φ(uk, vk) − φ(u⋆, v⋆) ≥ ε for k = O(c + 1/σ2 log(1/ε)), where c ≥ 0 does not depend on ε. This confirms the analysis matches the asymptotic behavior.

Local Suboptimality of SK:

The paper shows SK is locally suboptimal. For an L-smooth µ-strongly convex problem, the minimax optimal stepsize 2/(L+µ) > 1/L produces a better contraction. The authors show: "Plugging in L = 1 and µ = σ2, we see that the stepsize 2/(σ2+1) provides a local contraction of ((1−σ2)/(1+σ2))2 < (1−σ2)2. If σ2 → 1, this stepsize improves on the contraction of SK by nearly a factor of 4."

Accelerated Variants:

1. Global acceleration (Theorem 4.1): Under the same assumptions as Theorem 3.2, there exists an algorithm A1 that outputs an ε-approximate scaling in Kε = O(√(∥p∥1∥u⋆∥)/√ε) iterations, where each iteration has the same cost as SK. This improves on the O(∥p∥1D(2/3)/ε(2/3)) result from [2] in terms of ε.

2. Stochastic acceleration (Theorem 4.2): Under the same assumptions as Theorem 3.2, there exists an algorithm A2 that outputs û such that E[∥∇ζ(û)∥] ≤ ε with arithmetic complexity Õ(nnz(A) n(1/4) √(∥p∥∞ D)/√ε) in expectation.

3. Local online acceleration (Algorithm 2, OSMS): Theorem 4.3 (informal) states: Algorithm 2 outputs an ε-approximate scaling in O(1/σ2⋆ log(1/ε)) iterations.

4. Local Nesterov's acceleration (Algorithm 3, PAGD): Theorem 4.4 states: Suppose u1 = w1 satisfies the condition from Lemma 4.4, then Algorithm 3 outputs an ε-approximate scaling in O(2/√σ2 log(1/ε)) iterations.

Numerical Experiments:

The experiments compare SK, gradient descent (GD) on ζ, OSMS, and PAGD on entropy-regularized optimal transport instances. The results show: As our theory predicts, using a larger step size yields better local contraction than SK and Both PAGD and OSMS beat SK by several orders of magnitude within the budget.

Conclusion:

"This paper investigates the nonasymptotic local convergence of SK and provides new insights into its behavior. We also develop nonasymptotic accelerated variants through the semi-dual formulation. Our techniques extend to related problems, such as unbalanced optimal transport, and our analysis template applies to general two-block alternating minimization algorithms."

Improvements to AI Systems:

1. Provably Efficient Matrix Scaling for Optimal Transport (OT) Solvers:

Integrate the tight nonasymptotic local convergence bounds (Theorem 3.2) into AI systems that rely on Sinkhorn iterations for entropy-regularized OT (e.g., in generative models, domain adaptation, or attention mechanisms). The improved complexity guarantees (e.g., O(n cubed + (n/epsilon)) for doubly stochastic scaling) allow AI systems to certify convergence within a fixed iteration budget for large-scale OT problems, eliminating heuristic stopping criteria and reducing runtime unpredictability.

2. Accelerated Alternating Minimization for Two-Block Optimization:

Replace vanilla Sinkhorn-Knopp with the accelerated variants (OSMS and PAGD) in AI models that solve two-block convex-concave saddle point problems (e.g., Wasserstein GANs, adversarial training, or federated learning with dual variables). These variants achieve O(1/sqrt sigma 2 (1/epsilon)) local convergence (Theorem 4.4), enabling faster adaptation in online learning settings where the scaling matrix changes over time (e.g., dynamic transport plans in reinforcement learning or continual learning).

3. Adaptive Step-Size Control via Local Curvature Estimation:

Use the paper’s explicit computation of the local contraction rate 1 - sigma 2 (where sigma 2 is the second-largest singular value of a normalized matrix) to design self-tuning AI optimizers. The system can estimate sigma 2 from intermediate iterates and automatically switch to the optimal step size 2/(1+ sigma 2) (which improves contraction by up to 4x) for faster convergence, without manual hyperparameter tuning.

4. Polynomial-Time Guarantees for Sparse and Dense Matrix Scaling:

For AI systems handling sparse matrices (e.g., graph neural networks, sparse attention), the improved dense-matrix complexity O(n 7/4/sqrt epsilon) (vs. O(n/epsilon 2/3)) enables scalable inference on resource-constrained devices by reducing arithmetic operations. The system can now process larger sparse matrices (e.g., adjacency matrices in recommendation systems) with guaranteed error bounds, improving recommendation accuracy and fairness.

5. Stochastic Acceleration for Noisy Environments:

The stochastic acceleration (Theorem 4.2) provides an expected arithmetic complexity of (nnz(A) n 1/4/sqrt epsilon). This allows AI systems to handle noisy or approximate matrix inputs (e.g., from sensor data or partial observations) while maintaining convergence guarantees, making them robust in real-world deployment (e.g., robotics, autonomous navigation).

6. Template for Analyzing Other Alternating Algorithms:

The paper’s analysis template (for proving nonasymptotic local convergence of alternating minimization) can be applied to improve AI systems that use block coordinate descent (e.g., in non-negative matrix factorization, dictionary learning, or sparse coding). By deriving tight local rates for these algorithms, AI systems can predict and control their convergence behavior in deep learning pipelines, leading to more stable training.

We revisit the Sinkhorn-Knopp (SK) algorithm for the matrix scaling problem. Despite extensive literature on the global convergence of SK and its variants, its local linear convergence behavior remains less understood. We address this gap by providing the first nonasymptotic local analysis of SK that matches the rate obtained from existing asymptotic Jacobian-based arguments. We show that under certain connectivity conditions, SK is a polynomial-time algorithm for doubly stochastic matrix scaling. With the developed tools, we showcase the local suboptimality of SK and provide accelerated variants. Finally, for dense matrices, we improve the complexity of existing first-order matrix scaling algorithms from O(epsilon 2/3) to O(sqrt epsilon).