Interpolating between Optimal Transport and MMD using Sinkhorn Divergences

Introduction

The two main classes of losses $\text{{L}}(\alpha,\beta)$ which avoid these shortcomings are Optimal Transport distances and Maximum Mean Discrepancies: they are continuous with respect to the convergence in law and metrize its topology. That is, $\alpha_{n}\rightharpoonup\alpha\Leftrightarrow\text{{L}}(\alpha_{n},\alpha)\rightarrow 0$. The main purpose of this paper is to study the theoretical properties of a new class of geometric divergences which interpolates between these two families and thus offers an extra degree of freedom through a parameter $\varepsilon$ that can be cross-validated in typical learning scenarios.

1 Previous works

Out of this collection of methods, entropic regularization has recently emerged as a computationally efficient way of approximating OT costs. For $\varepsilon>0$, we define

MMD norms.

If $k$ is universal (Micchelli et al.,, 2006) (i.e. if the linear space spanned by functions $k(x,\cdot)$ is dense in $\mathcal{C}(\mathcal{X})$) we know that $\left\|\cdot\right\|_{k}$ metrizes the convergence in law. Such Euclidean norms, introduced for shape matching in (Glaunes et al.,, 2004), are often referred to as “Maximum Mean Discrepancies” (MMD) (Gretton et al.,, 2007). They have been extensively used for generative model (GANs) fitting in machine learning (Li et al.,, 2015; Dziugaite et al.,, 2015). MMD norms are cheaper to compute than OT and have a smaller sample complexity – i.e. approximation error when sampling a distribution.

2 Interpolating between OT and MMD using Sinkhorn divergences

Sinkhorn divergences. On the one hand, OT losses have appealing geometric properties; on the other hand, cheap MMD norms scales up to large batches with a low sample complexity. Why not interpolate between them to get the best of both worlds?

Following (Genevay et al.,, 2018) (see also (Ramdas et al.,, 2017; Salimans et al.,, 2018; Sanjabi et al.,, 2018)) we consider a new cost built from $\text{{OT}}_{\varepsilon}$ that we call a Sinkhorn divergence:

Such a formula satisfies $\text{{S}}_{\varepsilon}(\beta,\beta)=0$ and interpolates between OT and MMD (Ramdas et al.,, 2017):

In the literature, the formula (3) has been introduced more or less empirically to fix the entropic bias present in the $\text{{OT}}_{\varepsilon}$ cost: with a structure that mimicks that of a squared kernel norm (2), it was assumed or conjectured that $\text{{S}}_{\varepsilon}$ would define a positive definite loss function, suitable for applications in ML. This paper is all about proving that this is indeed what happens.

3 Contributions

The purpose of this paper is to show that the Sinkhorn divergences are convex, smooth, positive definite loss functions that metrize the convergence in law. Our main result is the theorem below, that ensures that one can indeed use $\text{{S}}_{\varepsilon}$ as a reliable loss function for ML applications – whichever value of $\varepsilon$ we pick.

Let $\mathcal{X}$ be a compact metric space with a Lipschitz cost function $\text{{C}}(x,y)$ that induces, for $\varepsilon>0$, a positive universal kernel $k_{\varepsilon}(x,y)\stackrel{{\scriptstyle\mbox{def.}}}{{=}}\exp(-\text{{C}}(x,y)/\varepsilon)$. Then, $\text{{S}}_{\varepsilon}$ defines a symmetric positive definite, smooth loss function that is convex in each of its input variables. It also metrizes the convergence in law: for all probability Radon measures $\alpha$ and $\beta\in\mathcal{M}_{1}^{+}(\mathcal{X})$,

This theorem legitimizes the use of the unbiased Sinkhorn divergences $\text{{S}}_{\varepsilon}$ instead of $\text{{OT}}_{\varepsilon}$ in model-fitting applications. Indeed, computing $\text{{S}}_{\varepsilon}$ is roughly as expensive as $\text{{OT}}_{\varepsilon}$ (the computation of the corrective factors being cheap, as detailed in Section 3) and the “debiasing” formula (3) allows us to guarantee that the unique minimizer of $\alpha\mapsto\text{{S}}_{\varepsilon}(\alpha,\beta)$ is the target distribution $\beta$ (see Figure 1). Section 3 details how to implement these divergences efficiently: our algorithms scale up to millions of samples thanks to freely available GPU routines. To conclude, we showcase in Section 4 the typical behavior of $\text{{S}}_{\varepsilon}$ compared with $\text{{OT}}_{\varepsilon}$ and standard MMD losses.

Proof of Theorem 1

We now give the proof of Theorem 1. Our argument relies on a new Bregman divergence derived from a weak∗ continuous entropy that we call the Sinkhorn entropy (see Section 2.2). We believe this (convex) entropy function to be of independent interest. Note that all this section is written under the assumptions of Theorem 1; the proof of some intermediate results can be found in the appendix.

First, let us recall some standard results of regularized OT theory (Peyré and Cuturi,, 2017). Thanks to the Fenchel-Rockafellar theorem, we can rewrite Cuturi’s loss (1) as

where $f\oplus g$ is the tensor sum $(x,y)\in\mathcal{X}^{2}\mapsto f(x)+g(y)$. The primal-dual relationship linking an optimal transport plan $\pi$ solving (1) to an optimal dual pair $(f,g)$ that solves (8) is

Crucially, the first order optimality conditions for the dual variables are equivalent to the primal’s marginal constraints $(\pi_{1}=\alpha,\pi_{2}=\beta)$ on (9). They read

where the “Sinkhorn mapping” $\text{{T}}:\mathcal{M}_{1}^{+}(\mathcal{X})\times\mathcal{C}(\mathcal{X})\rightarrow\mathcal{C}(\mathcal{X})$ is defined through

with a SoftMin operator of strength $\varepsilon$ defined through

The following proposition recalls some important properties of $\text{{OT}}_{\varepsilon}$ and the associated dual potentials. Its proof can be found in Section B.1.

The following proposition, whose proof is detailed in Section B.2, shows that the dual potentials are the gradients of $\text{{OT}}_{\varepsilon}$.

$\text{{OT}}_{\varepsilon}$ is weak* continuous and differentiable. Its gradient reads

where $(f,g)$ satisfies $f=\text{{T}}(\beta,g)$ and $g=\text{{T}}(\alpha,f)$ on the whole domain $\mathcal{X}$ and T is the Sinkhorn mapping (11).

Let us stress that even though the solutions of the dual problem (8) are defined $(\alpha,\beta)$-a.e., the gradient (15) is defined on the whole domain $\mathcal{X}$. Fortunately, an optimal dual pair $(f_{0},g_{0})$ defined $(\alpha,\beta)$-a.e. satisfies the optimality condition (10) and can be extended in a canonical way: to compute the “gradient” pair $(f,g)\in\mathcal{C}(\mathcal{X})^{2}$ associated to a pair of measures $(\alpha,\beta)$, using $f=\text{{T}}(\beta,g_{0})$ and $g=\text{{T}}(\alpha,f_{0})$ is enough.

2 Sinkhorn and Haussdorf divergences

Crucially, we now assume that $k_{\varepsilon}$ is a positive universal kernel on the space of signed Radon measures.

Under the assumptions above, we define the Sinkhorn negentropy of a probability Radon measure $\alpha\in\mathcal{M}_{1}^{+}(\mathcal{X})$ through

The following proposition is the cornerstone of our approach to prove the positivity of $\text{{S}}_{\varepsilon}$, providing an alternative expression of $\text{{F}}_{\varepsilon}$. Its proof relies on a change of variables $\mu=\exp(f/\varepsilon)\,\alpha$ in (8) that is detailed in the Section B.3 of the appendix.

The following proposition, whose proof can be found in the Section B.4 of the appendix, leverages the alternative expression (18) to ensure the convexity of $\text{{F}}_{\varepsilon}$.

Under the same hypotheses as Proposition 3, $\text{{F}}_{\varepsilon}$ is a strictly convex functional on $\mathcal{M}_{1}^{+}(\mathcal{X})$.

We now define an auxiliary “Hausdorff” divergence that can be interpreted as an $\text{{OT}}_{\varepsilon}$ loss with decoupled dual potentials.

Thanks to Proposition 2, the Sinkhorn negentropy $\text{{F}}_{\varepsilon}$ is differentiable in the sense of (14). For any probability measures $\alpha,\beta\in\mathcal{M}_{1}^{+}(\mathcal{X})$ and regularization strength $\varepsilon>0$, we can thus define

It is the symmetric Bregman divergence induced by the strictly convex functional $F_{\varepsilon}$ (Bregman,, 1967) and is therefore a positive definite quantity.

3 Proof of the Theorem

We are now ready to conclude. First, remark that the dual expression (8) of $\text{{OT}}_{\varepsilon}(\alpha,\beta)$ as a maximization of linear forms ensures that $\text{{OT}}_{\varepsilon}(\alpha,\beta)$ is convex with respect to $\alpha$ and with respect to $\beta$ (but not jointly convex if $\varepsilon>0$). $\text{{S}}_{\varepsilon}$ is thus convex with respect to both inputs $\alpha$ and $\beta$ as a sum of the functions $\text{{OT}}_{\varepsilon}$ and $\text{{F}}_{\varepsilon}$ – see Proposition 4.

Using (15) to get $\nabla_{2}\text{{OT}}_{\varepsilon}(\alpha,\alpha)=-\nabla\text{{F}}_{\varepsilon}(\alpha)$, $\nabla_{1}\text{{OT}}_{\varepsilon}(\beta,\beta)=-\nabla\text{{F}}_{\varepsilon}(\beta)$ and summing the above inequalities, we show that $\text{{H}}_{\varepsilon}\leqslant\text{{S}}_{\varepsilon}$, which implies (5).

To prove (6), note that $\text{{S}}_{\varepsilon}(\alpha,\beta)=0\Rightarrow\text{{H}}_{\varepsilon}(\alpha,\beta)=0$, which implies that $\alpha=\beta$ since $\text{{F}}_{\varepsilon}$ is a strictly convex functional.

Finally, we show that $\text{{S}}_{\varepsilon}$ metrizes the convergence in law (7) in the Section B.5 of the appendix.

Computational scheme

We have shown that Sinkhorn divergences (3) are positive definite, convex loss functions on the space of probability measures. Let us now detail their implementation on modern hardware.

1 The Sinkhorn algorithm(s)

Proposition 1 is key to the modern theory of regularized Optimal Transport: it allows us to compute the $\text{{OT}}_{\varepsilon}$ cost – and thus the Sinkhorn divergence $\text{{S}}_{\varepsilon}$, thanks to (3) – using dual variables that have the same memory footprint as the input measures: solving (8) in our discrete setting, we only need to store the sampled values of the dual potentials $f$ and $g$ on the measures’ supports.

denotes a (stabilized) log-sum-exp reduction.

If $(\bm{f},\bm{g})$ is an optimal pair of dual vectors that satisfies Equations (20-21), we deduce from (13) that

But how can we solve this coupled system of equations given $\bm{\alpha}$, $\bm{x}$, $\bm{\beta}$ and $\bm{y}$ as input data?

The Sinkhorn algorithm.

One simple answer: by enforcing (20) and (21) alternatively, updating the vectors $\bm{f}$ and $\bm{g}$ until convergence (Cuturi,, 2013). Starting from null potentials $\bm{f}_{i}=0=\bm{g}_{j}$, this numerical scheme is nothing but a block-coordinate ascent on the dual problem (8). One step after another, we are enforcing null derivatives on the dual cost with respect to the $\bm{f}_{i}$’s and the $\bm{g}_{j}$’s.

Convergence.

The “Sinkhorn loop” converges quickly towards its unique optimal value: it enjoys a linear convergence rate (Peyré and Cuturi,, 2017) that can be improved with some heuristics (Thibault et al.,, 2017). When computed through the dual expression (23), $\text{{OT}}_{\varepsilon}$ and its gradients (26-27) are robust to small perturbations of the values of $\bm{f}$ and $\bm{g}$: monitoring convergence through the $\text{{L}}^{1}$ norm of the updates on $\bm{f}$ and breaking the loop as we reach a set tolerance level is thus a sensible stopping criterion. In practice, if $\varepsilon$ is large enough – say, $\varepsilon\geqslant\texttt{.05}$ on the unit square with an Earth Mover’s cost $\text{{C}}(x,y)=\|x-y\|$ – waiting for 10 or 20 iterations is more than enough.

All in all, the baseline Sinkhorn loop provides an efficient way of solving the discrete problem $\text{{OT}}_{\varepsilon}(\alpha,\beta)$ for generic input measures. But in the specific case of the (symmetric) corrective terms $\text{{OT}}_{\varepsilon}(\alpha,\alpha)$ and $\text{{OT}}_{\varepsilon}(\beta,\beta)$ introduced in (3), we can do better.

The key here is to remark that if $\alpha=\beta$, the dual problem (8) becomes a concave maximization problem that is symmetric with respect to its two variables $f$ and $g$. Hence, there exists a (unique) optimal dual pair $(f,g=f)$ on the diagonal which is characterized in the discrete setting by the symmetric optimality condition: $\forall i\in[1,\text{{N}}],$

Fortunately, given $\bm{\alpha}$ and $\bm{x}$, the optimal vector $\bm{f}$ that solves this equation can be computed by iterating a well-conditioned fixed-point update:

This symmetric variant of the Sinkhorn algorithm can be shown to converge much faster than the standard loop applied to a pair $(\alpha,\beta=\alpha)$ on the diagonal, and three iterations are usually enough to compute accurately the optimal dual vector.

2 Computing the Sinkhorn divergence and its gradients

In this day and age, we could be tempted to rely on the automatic differentiation engines provided by modern libraries, which let us differentiate the result of twenty or so Sinkhorn iterations as a mere composition of elementary operations (Genevay et al.,, 2018). But beware: this loop has a lot more structure than a generic feed forward network. Taking advantage of it is key to a x2-x3 gain in performances, as we now describe.

Crucially, we must remember that the Sinkhorn loop is a fixed point iterative solver: at convergence, its solution satisfies an equation given by the implicit function theorem. Thanks to (15), using the very definition of gradients in the space of probability measures (14) and the intermediate variables in the computation of $\text{{S}}_{\varepsilon}(\alpha,\beta)$, we get that

Graph surgery with PyTorch.

Assuming convergence in the Sinkhorn loops, it is thus possible to compute the gradients of $\text{{S}}_{\varepsilon}$ without having to backprop through the twenty or so iterations of the Sinkhorn algorithm: we only have to differentiate the expression above with respect to $x$. But does it mean that we should differentiate C or the log-sum-exp operation by hand? Fortunately, no!

Modern libraries such as PyTorch (Paszke et al.,, 2017) are flexible enough to let us “hack” the naive autograd algorithm, and act as though the optimal dual vectors $\bm{f}_{i}$, $\bm{p}_{i}$, $\bm{g}_{j}$ and $\bm{q}_{j}$ did not depend on the input variables of $\text{{S}}_{\varepsilon}$. As documented in our reference code,

an appropriate use of the .detach() method in PyTorch is enough to get the best of both worlds: an automatic differentiation engine that computes our gradients using the formula at convergence instead of the baseline backpropagation algorithm. All in all, as evidenced by the benchmarks provided Figure 3, this trick allows us to divide by a factor 2-3 the time needed to compute a Sinkhorn divergence and its gradient with respect to the $\bm{x}_{i}$’s.

3 Scaling up to large datasets

The Sinkhorn iterations rely on a single non-trivial operation: the log-sum-exp reduction (22). In the ML literature, this SoftMax operator is often understood as a row- or column-wise reduction that acts on $[\text{{N}},\text{{M}}]$ matrices. But as we strive to implement the update rules (20-21) and (25) on the GPU, we can go further.

First, if the number of samples N and M in both measures is small enough, we can optimize the GPU usage by computing Sinkhorn divergences by batches of size B. In practice, this can be achieved by encoding the cost function C as a 3D tensor of size $[\text{{B}},\text{{N}},\text{{M}}]$ made up of stacked matrices $(\text{{C}}(\bm{x}_{i},\bm{y}_{j}))_{i,j}$, while $\bm{f}$ and $\bm{g}$ become $[\text{{B}},\text{{N}}]$ and $[\text{{B}},\text{{M}}]$ tensors, respectively. Thanks to the broadcasting syntax supported by modern libraries, we can then seamlessly compute, in parallel, loss values $\text{{S}}_{\varepsilon}(\alpha_{k},\beta_{k})$ for $k$ in $[1,\text{{B}}]$.

The KeOps library.

Unfortunately though, tensor-centric methods such as the one presented above cannot scale to measures sampled with large numbers N and M of Dirac atoms: as these numbers exceed 10,000, huge $[\text{{N}},\text{{M}}]$ matrices stop fitting into GPU memories. To alleviate this problem, we leveraged the KeOps library (Charlier et al.,, 2018) that provides online map-reduce routines on the GPU with full PyTorch integration. Performing online log-sum-exp reductions with a running maximum, the KeOps primitives allow us to compute Sinkhorn divergences with a linear memory footprint. As evidenced by the benchmarks of Figures 3-3, computing the gradient of a Sinkhorn loss with 100,000 samples per measure is then a matter of seconds.

Numerical illustration

In the previous sections, we have provided theoretical guarantees on top of a comprehensive implementation guide for the family of Sinkhorn divergences $\text{{S}}_{\varepsilon}$. Let us now describe the geometry induced by these new loss functions on the space of probability measures.

To compare MMD losses $\text{{L}}_{k}$ with Cuturi’s original cost $\text{{OT}}_{\varepsilon}$ and the de-biased Sinkhorn divergence $\text{{S}}_{\varepsilon}$, a simple yet relevant experiment is to let a model distribution $\alpha(t)$ flow with time $t$ along the “Wasserstein-2” gradient flow of a loss functional $\alpha\mapsto\text{{L}}(\alpha,\beta)$ that drives it towards a target distribution $\beta$ (Santambrogio,, 2015). This corresponds to the “non-parametric” version of the data fitting problem evoked in Section 1, where the parameter $\theta$ is nothing but the vector of positions $\bm{x}$ that encodes the support of a measure $\alpha=\tfrac{1}{\text{{N}}}\sum_{i=1}^{N}\delta_{\bm{x}_{i}}$. Understood as a “model free” idealization of fitting problems in machine learning, this experiment allows us to grasp the typical behavior of the loss function as we discover the deformations of the support that it favors.

with a Euler scheme and display the evolution of $\alpha(t)$ up to time $t=5$.

Interpretation.

In both figures, the fourth line highlights the entropic bias that is present in the $\text{{OT}}_{\varepsilon}$ loss: $\alpha(t)$ is driven towards a minimizer that is a “shrunk” version of $\beta$. As showed in Theorem 1, the de-biased loss $\text{{S}}_{\varepsilon}$ does not suffer from this issue: just like MMD norms, it can be used as a reliable, positive-definite divergence.

Going further, the dynamics induced by the Sinkhorn divergence interpolates between that of an MMD ($\varepsilon=+\infty$) and Optimal Transport ($\varepsilon=0$), as shown in (4). Here, $\text{{C}}(x,y)=\|x-y\|$ and we can indeed remark that the second and third lines bridge the gap between the flow of the energy distance $\text{{L}}_{-\|\cdot\|}$ (in the first line) and that of the Earth Mover’s cost $\text{{OT}}_{0}$ which moves particles according to an optimal transport plan.

Please note that in both experiments, the gradient of the energy distance with respect to the $\bm{x}_{i}$’s vanishes at the extreme points of $\alpha$’s support. Crucially, for small enough values of $\varepsilon$, $\text{{S}}_{\varepsilon}$ recovers the translation-aware geometry of OT and we observe a clean convergence of $\alpha(t)$ to $\beta$ as no sample lags behind.

Conclusion

Recently introduced in the ML literature, the Sinkhorn divergences were designed to interpolate between MMD and OT. We have now shown that they also come with a bunch of desirable properties: positivity, convexity, metrization of the convergence in law and scalability to large datasets.

To the best of our knowledge, it is the first time that a loss derived from the theory of entropic Optimal Transport is shown to stand on such a firm ground. As the foundations of this theory are progressively being settled, we now hope that researchers will be free to focus on one of the major open problems in the field: the interaction of geometric loss functions with concrete machine learning models.

Appendix A Standard results

Before detailing our proofs, we first recall some well-known results regarding the Kullback-Leibler divergence and the SoftMin operator defined in (12).

It can be rewritten as an $f$-divergence associated to

Dual formulation.

Lower bound on the sup. If $\alpha$ is not absolutely continuous with respect to $\beta$, there exists a Borel set $A$ such that $\alpha(A)>0$ and $\beta(A)=0$. Consequently, for $h=\lambda\,\textbf{1}_{A}$,

Going further, the density of continuous functions in the space of bounded measurable functions allows us to restrict the optimization domain:

Let $h=\sum_{i\in I}h_{i}\,\textbf{1}_{A_{i}}$ be a simple Borel function on $\mathcal{X}$, and let us choose some error margin $\delta>0$. Since $\alpha$ and $\beta$ are Radon measures, for any $i$ in the finite set of indices $I$, there exists a compact set $K_{i}$ and an open set $V_{i}$ such that $K_{i}\subset A_{i}\subset V_{i}$ and

Moreover, for any $i\in I$, there exists a continuous function $\varphi_{i}$ such that $\textbf{1}_{K_{i}}\leqslant\varphi_{i}\leqslant\textbf{1}_{V_{i}}$. The continuous function $g=\sum_{i\in I}h_{i}\varphi_{i}$ is then such that

We can then show that the Kullback-Leibler divergence is weakly lower semi-continuous:

If $\alpha_{n}\rightharpoonup\alpha$ and $\beta_{n}\rightharpoonup\beta$ are weakly converging sequences in $\mathcal{M}^{+}(\mathcal{X})$, we get

According to (31), the KL divergence is defined as a pointwise supremum of weakly continuous applications

A.2 SoftMin Operator

Under the assumptions of the definition (12), we get that

If $\varphi$ and $\psi$ are two continuous functions in $\mathcal{C}(\mathcal{X})$ such that $\varphi\leqslant\psi$,

Let $(\alpha_{n})$ be a sequence of probability measures converging weakly towards $\alpha$, and $(\varphi_{n})$ be a sequence of continuous functions that converges uniformly towards $\varphi$. Then, for $\varepsilon>0$, the SoftMin of the values of $\varphi_{n}$ on $\alpha_{n}$ converges towards the SoftMin of the values of $\varphi$ on $\alpha$, i.e.

Appendix B Proofs

The existence of an optimal pair $(f,g)$ of potentials that reaches the maximal value of the dual objective is proved using the contractance of the Sinkhorn map T, defined in (11), for the Hilbert projective metric (Franklin and Lorenz,, 1989).

While optimal potentials are only defined $(\alpha,\beta)$-a.e., as highlighted in Proposition 1, they are extended to the whole domain $\mathcal{X}$ by imposing, similarly to the classical theory of OT (Santambrogio,, 2015, Remark 1.13), that they satisfy

with $T$ defined in (11). We thus assume in the following that this condition holds. The following propositions studies the uniqueness and the smoothness (with respect to the spacial position and with respect to the input measures) of these functions $(f,g)$ defined on the whole space.

For $t\in$, let us define $f_{t}=f_{0}+t(f_{1}-f_{0})$, $g_{t}=g_{0}+t(g_{1}-g_{0})$ and

the value of the dual objective between the two optimal pairs. As $\varphi$ is a concave function bounded above by $\varphi(0)=\varphi(1)=\text{{OT}}_{\varepsilon}(\alpha,\beta)$, it is constant with respect to $t$. Hence, for all $t$ in $$,

This is only possible if, $\alpha\otimes\beta$-a.e. in $(x,y)$,

As we extend the potentials through (34), the SoftMin operator commutes with the addition of $K$ (33) and lets our result hold on the whole feature space. ∎

According to (34), $f$ is a SoftMin combination of $\kappa$-Lipschitz functions of the variable $x$; using the algebraic properties of the SoftMin operator detailed in (32-33), one can thus show that $f$ is a $\kappa$-Lipschitz function on the feature space. The same argument holds for $g$. ∎

Let $\alpha_{n}\rightharpoonup\alpha$ and $\beta_{n}\rightharpoonup\beta$ be weakly converging sequences of measures in $\mathcal{M}_{1}^{+}(\mathcal{X})$. Given some arbitrary anchor point $x_{o}\in\mathcal{X}$, let us denote by $(f_{n},g_{n})$ the (unique) sequence of optimal potentials for $\text{{OT}}_{\varepsilon}(\alpha_{n},\beta_{n})$ such that $f_{n}(x_{o})=0$.

Then, $f_{n}$ and $g_{n}$ converge uniformly towards the unique pair of optimal potentials $()$ for $\text{{OT}}_{\varepsilon}(\alpha,\beta)$ such that $f(x_{o})=0$. Up to the value at the anchor point $x_{o}$, we thus have that

Being equicontinuous and uniformly bounded on the compact set $\mathcal{X}$, the sequence $(f_{n},g_{n})_{n}$ satisfies the hypotheses of the Ascoli-Arzela theorem: there exists a subsequence $(f_{n_{k}},g_{n_{k}})_{k}$ that converges uniformly towards a pair $()$ of continuous functions. $k$ tend to infinity, we see that $f(x_{o})=0$ and, using the continuity of the SoftMin operator (Proposition 10) on the optimality equations (10), we show that $()$ is an optimal pair for $\text{{OT}}_{\varepsilon}(\alpha,\beta)$.

Now, according to Proposition 11, such a limit pair of optimal potentials $()$ is unique. $(f_{n},g_{n})_{n}$ is thus a compact sequence with a single possible adherence value: it has to converge, uniformly, towards $()$. ∎

B.2 Proof of Proposition 2

The proof is mainly inspired from (Santambrogio,, 2015, Proposition 7.17). Let us consider $\alpha$, $\delta\alpha$, $\beta$, $\delta\beta$ and times $t$ in a neighborhood of , as in the statement above. We define $\alpha_{t}=\alpha+t\delta\alpha$, $\beta_{t}=\beta+t\delta\beta$ and the variation ratio $\Delta_{t}$ given by

Using the very definition of $\text{{OT}}_{\varepsilon}$ and the continuity property of Proposition 13, we now provide lower and upper bounds on $\Delta_{t}$ as $t$ goes to .

As written in (13), $\text{{OT}}_{\varepsilon}(\alpha,\beta)$ can be computed through a straightforward, continuous expression that does not depend on the value of the optimal dual potentials $()$ at the anchor point $x_{o}$:

Combining this equation with Proposition 13 (that guarantees the uniform convergence of potentials for weakly converging sequences of probability measures) allows us to conclude.

Lower bound.

First, let us remark that $()$ is a suboptimal pair of dual potentials for $\text{{OT}}_{\varepsilon}(\alpha_{t},\beta_{t})$. Hence,

since $g$ and $f$ satisfy the optimality equations (10).

Upper bound.

Conversely, let us denote by $(g_{t},f_{t})$ the optimal pair of potentials for $\text{{OT}}_{\varepsilon}(\alpha_{t},\beta_{t})$ satisfying $g_{t}(x_{o})=0$ for some arbitrary anchor point $x_{o}\in\mathcal{X}$. As $(f_{t},g_{t})$ are suboptimal potentials for $\text{{OT}}_{\varepsilon}(\alpha,\beta)$, we get that

Conclusion.

Thanks to Proposition 13, we thus know that $f_{t}$ and $g_{t}$ converge uniformly towards $f$ and $g$. Combining the lower and upper bound, we get

since $\delta\alpha$ and $\delta\beta$ both have an overall mass that sums up to zero.

B.3 Proof of Proposition 3

The definition of $\text{{OT}}_{\varepsilon}(\alpha,\alpha)$ is that

Thanks to the symmetry of this concave problem with respect to the variables $f$ and $g$, we know that there exists a pair $(f,g=f)$ of optimal potentials on the diagonal, and

where $\star$ denotes the smoothing (convolution) operator defined through

for $k\in\mathcal{C}(\mathcal{X}\times\mathcal{X})$ and $\mu\in\mathcal{M}^{+}(\mathcal{X})$.

Optimizing on measures.

keeping in mind that $\alpha$ is a probability measure, we then get that

where we optimize on positive measures $\mu\in\mathcal{M}^{+}(\mathcal{X})$ such that $\alpha\ll\mu$ and $\mu\ll\alpha$.

Expansion of the problem.

As $k_{\varepsilon}(x,y)=\exp(-\text{{C}}(x,y)/\varepsilon)$ is positive for all $x$ and $y$ in $\mathcal{X}$, we can remove the $\mu\ll\alpha$ constraint from the optimization problem:

Existence of the optimal measure μ𝜇\mu.

In the expression above, the existence of an optimal $\mu$ is given as a consequence of the well-known fact from OT theory that optimal dual potentials $f$ and $g$ exist, so that the dual OT problem (8) is a max and not a mere supremum. Nevertheless, since this property of $\text{{F}}_{\varepsilon}$ is key to the metrization of the convergence in law by Sinkhorn divergences, let us endow it with a direct, alternate proof:

For any $\alpha\in\mathcal{M}_{1}^{+}(\mathcal{X})$, assuming that $\mathcal{X}$ is compact, there exists a unique $\mu_{\alpha}\in\mathcal{M}^{+}(\mathcal{X})$ such that

Moreover, $\alpha\ll\mu_{\alpha}\ll\alpha$.

Notice that for $(\alpha,\mu)\in\mathcal{M}_{1}^{+}(\mathcal{X})\times\mathcal{M}^{+}(\mathcal{X})$,

Since C is bounded on the compact set $\mathcal{X}\times\mathcal{X}$ and $\alpha$ is a probability measure, we can already say that

Upper bound on the mass of μ𝜇\mu.

Since $\mathcal{X}\times\mathcal{X}$ is compact and $k_{\varepsilon}(x,y)>0$, there exists $\eta>0$ such that $k(x,y)>\eta$ for all $x$ and $y$ in $\mathcal{X}$. We thus get

As we build a minimizing sequence $(\mu_{n})$ for $\text{{F}}_{\varepsilon}(\alpha)$, we can thus assume that $\langle\mu_{n},\,1\rangle$ is uniformly bounded by some constant $M>0$.

Weak continuity.

Crucially, the Banach-Alaoglu theorem asserts that

is weakly compact; we can thus extract a weakly converging subsequence $\mu_{n_{k}}\rightharpoonup\mu_{\infty}$ from the minimizing sequence $(\mu_{n})$. Using Proposition 8 and the fact that $k_{\varepsilon}$ is continuous on $\mathcal{X}\times\mathcal{X}$, we show that $\mu\mapsto\text{{E}}_{\varepsilon}(\alpha,\mu)$ is a weakly lower semi-continuous function: $\mu_{\infty}=\mu_{\alpha}$ realizes the minimum of $\text{{E}}_{\varepsilon}$ and we get our existence result.

Uniqueness.

We assumed that our kernel $k_{\varepsilon}$ is positive universal. The squared norm $\mu\mapsto\|\mu\|_{k_{\varepsilon}}^{2}$ is thus a strictly convex functional and using Proposition 6, we can show that $\mu\mapsto\text{{E}}_{\varepsilon}(\alpha,\mu)$ is strictly convex. This ensures that $\mu_{\alpha}$ is uniquely defined. ∎

B.4 Proof of Proposition 4

Let us take a pair of measures $\alpha_{0}\neq\alpha_{1}$ in $\mathcal{M}_{1}^{+}(\mathcal{X})$, and $t\in(0,1)$; according to Proposition 14, there exists a pair of measures $\mu_{0}$, $\mu_{1}$ in $\mathcal{M}^{+}(\mathcal{X})$ such that

which is enough to conclude. To show the strict inequality, let us remark that

B.5 Proof of the Metrization of the Convergence in Law

The regularized OT cost is weakly continuous, and the uniform convergence for dual potentials ensures that $\text{{H}}_{\varepsilon}$ and $\text{{S}}_{\varepsilon}$ are both continuous too. Paired with (6), this property guarantees the convergence towards of the Hausdorff and Sinkhorn divergences, as soon as $\alpha_{n}\rightharpoonup\alpha$.

Conversely, let us assume that $\text{{S}}_{\varepsilon}(\alpha_{n},\alpha)\rightarrow 0$ (resp. $\text{{H}}_{\varepsilon}(\alpha_{n},\alpha)$). Any weak limit $\alpha_{n_{\infty}}$ of a subsequence $(\alpha_{n_{k}})_{k}$ is equal to $\alpha$: since our divergence is weakly continuous, we have $\text{{S}}_{\varepsilon}(\alpha_{n_{\infty}},\alpha)=0$ (resp. $\text{{H}}_{\varepsilon}(\alpha_{n_{\infty}},\alpha)$), and positive definiteness holds through (6).

In the meantime, since $\mathcal{X}$ is compact, the set of probability Radon measures $\mathcal{M}_{1}^{+}(\mathcal{X})$ is sequentially compact for the weak-$\star$ topology. $\alpha_{n}$ is thus a compact sequence with a unique adherence value: it converges, towards $\alpha$.