Dear Reviewer,
Thank you for your positive and constructive feedback on our work. We are particularly pleased that you highlighted the convergence guarantees for decreasing regularization in OT and the use of online algorithm that is suitable for large-scale problems. We address the specific concerns below:

Weaknesses

1. We thank the reviewer for this insightful question. It is correct that our theoretical analysis assumes uniform weights for the target distribution, as stated in Assumption 1. This restriction was made primarily to simplify the convergence proof. However, we believe that the results should still hold when the target distribution has strictly positive but unequal weights, although extending the theory to this more general setting would require additional technical work.
   In the revised version, we will make this point clearer and explicitly state that Example 2 was included to illustrate that DRAG continues to perform well in practice even when the uniformity assumption is relaxed.

2. Since the result of Proposition 2 holds for any $q$, there is no negative trade-off: we can always make the decay arbitrarily fast by increasing $q$, even when $2a$ is close to $b$. Indeed, given any $\delta > 0$ and any constants $a, b$ such that $2a < b$, and for any $s > 0$, we can always choose $q$ large enough so that

We will clarify this point in the revised version and emphasize that there is no trade-off in the choice of $a$ and $b$ here, thanks to the fact that the result holds for all $q > 0$.

3. We agree that the number of datasets considered in both the generative modeling and Monge-Kantorovich quantiles experiments is currently limited, and that more examples would help strengthen our conclusions.
   To address this, we have already extended both experimental sections and will include the results in the appendix of the revised version. In particular, we now provide more comprehensive examples involving 16 different measures for both the generative modeling and MK quantile tasks. Moreover, we use discrete measures with up to $10^8$ points for MK quantiles, emphasizing that DRAG, as an SGD-based algorithm, is well suited for large-scale problems.
   While the rebuttal process does not allow us to include additional figures, we summarize the findings here:

• For the generative modeling task, DRAG consistently outperformed Adam-based optimization, especially on more complex synthetic measures, where standard SGD showed instability.
• For the MK quantile task, the general observations remained the same: Neural OT performed significantly worse than DRAG , both in terms of coverage and bias, despite being significantly slower. DRAG continued to outperform Online Sinkhorn in terms of accurately covering the support of the target distribution without entropic bias, which aligns with the theoretical expectation of convex MK regions.
  We will make sure to discuss these additional results in the final version to provide a more complete empirical validation.

Questions

1. We thank the reviewer for the suggestion. In the revised version, we will present results in general form with explicit cost dependence when possible, especially when introducing the background on (entropic) Optimal Transport in Section 2. We will then specialize to the quadratic cost case only when needed, to enhance clarity and avoid unnecessary notational overload.

2. The projection step in Algorithm 1 is indeed a thresholding operation. This step is crucial in the convergence analysis of the non-averaged iterates of DRAG, as established in Theorem 1. In practice, we also find this projection step helpful, particularly during the early iterations, where it can improve both the convergence speed and the stability of the algorithm.

3. The notational framework on lines 200–203 was intended to emphasize that the optimal potential $\mathbf{g}^{\*}$ is only unique up to transformations of the form $\mathbf{g}^{\*} + a \mathbf{1}_M$, with $a \in \mathbb{R}$. Therefore, we conduct all our convergence analysis in the orthogonal complement of $\mathrm{vect}(\mathbf{1}_M)$. We apologize for the lack of clarity and will provide a clearer explanation in the final version of the paper.

Dear Reviewer,
We would like to thank you for your time and thoughtful comments. Please find below our responses to the weaknesses and questions you raised.

Weaknesses

1. Originality of DRAG

With our algorithm DRAG, and in the semi-discrete setting, we provide, to the best of our knowledge, the first theoretical proof of accelerated convergence rates due to decreasing regularization for solving an OT problem in any setting. While the heuristic of $\varepsilon$-scaling / $\varepsilon$-annealing [1] in discrete OT has been observed in practice to enhance the convergence rate of the Sinkhorn algorithm, theoretical guarantees that a decreasing regularization scheme can indeed accelerate the convergence of an OT algorithm remained an open question prior to our work.

Whereas algorithms with fixed regularization typically suffer from adverse dependence on the regularization parameter, our use of a decreasing regularization sequence entirely avoids this issue.
Moreover, unlike discrete solvers that reuse all samples multiple times (multi-pass), DRAG is adaptive to the number of samples and operates in an online (one-pass) setting, without the need to store past samples. This enables us to solve the semi-discrete OT problem without incurring either discretization or regularization bias.

In the final version, we will further emphasize these distinctions, highlighting DRAG’s main contribution: it is the first OT algorithm to confirm, through theory, the heuristic (previously only observed in practice) that decreasing regularization leads to acceleration.

2. Scalability of the method (This reply is mostly the same for both Reviewer fbxE and Reviewer zxMk.)

We reran our synthetic data experiments (Examples 1, 2, and 3) in dimensions 100 and 1000, instead of 10. We observed the same convergence rate behavior across all experiments with respect to the number of samples, with similar final errors. The experiments were approximately 4 and 15 times slower to run, respectively, which is solely due to the increased time required to compute $c(x, y) = \frac{1}{2}||x - y||^2$ on our GPU. This similar behavior in higher dimensions is not surprising in the semi-discrete OT setting, which is known not to suffer from the curse of dimensionality, as shown in [2]. While changing the parameter $d$ does not affect the error in our setting, increasing the number of points $M$ does impact the convergence rate, as we are solving a stochastic optimization problem in $\mathbb{R}^M$ via the semi-dual formulation. Since our algorithm is SGD-based and has linear memory complexity $\mathcal{O}(M)$ and computational complexity $\mathcal{O}(dM)$, we believe it is highly scalable. We will include these comments in the final version and thank the reviewer for highlighting this point. Furthermore, we will replace our synthetic examples in Figures 1, 2, and 3 with their dimension-1000 versions to better illustrate this scalability.

Questions

1. The averaging of the gradient steps is indeed essential, both theoretically and numerically. Theoretically: Averaging the iterates leads to improved convergence guarantees, as shown in Theorem 2. For instance, for any $\alpha > 1/2$, choosing parameters $a = \frac{1}{3}^-$ and $b = 2/3$ results in a convergence rate of only $\mathcal{O}(1/t^b)$ for the non-averaged iterates, compared to the optimal rate of $\mathcal{O}(1/t)$ for the averaged iterates. Moreover, our convergence analysis requires $b < 1$, which prevents us from achieving the optimal $\mathcal{O}(1/t)$ rate for non-averaged iterates in theory. Note that acceleration by averagin is classic in the stochastic gradient descent litterature [2]. Numerically: This difference is also observed in practice. The non-averaged iterates converge at a rate of approximately $\mathcal{O}(1/t^b)$, while the averaged iterates converge faster and exhibit improved stability. To illustrate this, we will include the errors of the non-averaged iterates in Figure 1, showing that averaging indeed enhances convergence in practice.

2. We confirm that these are indeed typos, and we thank you for pointing them out.

We hope that we have addressed all of your concerns and that the clarifications and improvements above provide a clearer picture of the quality, originality, and significance of our work.

References

[1] B. Schmitzer. Stabilized sparse scaling algorithms for entropy regularized transport problems. SIAM Journal on Scientific Computing

[2] Polyak, B. T., & Juditsky, A. B. (1992). Acceleration of stochastic approximation by averaging. SIAM Journal on Control and Optimization, 30(4), 838–855.

Dear Reviewer,
Thank you for your positive and constructive feedback on our work. We are particularly pleased that you highlighted the convergence guarantees for decreasing regularization in OT. We also believe this is the first such result. We appreciate your insightful questions and helpful suggestions.

Questions

A) We apologize for the lack of clarity. What we intended to convey is that, when using entropic OT as in [45], decreasing the regularization with respect to the number of points is necessary to achieve the optimal convergence rate. Using a larger, fixed regularization would lead to a non-optimal bias. Of course, an estimator that does not rely on entropic OT could, in principle, also achieve the optimal rate. Although, to the best of our knowledge, no such estimator currently exists for the semi-discrete OT setting. Therefore, while decreasing regularization is necessary when using entropic OT, using entropic OT itself may not be strictly necessary. We will better explain our meaning in the revised version of the paper.

B) We believe this question may have arisen due to a typo in the comment following the statement of Theorem 2, on line 241, where it incorrectly states $\alpha \in (0,1]$ instead of the correct $\alpha \in (1/2, 1]$. This typo was corrected in the full version of the paper included in the supplementary material. We apologize for the inconvenience. Both theorems allow for $\alpha \in (0,1]$ and yield a convergence rate. There is indeed an improvement in the convergence rate for the averaged case: when $\alpha > 1/2$, choosing $a = 1/3^-$ and $b = 2/3$ leads to a convergence rate of $\mathcal{O}(1/t^{2/3})$ for the non-averaged iterates, compared to the optimal rate of $\mathcal{O}(1/t)$ for the averaged iterates. This acceleration due to averaging is well known in the context of SGD algorithms [1], where, without averaging, the iterates tend to oscillate around the minimizer due to stochastic noise, while averaging “smooths out” this noise by aggregating information from multiple iterates.

C) We acknowledge that restricted strong convexity is a key property for the success of the convergence analysis of DRAG. In our setting, it holds locally with a parameter that does not depend on $\varepsilon$, because, denoting by $\rho^\varepsilon_*$ the smallest eigenvalue of the Hessian of $H_\varepsilon$ on $\mathrm{vect}(\mathbf{1})^\bot$, we have a uniform lower bound $\rho_* = \inf_{\varepsilon \in (0,1]} \rho^\varepsilon_*$. In particular, we do not have $\rho^\varepsilon_* \to 0$ as $\varepsilon \to 0$. This is due to the fact that, under our assumptions, even the non-regularized semi-dual has $\rho^0_* > 0$. This provides a form of local strong convexity, stemming from the regularity of the quadratic cost and the uniqueness of the minimizer. Using Theorem 3.2 in [2], we have $\rho^\varepsilon_* \geq \frac{1}{2}\rho^{\*}$ for all $\varepsilon \in [0,1]$. However, while the local restricted convexity is common for all $\varepsilon$, we see that the area in which it is true depends on $\varepsilon$, which thus show that when we decrease the regularization, the radius around $\mathbf{g}_{\varepsilon}^*$ in which we have a favorible restrcited strong convexity is smaller.
We hope this explanation clarifies the local behavior, and we will incorporate it into the final version of the paper.

References

[1] Polyak, B. T., & Juditsky, A. B. (1992). Acceleration of stochastic approximation by averaging. SIAM Journal on Control and Optimization, 30(4), 838–855.

Dear Reviewer,
Thank you for your detailed feedback and constructive comments. We also thank you for pointing out some typos such as the one regarding equation (2), they will be corrected in the revised version. We adress each point below.

Weaknesses

1. We acknowledge that several of the proofs are technically involved, and we sincerely appreciate the time and attention you dedicated to reading through them. We would like to emphasize that, to the best of our knowledge, the proofs are both correct and original, and we have carefully verified their correctness. We appreciate the reviewer's concern and would welcome any specific suggestions that could help improve the accessibility of our proofs.

2. Our results indeed cover only the quadratic cost case, as we rely on Corollary 2 in [1] concerning the convergence of entropic potentials, which we stated in Proposition 1. This result applies only to the quadratic cost. Note that [2] also focuses exclusively on the quadratic cost, as their analysis is similarly based on [1].
   That said, the two key technical lemmas: Lemma 1 (projection step) and Lemma 2 (Restricted Strong Convexity) can be generalized to other cost functions. In particular, the lower bound on $\rho_*$ stated on line 213 is proven for a broader class of costs in Proposition 5.1 of [3].
   Proposition 1, which is restricted to the quadratic cost, is therefore the limiting factor. Generalizing it to other costs would directly enable us to extend our results without modifying the convergence algorithm or the structure of our proof. However, we believe this generalization is technically challenging and constitutes an interesting direction for future work.
   In the final version, we will further emphasize that Proposition 1 is the limiting factor and will explicitly state the generalization to other costs as a valuable research direction.

3. (This reply is mostly the same for both Reviewer fbxE and Reviewer zxMk) We reran our synthetic data experiments (Examples 1, 2, and 3) in dimensions 100 and 1000, instead of 10. We observed the same convergence rate behavior across all experiments with respect to the number of samples, with similar final errors. The experiments were approximately 4 and 15 times slower to run, respectively, which is solely due to the increased time required to compute $c(x, y) = \frac{1}{2}||x - y||^2$ on our GPU.
   This similar behavior in higher dimensions is not surprising in the semi-discrete OT setting, which is known not to suffer from the curse of dimensionality, as shown in [2]. While changing the parameter $d$ does not affect the error in our setting, increasing the number of points $M$ does impact the convergence rate, as we are solving a stochastic optimization problem in $\mathbb{R}^M$ via the semi-dual formulation.
   Since our algorithm is SGD-based and has linear memory complexity $\mathcal{O}(M)$ and computational complexity $\mathcal{O}(dM)$ per iteration, we believe it is highly scalable.
   We will include these comments in the final version and thank the reviewer for highlighting this point. Furthermore, we will replace our synthetic examples in Figures 1, 2, and 3 with their dimension-1000 versions to better illustrate this scalability.

4. For real-world applications, DRAG can be used in generative modeling tasks that rely on semi-discrete OT to improve performance, as in [4] with GANs, [5] with diffusion models, and in statistics for Monge–Kantorovich quantile estimation [7]. Besides the fact that DRAG has optimal convergence guarantees for semi-discrete OT, we also observed in practice its better performance compared to Adam, which is not tailored for the OT task. As mentioned in our reply to Reviewer Ueg5’s request, we have extended both the generative modeling and Monge–Kantorovich quantile experiments to include more complex and large-scale examples.
   Furthermore, unlike DRAG, Adam does not come with theoretical guarantees for convergence to the true OT potential and map, to the best of our knowledge.

Questions

1. a) (l. 153) You are correct, the assumption is unnecessary.
   b) The notational framework on lines 202–203 was intended to emphasize that the optimal potential $\mathbf{g}^*$ is only unique up to transformations of the form $\mathbf{g}^{\*} + a \mathbf{1}_{M}$, with $a \in \mathbb{R}$. Therefore, we conduct all our convergence analysis in the orthogonal complement of $\mathrm{vect}(\mathbf{1}_M)$.

2. While extending our results to the continuous setting is an interesting direction for future work, our method fundamentally relies on the fact that the semi-discrete setting leads to a finite-dimensional problem. This property does not hold in the fully continuous setting. Therefore, the semi-discrete nature of the problem is indeed a limitation of our approach. That said, we would like to emphasize that the semi-discrete setting is, in itself, a very interesting and challenging problem with numerous applications. For example, it arises in generative modeling [4,5] and statistics, such as with Monge–Kantorovich quantiles.

3. While our algorithm achieves the best convergence rate when $\mu$ is $\alpha$-Hölder with $\alpha > 1/2$, we still guarantee convergence for any $\alpha$, as stated in Theorem 2. The condition $\alpha > 1/2$ is somewhat conservative and arises from our convergence analysis.
   That said, we believe that at least continuity of the source measure is necessary, since without that the optimal potential $\mathbf{g}^{\*}$ may not be unique, even in the orthogonal complement of $\mathrm{vect}(\mathbf{1}\_{M})$. For example, consider the case where $\mu = \frac{1}{2}\mathcal{U}([-2, -1]) + \frac{1}{2}\mathcal{U}([1, 2])$ and $\nu = \frac{1}{2}\delta_{-1} + \frac{1}{2}\delta_{1}$. Then, for the quadratic cost, all vectors of the form $\mathbf{g}_{x} = (-x, x)$ with $x \in (-1, 1)$ are optimal for the semi-dual problem, even though they all lie in the orthogonal of $\mathrm{vect}(\mathbf{1}_M)$.

4. (Cost question) see weakness 2.

Comparison with [2]

While [2] provides statistical guarantees for the approximation of the semi-discrete OT map and establishes that an optimal rate of $\mathcal{O}(1/\sqrt{t})$ is achievable, this is done by discretizing the source measure using $t$ points and employing the discrete entropic OT map as an estimator, with entropic regularization $\varepsilon \simeq 1/\sqrt{t}$.

Their work offers strong statistical guarantees for OT map estimation, provided that the discrete OT map is accurately computed. However, it focuses on an offline setting where all $t$ samples from the source measure are stored, and a discrete entropic OT problem is solved. The analysis is limited to sample complexity, specifically the statistical bias between the semi-discrete and discrete entropic OT maps. The computational complexity of Sinkhorn required to accurately approximate the OT map is not discussed, even though it typically scales poorly with $\varepsilon^{-1}$. Moreover, in this offline or multi-pass setting, increasing the number of Sinkhorn iterations does not correspond to using more samples.

In contrast, our paper directly tackles the semi-discrete OT problem without introducing discretization bias. We propose an algorithm that operates in an online (single-pass) setting, does not require storing past samples from the source measure, and demonstrates that decreasing regularization can indeed accelerate convergence for OT algorithms. Their paper provides stronger statistical guarantees, achieving optimal rates for $\alpha \in (0,1)$, whereas our results apply to $\alpha > 1/2$. However, we directly solve the semi-discrete problem, avoid discretization bias, and offer a scalable algorithm with linear computational and memory complexity per iteration.

Overall contribution of the paper

We would like to emphasize that, while decreasing regularization is a popular heuristic in optimal transport, DRAG is, to the best of our knowledge, the first algorithm in OT where acceleration due to decreasing regularization is formally proven. Prior to our work, theoretical guarantees that a decreasing regularization scheme can indeed accelerate the convergence of an OT algorithm remained an open question.
Whereas a fixed regularization typically leads to a complexity that depends adversely on $\varepsilon^{-1}$ in OT algorithms, our approach achieves convergence rates without such adverse dependence, thanks to the careful design of our regularization schedule. This theoretical benefit is reflected in practice, as shown in Figure 2, where DRAG outperforms all fixed regularization baselines.

We hope that our response has clarified the originality and significance of our work. We appreciate the reviewer’s time and hope they find the revised explanations helpful.

References

[1] Delalande, A. (2022). Nearly tight convergence bounds for semi-discrete entropic optimal transport. International Conference on Artificial Intelligence and Statistics (pp. 1619-1642). PMLR

[2] A.-A. Pooladian, V. Divol, and J. Niles-Weed. Minimax estimation of discontinuous optimal transport maps: The semi-discrete case. In the International Conference on Machine Learning, 2023.

[3] Chizat, L., Delalande, A., & Vaškevičius, T. (2024). Sharper Exponential Convergence Rates for Sinkhorn's Algorithm in Continuous Settings. arXiv preprint arXiv:2407.01202.

[4] An, D., Guo, Y., Lei, N., Luo, Z., Yau, S. T., & Gu, X. (2019). AE-OT: A new generative model based on extended semi-discrete optimal transport. ICLR 2020.

[5] Li, Zezeng, et al. "DPM-OT: a new diffusion probabilistic model based on optimal transport." Proceedings of the ieee/cvf international conference on computer vision. 2023.

Dear Reviewer,
Thank you for your follow-up and for increasing your score based on our clarifications and additional experiments. We truly appreciate your engagement with our work.

We understand your concerns about scalability. While 1000 dimensions may not appear large by machine learning standards, it corresponds to the true intrinsic dimension of the data in our setting. Our goal was to demonstrate the absence of the curse of dimensionality. If the curse were present, we would expect a convergence rate of the form $\mathcal{O}(n^{-1/d})$, which deteriorates with increasing dimension. Instead, we observe a convergence rate of $\mathcal{O}(n^{-1/2})$, which is independent of the dimension and thus highlights that semi-discrete OT avoids this curse. We acknowledge that optimal transport can be applied in high-dimensional ambient spaces; however, in practice, the data often lies on a low-dimensional manifold, which enables favorable convergence behavior [1].

You are right that our theoretical results currently apply only to the quadratic cost, corresponding to the 2-Wasserstein distance. While this is indeed a restriction, we note that the quadratic cost plays a central role in optimal transport theory and has deep connections to PDEs, geometry, and many practical applications. Although DRAG can be used with other cost functions without modification, extending the theoretical analysis beyond the quadratic case is an interesting direction for future work.

As for the continuous setting, we agree that our focus is on semi-discrete OT. Still, the ability of DRAG to operate online, with only one distribution discretized and without storing past samples, contrary to discrete solvers, makes it a practical proxy for large-scale continuous problems, especially given its linear memory usage in the number of discretized points.

Reference: [1] Weed, J., & Bach, F. (2019). Sharp asymptotic and finite-sample rates of convergence of empirical measures in Wasserstein distance. Bernoulli, 25(4A), 2620-2648.

Dear Reviewer,
Thank you for your positive evaluation of our work and for your questions. We are also glad that you appreciate that our work is, to the best of our knowledge, the first to provide theoretical guarantees for acceleration through a decreasing regularization scheme in OT, whereas this approach had previously been used in the discrete setting only as a heuristic.

Please find below our replies to the weaknesses you pointed out and the questions you raised.

Weaknesses

1. You are right, it should be $\frac{d\pi}{d\mu \otimes d\nu}$. We will incorporate this correction in the final version.

2. While $\gamma_k = \frac{\gamma_1}{k^b}$ with $\gamma_1, b > 0$ is the standard decreasing step size (learning rate) sequence in the SGD algorithm, as in [1, Equation (1)], it is true that in our case, the decay parameter $b$ depends on $a$ and $\alpha$ to achieve the best convergence rate. This follows from the convergence analysis of DRAG in Theorem 2, where the optimal rate is always attained for $b = \frac{2}{3}$ as soon as $\alpha > 1/2$, as discussed in the paragraph "Dependence on $a$, $b$, and $\alpha$" (l. 224).
   In the final version, we will emphasize that the specific scaling of $\gamma_k$ is motivated by the convergence proof of DRAG.

3. The projection step in Algorithm 1 is indeed a thresholding operation. This step is crucial in the convergence analysis of the non-averaged iterates of DRAG, as established in Theorem 1. In practice, we also find this projection step helpful, particularly during the early iterations, where it can improve both the convergence speed and the stability of the algorithm.
   In the later iterations, when the estimator $\mathbf{g}_t$ is close to $\mathbf{g}^*$, the thresholding becomes unnecessary and typically does not occur. In the final version, we will further emphasize the usefulness of this projection step, both from a theoretical perspective and in practice, especially in the early stages of the algorithm.

Question on Proposition 1

$K_0$ depends not only on the target measure $\nu$, but also on the source measure $\mu$, through the radius $R$ and the constants $m_1, m_2 > 0$ such that $m_1 < d\mu < m_2$. We will update Proposition 1 to reflect this dependence and include a reference to Remark 2.1 in [2].

Additionally, we will add a clarification that our choice of $a$ in the decreasing scheme $\varepsilon_k = 1/k^a$ (used in the convergence rate of DRAG in Theorem 2) ensures that the terms involving $K_0$ become negligible as the number of iterations increases. More precisely, they are of the order $o(||\bar{\mathbf{g}}_t - \mathbf{g}^*||^2)$ and therefore do not affect the leading convergence rate.

References

[1] Moulines, E., & Bach, F. (2011). Non-asymptotic analysis of stochastic approximation algorithms for machine learning. Advances in neural information processing systems, 24.

[2] Delalande, A. (2022). Nearly tight convergence bounds for semi-discrete entropic optimal transport. International Conference on Artificial Intelligence and Statistics (pp. 1619-1642). PMLR.