When OT meets MoM: Robust estimation of Wasserstein Distance

Guillaume Staerman, Pierre Laforgue, Pavlo Mozharovskyi, Florence d'Alché-Buc

Introduction

Computing distances between probability distributions has become a central question in numerous modern Machine Learning applications, ranging from generative modeling to clustering. Optimal Transport (OT) offers an appealing and insightful tool to solve this problem, building upon the Wasserstein distance. Given two probability distributions, the latter is defined in terms of the solution to the Monge-Kantorovich optimal mass transportation problem. Interestingly, it relies on a ground distance between points to build a distance between probability distributions . For that reason, the Wasserstein distance stands out from the divergences usually exploited in generative modeling, like the f-divergences , by its ability to take into account the underlying geometry of the space, capturing the difference between probability distributions even when they have non-overlapping supports. This appealing property has been successfully exploited in Generative Adversarial Networks (GANs) , as well as in Variational Autoencoders (VAEs) , where the Wasserstein distance can advantageously replace an f-divergence as the loss function. Many other applications rely on the entropic-regularized approximations introduced by , which has considerably alleviated the inherent computational complexity of the Wasserstein distance in the discrete case, by drawing on the Sinkhorn-Knopp algorithm. A common feature to almost all these works is that the Wasserstein distance is estimated from finite samples. While this problem has long been theoretically studied under the i.i.d. assumption , it has never been tackled through the lens of robustness to outliers, a crucial issue in Reliable Machine Learning. Indeed, data is nowadays collected at a large scale in unmastered acquisition conditions, and through a large variety of devices and platforms. The resulting datasets often present undesirable influential observations, whether they are errors or rare observations. The presence of corrupted data may heavily damage the quality of estimators, calling for dedicated methods such as JS/TV-GANs in the particular case of robust shift-parameter estimation, Robust Divergences in variational inference , or more general tools from robust statistics .

The aim of this work is to propose outliers-robust estimators of the Wasserstein distance, and illustrate their application in generative modeling. To that end, we explore how to combine a Median-of-Means approach with Optimal Transport. The Median-of-Means (MoM) is a robust mean estimator firstly introduced in complexity theory during the 1980s . Following the seminal deviation study by , MoM has lately witnessed a surge of interest, mainly due to its attractive sub-gaussian behavior, under the sole assumption that the underlying distribution has finite variance . Originally devoted to scalar random variables, MoM has notably been extended to random vectors and UU-statistics . As a natural alternative to the empirical mean, MoM has become the cornerstone of several robust learning procedures in heavy-tailed situations, including bandits and MoM-tournaments . A more recent line of work now focuses on MoM’s ability to deal with outliers. Aside from concentration results in a contaminated context , it has yielded promising applications in robust mean embedding , and the more general MoM-minimization framework .

In this paper, we introduce and study outliers-robust estimators of the Wasserstein distance based on the MoM methodology. Our contribution is threefold:

Focusing on the Kantorovich-Rubinstein duality , we present three novel MoM-based estimators, leveraging in particular Medians of UU-statistics (MoU). In the realistic setting of contaminated data, we show their strong consistency, and provide non-asymptotic bounds as well.

We propose a dedicated algorithm to compute these three estimators in practice. Applied on a parametric family of Lipschitz functions, e.g. neural networks with clipped weights, it performs a MoM/MoU gradient descent algorithm. A sensitivity analysis of the unique parameter of these estimators is also provided throught numerical experiments on toy datasets.

We robustify WGANs (w.r.t. outliers) using a MoM-based estimator as loss function. We show the benefits of this approach through convincing numerical results on two contaminated well known benchmarks: CIFAR10 and Fashion MNIST.

Background and preliminaries

Given p∈[1,∞)p\in[1,\infty), the Wasserstein distance of order pp between two arbitrary measures μ\mu and ν\nu is defined through the resolution of the Monge-Kantorovitch mass transportation problem :

Of particular interest is the problem of estimating the Wasserstein distance between μ\mu and ν\nu given a finite number of observations. The usual assumption is to rely upon two samples X={X1,…,Xn}\mathbf{X}=\{X_{1},\ldots,X_{n}\} and Y={Y1,…,Ym}\mathbf{Y}=\{Y_{1},\ldots,Y_{m}\}, composed of i.i.d. realizations drawn respectively from μ\mu and ν\nu. The corresponding empirical distributions denoted by μ^n=(1/n)∑i=1nδXi\hat{\mu}_{n}=(1/n)\sum_{i=1}^{n}\delta_{X_{i}}, and ν^m=(1/m)∑j=1mδYj\hat{\nu}_{m}=(1/m)\sum_{j=1}^{m}\delta_{Y_{j}}. The natural questions are then: how to compute the estimator W(μ^n,ν^m)\mathcal{W}(\hat{\mu}_{n},\hat{\nu}_{m}), and does it converge towards W(μ,ν)\mathcal{W}(\mu,\nu)? In the dual formulation (2), computing W(μ^n,ν^m)\mathcal{W}(\hat{\mu}_{n},\hat{\nu}_{m}) is equivalent to replace the expectations with empirical means. The unit ball of Lipschitz functions can be replaced with a parameterized family of Lipschitz functions, more amenable for learning when W(μ^n,ν^m)\mathcal{W}(\hat{\mu}_{n},\hat{\nu}_{m}) is used as a loss function, see e.g. Wasserstein GANs . From the theoretical side, a substantial number of works have studied the convergence of W(μ^n,ν^m)\mathcal{W}(\hat{\mu}_{n},\hat{\nu}_{m}) under the i.i.d. setting described above. Statistical rates of convergence of the original OT problem are known to be slow rates with respect to the dimension dd of the input space, i.e. they are of order O(n−1/d)O(n^{-1/d}) .

2 Median-of-Means

When Wasserstein meets MoM

In this section, we investigate how MoM estimators can be leveraged to define and analyze new estimators of W(μ,ν)\mathcal{W}(\mu,\nu) that exhibit strong theoretical guarantees in presence of outliers. In order to assess robustness, we place ourselves in the realistic O∪I\mathcal{O}\cup\mathcal{I} framework, see e.g. , devoted to data contamination. In this setting, the i.i.d. assumption is relaxed, and the following assumption is instead adopted.

Starting from the expression of the dual expression (2), we observe that it can be considered with a two-fold perspective. The first one consists in considering the Wasserstein distance as the supremum of the difference between two expected values. The second one, obtained by linearity of the expectation, rather regards W(μ,ν)\mathcal{W}(\mu,\nu) as the supremum of single expected values, but taken with respect to the tuple (X,Y)(X,Y), and associated to the kernel: hϕ ⁣:(X,Y)↦ϕ(X)−ϕ(Y)h_{\phi}\colon(X,Y)\mapsto\phi(X)-\phi(Y).

Although quite elementary at first sight, this two-fold perspective gains complexity when applied to the empirical distributions μ^n\hat{\mu}_{n} and ν^m\hat{\nu}_{m}. Indeed, following the first perspective, the natural estimator obtained is the supremum of the differences between two empirical averages, while the second one leads to the supremum of 22-samples UU-statistics of degrees (1,1)(1,1) and kernels hϕh_{\phi}. So far, both points of view are strictly equivalent by linearity of the expectation and the empirical mean. However, this equivalence breaks down as soon as non-linearities are introduced, through MoM-like estimators for instance. We therefore introduce three distinct estimators of W(μ,ν)\mathcal{W}(\mu,\nu), that differ upon which estimator of Section 2.2 is used.

We define the Median-of-Means and the Median-of-UU-statistics estimators of the 11-Wasserstein distance as follows:

While WMoM\mathcal{W}_{\text{MoM}} relies on the difference between individual median blocks, WMoU-diag\mathcal{W}_{\text{MoU-diag}} considers the median over all possible combinations of blocks between X\mathbf{X} and Y\mathbf{Y}. As an intermediate step, WMoU-diag\mathcal{W}_{\text{MoU-diag}} looks after diagonal blocks only. The latter formulation is used in to derive robust mean embedding and Maximum Mean Discrepancy estimators. The theoretical analysis is made simpler by the independence between the blocks, but the estimator suffers from an increased variance due to the important loss of information, see Figure 1(c) and . It should be noticed however that WMoU-diag\mathcal{W}_{\text{MoU-diag}} enjoys a much lower computational cost in practice.

Another important question to be addressed is: how to handle the non-differentiability introduced by the median operator? Indeed, the Wasserstein distance often acts as a loss function, e.g. in generative modeling (VAEs, GANs), and optimizing through a MoM/MoU-based criterion then becomes crucial. One answer is to uses a MoM-gradient descent algorithm . It consists in performing a mini-batch gradient step based on the median block. In order to avoid local minima, authors propose shuffle the partition at each step of the descent, leading to the minimization of an expected MoM loss (w.r.t. the shuffling) that is more stable. Notice that this method goes beyond random partitions, and easily adapts to the randomized extensions discussed above.

2 Theoretical guarantees

There exist CO≥1C_{\mathcal{O}}\geq 1 and 0≤αO<10\leq\alpha_{\mathcal{O}}<1 such that nO≤CO2 nαOn_{\mathcal{O}}\leq C_{\mathcal{O}}^{2}~{}n^{\alpha_{\mathcal{O}}} and mO≤CO2 mαOm_{\mathcal{O}}\leq C_{\mathcal{O}}^{2}~{}m^{\alpha_{\mathcal{O}}}.

We start by an asymptotic result establishing the strong consistency of estimators in Definition 2. It highlights the different outlier configurations allowed through conditions on the proportions of outliers τX\tau_{\mathbf{X}} and τY\tau_{\mathbf{Y}}.

Suppose that samples X\mathbf{X} and Y\mathbf{Y} satisfy Assumptions 1 and 3. Then, choosing KX=⌈2τX n⌉K_{\mathbf{X}}=\lceil\sqrt{2\tau_{\mathbf{X}}}~{}n\rceil, it holds:

If finally τX+τY<1/2\tau_{\mathbf{X}}+\tau_{\mathbf{Y}}<1/2 and n=mn=m, then choosing KX=KY=⌈2(τX+τY) n⌉K_{\mathbf{X}}=K_{\mathbf{Y}}=\lceil\sqrt{2(\tau_{\mathbf{X}}+\tau_{\mathbf{Y}})}~{}n\rceil, it holds:

Suppose that samples X\mathbf{X} and Y\mathbf{Y} satisfy Assumption 1, and define Γ ⁣:τ↦1+2τ/1−2τ\Gamma\colon\tau\mapsto\sqrt{1+\sqrt{2\tau}}/\sqrt{1-2\tau}. Then, for all δ∈]0,exp⁡(−4n2τX)]\delta\in]0,\exp(-4n\sqrt{2\tau_{\mathbf{X}}})], choosing KX=⌈2τX n⌉K_{\mathbf{X}}=\lceil\sqrt{2\tau_{\mathbf{X}}}~{}n\rceil, it holds with probability at least 1−δ1-\delta:

If furthermore τX+τY<1/2\tau_{\mathbf{X}}+\tau_{\mathbf{Y}}<1/2 and n=mn=m, then for all δ∈]0,exp⁡(−4n2(τX+τY))]\delta\in]0,\exp(-4n\sqrt{2(\tau_{\mathbf{X}}+\tau_{\mathbf{Y}})})], choosing KX=KY=⌈2(τX+τY) n⌉K_{\mathbf{X}}=K_{\mathbf{Y}}=\lceil\sqrt{2(\tau_{\mathbf{X}}+\tau_{\mathbf{Y}})}~{}n\rceil, it holds with probability at least 1−δ1-\delta:

The proof derives from concentration results established in , combined with a generic chaining argument. It should be noticed that constant C2(τX)C_{2}(\tau_{\mathbf{X}}) explodes as τX\tau_{\mathbf{X}} goes to 1/21/2, which is expected: the more outliers, the more difficult it is to estimate W(μ,ν)\mathcal{W}(\mu,\nu). We also stress that the dependence in 1/1−2τX1/\sqrt{1-2\tau_{\mathbf{X}}} is better than the 1/(1−2τX)3/21/(1-2\tau_{\mathbf{X}})^{3/2} term exhibited in . Integrating the deviation probabilities of Proposition 5 and using Assumption 3, we finally obtain our main theorem, that provides a nonasymptotic control on the expected value of our estimators deviations from W(μ,ν)\mathcal{W}(\mu,\nu).

Suppose that samples X\mathbf{X} and Y\mathbf{Y} satisfy Assumptions 1 and 3, and recall the notation used in Proposition 5. Let β∈\beta\in, then for all nn such that n1d+2+1−β2≥C1(τX)/(2C2(τX)(2τX)14)n^{\frac{1}{d+2}+\frac{1-\beta}{2}}\geq C_{1}(\tau_{\mathbf{X}})/(2C_{2}(\tau_{\mathbf{X}})(2\tau_{\mathbf{X}})^{\frac{1}{4}}), it holds:

with κ1(τ)=C1(τ)\kappa_{1}(\tau)=C_{1}(\tau), κ2(τ)=2COC2(τ)(2/τ)1/4\kappa_{2}(\tau)=2C_{\mathcal{O}}C_{2}(\tau)(2/\tau)^{1/4}, and κ3(τ)=πC2(τ)/2\kappa_{3}(\tau)=\sqrt{\pi}C_{2}(\tau)/2.

Of course, the above bound only makes sense if β>αO\beta>\alpha_{\mathcal{O}}. In particular, if αO≤d/(d+2)\alpha_{\mathcal{O}}\leq d/(d+2), setting β=1\beta=1 gives that for all nn such that n1d+2≥C1(τX)/(2C2(τX)(2τX)14)n^{\frac{1}{d+2}}\geq C_{1}(\tau_{\mathbf{X}})/(2C_{2}(\tau_{\mathbf{X}})(2\tau_{\mathbf{X}})^{\frac{1}{4}}), with the notation κ=κ1+κ2+κ3\kappa=\kappa_{1}+\kappa_{2}+\kappa_{3}, it holds:

If furthermore τX+τY<1/2\tau_{\mathbf{X}}+\tau_{\mathbf{Y}}<1/2 and n=mn=m, then for all nn s.t. n1d+2≥C1(τX+τY)/(2C2(τX+τY)(2(τX+τY))14)n^{\frac{1}{d+2}}\geq C_{1}(\tau_{\mathbf{X}}+\tau_{\mathbf{Y}})/(2C_{2}(\tau_{\mathbf{X}}+\tau_{\mathbf{Y}})(2(\tau_{\mathbf{X}}+\tau_{\mathbf{Y}}))^{\frac{1}{4}}), with the notation κ′=2κ1+22κ2+2κ3\kappa^{\prime}=2\kappa_{1}+2\sqrt{2}\kappa_{2}+2\kappa_{3}, it holds:

The unique property of the Wasserstein distance we used, compared to other Integral Probability Metrics (IPMs) , is the way to bound the entropy of the unit ball of Lipschitz functions. The present analysis can thus be extended in a direct fashion to any other IPM that has finite entropy.

MoM-based estimators in practice

In this section, we first propose a novel algorithm to approximate the MoM/MoU-based estimators using neural networks and provide an empirical study of its behaviour on two toy datasets. Then, we show how to robustify Wasserstein-GANs and present MoMWGAN, a MoM-based variant of GAN, which is evaluated on two well-known image benchmarks.

As show in Section 3, MoM/MoU-based estimation of the Wasserstein distance offers a robust alternative to the classical empirical estimator of W\mathcal{W}. Indeed, the empirical estimator of W\mathcal{W} would not converge towards the target in the O∪I\mathcal{O}\cup\mathcal{I} framework. The proposed estimators are consistent and have convergence rates of order O(n−1/(d+2))O(n^{-1/(d+2)}) with the O∪I\mathcal{O}\cup\mathcal{I} framework.These convergence rates are similar, when dd is not too small, to those of the empirical estimator of W\mathcal{W} in a non-contaminated setting. Nevertheless, the question of computing these estimators raises two major difficulties: (i) the optimization over the unit ball of Lipschitz functions is intractable, which is a difficulty common to the approximation of the standard Wasserstein distance, and (ii) the non-differentiability of the median-based loss. The first issue is well known of the practioners of the Wasserstein distance who usually prefer to rely on its primal definition with an entropy-based regularization . However, learning algorithms devoted to Wasserstein GANs overcome this by weight clipping or gradient penalization to impose to the GAN a Lipchitz constraint. Similarly we propose here to limit Φ\Phi to be a neural network with similar constraints on weights to ensure its MM-Lipschitzianity. This enables to approximate the Wasserstein distance up to a (unknown) multiplicative coefficient MM.To overpass (ii), one can adopt MoM/MoU gradient descent. Exploited in the context of robust classification , using MoM/MoU gradient descent has been proved to be equivalent to minimize the expectation over the sampling strategy of blocks of WMoM,WMoU-diag\mathcal{W}_{\text{MoM}},\mathcal{W}_{\text{MoU-diag}} and WMoU\mathcal{W}_{\text{MoU}}. Combining these techniques, we design novel algorithms to compute approximations of the proposed estimators: W~MoM\widetilde{\mathcal{W}}_{\text{MoM}} (see Algorithm 1), W~MoU-diag\widetilde{\mathcal{W}}_{\text{MoU-diag}} and W~MoU\widetilde{\mathcal{W}}_{\text{MoU}} (see the Supplementary Material).

2 Empirical study

Two simulated datasets in 2D space with different kinds of anomalies are used in the experiments. The random vectors X1X_{1} and X2X_{2} are chosen to be distributed according a mixture of a standard Gaussian distribution and an "anomaly" distribution, respectively A1\mathcal{A}_{1} and A2\mathcal{A}_{2} defined as follows. A1\mathcal{A}_{1} is the uniform distribution U[−50,50]\mathcal{U}[\mathbf{-50},\mathbf{50}] that mimics isolated outliers while A2\mathcal{A}_{2} is the standard Cauchy distribution shifted by 25, defined to mimic aggregate outliers (see e.g. ). The random vector YY is Gaussian with Y∼N(5,I2)Y\sim\mathcal{N}(\mathbf{5},I_{2}), Datasets D1=(X1,Y)\mathcal{D}_{1}=(\mathbf{X}_{1},\mathbf{Y}) and D2=(X2,Y)\mathcal{D}_{2}=(\mathbf{X}_{2},\mathbf{Y}) contain 500 independent and identical copies of (X1,Y)(X_{1},Y), (X2,Y)(X_{2},Y) respectively, with the same proportion of outliers τX\tau_{X}.

Evaluation metrics.

The Lipchitz constant MM being unknown and highly depending of the clipping parameter choice, it wouldn’t be appropriate to compare the true 1-Wasserstein value, equal to 50\sqrt{50}, with W~MoM,W~MoU-diag\widetilde{\mathcal{W}}_{\text{MoM}},\widetilde{\mathcal{W}}_{\text{MoU-diag}} and W~MoU\widetilde{\mathcal{W}}_{\text{MoU}}. Therefore, we propose to compare W~MoM\widetilde{\mathcal{W}}_{\text{MoM}}, W~MoU-diag\widetilde{\mathcal{W}}_{\text{MoU-diag}} and W~MoU\widetilde{\mathcal{W}}_{\text{MoU}} to W~\widetilde{\mathcal{W}}, the 1-Wasserstein distance approximated by Algorithm 1, when MoM is not used, e.g. KX=KY=1K_{\mathbf{X}}=K_{\mathbf{Y}}=1, by measuring the absolute error between them.

The numbers of blocks, KXK_{\mathbf{X}} and KYK_{\mathbf{Y}}, are crucial parameters for computation. They define the trade-off between the robustness of the estimator and computational burden. However the theory does not give enough insights about their value: the necessary assumption for the consistency is only that they should be greater than 2τXn2\tau_{X}n (see Section 3.2). An empirical study of the influence of their values on the behavior of the approximations of WMoM,WMoU-diag\mathcal{W}_{\text{MoM}},\mathcal{W}_{\text{MoU-diag}} and WMoU\mathcal{W}_{\text{MoU}} is therefore much useful. For sake of simplicity, we set KX=KYK_{\mathbf{X}}=K_{\mathbf{Y}} in the subsequent experiments.

In a first experiment, we explore the ability of algorithm 1 and variants described in the supplements to override outliers according to the values of KXK_{\mathbf{X}} and with different rates of outliers τX\tau_{\mathbf{X}}. The approximations W~MoM\widetilde{\mathcal{W}}_{\text{MoM}},W~MoU-diag\widetilde{\mathcal{W}}_{\text{MoU-diag}} and W~MoU\widetilde{\mathcal{W}}_{\text{MoU}} are computed using a simple multi-layer perceptrons with one hidden layer and MoM gradient descent over several τX\tau_{\mathbf{X}} and KXK_{\mathbf{X}} on both datasets. The experiment is repeated 20 times with various seeds. Mean results are displayed. Figure 2 represents absolute deviations between the 1-Wasserstein distance approximated with a MLP when τX=0\tau_{\mathbf{X}}=0 and W~MoU-diag\widetilde{\mathcal{W}}_{\text{MoU-diag}} with various anomalies settings and different values of KXK_{X}. The reader is invited to refer to Section B of the supplements to see similar results for W~MoU\widetilde{\mathcal{W}}_{\text{MoU}} and W~MoM\widetilde{\mathcal{W}}_{\text{MoM}}). Shaded areas, in Figure 2 represent 25%-75% quantiles over the 20 repetitions. On both datasets, we observe that the approximation algorithm succeeds to provide an estimation of WMoU-diag\mathcal{W}_{\text{MoU-diag}}, able to override outliers with different τX\tau_{\mathbf{X}} while KXK_{\mathbf{X}} is high enough. From Section 3.2, we know that KXK_{\mathbf{X}} needs to be higher than 2τXn2\tau_{\mathbf{X}}n to have theoretical guarantees. Experiments show that in practice, this condition is not necessary in every situations. For example, when τX=0.1\tau_{\mathbf{X}}=0.1 (i.e. 50 anomalies) in Figure 2 (left), only 70 blocks are needed to override outliers. The reason is that hypothesis makes things work in the worst case, i.e., when each outlier is isolated in one block which lead to have τXn\tau_{\mathbf{X}}n contaminated blocks. This is rarely the case in practice, several blocks can be contaminated by many outliers and this is why fewer blocks are needed.

In a second experiment illustrated by Figure 3, we study the convergence of the approximation algorithm with and without anomalies for different values of KXK_{\mathbf{X}} on D1\mathcal{D}_{1}. To get a fair comparison between the different settings of the algorithm, we compare the predicted values across the "learning" epoch. Here during one epoch, the algorithm has made a gradient pass over the whole dataset, which means that one epoch corresponds one iteration of the approximation algorithm if KX=1K_{\mathbf{X}}=1 (no MoM estimation), and to KXK_{\mathbf{X}} iterations, in the other cases. In both cases (with or without anomalies), the higher KXK_{\mathbf{X}} is, the faster the approximation algorithm converges. Without surprise, the MoM approach benefits from the same properties than a mini-batch approach. When there is no anomalies, the distance values reached after convergence are close to the "true" value (obtained with the plain estimator when KX=1K_{\mathbf{X}}=1), especially when KXK_{\mathbf{X}} is lower. This means that the MoM-based algorithm can be used routinely instead of the plain estimator. With 5% of anomalies, one can see that distance values reached after convergence get closer to the target as KXK_{\mathbf{X}} grows.

3 Application to robust Wasserstein GANs

In this part, we introduce a robust modification of WGANs, named MoMWGAN, using one of the three proposed estimators in Section 3.

The behaviour of likelihood-free generative modeling such as Generative Adversarial Networks in the presence of outliers, i.e., with heavy-tails distributions or contaminated data, has been poorly investigated up to very recently. At our knowledge, the unique reference is . In particular, Gao et al. have studied theoretically and empirically the robustness of f-GAN in the special case of mean estimation for elliptical distributions. In contrast, we illustrate here the theoretical results of section 3 by applying a MoM approach to robustify WassersteinGAN and show on two real-world image benchmarks how this new variant of GAN behaves when learned with contaminated data.

Reminder on GAN: Let us briefly recall the GAN principle. A GAN learns a function gθ:Z→Xg_{\theta}:\mathcal{Z}\rightarrow\mathcal{X} such that samples generate by gθ(z)∼Pθg_{\theta}(z)\sim P_{\theta}, taking as input a sample zz (from some reference measure ξ\xi, often Gaussian) in a latent space Z\mathcal{Z}, are close to those of the true distribution PrP_{r} of data. Wasserstein GANs use the 1-Wasserstein Distance under its Kantorovich-Rubinstein dual formula as the loss function. Instead of maximizing over the unit ball of Lipschitz functions, one uses a parametric family of M-Lipschitz functions under the form of neural net with clipped weights ww . Following up the theoretical analysis of Section 3, we introduce a MoM-based WGAN (MoMWGAN) model, combining the WMoM\mathcal{W}_{\text{MoM}} estimator studied in 3 and WGAN’s framework. Following the weight clipping approach, MoMWGAN boils down to the problem:

Note that the MoM procedure is chosen to be only applied on the observed contaminated sample. It is not clear in which way the sample drawn from the currently learned density is polluted and thus defining the number of blocks would be an issue. Optimization in WGAN is usually performed by taking mini-batches to reduce the computational load. In the same spirit, we apply MoM inside contaminated mini-batches as described in Algorithm 4. To get the outliers-robust property observed in the numerical experiments, we pay the price of finding the median block at each step by evaluating the loss which significantly increases the computational complexity.

To test the robustness of MoMWGAN we contaminated two well-known image datasets, CIFAR10 and Fashion MNIST, with two anomalies settings. Noise based-anomalies are added to CIFAR10, i.e., images with random intensity pixels drawn from a uniform law. For Fashion MNIST, the five first classes are considered as "informative data" while the sixth (Sandal) contains the anomalies. In both settings, WGAN and MoMWGAN are trained on the training samples contaminated in a uniform fashion with a proportion of 1.5% of outliers in both datasets. Both models use standard parameters of WGAN. KX=4K_{\mathbf{X}}=4 blocks have been used by MoMWGAN in both experiments. To assess performance of the resulting GANs, we generated 50000 generated images using each model (WGAN and MoMGAN) and measured the Fréchet Inception Distance (FID) between the generated examples in both cases and the (real) test sample. Table 1 shows that MoMWGAN improves upon WGAN in terms of outliers-robustness. Furthermore, some generated images are represented in Figure 4. One can see that outliers do not affect MoMWGAN generated samples while WGAN reproduce noise on contaminated CIFAR10 dataset. For Fashion MNIST, one may see that fewer images are degraded with MoMWGAN generator.

Conclusion and perspectives

In this paper, we have introduced three robust estimators of the Wasserstein distance based on MoM methodology. We have shown asymptotic and non-asymptotic results in the context of polluted data, i.e. the O∪I\mathcal{O}\cup\mathcal{I} framework. Surpassing computational issues, we have designed an algorithm to compute, in a efficient way, these estimators. Numerical experiments have highlighted the behavior of these estimators over their unique tuning parameter. Finally, we proposed to robustify WGANs using one of the introduced estimators and have shown its benefits on convincing numerical results. The theoretically well-founded MoM approaches to robustify the Wasserstein distance open the door to numerous applications beyond WGAN, including variational generative modeling. The promising MoMGAN deserves more attention and future work will concern the analysis of the estimator it provides.

Acknowlegments

The authors thank Pierre Colombo for his helpful remarks. This work has been funded by BPI France in the context of the PSPC Project Expresso (2017-2021).

References

A Technical Proofs

In this section are detailed the proofs of the theoretical claims stated in the core article. We first recall a simple lemma on the difference between two median vectors.

Thus, for all b\bm{b} within the infinite ball of center a\bm{a} and radius ϵ\epsilon it holds:

where BmedX\mathcal{B}^{\mathbf{X}}_{\text{med}} and BmedY\mathcal{B}^{\mathbf{Y}}_{\text{med}} are the median blocks of ϕ‾X,k−ϕ‾Y,l\overline{\phi}_{\mathbf{X},k}-\overline{\phi}_{\mathbf{Y},l} for 1≤k≤KX1\leq k\leq K_{\mathbf{X}} and 1≤l≤KY1\leq l\leq K_{\mathbf{Y}}. From Sections A.1 and A.1, we deduce that:

where we have used the fact that IX×IY\mathcal{I}_{\mathbf{X}}\times\mathcal{I}_{\mathbf{Y}} represents a majority of blocks, and the subadditivity of the supremum. By independence between samples X and Y, and between the blocks, it holds:

Now, the arguments to get the right-hand side equal to 11 are similar to those used in Lemma 3.1 and Proposition 3.2 in . We expose them explicitly for the sake of clarity.

Therefore F(x)F(x) is finite, and following Lemma 3.1. in we have

Since H(ε,BL,L1(μ^n))≤H(ε,BL,∥⋅∥∞)\mathcal{H}(\varepsilon,\mathcal{B}_{L},L^{1}(\hat{\mu}_{n}))\leq\mathcal{H}(\varepsilon,\mathcal{B}_{L},\|\cdot\|_{\infty}) and H(ε,BL,L1(ν^m))≤H(ε,BL,∥⋅∥∞)\mathcal{H}(\varepsilon,\mathcal{B}_{L},L^{1}(\hat{\nu}_{m}))\leq\mathcal{H}(\varepsilon,\mathcal{B}_{L},\|\cdot\|_{\infty}) then when, respectively, nn and mm go to infinity, we have

It is then direct to adapt the reasoning from Equation 6. ∎

A.2 Proof of Proposition 5

Let ψ∈BL\psi\in\mathcal{B}_{L}. From Equation 7, we know that −diam(K)≤ψ(X)≤diam(K)-\text{diam}(\mathcal{K})\leq\psi(X)\leq\text{diam}(\mathcal{K}), so that ψ(X)\psi(X) is in particular sub-Gaussian with parameter λ=diam(K)\lambda=\text{diam}(\mathcal{K}). A direct application of Proposition 1 in then gives that for all δ∈]0,e−4n2τX]\delta\in]0,e^{-4n\sqrt{2\tau_{\mathbf{X}}}}] and KX=⌈2τXn⌉K_{\mathbf{X}}=\lceil\sqrt{2\tau_{\mathbf{X}}}n\rceil , it holds with probability at least 1−δ1-\delta:

with Γ ⁣:τX↦1+2τX/1−2τX\Gamma\colon\tau_{\mathbf{X}}\mapsto\sqrt{1+\sqrt{2\tau_{\mathbf{X}}}}/\sqrt{1-2\tau_{\mathbf{X}}}. Using 8, observe also that ∀(ϕ,ψ)∈BL2\forall(\phi,\psi)\in\mathcal{B}_{L}^{2} it holds:

Now, let ζ>0\zeta>0, and ψ1,…,ψN(ζ,BL,∥⋅∥∞)\psi_{1},\ldots,\psi_{\mathcal{N}(\zeta,\mathcal{B}_{L},\|\cdot\|_{\infty})} be a ζ\zeta-coverage of BL\mathcal{B}_{L} with respect to ∥⋅∥∞\|\cdot\|_{\infty}. We know from that there exists CL>0C_{L}>0 such that for all ζ>0\zeta>0 it holds:

From now on, we use N=N(ζ,BL,∥⋅∥∞)\mathcal{N}=\mathcal{N}(\zeta,\mathcal{B}_{L},\|\cdot\|_{\infty}) for notation simplicity. Let ϕ\phi be an arbitrary element of BL\mathcal{B}_{L}. By definition, there exists i≤Ni\leq\mathcal{N} such that ∥ϕ−ψi∥∞≤ζ\|\phi-\psi_{i}\|_{\infty}\leq\zeta. Section A.2 then gives:

Applying Equation 8 to every ψi\psi_{i}, the union bound gives that with probability at least 1−δ1-\delta it holds:

Taking the supremum in both sides of Equation 11, it holds with probability at least 1−δ1-\delta:

Choosing ζ∼1/n1/(d+2)\zeta\sim 1/n^{1/(d+2)} and breaking the square root finally gives that it holds with probability at least 1−δ1-\delta:

with C1(τX)=2+CLC2(τX)C_{1}(\tau_{\mathbf{X}})=2+C_{L}C_{2}(\tau_{\mathbf{X}}), and C2(τX)=4 diam(K) Γ(τX)C_{2}(\tau_{\mathbf{X}})=4~{}\text{diam}(\mathcal{K})~{}\Gamma(\tau_{\mathbf{X}}).

Adaptation to MoU. From Equation 7, we get that the kernel hϕ ⁣:(X,Y)↦ϕ(X)−ϕ(Y)h_{\phi}\colon(X,Y)\mapsto\phi(X)-\phi(Y) has finite essential supremum ∥hϕ(X,Y)∥∞≤diam(K)\|h_{\phi}(X,Y)\|_{\infty}\leq\text{diam}(\mathcal{K}). Using Proposition 4 in with the same reasoning as above leads to the desired result, multiplying constants by a 22 factor. ∎

A.3 Proof of Theorem 7

Since n1d+2+1−β2≥C1(τX)/(2C2(τX)(2τX)14)n^{\frac{1}{d+2}+\frac{1-\beta}{2}}\geq C_{1}(\tau_{\mathbf{X}})/(2C_{2}(\tau_{\mathbf{X}})(2\tau_{\mathbf{X}})^{\frac{1}{4}}), then for all δ∈]0,e−4n2τX]\delta\in]0,e^{-4n\sqrt{2\tau_{\mathbf{X}}}}], it holds:

Combining with the first results of Proposition 4, for all δ∈]0,e−4n2τX]\delta\in]0,e^{-4n\sqrt{2\tau_{\mathbf{X}}}}], it holds with probability at least 1−δ1-\delta:

Reverting the inequation gives that it holds

One may finally use that for a nonnegative random variable it holds:

Where the second line holds thanks to Assumption 6.

Adaptation to MoU. The adaptation is straightforward, up to Equation 13, that now writes:

Using Assumption 6 on both samples X and Y, it leads to the desired results. ∎

B Additional material of the numerical part

In this part, we introduce algorithms and additional experiments that could not be in the paper for lack of space.

Here, algorithms to compute WMoU-diag(μn,νn)\mathcal{W}_{\text{MoU-diag}}(\mu_{n},\nu_{n}) and WMoU(μn,νn)\mathcal{W}_{\text{MoU}}(\mu_{n},\nu_{n}) are displayed.

B.2 Additional experiments

In this part, numerical results for W~MoU\widetilde{\mathcal{W}}_{\text{MoU}} and W~MoM\widetilde{\mathcal{W}}_{\text{MoM}}, related to the Section 4.2 of the paper, are displayed. Results of both experiments, depicted in Figure 5 and 6, are quite similar due to the simplicity of the problem.