A Provably Convergent and Practical Algorithm for Gromov–Wasserstein Optimal Transport¶
会议: NeurIPS 2026
arXiv: 2605.04175
OpenReview: https://openreview.net/forum?id=umKsf7P52m
领域: 优化/理论
关键词: 格罗莫夫-瓦瑟斯坦最优传输, 投影梯度法, 不精确投影, 原始可行性残差, 收敛性分析
一句话总结¶
针对非凸格罗莫夫-瓦瑟斯坦最优传输中精确投影代价高昂导致理论与工程实现脱节的问题，本文提出了一种基于易验证原始可行性残差的不精确投影梯度（iPG）框架，不仅在理论上严格证明了全序列收敛到驻点，而且在图对齐任务上兼顾了高稀疏度、近可行性与优异的目标函数解质量。
研究背景与动机¶
格罗莫夫-瓦瑟斯坦最优传输（GWOT）通过比对度量测度空间内的相对几何结构实现跨空间对齐，因而不依赖于跨域显式距离，在非同构空间匹配、跨模态单细胞对齐和异构图对齐中成为核心基础工具。然而，GWOT的目标函数包含双线性耦合项，整体具有非凸二次规划特性，其约束集合为具有非负与双边边际分布约束的输运多面体。在无正则化的纯非凸优化中，欧氏投影梯度法（PG）因单步形式直观且理论上能保持精确最优传输解固有的极高稀疏性，原本是最理想的第一阶基准方案。
但大规模场景下将梯度步投影回输运多面体的精确欧氏投影本身是一个带约束强凸二次规划，计算开销极其沉重；工程实践中往往仅能调用内部迭代求解器（如加速梯度法或牛顿类方法）进行截断逼近。现有的一阶投影梯度理论大多强假设子问题能求解至精确机器精度，缺乏面向不精确子求解的严格收敛机制，从而迫使研究界转向熵正则化平滑（如 EPGD、PPA）或无内部投影的双线性松弛单环算法（如 BAPG）。然而，熵平滑不可避免地会引入平滑偏差，破坏传输方案的稀疏结构，导致匹配解完全稠密且对正则化参数极度敏感；而完全依赖线搜索的条件梯度法（CG）则收敛极其迟缓且极易陷入较劣的局部极值。
面对实践逼近与严格理论之间的鸿沟，关键突破在于能否在不依赖未知精确投影点的前提下，给出一种计算代价极低且可验证的不精确停机准则。本文的核心 idea 是通过构建基于对偶迭代点归一化的可计算原始可行性残差准则，提出一种不精确投影梯度算法（iPG），仅凭单步内部求解的可行性残差即能统一子步误差与下降充分性，建立非凸 GWOT 序列的全局驻点收敛性与全序列收敛性。
方法详解¶
整体框架¶
iPG 针对二次损失 GWOT 展开：给定两度量空间的内部相异度矩阵 \(C_{\mathcal{X}}\in\mathbb{R}^{n\times n}\) 与 \(C_{\mathcal{Y}}\in\mathbb{R}^{m\times m}\) 及离散测度分布 \(a\in\Delta_n, b\in\Delta_m\)，将目标转化为在输运多面体 \(\mathcal{U}(a,b)=\{\Pi\ge 0: \Pi 1_m=a, \Pi^\top 1_n=b\}\) 上极小化平滑非凸目标函数 \(f(\Pi):=-\langle C_{\mathcal{X}}\Pi C_{\mathcal{Y}},\Pi\rangle\)。算法采用外层投影梯度迭代与内层对偶凸平滑优化嵌套的双层机制：外层计算当前耦合矩阵处的平滑梯度 \(\nabla f(\Pi^k)\) 并构造欧氏近邻目标；内层通过加速一阶算法在对偶空间求解无约束凸函数，一旦对偶变量导出的非负矩阵满足基于自身范数加权的原始边际残差容差条件，即截断内层并更新传输矩阵。
%%{init: {'flowchart': {'rankSpacing': 24, 'nodeSpacing': 28, 'padding': 6, 'wrappingWidth': 400}}}%%
flowchart TD
A["输入：相异度矩阵 C_X, C_Y 与测度 a, b"] --> B["外层梯度步：非凸GW梯度计算<br/>∇f(Π^k) = -2 C_X Π^k C_Y"]
B --> C["对偶子问题重构：无约束平滑极小化<br/>min g_k(y) 建立梯度与仿射约束对偶映射"]
C --> D["内层可验证停机判定：归一化残差控制<br/>||A Π^(k+1) - r|| ≤ ε_k / (1 + ||y^(k+1)||)"]
D --> E["双重影子序列桥接：Hoffman误差界与KL收敛<br/>控制与精确投影点偏离，完成全序列收敛"]
E --> F["输出：近乎完全可行的高稀疏匹配方案 Π^k"]关键设计¶
1. 对偶子问题重构：将输运投影转化为平滑无约束凸优化
欧氏投影子问题要求将矩阵 \(Z^k = \Pi^k - \frac{1}{\gamma}\nabla f(\Pi^k)\) 投影到输运多面体 \(\mathcal{U}(a,b)\)，其标准形式为在仿射约束 \(\mathcal{A}\Pi = r\)（其中 \(\mathcal{A}\Pi = [\Pi 1_m; \Pi^\top 1_n]\)，\(r=[a; b]\)）与非负锥 \(\Pi\ge 0\) 下极小化 \(\frac{1}{2}\|\Pi - Z^k\|_F^2\)。若直接处理原问题约束，内外层边界难以协调。iPG 将等式约束引入拉格朗日乘子向量 \(y\in\mathbb{R}^{n+m}\)，利用非负锥投影具有简单显式截断 \(\mathrm{Proj}_{\mathbb{R}_+^{n\times m}}(M) = \max(M, 0)\) 的特性，消去原变量得到其对偶无约束凸优化目标： $\(g_k(y) := \frac{1}{2}\left\|\mathrm{Proj}_{\mathbb{R}_+^{n\times m}}\left(\mathcal{A}^* y + \Pi^k - \frac{1}{\gamma}\nabla f(\Pi^k)\right)\right\|_F^2 - \langle r, y\rangle\)$ 该对偶函数 \(g_k(y)\) 处处连续可微，且其梯度恰好为原变量在截断点处的仿射约束违差：\(\nabla g_k(y) = \mathcal{A}\,\mathrm{Proj}_{\mathbb{R}_+^{n\times m}}(\mathcal{A}^* y + Z^k) - r\)。这使得内层子问题完全变为无约束平滑凸优化，可直接调用 Nesterov 加速梯度法或半光滑牛顿法极速求解，避免了复杂的交替方向分裂或高维锥求解。
2. 内层可验证停机判定：对偶范数加权的原始可行性残差控制
在传统不精确一阶算法理论中，停止准则常依赖于当前迭代点与未知精确解的距离 \(\|\Pi^{k+1} - \Pi^{k,*}\|_F\)，在工程中无法直接计算。iPG 创新地提出了只需计算当前对偶输出 \(y^{k+1}\) 及其诱导的非负矩阵 \(\Pi^{k+1}\) 的显式停机条件： $\(\|\mathcal{A}\Pi^{k+1} - r\|_2 \le \frac{\varepsilon_k}{1 + \|y^{k+1}\|_2}\)$ 其中 \(\{\varepsilon_k\}\) 为预先设定的正可和衰减序列。该设计的精妙之处在于分母中的归一化因子 \((1+\|y^{k+1}\|_2)\)：在无严格对偶有界性假设的通用设置下，对偶乘子可能随外层演化而发生范数漂移。由于一阶对偶间隙项涉及对偶乘子与原始可行性违差的内积 \(\langle y^{k+1}, \mathcal{A}\widetilde{\Pi}^{k+1} - \mathcal{A}\Pi\rangle\)，除以该对偶范数项能严格抹平潜在的乘子放大效应，将次微分摄动 \(\widehat{\varepsilon}_k\) 稳固地绑定在常数倍的 \(\varepsilon_k\) 尺度上，构造出单调可达且易于在内层循环每步监控的判定屏障。
3. 双重影子序列桥接：基于 Hoffman 误差界与 KL 几何的严格收敛保障
为了跨越不精确非可行迭代点与非凸收敛性之间的理论鸿沟，iPG 在证明中引入了双重“影子序列（Shadow Sequences）”。首先构造显式满足边缘约束的紧邻可行序列 \(\{\widetilde{\Pi}^k\}\subset \mathcal{U}(a,b)\)，其偏离度满足 \(\|\widetilde{\Pi}^k - \Pi^k\|_F \le 2\sqrt{m+n}\|\mathcal{A}\Pi^k - r\|_2\)，利用外层函数 Lipschitz 光滑性质与弱次微分性质证明了近似充分下降不等式，由此得出子序列必然收敛到 GWOT 原问题的驻点。
进一步，为了证明全序列收敛而非仅子序列收敛，算法引入第二影子序列——精确投影点序列 \(\Pi^{k,*} = \mathrm{Proj}_{\mathcal{U}(a,b)}(Z^k)\)。由于多面体互补约束系统在仿射变换下属于有限个多面体系统的并集，利用一致 Hoffman 误差界常数 \(\omega>0\)，证明了近似解误差由原始残差严格线性控制：\(\|\Pi^{k+1} - \Pi^{k,*}\|_F \le \omega\|\mathcal{A}\Pi^{k+1} - r\|_2 \le \omega\varepsilon_k\)。随后引入增广半代数功函数 \(\Psi_\tau(\Pi, t) := f(\Pi) + \delta_{\mathcal{U}(a,b)}(\Pi) + \frac{t^\tau}{\tau}\)，当容差序列满足 \(\sum_{k=0}^\infty \varepsilon_k < \infty\) 及 \(\sum_{k=0}^\infty (\sum_{i=k}^\infty \varepsilon_i^2)^{1-1/\tau} < \infty\)（如选取 \(\varepsilon_k = (1+k)^{-\alpha}, \alpha > 1\)）时，借助 Kurdyka-Łojasiewicz（KL）不等式证明了 \(\{\Pi^{k,*}\}\) 为柯西序列，最终确立了全序列向局部最优驻点的确定性收敛。
实验关键数据¶
主实验¶
在合成图对齐任务上进行评测：两图节点数 \(n\in\{100, 200, 300, 400, 500\}\)，基于 Erdős–Rényi 随机图（连边概率 \(p_{\mathrm{edge}}=0.2\)）添加置换与随机边翻转扰动（\(\eta=0.1\)），计算全源最短路径相异度矩阵。比较方法包括条件梯度法（CG）、熵正则化投影梯度（EPGD）、熵近端点算法（PPA）、Bregman交替投影梯度（BAPG）以及本文方法（iPG）。各方法均采用统一均匀边际分布与均匀初值，评测指标包含 GW 损失（越低越好）、稀疏度（0元素比例）、可行性残差 \(\|\Pi 1 - p\|_2 + \|\Pi^\top 1 - q\|_2\)、Hungarian对齐匹配精度（Accuracy）以及耗时（秒）。
| 问题规模 \(n\) | 算法配置 | GW Loss \(\downarrow\) | 稀疏度 (Sparsity) \(\uparrow\) | 可行性残差 (Feasibility) \(\downarrow\) | 对齐精度 (Accuracy) \(\uparrow\) | 耗时 (s) \(\downarrow\) |
|---|---|---|---|---|---|---|
| \(n=100\) | CG | 3.60e-02 | 0.99 | 0.00e+00 | 0.02 | 0.00 |
| \(n=100\) | EPGD (\(\epsilon=10^{-1}\)) | 5.36e-02 | 0.00 | 5.83e-17 | 0.04 | 0.00 |
| \(n=100\) | PPA (\(\epsilon=10^{-2}\)) | 3.08e-02 | 0.43 | 1.04e-05 | 0.25 | 2.61 |
| \(n=100\) | BAPG (\(\epsilon=10^{-1}\)) | 3.10e-02 | 0.48 | 2.27e-02 | 0.12 | 0.50 |
| \(n=100\) | iPG (Ours) | 3.21e-02 | 0.99 | 7.63e-09 | 0.17 | 0.92 |
| \(n=300\) | CG | 8.34e-02 | 1.00 | 0.00e+00 | 0.01 | 0.05 |
| \(n=300\) | EPGD (\(\epsilon=10^{-3}\)) | 6.07e-02 | 0.00 | 1.85e-05 | 0.51 | 6.18 |
| \(n=300\) | PPA (\(\epsilon=10^{-2}\)) | 6.10e-02 | 0.49 | 6.25e-06 | 0.52 | 13.11 |
| \(n=300\) | BAPG (\(\epsilon=10^{-1}\)) | 6.03e-02 | 0.29 | 6.15e-03 | 0.51 | 3.68 |
| \(n=300\) | iPG (Ours) | 5.85e-02 | 1.00 | 1.15e-09 | 0.57 | 15.66 |
| \(n=500\) | CG | 7.52e-02 | 1.00 | 0.00e+00 | 0.00 | 0.21 |
| \(n=500\) | EPGD (\(\epsilon=10^{-3}\)) | 5.66e-02 | 0.00 | 1.59e-05 | 0.36 | 2.24 |
| \(n=500\) | PPA (\(\epsilon=10^{-2}\)) | 4.76e-02 | 0.45 | 3.78e-06 | 0.57 | 30.32 |
| \(n=500\) | BAPG (\(\epsilon=10^{-1}\)) | 4.71e-02 | 0.20 | 3.62e-03 | 0.57 | 12.33 |
| \(n=500\) | iPG (Ours) | 3.25e-02 | 1.00 | 9.26e-10 | 0.85 | 38.41 |
消融与敏感性分析¶
在 \(n=500\) 规模下对比熵正则化系数 \(\epsilon \in \{10^{-3}, 10^{-2}, 10^{-1}, 10^0\}\) 对代表性基线模型的影响，并与无平滑参数依赖的 iPG 进行对比：
| 求解器 | 正则化参数 \(\epsilon\) | GW Loss \(\downarrow\) | 稀疏度 (Sparsity) \(\uparrow\) | 可行性残差 \(\downarrow\) | 对齐精度 (Accuracy) \(\uparrow\) | 运行时间 (s) \(\downarrow\) |
|---|---|---|---|---|---|---|
| iPG (Ours) | 无正则化参数 | 3.25e-02 | 1.00 | 9.26e-10 | 0.85 | 38.41 |
| PPA | \(10^{-3}\) | 4.52e-02 | 0.76 | 3.68e-05 | 0.62 | 10.46 |
| PPA | \(10^{-2}\) | 4.76e-02 | 0.45 | 3.78e-06 | 0.57 | 30.32 |
| PPA | \(10^{-1}\) | 5.00e-02 | 0.01 | 6.57e-07 | 0.52 | 71.51 |
| PPA | \(10^{0}\) | 5.31e-02 | 0.00 | 8.47e-08 | 0.48 | 100.31 |
| EPGD | \(10^{-3}\) | 5.66e-02 | 0.00 | 1.59e-05 | 0.36 | 2.24 |
| EPGD | \(10^{-2}\) | 9.19e-02 | 0.00 | 4.04e-17 | 0.01 | 0.05 |
| EPGD | \(10^{-1}\) | 9.19e-02 | 0.00 | 4.22e-17 | 0.01 | 0.05 |
| BAPG | \(10^{-3}\) | 异常崩溃 | - | - | - | - |
| BAPG | \(10^{-2}\) | 9.09e-01 | 0.99 | 5.86e-01 | 0.00 | 7.45 |
| BAPG | \(10^{-1}\) | 4.71e-02 | 0.20 | 3.62e-03 | 0.57 | 12.33 |
| BAPG | \(10^{0}\) | 5.01e-02 | 0.00 | 4.38e-04 | 0.52 | 18.35 |
关键发现¶
- 结构稀疏度决定下游离散化质量：经典 CG 虽然输出完全稀疏且可行，但陷入了较差的驻点，导致 \(n=500\) 时精度仅为 \(0.00\)；熵基线（EPGD 与 PPA）当 \(\epsilon\) 增大时耦合方案完全稠密化（稀疏度降为 \(0.00\)），将其转为硬对齐时精度出现明显退化。iPG 保持了 \(1.00\) 的极高稀疏性，且在 \(n=500\) 获得了高达 \(0.85\) 的对齐精度。
- 算法对超参数的稳健性优势显著：BAPG 在较小正则化 \(\epsilon=10^{-3}\) 下数值溢出直接崩溃，在 \(\epsilon=10^{-2}\) 时可行性残差高达 \(0.586\) 彻底失效；PPA 随着 \(\epsilon\) 增加计算耗时呈指数攀升（从 10s 激增至 100s）。iPG 采用多项式阶容差衰减 \(\varepsilon_k=(1+k)^{-3}\)，无需试凑熵参数，表现出极强的数值鲁棒性与高效性。
- 目标函数解质量与残差收敛同步占优：在各节点规模下，iPG 的最终 GW 损失均显著低于其他基线（\(n=500\) 时 iPG 为 3.25e-02，明显优于 PPA 的 4.76e-02 及 CG 的 7.52e-02），且可行性残差稳定维持在 \(10^{-9}\) 极低量级。
亮点与洞察¶
- 真正可监控的可行性停机准则：巧妙利用对偶问题梯度即为原始约束残差的解析性质，将不可知的空间投影距离转化为内层随时可算的对偶加权残差，彻底消除了算法实现与理论证明之间的鸿沟。
- 结合 Hoffman 误差界的非凸 KL 收敛分析范式：非凸复合优化通常难以在不精确投影下建立全序列收敛性。该文证明多面体投影的误差界常数 \(\omega\) 独立于外层位移，从而成功把不精确残差嵌入半代数势函数与 KL 框架，为非凸有约束优化提供了极具启发性的证明范式。
- 非熵欧氏投影在组合匹配上的本质优势：证明了“保留原始欧氏非凸模型 + 快速不精确投影”在离散几何匹配任务上不仅可行，而且比强制平滑的 Sinkhorn 类方法能更好地保留局部几何刚性。
局限与展望¶
- 作者承认的局限：目前的实证研究集中在合成图结构对齐上，尚未拓展到大规模跨模态生物单细胞或海量真实世界网络；同时算法目前主要依赖单核串行加速梯度，尚未针对现代 GPU/分布式异构架构进行显存与并行性定制。
- 潜在的研究局限：外层步长需满足反学习率 \(\gamma > L_f = 2\|C_{\mathcal{X}}\|_2\|C_{\mathcal{Y}}\|_2\)，对于超大规模相异度矩阵，严格计算谱范数上界可能偏向保守，带来较多的外层迭代次数；缺乏自适应回溯线搜索机制。
- 后续改进方向：可进一步探索结合方差缩减（VR）的随机或小批量（mini-batch）GWOT 扩展版本，并在内层引入稀疏牛顿/预条件共轭梯度以应对百万级节点对齐。
相关工作与启发¶
- vs 条件梯度法 (CG / Frank-Wolfe)：CG 每轮求解线性规划子问题，虽然迭代轻量且保证可行，但收敛速率受限为 \(\mathcal{O}(1/k)\) 且极易陷入较差的非凸驻点。iPG 采用带二次邻近项的投影梯度，配合加速一阶内层求解，能跳出劣质局部解并显著降低 GW 损失。
- vs 熵投影梯度与近端点法 (EPGD / PPA)：EPGD 与 PPA 引入熵正则化后借助 Sinkhorn 算法实现快速投影，但平滑项导致传输矩阵几乎完全失去稀疏性，且对超参敏感。iPG 直接作用于无正则化原问题，既保持高稀疏结构又无需调参。
- vs Bregman 交替投影梯度法 (BAPG)：BAPG 试图构建单环无内层投影的一阶格式，但依赖于松弛的 GW 模型，且在实际数值实现中缺乏对微小正则化的容忍度，极易遭遇数值发散。iPG 采用有原则的双层嵌套与严格残差控制，数值稳定性远胜 BAPG。
评分¶
- 新颖性: ⭐⭐⭐⭐ [提出了基于对偶范数加权残差的无未知量停机准则，设计了新颖的非凸GWOT不精确分析范式]
- 实验充分度: ⭐⭐⭐⭐ [基线覆盖全面，细致评测了目标值、稀疏度、残差、下游匹配精度和运行时间，统计量详实]
- 写作质量: ⭐⭐⭐⭐⭐ [数学推导清晰严密，从 Lemma A.1、Hoffman 误差界到 KL 收敛性的链条极其完整自洽]
- 价值: ⭐⭐⭐⭐ [为大规模非凸最优传输提供了兼顾数学严格性与计算实用性的非熵第一阶解决路径]