1 Introduction↩︎

Let \(f,g\) be two probability densities over \(\mathbb{R}^d\), supported on the same set. For any measurable function \(\varphi\), one can write \[\label{eq:is} \mathbb{E}_{Y \sim g}[\varphi(Y)] = \mathbb{E}_{X \sim f}\left[\varphi(X) \frac{g(X)}{f(X)}\right].\tag{1}\] In particular, for any \(\varphi\), if we have samples \(X_i\) from \(f\), we can approximate the expectation \(\mathbb{E}_{Y \sim g}[\varphi(Y)]\) by the empirical average \[\label{eq:is950} \frac{1}{n}\sum_{i=1}^n \varphi(X_i) \frac{g(X_i)}{f(X_i)}.\tag{2}\] This random variable has mean \(\mathbb{E}_{Y \sim g}[\varphi(Y)]\), and if the samples are independent, by the law of large numbers, it converges to this expectation as \(n\to \infty\), which gives an estimator of the integral \(\mathbb{E}_{Y \sim g}[\varphi(Y)]\) using only samples from \(f\). However, 2 implicitly assumes that we can evaluate the density \(g\) itself. In many applications, the target is only available up to a multiplicative constant, and the normalization constant \(\int g\) could be impossible to compute. This is the usual situation for unnormalized measures, e.g. in Bayesian inference where the posterior is proportional to likelihood times prior. In that case, the weights \(g(X_i)/f(X_i)\) in 2 cannot be computed.

Self-normalized importance sampling resolves this issue by replacing 2 with a ratio in which any unknown multiplicative constant cancels, \[\label{eq:is95sn} \frac{\sum_{i=1}^n \varphi(X_i) g(X_i)/f(X_i)}{\sum_{i=1}^n g(X_i)/f(X_i)}.\tag{3}\] While this estimator is in general biased, the law of large numbers for 2 still ensures that the numerator converges to \(\int \varphi g\) and the denominator converges to \(\int g\), hence the estimator converges to \(\mathbb{E}_{X \sim g}[\varphi(X)]\) even if \(g\) is only known up to a normalization constant.

The estimator 3 can be seen as the integral of \(\varphi\) against the self-normalized importance sampling measure \(\hat{g}_n\), which motivates the following definition. From now on, we will restrict to the cube \((0,1)^d\).

Importance sampling and its variants are used in a large variety of contexts, and it is crucial to understand the accuracy of these methods and how it behaves depending on the sample size \(n\). Many works were devoted to this topic, as we will review in the next section. It is generally acknowledged that the accuracy of 3 depends on the distribution and concentration of the importance weights \(W_i = g(X_i)/f(X_i)\): if the sampling density \(f\) and the target density \(g\) are equal, these weights are all equal to 1. Conversely, any significant discrepancy between \(f\) and \(g\) will lead to a poor concentration of the importance weights (they could have large variance, or be heavy-tailed), resulting in a very poor accuracy of 3 .

In 2018, a landmark paper by Chatterjee and Diaconis [1] studied how \(\int \varphi d\hat{g}_n\) is close to \(\int \varphi dg\) for a fixed function \(\varphi\) with unit square-norm. They essentially showed that there is a cutoff at the sample size \(n = e^{\mathrm{d}_{\mathrm{KL}}(g\mid f)}\), where \(\mathrm{d}_{\mathrm{KL}}(g|f)\) is the Kullback-Leibler divergence between \(f\) and \(g\), in the sense that if \(n \gg e^{\mathrm{d}_{\mathrm{KL}}(g\mid f)}\), then \[\mathbb{P}\left(\left|\int \varphi d\hat{g}_n - \int \varphi dg\right|>\varepsilon\right) = O(\varepsilon),\] while if \(n \ll e^{\mathrm{d}_{\mathrm{KL}}(f|g)}\), then one can find (at least) one square-integrable function \(\varphi\) such that \(\int \varphi d\hat{g}_n = 1\) with overwhelming probability, but \(\int \varphi dg \leqslant \varepsilon\). This is often interpreted as a curse of dimensionality phenomenon, because in high dimension, the divergence \(\mathrm{d}_{\mathrm{KL}}(g\mid f)\) can grow linearly with \(d\) in many contexts (see [2] for a prototypical example).

In this work, we are interested in a more quantitative study of the estimation error using the self-normalized importance sampling measure. Our results are global in the sense that we directly measure the \(p\)-Wasserstein distance \(\mathscr{W}_p(\hat{g}_n, g)\) for \(p\geqslant 1\), and prove that it has order \(n^{-1/d}\).

2 Main result↩︎

2.1 Asymptotic behavior of the Wasserstein distance↩︎

The \(p\)-Wasserstein distance is defined for pairs of Borel measures \(\mu\) and \(\nu\) on \((0,1)^d\) by the formula \[\label{def:Wp}\mathscr{W}_p^p(\mu, \nu)= \inf_{\Pi \in \mathscr{C}(\mu, \nu)}\int |x - y|^p \mathrm{d}\Pi(x, y),\tag{4}\] where \(\mathscr{C}(\mu, \nu)\) is the set of couplings between \(\mu\) and \(\nu\). A coupling between \(\mu\) and \(\nu\) is a measure \(\Pi\) on the product space with marginals \(\mu\) and \(\nu\). Note that we did not impose any normalization on \(\mu\) and \(\nu\), which may or may not be probability measures—but they must have the same mass.

2.2 Remarks↩︎

 

2.2.0.1 A

The constants \(\beta_{p,d}\) and \(\beta^{\mathrm{low}}_{p,d}\) are respectively defined in 14 and 29 . It is conjectured that they are equal. This was only proved in the case where \(p=2\) and \(d\geqslant 3\) in [3].

2.2.0.2 B

Our proof has the decisive advantage of being valid for the whole range of \(p\geqslant 1\) and \(d\geqslant 3\). Other simpler proofs (like the one in [4]) could work, but only for tighter ranges like \(p\leqslant d/2\).

2.2.0.3 C

The hypothesis that the common support of \(f\) and \(g\) is the cube \((0,1)^d\) could be relaxed to cover compact sets with Lipschitz boundary at the price of some technicalities, see [5].

2.2.0.4 D

When \(g=f\), this theorem reduces to [5], where it is proved that \(\mathscr{W}^p_p(\hat{g}_n, g)\) is of order \(n^{-p/d}\beta_{p,d}\int g^{1 - p/d}\). The extra term \(\int g^{1 - p/d}\) is equal to 1 if \(g\equiv 1\), hence this term accounts for the non-uniformity of \(g\). Since \(x\mapsto x^{1 - p/d}\) is convex if \(p\geqslant d\) and concave if \(p < d\), Jensen’s inequality yields \[\label{eq:g941-p47d} \int_{(0,1)^d} g(x)^{1 - p/d} \leqslant 1 \quad \text{ iff } p < d.\tag{5}\] The interpretation of 5 is that, when \(p < d\) and \(g\) is nonuniform, the IS error is smaller than if \(g\) were uniform. The regime \(p< d\) is the typical one in high dimension; in practice, one often uses \(p=1\) or \(p=2\).

2.3 Asymptotically optimal sampling distribution↩︎

It might sound natural that, when approximating \(g\) by IS, using reweighted samples from a distribution \(f\neq g\), the error should be higher, i.e. the extra term \(\int g f^{-p/d}\) should be greater than \(\int g^{1 - p/d}\). Again, this is not true in general. Suppose for example that the sampling distribution is uniform, ie \(f=1\). Then, 5 shows that \(\int g^{1 - p/d} > 1=\int g f^{-p/d}\) if \(p \geqslant d\). In other words, the IS error using real samples from \(g\) is greater than the IS error using uniform samples from \(g\). So it is only natural to ask, for a fixed \(g\), what is the sampling distribution \(f\) that minimizes the error \(\mathscr{W}_p(\hat{g}_n, g)^p\). Our results do not provide the exact asymptotic of \(\mathscr{W}_p(\hat{g}_n, g)^p\) as \(n \to \infty\), because the two constants \(\beta_{p,d}\) and \(\beta^{\mathrm{low}}_{p,d}\) are not known to be equal — but they are generally conjectured to be equal. Under this conjecture, we would have \[\mathscr{W}_p^p(\hat{g}_n, g) \sim n^{-p/d}\beta_{p,d}\int_{(0,1)^d} g(x) f(x)^{-p/d}\mathrm{d} x\] and it is natural to minimize this asymptotic equivalent with respect to \(f\). The following proposition gives the minimizer.

It is noticeable that if \(g\) is not uniform, then \(f^* \neq g\). Moreover, since \(g\) is admissible, we necessarily have \(J_g(f^*) \leqslant J_g(g)\). This result might feel counter-intuitive: when using importance sampling to approximate \(g\), the best sampling distribution is not \(g\) itself if we measure the error by the transport cost with \(\mathscr{W}^p_p\). For example, in 3d problems where the error is measured in the euclidean distance \((p=2)\), the best sampling distribution is (proportional to) \(g^{3/5}\), which is very different from \(g\).

2.4 Proof outline↩︎

The proof of this theorem spans Sections 5 and 6, and adapts the general strategy developed in a series of papers studying optimal matchings between point sets and measures, [8],[9], [5]. Appendices 8 and 9 contain classical results on concentration inequalities and functional inequalities used in the proof.

Acknowledgements↩︎

The authors thank Dario Trevisan for useful discussions and suggestions. Part of this research was supported by the ANR project STOIQUES. This work was supported by a public grant from the Fondation Mathématique Jacques Hadamard.

3 Related work and perspectives↩︎

3.1 Literature review↩︎

Importance sampling (IS) is a classical Monte Carlo device, tracing back at least to the early work of Kahn [10], and now ubiquitous in statistics and applied probability. We refer to the surveys by Agapiou et al. [11] or Elvira et al. [7] for a broad overview.

A large body of work on IS aims at stabilizing weights or diagnosing when stabilization is needed. Truncated importance sampling [12] and Pareto-smoothed importance sampling [13] are widely used practical fixes, and recent work studies regimes with unbounded weight functions, like [14]. IS has become a core tool in variational inference, see for example the importance-weighted VAE [15] or the recent work [16]. IS is also used in the training or sampling of generative models in machine learning, see [2], [17]–[20] for recent examples.

These works typically measure accuracy through function-by-function errors (variance, CLT, concentration for a fixed test function) or through scalar summaries of the weights (ESS and its variants), which are indispensable in practice but do not directly quantify how close the random measure \(\hat{g}_n\) is to the target distribution as an object. In general, there are few papers on a rigorous study of the accuracy of IS beyond the landmark paper [1], and more recently [21] but with a focus on the normalization constant estimation in a dynamic setting (like AIS or using Jarzynski’s inequality). By treating self-normalized IS as a weighted empirical measure and analyzing it in \(p\)-Wasserstein distance, we obtain sharp scaling laws for \(\mathscr{W}_p(\hat{g}_n,g)\). This provides a genuinely different guarantee: it controls the full output of IS as a probability measure (for example, if \(p=1\), all Lipschitz test functions simultaneously), and it makes explicit how the proposal \(f\) enters through a transport-relevant functional of \(f\) and \(g\).

Our contribution is to bridge IS with the geometric theory of empirical measures in optimal transport. Rates for Wasserstein distances between empirical measures and their limits are connected to optimal matching problems, starting from Ajtai–Komlós–Tusnády [22] and further developed in, e.g., [4], [5], [5], [9], [23]–[27]. Another line of research focused on general metric (Polish) spaces, see [28]–[31] ; but we are not aware of any work on importance sampling or weighted empirical measures in Wasserstein distance.

3.2 Future work↩︎

3.2.1 Unbounded weights↩︎

The main limitation of our work lies in the restriction that \(f\) and \(g\) have full support on \((0,1)^d\) and are bounded from above and below. This is a strong assumption: in this case the weights \(W = g(X)/f(X)\) are bounded and cannot be heavy-tailed, for example. But their distribution is crucial to assess the performance of an importance sampling strategy. In fact, \(\mathbb{E}_{X \sim f}[W] = 1\), but the variance is given by \[\mathbb{E}_{X \sim f}[W^2] - 1 = \int \frac{g(x)^2}{f(x)}dx -1 = \mathbb{E}_{X \sim g}[W] - 1.\] This quantity can be infinite, leading to a high variability of 2 , which would be dominated by a few large weights. Unfortunately, our method does not apply to this case, and extending it to the case where \(f,g\) are allowed to have unbounded or zero values is hard. See however [3], [9] for some results in this direction.

3.2.2 Annealed importance sampling↩︎

Annealed Importance Sampling [32] is a popular refinement of IS. The end goal is to estimate means of functionals with respect to a density \(g=g_0\), starting from samples \(X_i\) from a density, say, \(g_n\). Given a chain of densities \(g_0, g_1, \ldots, g_n\), typically \(g_k\) is proportional to \(g_0^{\beta_k}g_n^{1-\beta_k}\) for some \(1=\beta_1 > \beta_2 > \ldots > \beta_n=0\), the annealed importance sampling estimator is defined as follows: starting at \(X_i^{(n)}=X_i\), one samples \(X_i^{(n-1)}\). This is typically done by running a Markov chain starting from \(X_i^{(n)}\), whose invariant measure is \(g_{n-1}\). The chain is run for several steps; the resulting sample \(X_i^{(n-1)}\) might not be exactly distributed according to \(g_{n-1}\). Then one samples \(X_i^{(n-2)}\) by running a Markov chain starting from \(X_i^{(n-1)}\), whose invariant measure is \(g_{n-2}\), and so on, until one reaches \(X_i^{(0)}\). The annealed importance sampling measure is then defined as \[\mu = \frac{1}{n}\sum_{i=1}^n W_i \delta_{X_i^{(0)}}\] where \[W_i = \frac{g_{n-1}(X_i^{(n-1)})}{g_n(X_i^{(n-1)})}\times \frac{g_{n-2}(X_i^{(n-2)})}{g_{n-1}(X_i^{(n-2)})}\times \cdots \times \frac{g_0(X_i^{(0)})}{g_1(X_i^{(0)})}.\] One difference with classical importance sampling is that the importance weights do not only depend on \(X_i^{(0)}\), the points effectively sampled at the end of the chain. To study AIS, we would need a version of our main result which is suitable for random weights \(W\) which only satisfy \[\mathbb{E}[W \mid X] = \frac{g(X)}{f(X)}.\] Our preliminary work on this general problem suggests that the constant \(\beta_{p,d}\) in our main result would now depend on the whole distribution of \(W\) in a non-trivial way.

3.2.3 Girsanov change-of-measure↩︎

Up to now, we have mostly been interested in the high but finite-dimensional setting; but the problem of finding estimates for infinite-dimensional quantities is also of interest. Let \((B_n)\) be a sequence of independent Brownian motions from \([0,1]\) to \(\mathbb{R}^d\). We equip \(C=\mathscr{C}([0,1],\mathbb{R}^d)\) with the topology of uniform convergence. The Wasserstein distance between Gaussian measures on \(C\) and the empirical measure \(\mu_n = n^{-1}\sum_{i=1}^n \delta_{B_i}\) was studied in¹ [29] (see also references therein), although the authors only seem to prove upper bounds. But in applied stochastic calculus, many problems revolve around importance sampling through the Girsanov identity. In its simplest form, this theorem states that if \(B\) is a Brownian motion and \((h_t)_{t \in [0,1]}\) is a process adapted to the filtration generated by \(B\), then for any functional \(\varphi\), \[\mathbb{E}[\varphi(B)Z_T] = \mathbb{E}[\varphi(X)]\] where \(Z_T\) is the terminal value of the “change of measure" process \[Z_t = \exp\!\left(-\int_0^t h_s\,dB_s - \frac{1}{2}\int_0^t h_s^2\,ds\right),\] and \(X = (X_t)_{t \in [0,1]}\) is the diffusion process defined by \(dX_t = dB_t + h_t dt\). This allows one to compute observables from \(X\) using only samples from \(B\), thus falling exactly in the framework of importance sampling. We plan on extending our results to this setting in the near future.

4 Background and notation↩︎

Notation↩︎

For any Borel set \(S\), we denote by \(|S|\) its volume and \(\mathbf{1}_S\) its indicator function. For any real quantities \(A,B\), the notation \(A \lesssim B\) or equivalently \(A = O(B)\) means that there exists a constant \(C>0\) (possibly depending on \(p\) and \(d\) only) such that \(A \leqslant C B\). We write \(A \sim B\) if \(A \lesssim B\) and \(B \lesssim A\).

For a measure \(\mu\) and a set \(S\), we denote by \(\mu{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}S\) the restriction of \(\mu\) to \(S\).

Elementary inequalities↩︎

Throughout the proof, we will repeatedly use the two elementary inequalities: \[\label{ineq:elementary}\forall a, b\geqslant 0, \quad\forall p\geqslant 1, \qquad (a+b)^p \leqslant (1+\varepsilon)a^p + c_p\varepsilon^{1-p}b^p\tag{7}\] for some constant \(c_p\); and for all \(a_1, \ldots, a_k \geqslant 0\), \[\label{ineq:elementary95sum}\forall p\geqslant 1, \qquad (a_1 + \dots + a_k)^p \leqslant k^p (a_1^p + \dots + a_k^p).\tag{8}\]

The Wasserstein distance and its properties↩︎

For two measures \(\mu\) and \(\nu\) on the same measurable space \(S\) and with the same mass, the Wasserstein distance of order \(p\) is defined by \[\mathscr{W}_p(\mu, \nu)= \left(\inf_{\Pi \in \mathscr{C}(\mu, \nu)}\int_{S \times S} |x - y|^p \mathrm{d}\Pi(x, y)\right)^{1/p},\] where \(\mathscr{C}(\mu, \nu)\) is the set of couplings between \(\mu\) and \(\nu\). A coupling between \(\mu\) and \(\nu\) is a measure \(\Pi\) on the product space \(S \times S\) with marginals \(\mu\) and \(\nu\). Note that we did not impose any normalization on \(\mu\) and \(\nu\), which may or may not be probability measures—but they must have the same mass. The quantity \(\mathscr{W}_p^p(\mu, \nu)\) is called the \(p\)-Wasserstein cost between \(\mu\) and \(\nu\), and it is linear in the mass of \(\mu\) and \(\nu\) in the sense that \[\label{wasserstein95linearity} \mathscr{W}_p^p(\alpha \mu, \alpha \nu) = \alpha \mathscr{W}_p^p(\mu, \nu).\tag{9}\]

If \(S\) is bounded, then clearly the term \(|x-y|^p\) is bounded by \(\mathrm{diam}(S)^p\), so that \(\int_{S \times S} |x - y|^p \mathrm{d}\Pi(x, y) \leqslant \mathrm{diam}(S)^p \Pi(S \times S)\), and since \(\Pi\) is a coupling between \(\mu\) and \(\nu\), we have \(\Pi(S \times S) = \mu(S) = \nu(S)\). Consequently, \[\label{wasserstein95bound95bounded95domain} \mathscr{W}_p^p(\mu, \nu) \leqslant \mathrm{diam}(S)^p \mu(S).\tag{10}\]

If \(R\) is a subset of \(S\) such that \(\mu(R) = \nu(R)\), we denote by \(\mathscr{W}_{p,R}(\mu, \nu)\) the Wasserstein distance of order \(p\) between \(\mu{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}R\) and \(\nu{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}R\). If \(\mathscr{R}\) is a collection of disjoint subsets of \(S\) such that \(\cup_{R \in \mathscr{R}} R = S\) and such that \(\mu(R) = \nu(R)\) for all \(R \in \mathscr{R}\), then for any collection of couplings \(\Pi_R\) between \(\mu{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}R\) and \(\nu{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}R\) for \(R \in \mathscr{R}\), the measure \(\Pi = \sum_{R \in \mathscr{R}} \Pi_R\) is a coupling between \(\mu\) and \(\nu\). Consequently, \[\mathscr{W}_{p,S}^p(\mu, \nu)\leqslant \sum_{R \in \mathscr{R}} \mathscr{W}_{p,R}^p(\left.\mu\right|_R, \left.\nu\right|_R).\] This property is called subadditivity of the Wasserstein distance.

5 Proof of the upper bound↩︎

5.1 The setup↩︎

Let \(g\) be a probability density on \((0,1)^d\), and let \(X\) be a random variable on \((0,1)^d\) with probability density \(f\). We set \[\label{def:W}W = \frac{g(X)}{f(X)}.\tag{11}\] Then, for any test function \(\varphi\), \[\mathbb{E}_{X \sim g}[\varphi(X)] = \mathbb{E}_{X \sim f}[\varphi(X) W].\] The density \(f\) of \(X\) will often be called the sampling density while \(g\) is called the target density.

Several of our estimates are formulated in terms of the relative second moment of \(W\), denoted \[\label{def:v}v_A = \frac{\mathrm{Var}(W \mathbf{1}_{X \in A})}{\mathbb{E}[W \mathbf{1}_{X \in A}]^2}\qquad \forall A \subset (0, 1)^d.\tag{12}\] In this paper, we only use very crude bounds on \(v_A\). Indeed, since we assumed in Hypothesis 1 that \(W \leqslant c\), it is easily seen that \(v_A \leqslant \frac{c}{\int_A g}\); and since \(g \geqslant c^{-1}\), this is smaller than \(c^2 |A|^{-1}\): \[\label{bound95v}\forall A, \qquad v_A \lesssim|A|^{-1}.\tag{13}\]

Let \(n\) be an integer and \((X_i)_{i \in [n]}\) be a sequence of independent random variables with density \(f\), with their weights \(W_i = g(X_i)/f(X_i)\). We set \[\mu_n = \sum_{i=1}^n W_i \delta_{X_i}\] the empirical measure of the \(X_i\), weighted by the \(W_i\), and \[Z_n = \sum_{i=1}^n W_i\] the normalization constant, so that \(Z_n^{-1}\mu_n\) is a probability measure. When \(g=f\), the weights \(W_i\) are all equal to 1. Under the additional hypothesis that \(g\equiv 1\), the main result of [8] is the existence of the following limit, which will be our definition of \(\beta_{p,d}\): \[\label{def:beta}\limsup_{n \to \infty} n^{p/d} \mathbb{E}\mathscr{W}_p^p\left(\frac{1}{n}\sum_{i=1}^n \delta_{X_i}, 1\right) = \beta_{p,d}.\tag{14}\] Our goal is to prove that, under Hypothesis 1, \[\label{eq:main950} \limsup_{n \to \infty} n^{p/d} \mathbb{E}\mathscr{W}_p^p\left({Z_n}^{-1}\mu_n, g\right) \leqslant \beta_{p,d}\int_{(0,1)^d} g(x) f(x)^{-p/d} \mathrm{d} x.\tag{15}\] We will first prove a non-normalized version of 15 : \[\label{eq:main951} \limsup_{n \to \infty} n^{p/d-1} \mathbb{E}\mathscr{W}_p^p\left(\mu_n, Z_n g\right) \leqslant \beta_{p,d}\int_{(0,1)^d} g(x) f(x)^{-p/d} \mathrm{d} x.\tag{16}\] We will then deduce 15 from 16 .

5.2 Dyadic decomposition↩︎

We divide the cube \(Q^0=(0,1)^d\) dyadically into \(2^{dk}\) dyadic cubes \(Q^k_i\) with side length \[\label{def:l}\ell_k=2^{-k}\tag{17}\] and we note \(\mathscr{D}_k = \{Q^k_i : i \in [2^{dk}]\}\) the set of these cubes. Note that the diameter of any \(Q \in \mathscr{D}_k\) is \[\label{diam} \mathrm{diam}(Q) = \sqrt{d}\,2^{-k} \sim \ell_k.\tag{18}\] We will study how the points \(X_i\) are distributed in these cubes and how this distribution changes when we change the scale \(k\). For this, we note \[\label{def:nu}\nu_n = \sum_{i=1}^n \delta_{X_i}\tag{19}\] the unweighted empirical measure of the \(X_i\). For every measurable set \(Q\subset (0,1)^d\), we define \[\begin{align} \label{def:kappa95gamma95delta} &\kappa_Q=\frac{\mu_n(Q)}{\int_Q g}, &&\gamma_Q = \frac{\mu_n(Q)}{\nu_n(Q)} &&\delta_Q = \frac{\mu_n(Q)}{\int_Q f}. \end{align}\tag{20}\]

In the case where \(\nu_n(Q)=0\), then we also have \(\mu_n(Q)=0\) and we set \(\gamma_Q=0\). For any scale \(k\), we define the measure \[\lambda^k=\sum_{Q \in \mathscr{D}_k} \kappa_Q \mathbf{1}_Q g.\] Note that, since \(\kappa_{Q^0} = Z_n\), we clearly have \(\lambda^0 = Z_n g\). We will go down to a maximal scale \(K\) defined by \[\label{def:K} K = \lfloor d^{-1} \log_2 n \rfloor - \omega \quad \qquad \omega = \lfloor\log_2 ((\log n)^\alpha)\rfloor\tag{21}\] for some \(\alpha > 0\) — it will be clear later on that any \(\alpha \geqslant 3\) works. For simplicity, we could choose \(\alpha=3\). The side length \(\ell_K\) of any cube \(Q \in \mathscr{D}_K\) is thus equal to \[\label{eq:ell95K} \ell_K = 2^{-K} \sim \frac{2^{\omega}}{n^{1/d}} = \frac{(\log(n))^{\alpha}}{n^{1/d}}\tag{22}\] and the volume of any cube \(Q \in \mathscr{D}_K\) is equal to \[\label{eq:vol95DK} |Q| = \ell_K^d \sim \frac{(\log(n))^{\alpha d}}{n}.\tag{23}\]

5.3 Reduction to the dyadic decomposition at the maximal scale↩︎

This lemma, whose proof is postponed to Section 7, reduces the problem to the study of local expectations like \(\mathbb{E}\left[\mathscr{W}_{p,Q}^p(\mu_n , \kappa_{Q} g)\right]\) for any fixed dyadic cube \(Q\), which is what we do now. We will simply write \(\gamma, \delta, \kappa\) instead of \(\gamma_Q, \delta_Q, \kappa_Q\). On the dyadic cube \(Q\), we will compare \(\mu_n\) with the measures \(\gamma \nu_n, \delta f\) and \(\kappa g\) which all have the same mass \(\mu_n(Q)\) on \(Q\). By the triangle inequality, \[\begin{align} \mathscr{W}_{p,Q}^p(\mu_n , \kappa g)\le \left( \mathscr{W}_{p,Q}(\mu_n , \gamma \nu_n)+ \mathscr{W}_{p,Q}(\gamma \nu_n, \delta f)+ \mathscr{W}_{p,Q}( \delta f,\kappa g)\right)^p. \end{align}\] It will appear thereafter that the leading contribution is the ‘flattened’ term \(\mathscr{W}_{p,Q}^p(\gamma \nu_n, \delta f)\), which transports a point measure with all weights equal to \(\gamma\) towards the reference measure \(\delta f\). Using 7 and averaging, we have \[\begin{align} \nonumber \mathbb{E}[\mathscr{W}_{p,Q}^p(\mu_n , \kappa g)] &\leqslant (1+\varepsilon) \mathbb{E}[\mathscr{W}_{p,Q}^p(\gamma \nu_n, \delta f)]+ \frac{c_p}{\varepsilon^{p-1}}\left( \mathbb{E}[\mathscr{W}_{p,Q}^p(\mu_n , \gamma \nu_n)]+ \mathbb{E}[\mathscr{W}_{p,Q}^p( \delta f,\kappa g)]\right)\\ &=: (1+\varepsilon)(\mathrm{I}) + \frac{c_p}{\varepsilon^{p-1}}[(\mathrm{II}) + (\mathrm{III})].\label{threeterms} \end{align}\tag{24}\] We bound each of these three terms separately, and uniformly in \(Q \in \mathscr{D}_K\). The proofs are again postponed to Section 7.

Due to Hypothesis 1, the term \(\int_Q gf^{-p/d}\) has order \(|Q|=(\log n)^{\alpha d}/n\), so the right-hand side is of the order of \(\gtrsim n^{ - p/d}\log(n)^{\alpha d}\). The next lemmas show that \((\mathrm{II})\) and \((\mathrm{III})\) are negligible at this scale.

5.4 Proof of the main theorem↩︎

These lemmas together with 24 imply that for \(n\) large enough, \[\forall Q \in \mathscr{D}_K, \quad \mathbb{E}[\mathscr{W}_{p,Q}^p(\mu_n , \kappa g)] \leqslant (1+\varepsilon)^{4+p} \beta_{p,d}n^{1 - p/d} \int_{Q} g f^{-p/d}.\] Summing over all \(Q \in \mathscr{D}_K\) as in the right-hand side of ?? then gives \[\limsup_{n \to \infty} \sum_{Q \in \mathscr{D}_K} n^{p/d - 1}\mathbb{E}[\mathscr{W}_{p,Q}^p(\mu_n , \kappa_Q g)] \leqslant (1+\varepsilon)^{4+p} \beta_{p,d}\int_{(0,1)^d} g f^{-p/d}.\] This is valid for all \(\varepsilon > 0\), so that we obtain 16 . To obtain 15 , we note that \[\mathscr{W}_{p}^p(\hat{g}_n, g) = Z_n^{-1} \mathscr{W}_{p}^p(\mu_n, Z_n g) = (nZ_n^{-1}) n^{-1}\mathscr{W}_{p}^p(\mu_n, Z_n g).\]

Lemma 10 with \(Q = (0,1)^d\) shows that \(\mathbb{P}(|Z_n / n - 1| > \varepsilon) \leqslant 2e^{-nC}\) for some constant \(C\) depending on \(\varepsilon\). Therefore, using \(\mathscr{W}_{p,Q}(\mu_n, Z_n g)^p \lesssim Z_n\lesssim n\), we get \[\begin{align} \mathbb{E}[\mathscr{W}_{p}^p(\hat{g}_n, g)] &\leqslant (1+\varepsilon)n^{-1}\mathbb{E}[\mathscr{W}_{p}^p(\mu_n, Z_n g)] + C\mathbb{E}[n\mathbf{1}_{|Z_n / n - 1| > \varepsilon}]\\ &\leqslant (1+\varepsilon)n^{-1}\mathbb{E}[\mathscr{W}_{p}^p(\mu_n, Z_n g)] + Ce^{-nC}. \end{align}\] The final conclusion 15 follows from 16 .

6 Proof of the lower bound↩︎

6.1 The Boundary Wasserstein distance↩︎

For a bounded convex set \(\Omega\subset \mathbb{R}^d\) and two positive measures \(\mu, \nu\) we define the boundary Wasserstein distance as \[\label{defWb} \mathscr{W}b_{p,\Omega}(\mu,\nu)=\left(\inf_{\Pi\in \mathscr{C}_\Omega(\mu, \nu)}\int_{\overline{\Omega}\times\overline{\Omega}}|x-y|^p d\Pi(x,y)\right)^\frac{1}{p}\tag{25}\] where \(\mathscr{C}_\Omega(\mu, \nu)\) is the set of couplings \(\Pi\) such that \(\Pi_1{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}\Omega=\mu{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}\Omega\) and \(\Pi_2{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}\Omega=\nu{{\mathchoice {\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\>\vrule width.2pt \vbox to7pt{\vfill \hrule width7pt height.2pt}\>}{\,\vrule width.2pt \vbox to5pt{\vfill \hrule width5pt height.2pt}\,}{\,\vrule width.2pt \vbox to3.35pt{\vfill \hrule width3.35pt height.2pt}\,}}}\Omega\), see [33]. Let us point out that as opposed to \(\mathscr{W}_p\), we do not need to require \(\mu(\Omega)=\nu(\Omega)\) for \(\mathscr{W}b_{p,\Omega}(\mu,\nu)\) to be well-defined. However, if both measures have the same mass, we always have

\[\label{ineq:wb-w} \mathscr{W}b_{p,\Omega}^p(\mu,\nu)\leqslant \mathscr{W}_p(\mu,\nu)^p.\tag{26}\] As shown in, e.g., [3], [5], [9], the boundary Wasserstein distance is super-additive: for any collection of disjoint bounded convex sets \(A_1, \ldots, A_k \subset \Omega\), \[\label{eq:superad} \mathscr{W}b_{p,\Omega}^p(\mu,\nu)\geqslant \sum_{i=1}^k \mathscr{W}b_{p,A_i}^p(\mu,\nu).\tag{27}\]

6.2 The lower constant↩︎

We defined \(\beta_{p,d}\) in 14 . We may define the constant \(\beta^{\mathrm{low}}_{p,d}\) in a similar way, using the boundary Wasserstein distance. For an integer \(n\) and a family of iid random variables \(X_i\) uniformly distributed in \(Q_0=(0,1)^d\) we set \[\label{defFb} Fb_p(n)=n^{p/d}\mathbb{E}\left[ \mathscr{W}b_{p,Q_0}^p\left(\frac{1}{n}\sum_{i=1}^n \delta_{X_i},1\right) \right].\tag{28}\] Then, from [9], the following limit exists and belongs to \(]0, \infty[\): \[\label{def:betalow} \beta^{\mathrm{low}}_{p,d}=\lim_{n\to \infty} Fb_p(n).\tag{29}\]

6.3 The lower bound↩︎

We aim to prove that under Hypothesis 1 on \(f\) and \(g\), we have \[\label{mainclaimlower} \liminf_{n\to \infty} n^{p/d-1} \mathbb{E}\left[ \mathscr{W}_{p}^p(\mu_n,Z_ng) \right]\geqslant \beta^{\mathrm{low}}_{p,d}\int_{(0,1)^d} g f^{-p/d}.\tag{30}\] While we could in principle argue along the same lines as for the upper bound 16 , as shown in [9] the lower bound 30 is actually simpler. Fix \(k\in \mathbb{N}\) with \(k\gg1\) and set \(L=2^k\). For \(n\in \mathbb{N}\) we set \[K = \lfloor d^{-1} \log_2 n \rfloor\] which is not exactly the same as in 21 , and we define \[\ell= L 2^{-K}=2^{k-K}\sim L n^{-1/d}.\] We then consider the partition \(\mathscr{D}\) of dyadic cubes with side length \(\ell\). Using 26 and the super-additivity of \(\mathscr{W}b_{p,Q_0}\), we have \[\label{superad} \mathbb{E}\left[ \mathscr{W}_{p}^p(\mu_n,Z_ng) \right]\geqslant \sum_{Q\in \mathscr{D}} \mathbb{E}\left[ \mathscr{W}b_{p,Q}^p(\mu_n,Z_ng) \right].\tag{31}\] We now bound from below \(\mathbb{E}\left[ \mathscr{W}b_{p,Q}(\mu_n,Z_ng)^p \right]\) for every fixed \(Q\). We recall that \(\nu_n=\sum_{i=1}^n \delta_{X_i}\) and that \(\kappa=\kappa_Q, \gamma=\gamma_Q, \delta=\delta_Q\) were defined in 20 . We also introduce \[\lambda=Z_n\frac{\delta}{\kappa}=Z_n \frac{\int_Q g}{\int_Q f}.\] By Young’s inequality and \(\mathscr{W}b_{p,Q}\le \mathscr{W}_{p,Q}\) for measures with the same mass, we have for \(\varepsilon\in(0,1)\), \[\label{splitDir} \mathscr{W}b_{p,Q}^p(\mu_n,Z_ng)\geqslant (1-\varepsilon)\mathscr{W}b_{p,Q}^p(\gamma\nu_n,\lambda f) -\frac{C}{\varepsilon^{p-1}}\left(\mathscr{W}_{p,Q}^p(\mu_n,\gamma \nu_n) + \mathscr{W}_{p,Q}^p(\lambda f,Z_ng)\right).\tag{32}\] Notice that the main difference with 24 is that here \(\gamma \nu_n\) and \(\lambda f\) do not have the same mass.
We first state the estimates for the last two terms, see Section 7 for the proofs.

We now state the estimate for the first and main term, see again Section 7 for the proof.

Combining the results of both lemmas with 32 and summing over \(Q\) yields \[n^{p/d-1}\mathbb{E}\left[ \mathscr{W}_{p}^p(\mu_n,Z_ng) \right]\geqslant (1-\varepsilon)(1-o_{L}(1)) \beta^{\mathrm{low}}_{p,d}\int_{Q_0} g f^{-p/d} -C_\varepsilon\left( o_L(1) + C_L(n^{-p/d}+n^{-1/d} )\right).\] Sending first \(n\to \infty\) we get \[\liminf_{n\to \infty} n^{p/d-1}\mathbb{E}\left[ \mathscr{W}_{p}^p(\mu_n,Z_ng) \right]\geqslant (1-\varepsilon)(1-o_{L}(1)) \beta^{\mathrm{low}}_{p,d}\int_{Q_0} g f^{-p/d} -C_\varepsilon o_L(1).\] Sending then \(L\to \infty\) and finally \(\varepsilon\to 0\) concludes the proof of 30 . Once 30 is established, we can deduce that \[\label{zvocfqip} \limsup_{n \to \infty} n^{p/d} \mathbb{E}\mathscr{W}_p^p\left({Z_n}^{-1}\mu_n, g\right) \leqslant \beta_{p,d}\int_{(0,1)^d} g(x) f(x)^{-p/d} \mathrm{d} x\tag{34}\] using the same arguments as in the proof of 15 in Section 5.4.

7 Proofs of the technical lemmas↩︎

7.1 Proof of Lemma 1↩︎

We start by writing a multi-scale decomposition of \(\mathscr{W}_p^p(\mu_n, Z_n g)\).

To prove Lemma 1, we only need to prove that the error term \(E\) in ?? is negligible, that is \[\label{eq:E} \limsup_{n \to \infty} n^{p/d - 1}E = 0.\tag{35}\]

The proof relies on the following estimate.

7.2 Proof of Lemma 2↩︎

By 9 , \((\mathrm{I})\) is equal to \(\mathbb{E}[\gamma \mathscr{W}_{p,Q} ( \nu_n, (\delta/\gamma) f)^p]\). We now define \[\begin{align} &\mathscr{A}_Q = \left\lbrace \gamma_Q > (1 + \varepsilon) \frac{\int_Q g}{\int_Q f} \right\rbrace &&\mathscr{A} = \bigcup_{Q \in \mathscr{D}_K} \mathscr{A}_Q \\ &\mathscr{B}_Q = \left\lbrace \left|\frac{\nu_n(Q)}{n\int_Q f} - 1\right| > \varepsilon \right\rbrace &&\mathscr{B} = \bigcup_{Q \in \mathscr{D}_K} \mathscr{B}_Q \end{align}\] and finally \(\mathscr{E} = \mathscr{A} \cup \mathscr{B}\). Then, \[\begin{align} \mathbb{E}[\gamma \mathscr{W}_{p,Q}^p ( \nu_n, (\delta/\gamma) f)] &= \mathbb{E}[\mathbf{1}_{\mathscr{E}^c} \gamma \mathscr{W}_{p,Q}^p ( \nu_n, (\delta/\gamma) f)] + \mathbb{E}[\mathbf{1}_{\mathscr{E}} \gamma \mathscr{W}_{p,Q}^p ( \nu_n, (\delta/\gamma) f)] \nonumber \\ &\leqslant (1+\varepsilon)\frac{\int_Q g}{\int_Q f} \mathbb{E}[\mathbf{1}_{\mathscr{E}^c}\mathscr{W}_{p,Q}^p ( \nu_n, (\delta/\gamma) f)] + \mathbb{E}[\mathbf{1}_{\mathscr{E}} \gamma \mathscr{W}_{p,Q}^p ( \nu_n, (\delta/\gamma) f)].\label{eq:I95bound} \end{align}\tag{40}\]

7.2.1 Event \(\mathscr{E}\) has small probability↩︎

By the union bound and the fact that \(|\mathscr{D}_K|=2^{dK}\leqslant n\), we have \(\mathbb{P}(\mathscr{E}) \leqslant n \sup_{Q \in \mathscr{D}_K}(\mathbb{P}(\mathscr{A}_Q)+\mathbb{P}(\mathscr{B}_Q))\). Lemma 11 ensures that \(\mathbb{P}(\mathscr{A}_Q) \leqslant e^{-nC|Q|}\) for some constant \(C\) depending on \(\varepsilon\) only. Since \(n|Q|\sim \log(n)^{\alpha d}\) (see 23 ) we see that \(\mathbb{P}(\mathscr{A}_Q) \leqslant e^{-C \log(n)^{\alpha d}}\) which is itself smaller than \(n^{-p/d-10}\) for \(n\) large enough. Proposition 9 also ensures that \(\mathbb{P}(\mathscr{B}_Q) \leqslant n^{-p/d-10}\) for \(n\) large enough. By the union bound, we finally obtain \[\label{eq:prob95bad95event} \mathbb{P}(\mathscr{E}) \lesssim n^{-p/d-9}.\tag{41}\]

7.2.2 Bounding the second term of 40↩︎

By 10 and the fact that \(\mathrm{diam}(Q)\leqslant 1\) we see that \(\mathscr{W}_{p,Q} ( \nu_n, (\delta/\gamma) f)^p \leqslant \nu_n(Q) \leqslant n\). Consequently, by the Cauchy-Schwarz inequality, \[\begin{align} \mathbb{E}[\mathbf{1}_{\mathscr{E}} \gamma \mathscr{W}_{p,Q} ( \nu_n, (\delta/\gamma) f)^p] &\leqslant n \mathbb{E}[\mathbf{1}_{\mathscr{E}} \gamma] \leqslant n \sqrt{\mathbb{P}(\mathscr{E})}\sqrt{\mathbb{E}[\gamma^2]}. \end{align}\]

Lemma 12 and 41 imply that when \(n\) is large enough, this is smaller than \(n^{-p/d-7}\).

7.2.3 A uniform bound for unweighted Wasserstein cost↩︎

We define \[F_p(k, \varrho) = \mathbb{E}[\mathscr{W}_{Q_0}^p(k\varrho, \eta_k)],\] where \(\varrho\) is any probability density on \(Q_0\), \(k\) is an integer and \(\eta_k\) is a point process with \(k\) independent points whose common probability density is \(\varrho\). We will then use the main result of [8].

In particular, there is an integer \(N\) such that for all \(k \geqslant N\), \[\label{eq:F95p95bound95uniform} k^{p/d - 1} F_p(k, 1) \leqslant (1+\varepsilon) \beta_{p,d}.\tag{42}\]

7.2.4 Bounding the first term of 40↩︎

The proof is analogous to [26], [34]. We note that \(\delta / \gamma = \nu_n(Q) / \int_Q f\). By scaling, we have \[\mathbb{E}[\mathscr{W}_{p,Q}^p( \delta/\gamma f, \nu_n) \mid \nu_n(Q) = k] = \ell^{p}F_p(k, f_Q)\] where \(f_Q\) is the restriction of \(f\) on \(Q\), rescaled to \(Q_0\): \[\label{eq:f95Q} \forall x \in Q_0, \quad f_Q(x) = \frac{\ell^d f(z + x\ell)}{\int_Q f}\tag{43}\] where \(z\) is the corner of \(Q\) with the smallest coordinates (in 2d, that would be the bottom left corner). The gradient of \(f_Q\) is \(\lesssim\ell\), and there is at least one \(x_0 \in (0,1)^d\) such that \(f_Q(x_0)=1\), hence for any \(x \in Q_0\), we have \(|f_Q(x) -1| \lesssim\ell |x - x_0| \lesssim\ell\). The density \(f_Q\) thus satisfies \(|f_Q - 1|_\alpha \leqslant c \ell\) when \(n\) is large enough (see 51 for the definition of the \(\alpha\)-Hölder norm). Lemma 13 ensures that there is a map \(T : (0,1)^d \to (0,1)^d\) preserving the boundary of \((0,1)^d\) and such that \(T_\# f_Q = 1\) and \(\mathrm{Lip}(T), \mathrm{Lip}(T^{-1}) \leqslant 1+ c_{\alpha,d}|f_Q - 1|_\alpha =: L\).

For any \((X,Y)\) a coupling of \(k f_Q\) and \(\eta_k\), we thus see that \((T(X), T(Y))\) is a coupling of \(k\times 1\) and a point process with \(k\) independent points whose common probability density is \(1\). Consequently, \[\mathbb{E}[|X-Y|^p] = \mathbb{E}[|T^{-1}T(X)-T^{-1}T(Y)|^p] \leqslant L^p \mathbb{E}[|T(X)-T(Y)|^p].\] We conclude that \[\label{eq:W95p95bound95uniform} \mathbb{E}[\mathscr{W}_{p,Q}^p(k f_Q, \eta_k) \mid \nu_n(Q) = k] \leqslant L^p F_p(k, 1).\tag{44}\]

On the event \(\mathscr{B}^c\), we know that for all \(Q \in \mathscr{D}_K\), \(\nu_n(Q)\geqslant (1-\varepsilon)n\int_Q f\), which by Hypothesis 1 is \(\gtrsim\log(n)^{\alpha d}\). In particular, when \(n\) is large enough, on the event \(\mathscr{B}^c\), we have \(\nu_n(Q)\geqslant N\) with \(N\) defined for 42 . Consequently, uniformly over all \(Q \in \mathscr{D}_K\), we have \[\begin{align} \mathbb{E}[\mathbf{1}_{\mathscr{E}^c}\mathscr{W}_{p,Q}^p( \delta/\gamma f, \nu_n)] &\leqslant \mathbb{E}[\mathbf{1}_{\mathscr{B}^c} \mathbb{E}[\mathscr{W}_{p,Q}^p( \delta/\gamma f, \nu_n) \mid \nu_n(Q)]] \\ & \leqslant \ell^p L^p \mathbb{E}[\mathbf{1}_{\mathscr{B}^c} F_p(\nu_n(Q), 1)] \\ & \leqslant \ell^p L^p \mathbb{E}\left[\mathbf{1}_{\mathscr{B}^c} (1+\varepsilon)\nu_n(Q)^{1 - p/d}\beta_{p,d}\right] \\ &\leqslant (1+\varepsilon)^{2 - p/d}\ell^p L^p \left(n \int_Q f\right)^{1 - p/d}\beta_{p,d}\\ & = (1+\varepsilon)^{p+1}\ell^p n^{1 - p/d}\beta_{p,d}\left(\int_Q f\right)^{1 - p/d} \end{align}\] where in the last line, we used that \(L-1 \lesssim\ell\), hence \(L \leqslant 1+\varepsilon\) when \(n\) is large enough.

7.2.5 Gathering everything↩︎

Let us go back to 40 . If \(n\) is large enough, then for all \(Q \in \mathscr{D}_K\), we have \[\begin{align} \mathbb{E}[\mathscr{W}_{p,Q}^p ( \gamma \nu_n, \delta f)] &\leqslant (1+\varepsilon)\frac{\int_Q g}{\int_Q f} \mathbb{E}[\mathbf{1}_{\mathscr{E}^c}\mathscr{W}_{p,Q}^p ( \nu_n, (\delta/\gamma) f)] + n^{-p/d-7} \\ &\leqslant (1+\varepsilon)^{1+p} n^{1 - p/d} \beta_{p,d} \frac{\int_Q g}{\left(\int_Q f \right)^{p/d}} + n^{-p/d-7} \end{align}\] By the uniform continuity of \(f\) over \(Q_0\), for all \(x \in Q\) we have \[\frac{\int_Q g}{\left(\int_Q f\right)^{p/d}} \leqslant (1+\varepsilon)g(x)f(x)^{-p/d}.\] Plugging this into the preceding equation, we obtain that (uniformly over all \(Q \in \mathscr{D}_K\)), when \(n\) is large enough, we have \[\begin{align} \mathbb{E}[\mathscr{W}_{p,Q}^p ( \gamma \nu_n, \delta f)] &\leqslant (1+\varepsilon)^{2+p} n^{1 - p/d} \beta_{p,d} \int_Q gf^{-p/d} + n^{-p/d-7}\\ &\leqslant (1+\varepsilon)^{3+p} n^{1 - p/d} \beta_{p,d} \int_Q gf^{-p/d} \end{align}\] where in the last line we used the fact that \(\int_Q gf^{-p/d} \gtrsim\int_Q g \gtrsim n^{-1}\), and thus \(n^{-p/d-7}\) is negligible with respect to \(n^{1 - p/d} \int_Q gf^{-p/d}\).

7.3 Proof of Lemma 3 and the first bound of Lemma 5↩︎

We prove the following general bound: for any \(Q \subset (0,1)^d\), \[\label{eq:general95bound95II} \mathbb{E}[\mathscr{W}_{p,Q}^p(\delta f, \kappa g)] \lesssim n \ell^{2p} |Q|.\tag{45}\] Applying this first to \(\ell=\ell_K\) and observing that \(\ell_K = (\log n)^{\alpha} n^{-1/d}\), this would prove Lemma 3. Moreover, this proves also the first bound in Lemma 5.

7.4 Proof of Lemma 4 and second bound of Lemma 5↩︎

We prove the more general bound: for any \(Q \subset (0,1)^d\), \[\label{eq:general95bound95III} \mathbb{E}[\mathscr{W}_{p,Q}^p(\mu_n, \gamma \nu_n)] \lesssim n \ell^{p+1} |Q|.\tag{46}\] Applying this first to \(\ell=\ell_K\) this would conclude the proof of Lemma 4. Moreover, this proves also the second bound in Lemma 5.
We start with the simple case where \(\gamma=1/n\) so that both \(\mu_n\) and \(\gamma \nu_n\) are probability measures. For any coupling \((X,Y)\) of \(\mu_n\) and \(\gamma\nu_n\), the cost of this coupling is equal to \[\sum_{x, y}\mathbb{P}(X=x, Y=y) |x-y|^p \leqslant \mathrm{diam}(Q)^p \mathbb{P}(X\neq Y)\] which proves that \(\mathscr{W}_{p,Q}^p(\mu_n , \gamma \nu_n) \leqslant \mathrm{diam}(Q)^p \mathrm{d}_{\mathrm{TV}}(\mu_n, \gamma \nu_n)\). This remains true for any \(\gamma\) by scaling. Now, the measures \(\mu_n\) and \(\gamma \nu_n\) are supported on the same discrete set, so their total variation distance is equal to \[\mathrm{d}_{\mathrm{TV}}(\mu_n, \gamma \nu_n) = \frac{1}{2}\sum_{X_i \in Q}|W_i- \gamma|.\]

We also have \(W_i - \gamma = \nu_n(Q)^{-1}\sum_{X_j \in Q} (W_i - W_j)\), so that \[\begin{align} \mathrm{d}_{\mathrm{TV}}(\mu_n, \gamma \nu_n) &\leqslant \frac{\nu_n(Q)^{-1}}{2}\sum_{X_i, X_j \in Q} |W_i - W_j|. \end{align}\] Remember that \(W_i = g(X_i)/f(X_i)\); thanks to Hypothesis 1, the function \(g/f\) is Lipschitz. Indeed, \[\left|\nabla \frac{g}{f}\right| = \left|\frac{\nabla g}{f} - \frac{g}{f}\frac{\nabla f}{f}\right| \leqslant c^2 + c^4 \leqslant 2c^4.\] We thus have \[\begin{align} \mathbb{E}\mathrm{d}_{\mathrm{TV}}(\mu_n, \gamma \nu_n) &\lesssim\mathbb{E}\frac{\nu_n(Q)^{-1}}{2}\sum_{X_i, X_j \in Q} |X_i - X_j| \\ &\lesssim\mathbb{E}\nu_n(Q)\mathrm{diam}(Q)\\ &\lesssim n\ell |Q|. \end{align}\] This concludes the proof of 46 .

7.5 Proof of Lemma 6↩︎

Lemma 13 gives a map \(T\) which is a \((1+\varepsilon)\)-bi-Lipschitz map on \(Q\) and transports \(f/\int_{Q} f\) to \(1/|Q|\). Let \(N=\nu_n(Q)\) and let \(\hat{\nu}_N=T\sharp \nu_n\) so that \[\hat{\nu}_N=\sum_{i=1}^N \delta_{Y_i}\] with \(Y_i=T(X_i)\) which are iid uniformly distributed in \(Q\). We then have \[\mathscr{W}b_{p,Q}^p(\gamma\nu_n,\lambda f)\geqslant (1-\varepsilon)^p \mathscr{W}b_{p,Q}^p\left(\gamma\hat{\nu}_N,\frac{Z_n}{ |Q|} \int_Q g\right).\] We set \(M=Z_n \gamma^{-1} \int_Q g\sim Z_n |Q|\). Using the homogeneity of \(\mathscr{W}b_{p,Q}\), the triangle inequality, and the Young inequality, we see that \[\begin{gather} \label{splitmainterm} \mathscr{W}b_{p,Q}^p\left(\gamma\hat{\nu}_N, \frac{Z_n}{ |Q|} \int_Q g\right)=\gamma \mathscr{W}b_{p,Q}^p\left(\hat{\nu}_N,M/|Q|\right)\\ \geqslant \gamma\left[(1-\varepsilon) \mathscr{W}b_{p,Q}^p(\hat{\nu}_N,N/|Q|)-\frac{C}{\varepsilon^{p-1}}\mathscr{W}b_{p,Q}^p(N/|Q|, M/|Q|)\right]. \end{gather}\tag{47}\] To estimate the last right-hand side term we appeal to [9] to obtain \[\mathscr{W}b_{p,Q}^p(N/|Q|, M/|Q|)\lesssim\ell^{p+d} \frac{|N/|Q|- M/|Q||^p}{(M/|Q|)^{p-1}}=\ell^{p} \frac{|N-M|^p}{M^{p-1}}.\] We now write that \[\begin{gather} |N-M|= \frac{N}{\mu_n(Q)}|\mu_n(Q)-Z_n\int_Q g|\le \frac{N}{\mu_n(Q)}|\mu_n(Q)-n\int_Q g| + \frac{N}{\mu_n(Q)} |n-Z_n|\int_Q g\\ \lesssim|\mu_n(Q)-n\int_Q g| + |n-Z_n||Q|=|\mu_n(Q)-\mathbb{E}\left[ \mu_n(Q) \right]| + |Z_n-\mathbb{E}\left[ Z_n \right]||Q| \end{gather}\] to obtain by Cauchy-Schwarz, \[\begin{gather} \mathbb{E}\left[ \frac{|N-M|^p}{M^{p-1}} \right]\\ \lesssim|Q|^{1-p} \mathbb{E}\left[ Z_n^{2(1-p)} \right]^{1/2}\left(\mathbb{E}\left[ |\mu_n(Q)-\mathbb{E}\left[ \mu_n(Q) \right]|^{2p} \right]^{1/2}+ |Q|^p\mathbb{E}\left[ |Z_n-\mathbb{E}\left[ Z_n \right]|^{2p} \right]^{1/2}\right)\\ \lesssim|Q|^{1-p} n^{1-p}\left( (n|Q|)^{p/2}+ |Q|^pn^{p/2}\right) \lesssim(|Q| n)^{1-p/2}. \end{gather}\] In conclusion, recalling that by definition \(n|Q|=L^d\), \[\begin{align} \mathbb{E}\left[ \mathscr{W}b_{p,Q}^p(N/|Q|, M/|Q|) \right]&\lesssim n^{1-p/d} (|Q|n)^{p(\frac{1}{d}-\frac{1}{2})} |Q|\label{estimchangemass}\\ &=n^{1-p/d} L^{-p(\frac{d}{2}-1)} |Q|\\ &=n^{1-p/d}o_L(1) |Q|. \end{align}\tag{48}\]

We now estimate the first right-hand side term of 47 . By the uniform continuity of \(f\) and \(g\), when \(n\) is large enough we have \(\gamma\geqslant (1-\varepsilon) g(x)/f(x)\). Recalling the definition 28 of \(Fb_p\), we thus have by scaling \[\gamma \mathbb{E}\left[ \mathscr{W}b_{p,Q}^p(\hat{\nu}_N,N/|Q|) \right]\geqslant (1-\varepsilon) \frac{g(x)}{ f(x)} \ell^p \mathbb{E}\left[ N^{1-p/d}Fb_p(N) \right].\] Since \(\mathbb{E}\left[ N \right]=n\int_Q f\geqslant (1-\varepsilon) f(x) n|Q|\) which tends to infinity as \(L\to \infty\) we have by definition of \(\beta^{\mathrm{low}}_{p,d}\) and the good concentration properties of \(N\) that \[\mathbb{E}\left[ N^{1-p/d}Fb_p(N) \right]\geqslant (1- o_L(1))\beta^{\mathrm{low}}_{p,d}\left(n\int_Q f\right)^{1-p/d}\] and thus \[\begin{align} \gamma \mathbb{E}\left[ \mathscr{W}b_{p,Q}^p(\hat{\nu}_N,N/|Q|) \right]&\geqslant (1-\varepsilon) (1- o_L(1))\beta^{\mathrm{low}}_{p,d}\frac{g(x)}{ f(x)} \ell^p (n|Q| f(x))^{1-p/d}\\ &=n^{1-p/d}(1-\varepsilon) (1- o_L(1))\beta^{\mathrm{low}}_{p,d}g(x) f(x)^{-p/d} |Q|. \end{align}\] Averaging over \(x\) we find \[\label{maintermlocDir} \gamma \mathbb{E}\left[ \mathscr{W}b_{p,Q}^p(\hat{\nu}_N,N/|Q|) \right]\geqslant n^{1-p/d}(1-\varepsilon) (1- o_L(1))\beta^{\mathrm{low}}_{p,d}\int_Q g f^{-p/d}.\tag{49}\] Combining 48 and 49 together yields ?? .

8 Concentration inequalities↩︎

We need concentration results for sums of random variables, typically of the form \[\sum_{i=1}^n W_i \mathbf{1}_{X_i \in Q}\] where \(Q\) is a very small measurable set — in particular, most of the summands are zero, and the variance of the summands is small. In this context, Bennett-type inequalities are tighter than Hoeffding’s. For generality, we will formulate our results for any pair \((X,W)\) of random variables, even if in the main text we will only use the case where \(X\) is a random variable with density \(f\) and \(W = g(X)/f(X)\).

8.1 Bennett’s inequality↩︎

Let \(Z_1, \ldots, Z_n\) be a sequence of iid random variables and let \(Z = \sum_{i=1}^n Z_i\). Then, Bennett’s inequality states that for any \(t>0\), \[\mathbb{P}(Z > \mathbb{E}[Z] + t) \leqslant \exp\left(-n (\sigma/m)^2 \times h\left(\frac{mt}{n\sigma^2}\right)\right),\] where \(m = |Z_i|_\infty\) and \(\sigma\) is the standard deviation of \(Z_i\). The function \(h\) is defined by \[h(x) = (1+x)\ln(1+x) - x.\] For future reference, we state the elementary properties of \(h\) as a lemma.

8.2 Concentration for linear functionals↩︎

Let \((X_i)\) be iid random variables with density \(f\) and let \(W_i = g(X_i)/f(X_i)\). For any bounded measurable function \(\varphi\), we are interested in the concentration of the functionals \[I(\varphi) = \sum_{i=1}^n W_i \varphi(X_i).\] The expectation of \(I\) is given by \(\mathbb{E}[I(\varphi)] = n\mathbb{E}_{X\sim g}[\varphi(X)]= n\int g\varphi\). Directly applying Bennett’s inequality, we get the following concentration result.

8.3 Concentration for \(\kappa\)↩︎

We will use the precedings bounds with \(\varphi = \mathbf{1}_{X \in Q}\) to get concentration bounds for \(\kappa_Q\). Recall that \[\kappa_Q = \frac{\sum W_i \mathbf{1}_{X_i \in Q}}{\int_Q g}\] so that \(\mathbb{E}[\kappa_Q] = n\). For simplicity we will note \(m,\sigma^2\) instead of \(m_Q,\sigma^2_Q\) and we will often use the notation \[b = \frac{m\mathbb{E}[W \mathbf{1}_{X \in Q}]}{\sigma^2}.\] The preceding concentration results thus yield that \[\label{i:kappa} \mathbb{P}\left(\left| \frac{\kappa_Q}{n} - 1 \right| > \delta\right) \leqslant 2e^{-n \sigma^2/m^2 h(\delta b)}\tag{50}\] and the same inequality holds for the lower deviation.

8.4 Concentration for \(\gamma\)↩︎

We recall that \[\gamma_Q = \frac{\sum_{i=1}^n W_i \mathbf{1}_{X_i \in Q}}{\sum_{i=1}^n \mathbf{1}_{X_i \in Q}} .\]

We will also need the following very crude bound on \(\gamma\).

8.5 Concentration for the number of points↩︎

Finally, we will need a uniform control of the number of points \(X_i\) falling in each cube \(Q \in \mathscr{D}_K\). We recall that \(\nu_n(A) = \sum_{i=1}^n \delta_{X_i}(A)\) is the unweighted empirical measure of the \(X_i\). Of course, the \(n^{-p/d-10}\) here is somehow arbitrary; any term which is sufficiently small with respect to \(n^{-p/d}\) would work.

9 Functional inequalities and transport↩︎

9.1 Functional inequalities↩︎

We use the Poincaré-Wirtinger inequality in the following form.

We also made repeated use of the following Proposition, which is proven in [8] (Lemma 3.4).

9.2 Transport map between almost uniform measures↩︎

We end this section by stating Proposition 2.4 from [26], a result on the existence of a (non-optimal) transport map between almost uniform measures. The \(\alpha\)-Hölder norm used in the lemma is defined by \[\label{def:alpha95holder95norm} |h|_\alpha = \sup_{x,y \in (0,1)^d, x\neq y} \frac{|h(x) - h(y)|}{|x-y|^\alpha} + \sup_{x \in (0,1)^d} |h(x)|\tag{51}\] and the Lipschitz constant of a function \(h\) is defined by \[\mathrm{Lip}(h) = \sup_{x,y \in (0,1)^d, x\neq y} \frac{|h(x) - h(y)|}{|x-y|}.\]

References↩︎

---

1. Their results apply to the Wiener measure and also to the fractional Brownian measure and the Brownian sheet.↩︎

2. \(v_R\) was defined in 12 .↩︎