Semidefinite Relaxations of the Gromov-Wasserstein Distance

Detailed discussion on related work

Could we check if this comment was intended as a weakness? It does not sound like one. We assume this comment was misplaced. Could the Reviewer clarify?

High run-time

We agree with the Reviewer’s concerns about run-time. We address some of these concerns in the global response. We have been working on developing faster algorithms for solving the GW-SDP. Unfortunately, this is a difficult task. There are many promising directions to develop more scalable algorithms for solving the GW-SDP problem. Some of these were discussed in the manuscript. We discuss first-order methods in the response to another Reviewer. These directions are promising but the scope of the work goes comfortably the scope of a single paper.

As a note, for the closely related semidefinite relaxation of the QAP, state-of-the-art methods are only able to solve problem instances where the dimension is around 40 [OWX:18]. The largest dimension in our experiments is 32, which is not far off the state of the art.

References

[OWX:18] D. E. Oliveira, H. Wolkowicz, & Y. Xu, (2018). ADMM for the SDP relaxation of the QAP. Mathematical Programming Computation, 10(4), 631-658.

When we wrote that "... our work is the only work we are aware of that addresses global optimality ... " and "The first contribution is that we compute globally optimal solutions to the GW problem. We can prove that our solutions are globally optimal.", it was to be understood in the context that global optimality was attained if the computed ratio was one.

When writing the rebuttal, we believed that all Reviewers understood the conditions under which global optimality held. Hence in writing a rebuttal, when summarizing our contributions, we did not think it was necessary to repeat the conditions under which global optimality held.

Unfortunately, this summary has given the impression that we claimed global optimality in all instances. We do not claim such a thing. We hope it was clear from the original paper we do not claim global optimality in all instances. We apologize for the confusion. This was not our intention. We hope the paper was very clear about its claims and limitations.

Given just a solution to GW problem, the paper's approximation ratio cannot be used to check whether the given solution is optimal or not.

To check if a solution is optimal, one typically needs to compute the dual problem (of some suitable form). Computing the dual is (sort of) unavoidable if you want a proof of optimality. The GW-SDP does this (in a way) because it is a convex program.

Also, given the solution to the GW-SDP formulation, one can certify the optimality, or the gap to optimality, of any proposed solution to the GW problem. Compute the approximation ratio. If it is equal to one, the proposed solution is optimal. If it is close to one, we know it is quite good.

In the (limited) additional experiments during the rebuttal phase, the paper has also shown empirically that in the experiment that were performed, their approach reach optimal solution when m is a multiple of n. However, no theoretical justification is provided.

We do not provide theoretical justification. In fact, we should not expect the SDP to provide exact solutions to all instances because the GW problem is not known to be tractable to solve. There must be some inputs for which the SDP is not exact. The surprising aspect here is that the SDP is exact for many inputs, and that is the point our experiments make.

(As to providing theoretical justification, such as proving such a phenomenon holds for random instances, is a difficult research problem. We don't have an answer for this at the moment.)

The Reviewer posed some further questions about outliers. We need a bit more time to respond to this and will provide a reply later.

Limited novelty

We emphasize, the paper has multiple contributions. The first contribution is that we compute globally optimal solutions to the GW problem. We can prove that our solutions are globally optimal. There is no existing work (as far as we are aware of), published or unpublished, that achieves this. This contribution alone, we believe, is a significant contribution.

The second contribution is to point out that many existing methods do not compute globally optimal solutions. This point, we believe, is a very concerning development, and researchers working on GW need to be aware.

The Reviewer is concerned with the limited novelty because the relaxation is a variant of an existing relaxation for the QAP. The point we make here is that this particular variant we propose is the simplest variant that answers the previous two questions.

We hope the Reviewer appreciates the conceptual aspects of our paper, namely that global optimality and certificate of optimality, are fundamental aspects of optimization.

Sparsity

When $m=n$, the feasible set of the GW problem is the set of doubly stochastic matrices. The extreme points of the set of doubly stochastic matrices is the collection of permutation matrices – this is known as the Birkhoff von-Neumann theorem. When $m \neq n$, the feasible set is known as the transportation polytope. The extreme points correspond to bipartite graphs (the set of vertices are $\alpha$ and $\beta$) that have no cycles. These matrices are also sparse.

The optimal solution of a convex program whose feasible region is the set of doubly stochastic matrices or more generally the transportation polytope tends to be sparse. This is because optimal solutions tend to lie on low-dimensional faces, which are sparse (this is a known phenomenon in convex geometry).

In short, the optimal solution to the GW problem should generally be sparse. The solutions to GW-SDP are sparse, and this is a good sign. The solutions obtained by Conjugate Gradient tend to be sparse, which is also a good sign. On the other hand, if the solution is non-zero everywhere, we can be certain the solution is sub-optimal. This is the case for e-GW solvers, and the reason is that e-GW adds a bias to the objective; that is, they solve a slightly different problem.

Heuristic solver

We will place the heuristic solver in the appendix, as recommended.

Unless this issue is resolved, the overall contribution of the paper must be considered limited.

We hope the Reviewer will take into consideration that developing scalable algorithms for solving large-scale SDPs is an active research area and an open challenge. In fact, for the closely related semidefinite relaxation of the QAP, state-of-the-art methods are only able to solve problem instances where the dimension is around 40 [OWX:18]. The largest dimension in our experiments is 32, which is not far off the state of the art.

There are many promising directions to develop more scalable algorithms for solving the GW-SDP problem. Some of these were discussed in the manuscript. We discuss first-order methods in the response to another Reviewer. These directions are promising but the scope of the work goes comfortably the scope of a single paper. We urge the Reviewer to moderate expectations as to what is technically feasible within the scope of a single paper.

Again, we reiterate an earlier point, that there is no other work we are aware of that addresses the global optimality of the GW problem. Addressing global optimality is one of the key contributions of this paper.

References

[OWX:18] D. E. Oliveira, H. Wolkowicz, & Y. Xu, (2018). ADMM for the SDP relaxation of the QAP. Mathematical Programming Computation, 10(4), 631-658.

We thank the reviewer for their quick response to our rebuttals.

Originality: The statement in the rebuttal, "The first contribution is that we compute globally optimal solutions to the GW problem," is inaccurate. To be precise, a solution is guaranteed to be globally optimal only when the right-hand side of equation (4) equals 1.

In our paper, we do make it clear that we only claim global optimality when the RHS of equation (4) equals one.

Sparsity: My concerns regarding sparsity have not been fully addressed. While it is correct that an optimal solution should be sparse, the argument that a sparse solution is close to the optimal one is not valid. Therefore, I find little merit in discussing the sparsity of the solutions.

The Reviewer asked why do solutions tend to be sparse. The initial comments do not appear to be a concern and we answered the question as it is. Our response simply states that solutions that are not sparse are certainly not optimal. Our proof of global optimality comes from the RHS of equation (4), not sparsity.

Real data

We use a publicly available database of triangular meshes (Sumner et al. 2004). We obtain 18 points and compute distance matrices using Dijkstra's algorithm. Each object's probability measure is chosen to be uniform. We apply (GW-SDP) to the corresponding metric-measure spaces to determine the correspondence between the selected vertices across different objects. Two representative examples are given in Figure 1 in the attached pdf file.

Following are the results when we do matching of distances matrices of different objects. In general, we expect shapes of the same animals to have a smaller GW distance than shapes of different animals, which is indeed the case for the three GW formulations. We still notice that GW-SDP consistently returns the smallest value when performing the same matching task.

Frank Wolfe with more samples

This is our interpretation of the Reviewer's comment. Fix both sets of distributions. Increase the number of samples for Frank-Wolfe but fix the number of samples for GW-SDP. What happens?

We conducted this experiment and we noticed that the objective value for Frank-Wolfe (GW-CG in table) decreases as we increase the number of samples and approach the objective value obtained by the SDP relaxation, which is obtained using few samples.

n GW-SDP GW-SDP Runtime (s) GW-CG GW-CG runtime (s) 10 0.4577 6.3753 1.135940 0.000389 100 0.629425 0.007571 1000 0.540984 2.520011 10000 0.496796 138.358954

We noticed that the objective value from Frank-Wolfe is greater than the objective value from GW-SDP. We do not know if this is because the solution from Frank Wolfe is sub-optimal, or if GW-SDP returns a smaller objective than the GW distance between the continuous distributions because of discretization error (from finite samples). This is an interesting question to investigate.

In any case, the experiment suggests it may be possible to mitigate the (lack of) global optimality issues we raise about Frank Wolfe simply by increasing the number of sample points. This is indeed a cheap way. Nevertheless, without the relaxation we propose, it is not possible to tell if the improved objective value obtained by increasing the number of samples is indeed globally optimal. Also, this experiment suggests that the SDP relaxation is able to obtain a far superior objective value with fewer data points compared to Frank Wolfe. One can further investigate trade-offs between these methods (as future work).

We caution against generalizing this approach of increasing samples to attain a better minima too broadly. For Frank Wolfe to converge to the global optimal solution, we need the optimization landscape to be favorable. This might not happen for more complicated data distributions or perhaps graphs. Again, without a relaxation (or a suitable dual problem), one cannot be certain if Frank Wolfe obtains the global solution.

Approximation gap

We checked this. The value is numerically equal to one. The plot appears to be less than one due to an issue with the graph plotting software.

Readability

Thanks for the suggestion. We will make these changes in the revision.

$m \neq n$

We performed an additional experiment where the number of samples in one distribution is fixed ($n=8$) and we vary the number of samples $m$ in the other distribution; see Figure 2 of the attached PDF. We notice that the relaxation is exact whenever $m$ is a multiple of $n$. On the other hand, the relaxation fails to be exact if $m$ is not a multiple of $n$.

Heuristic solver

Note that the heuristic solver will not be in the main body of the revised manuscript, based on a suggestion by another Reviewer.

The parameters are chosen heuristically. The parameters $\gamma_2​$ and $\gamma_3$ control the step size in each iteration, while $\gamma_1$​ is used for backtracking to prevent $\gamma_2$​ and $\gamma_3​$ from becoming too large. Through experimentation, we found that the algorithm performs better for relatively small choices.

Behavior of SDP relaxation when dimensions of the Gaussians change

We conducted experiments using Gaussian distributions with varying dimensions and found that the performance of the GW-SDP remained consistent across all tested dimensions. As such, we report results for a single representative dimension.

Barycenter algorithm

Could the Reviewer clarify which papers [1] and [2] refer to? Thanks.

Extensions

It is generally possible to extend the semidefinite relaxation to variants of the GW problem. We briefly describe the semidefinite relaxation to the outlier-robust GW problem by Kong et. al [KLTS:23]. Here, $(X,d_{X})$ and $(Y,d_{Y})$ are two metric spaces with accompanying measures $\mu$ and $\nu$. The distance between $\mu$ and $\nu$ is

$\min~~\langle L, P\rangle+\tau_1 d_{KL} (\pi 1,\alpha) + \tau_2 d_{KL} (\pi^T 1,\beta)$

$\mathrm{s.t.}\left( \begin{array}{cc} P & \mathrm{vec}(\pi)^T \\\ \mathrm{vec}(\pi) & 1 \end{array} \right) \succeq 0$

$\qquad\sum_{i}P_{(i,j),(k,l)} = f_{j}^{k,l}$, $\Sigma_{j} P_{(i,j),(k,l)} = g_i^{k,l}$

$\qquad P\geq 0$

$\qquad d_{KL}(\mu,\alpha)\leq\rho_1,d_{KL}(\nu,\beta)\leq\rho_2$

Note: The resulting formulation is convex but not an SDP because of the KL divergence. We discuss further extensions in the revised manuscript.

References

[KLTS:23] L Kong, J Li, J Tang, & A M-C So. (2023). Outlier-robust Gromov-Wasserstein for Graph Data. In Proceedings of the 37th International Conference on Neural Information Processing Systems.

[SP:04] RW Sumner and J Popovíc. (2004). Deformation transfer for triangle meshes. ACM Transactions on Graphics.

Run time comparisons

We wish to clarify that the current manuscript compares run times of the Conditional Gradient (CG) method, the entropic Gromov-Wasserstein method and our method. The comparisons are found in Page 5 of the submitted manuscript. For a fixed number of samples, the SDP runtime is far more costly than the other two algorithms.

We performed additional experiments to compare the runtime with algorithms suggested by the Reviewer. These results will be included in the revised manuscript.

We also conducted additional experiments where we compared the performance of GW-SDP with GW-CG and eGW-PPA. (These are described in the global comment and in the attached PDF.) Interestingly, in the experiment with real data, the objective values obtained from eGW-PPA are smaller than the values obtained from GW-CG, though still larger than those obtained by GW-SDP. (A smaller objective value is better in this instance.) This suggests the GW-SDP is the best algorithm for finding optimal solutions (in terms of the smallest objective value), followed by the eGW-PPA and then GW-CG.

First-order methods

First-order methods should improve the run-time in general. We investigated this approach but encountered numerous difficulties.

It is difficult to apply most first-order methods for solving SDPs to our relaxation. The main reason is that our formulation contains many constraints whereas most numerical techniques typically rely on having few affine constraints.

There is a line of work that uses proximal methods in combination with numerical schemes that allows families of constraints to be decoupled such as the ADMM method. Unfortunately, these ideas also do not translate easily. Our formulation has many different types of constraints: PSD, doubly stochasticity, non-negativity, and marginal sums. It is not clear how one designs proximal operators that project onto the intersection of these constraints, or design ADMM schemes that combine several different proximal operators.

There are some prior works that develop numerical schemes for solving semidefinite relaxations for the closely related Quadratic Assignment Problem (QAP) [KKBL:15,DML:17,OWX:18]. These works do point out that these semidefinite relaxations, while powerful, are not easy to solve. This suggests that it is equally difficult to develop numerical schemes for our semidefinite relaxations.

In fact, the scale of the problems we work on are not far from state-of-the-art for semidefinite relaxations of the QAP. In a work by Oliveira et al. -- a work purely about numerical optimization -- the authors propose an ADMM scheme to solve QAP instances where the dimension is 40 while the largest GW instance we solve has dimension m=n=32. This means that our problem size is not far from the state of the art.

We hope that the Reviewer accepts and values the conceptual and theoretical contributions already present in this paper, and recognizes that the scope of developing scalable numerical schemes, while important, goes comfortably beyond the scope of a single paper.

References

[KKBL:15] I. Kezurer, S. Z. Kovalsky, R. Basri, & Y. Lipman. (2015). Tight relaxation of quadratic match. In Computer Graphics Forum, volume 34, pages 115–128. Wiley Online Library.

[DML:17] N. Dym, H. Maron, & Y. Lipman. (2017). DS++: A Flexible, Scalable and Provably Tight Relaxation for Matching Problems. ACM Transactions on Graphics, 36(184):1–14.

[OWX:18] D. E. Oliveira, H. Wolkowicz, & Y. Xu, (2018). ADMM for the SDP relaxation of the QAP. Mathematical Programming Computation, 10(4), 631-658.