summaries: This article proposes a study on the rates of estimation of Brenier maps using the plug in "barycentric mapping" estimator. As a consequence of this study authors obtain rates of convergence for $W_2^{2}$ in various context (depending on the regularity of the measures and on the setting: discrete-discrete of semi-discrete). They apply their results for Wasserstein barycenter estimation and Nonparametric independence testing.  This paper analyze the statistical behavior of the general plug-in estimators of OT defined by barycentric projections. The authors provide a thorough analysis of the rate of convergence for the transport cost and the transport map. They also consider kernel smoothed plug-in estimators and relate its rate of convergence to the smoothness of the densities, which alleviates the curse of dimensionality suffered by the plug-in estimator. Summary =======  Context: --------  This paper tackles the problem of estimating the rate of convergence of empirical Monge maps (through its Kantorovich barycentric map approximator) to the true Monge map for  absolutely continuous Lebesgue measures.    Contribution ------------  - The authors provide new rates of convergence of the barycentric map using standard  empirical plug-in estimators (O(n^{-2/d})) - If the measures are regular enough (Holder smooth), the dependence on the dimension completely vanishes (which was established for Gaussian distributions https://arxiv.org/pdf/1905.10155.pdf) - The authors propose to take advantage of this smoothness result to sample from smoothened distributions (using a Kernel density estimator) and provide the sufficient conditions to obtain a desired rate. Higher smoothness however requires more samples: a tradeoff between statistical accuracy and computational cost cannot be avoided.   In standard OT problem, marginal distributions of the data are unknown. The current work considers plug-in estimation using Barycentric projection and derived rate of estimation. The idea of Barycentric projection is to plug the estimated marginal distributions of the data in the OT problem. Conditional mean of given us estimated and rate of estimation is derived. Various plug-in estimators of density of are considered such as empirical CDF and kernel density estimator. This paper studies the problem of estimating an optimal transport map between two probability distributions on Euclidean space, given IID samples from each of the distributions. The main result (Theorem 2.1) gives a general upper bound, in terms of the dual representation of the 2-Wasserstein distance. This is used to give upper bounds, in probability, on the rate of decay of the estimation risk both for the plugging in empirical distributions in the absence of smoothness assumptions and for kernel estimates under smoothness assumptions. Finally, implications for two applications (Wasserstein barycenter estimation and nonparametric independence testing) are discussed.
abstractText: Optimal transport maps between two probability distributions μ and ν on R have found extensive applications in both machine learning and statistics. In practice, these maps need to be estimated from data sampled according to μ and ν. Plugin estimators are perhaps most popular in estimating transport maps in the field of computational optimal transport. In this paper, we provide a comprehensive analysis of the rates of convergences for general plug-in estimators defined via barycentric projections. Our main contribution is a new stability estimate for barycentric projections which proceeds under minimal smoothness assumptions and can be used to analyze general plug-in estimators. We illustrate the usefulness of this stability estimate by first providing rates of convergence for the natural discretediscrete and semi-discrete estimators of optimal transport maps. We then use the same stability estimate to show that, under additional smoothness assumptions of Sobolev type or Besov type, kernel smoothed or wavelet based plug-in estimators respectively speed up the rates of convergence and significantly mitigate the curse of dimensionality suffered by the natural discrete-discrete/semi-discrete estimators. As a by-product of our analysis, we also obtain faster rates of convergence for plug-in estimators of W2(μ, ν), the Wasserstein distance between μ and ν, under the aforementioned smoothness assumptions, thereby complementing recent results in Chizat et al. (2020). Finally, we illustrate the applicability of our results in obtaining rates of convergence for Wasserstein barycenter between two probability distributions and obtaining asymptotic detection thresholds for some recent optimaltransport based tests of independence.
authors.name: Nabarun Deb
references.authors: Pankaj K Agarwal R Sharathkumar
references.title: Approximation algorithms for bipartite matching with metric and geometric costs
references.venue: In Proceedings of the forty-sixth annual ACM symposium on Theory of computing,
references.authors: Martial Agueh Guillaume Carlier
references.title: Barycenters in the wasserstein space
references.venue: SIAM Journal on Mathematical Analysis,
references.authors: Luigi Ambrosio Nicola Gigli Giuseppe Savaré
references.title: Gradient flows in metric spaces and in the space of probability measures
references.venue: Lectures in Mathematics ETH Zürich. Birkhäuser Verlag, Basel, second edition,
references.authors: Ery Arias-Castro David Mason Bruno Pelletier
references.title: On the estimation of the gradient lines of a density and the consistency of the mean-shift algorithm
references.venue: The Journal of Machine Learning Research,
references.authors: Arnab Auddy Nabarun Deb Sagnik Nandy
references.title: Exact detection thresholds for chatterjee’s correlation
references.venue: arXiv preprint arXiv:2104.15140,
references.authors: Franck Barthe Charles Bordenave
references.title: Combinatorial optimization over two random point sets. In Séminaire de Probabilités XLV, volume 2078 of Lecture Notes in Math., pages 483–535
references.authors: Jean-David Benamou Yann Brenier
references.title: A computational fluid mechanics solution to the monge-kantorovich mass transfer problem
references.venue: Numerische Mathematik,
references.authors: Jean-David Benamou Brittany D Froese Adam M Oberman
references.title: Numerical solution of the optimal transportation problem using the monge–ampère equation
references.venue: Journal of Computational Physics,
references.authors: Robert J Berman
references.title: Convergence rates for discretized monge–ampère equations and quantitative stability of optimal transport
references.venue: Foundations of Computational Mathematics,
references.authors: Espen Bernton Pierre E Jacob Mathieu Gerber Christian P Robert
references.title: Inference in generative models using the Wasserstein distance
references.venue: arXiv preprint arXiv:1701.05146,
references.authors: T.B. Berrett R.J. Samworth
references.title: Nonparametric independence testing via mutual information
references.authors: Jérémie Bigot Raúl Gouet Thierry Klein Alfredo López
references.title: Geodesic PCA in the Wasserstein space by convex PCA
references.venue: In Annales de l’Institut Henri Poincaré, Probabilités et Statistiques,
references.authors: Jérémie Bigot Raúl Gouet Thierry Klein Alfredo Lopez
references.title: Upper and lower risk bounds for estimating the wasserstein barycenter of random measures on the real line
references.venue: Electronic journal of statistics,
references.authors: Jérémie Bigot Thierry Klein
references.title: Characterization of barycenters in the wasserstein space by averaging optimal transport
references.venue: maps. ESAIM: Probability and Statistics,
references.authors: Adrien Blanchet Guillaume Carlier
references.title: Optimal transport and cournot-nash equilibria
references.venue: Mathematics of Operations Research,
references.authors: Melf Boeckel Vladimir Spokoiny Alexandra Suvorikova
references.title: Multivariate brenier cumulative distribution functions and their application to non-parametric testing
references.venue: arXiv preprint arXiv:1809.04090,
references.authors: Emmanuel Boissard Thibaut Le Gouic Jean-Michel Loubes
references.title: Distribution’s template estimate with wasserstein
references.venue: metrics. Bernoulli,
references.authors: Nicolas Bonneel Michiel Van De Panne Sylvain Paris Wolfgang Heidrich
references.title: Displacement interpolation using lagrangian mass transport
references.venue: In Proceedings of the 2011 SIGGRAPH Asia Conference,
references.authors: Yann Brenier
references.title: Polar factorization and monotone rearrangement of vector-valued functions
references.venue: Comm. Pure Appl. Math.,
references.authors: Luis A. Caffarelli
references.title: Boundary regularity of maps with convex potentials
references.venue: Comm. Pure Appl. Math.,
references.authors: Luis A. Caffarelli
references.title: The regularity of mappings with a convex potential
references.venue: J. Amer. Math. Soc.,
references.authors: Luis A. Caffarelli
references.title: Boundary regularity of maps with convex potentials
references.venue: II. Ann. of Math. (2),
references.authors: Guillaume Carlier Ivar Ekeland
references.title: Matching for teams
references.venue: Economic theory,
references.authors: Guillaume Carlier Adam Oberman Edouard Oudet
references.title: Numerical methods for matching for teams and wasserstein barycenters
references.venue: ESAIM: Mathematical Modelling and Numerical Analysis,
references.authors: Rick Chartrand Brendt Wohlberg Kevin Vixie Erik Bollt
references.title: A gradient descent solution to the monge-kantorovich problem
references.venue: Applied Mathematical Sciences,
references.authors: Victor Chernozhukov Alfred Galichon Marc Hallin Marc Henry
references.title: Monge-Kantorovich depth, quantiles, ranks and signs
references.venue: Ann. Statist.,
references.authors: Sinho Chewi Tyler Maunu Philippe Rigollet Austin J Stromme
references.title: Gradient descent algorithms for bures-wasserstein barycenters
references.venue: In Conference on Learning Theory,
references.authors: Pierre-André Chiappori Robert J McCann Lars P Nesheim
references.title: Hedonic price equilibria, stable matching, and optimal transport: equivalence, topology, and uniqueness
references.venue: Economic Theory,
references.authors: Lenaic Chizat Gabriel Peyré Bernhard Schmitzer François-Xavier Vialard
references.title: Scaling algorithms for unbalanced optimal transport problems
references.venue: Mathematics of Computation,
references.authors: Lenaic Chizat Pierre Roussillon Flavien Léger François-Xavier Vialard Gabriel Peyré
references.title: Faster Wasserstein distance estimation with the sinkhorn divergence
references.venue: Advances in Neural Information Processing Systems,
references.authors: Evgenii Chzhen Christophe Denis Mohamed Hebiri Luca Oneto Massimiliano Pontil
references.title: Leveraging Labeled and Unlabeled Data for Consistent Fair Binary Classification
references.venue: NeurIPS
references.authors: Sebastian Claici Edward Chien Justin Solomon
references.title: Stochastic wasserstein barycenters
references.venue: In International Conference on Machine Learning,
references.authors: Albert Cohen
references.title: Numerical analysis of wavelet methods, volume 32 of Studies in Mathematics and its Applications
references.authors: Nicolas Courty Rémi Flamary Amaury Habrard Alain Rakotomamonjy
references.title: Joint distribution optimal transportation for domain adaptation
references.venue: In Proceedings of the 31st International Conference on Neural Information Processing Systems,
references.authors: Nicolas Courty Rémi Flamary Devis Tuia Alain Rakotomamonjy
references.title: Optimal transport for domain adaptation
references.venue: IEEE transactions on pattern analysis and machine intelligence,
references.authors: Marco Cuturi
references.title: Sinkhorn distances: Lightspeed computation of optimal transport
references.venue: Advances in neural information processing systems,
references.authors: Marco Cuturi Arnaud Doucet
references.title: Fast computation of wasserstein barycenters
references.venue: In International conference on machine learning,
references.authors: Bharath Bhushan Damodaran Benjamin Kellenberger Rémi Flamary Devis Tuia Nicolas Courty
references.title: Deepjdot: Deep joint distribution optimal transport for unsupervised domain adaptation
references.venue: In Proceedings of the European Conference on Computer Vision (ECCV),
references.authors: Nabarun Deb Bhaswar B Bhattacharya Bodhisattva Sen
references.title: Efficiency lower bounds for distribution-free hotelling-type two-sample tests based on optimal transport
references.venue: arXiv preprint arXiv:2104.01986,
references.authors: Nabarun Deb Promit Ghosal Bodhisattva Sen
references.title: Measuring association on topological spaces using kernels and geometric graphs
references.venue: arXiv preprint arXiv:2010.01768,
references.authors: Nabarun Deb Bodhisattva Sen
references.title: Multivariate rank-based distribution-free nonparametric testing using measure transportation
references.venue: Journal of the American Statistical Association,
references.authors: Steffen Dereich Michael Scheutzow Reik Schottstedt
references.title: Constructive quantization: approximation by empirical measures
references.venue: Ann. Inst. Henri Poincaré Probab. Stat.,
references.authors: David L. Donoho Iain M. Johnstone Gérard Kerkyacharian Dominique Picard
references.title: Density estimation by wavelet thresholding
references.venue: Ann. Statist.,
references.authors: Dominique Drouet Mari Samuel Kotz
references.title: Correlation and dependence. Imperial College Press, London; distributed by World Scientific Publishing Co., Inc., River
references.authors: Ivar Ekeland Alfred Galichon Marc Henry
references.title: Optimal transportation and the falsifiability of incompletely specified economic models
references.venue: Economic Theory,
references.authors: Tarek A El Moselhy Youssef M Marzouk
references.title: Bayesian inference with optimal maps
references.venue: Journal of Computational Physics,
references.authors: Sira Ferradans Nicolas Papadakis Gabriel Peyré Jean-François Aujol
references.title: Regularized discrete optimal transport
references.venue: SIAM Journal on Imaging Sciences,
references.authors: Aden Forrow Jan-Christian Hütter Mor Nitzan Philippe Rigollet Geoffrey Schiebinger Jonathan Weed
references.title: Statistical optimal transport via factored couplings
references.venue: In The 22nd International Conference on Artificial Intelligence and Statistics,
references.authors: Nicolas Fournier Arnaud Guillin
references.title: On the rate of convergence in Wasserstein distance of the empirical measure
references.venue: Probab. Theory Related Fields,
references.authors: Terry L Friesz J Enrique Fernandez
references.title: A model of optimal transport maintenance with demand responsiveness
references.venue: Transportation Research Part B: Methodological,
references.authors: Charlie Frogner Chiyuan Zhang Hossein Mobahi Mauricio Araya-Polo Tomaso Poggio
references.title: Learning with a Wasserstein loss
references.venue: In Proceedings of the 28th International Conference on Neural Information Processing Systems - Volume 2,
references.authors: Kenji Fukumizu Bharath K Sriperumbudur Arthur Gretton Bernhard Schölkopf
references.title: Characteristic kernels on groups and semigroups
references.venue: In NIPS,
references.authors: Harold N Gabow Robert E Tarjan
references.title: Faster scaling algorithms for network problems
references.venue: SIAM Journal on Computing,
references.authors: Alfred Galichon
references.title: Optimal transport methods in economics
references.authors: Promit Ghosal Bodhisattva Sen
references.title: Multivariate ranks and quantiles using optimal transportation and applications to goodness-of-fit testing
references.venue: arXiv preprint arXiv:1905.05340,
references.authors: Nicola Gigli
references.title: On Hölder continuity-in-time of the optimal transport map towards measures along a curve
references.venue: Proc. Edinb. Math. Soc. (2),
references.authors: Evarist Giné Richard Nickl
references.title: Uniform central limit theorems for kernel density estimators
references.venue: Probab. Theory Related Fields,
references.authors: Evarist Giné Richard Nickl
references.title: Mathematical foundations of infinite-dimensional statistical models
references.venue: Cambridge Series in Statistical and Probabilistic Mathematics,
references.authors: Joan Glaunes Alain Trouvé Laurent Younes
references.title: Diffeomorphic matching of distributions: A new approach for unlabelled point-sets and sub-manifolds matching
references.venue: In Proceedings of the 2004 IEEE Computer Society Conference on Computer Vision and Pattern Recognition,
references.authors: Ian J Goodfellow Jean Pouget-Abadie Mehdi Mirza Bing Xu David Warde-Farley Sherjil Ozair Aaron Courville Yoshua Bengio
references.title: Generative adversarial networks
references.venue: Advances in Neural Information Processing Systems
references.authors: Raghuraman Gopalan Ruonan Li Rama Chellappa
references.title: Domain adaptation for object recognition: An unsupervised approach
references.venue: In 2011 International Conference on Computer Vision,
references.authors: Paula Gordaliza Eustasio Del Barrio Gamboa Fabrice Jean-Michel Loubes
references.title: Obtaining fairness using optimal transport theory
references.venue: In International Conference on Machine Learning,
references.authors: A. Gretton R. Herbrich A. Smola O. Bousquet B. Schölkopf
references.title: Kernel methods for measuring independence
references.venue: J. Mach. Learn. Res.,
references.authors: Arthur Gretton Olivier Bousquet Alex Smola Bernhard Schölkopf
references.title: Measuring statistical dependence with hilbert-schmidt norms
references.venue: In International conference on algorithmic learning theory,
references.authors: Arthur Gretton Kenji Fukumizu Choon Hui Teo Le Song Bernhard Schölkopf Alexander J Smola
references.title: A kernel statistical test of independence
references.venue: In Nips,
references.authors: Florian Gunsilius Susanne M Schennach
references.title: Independent nonlinear component analysis
references.venue: Technical report, cemmap working paper,
references.authors: Florian F Gunsilius
references.title: On the convergence rate of potentials of brenier maps
references.venue: Econometric Theory,
references.authors: Marc Hallin Eustasio Del Barrio Juan Cuesta-Albertos Carlos Matrán
references.title: Distribution and quantile functions, ranks and signs in dimension d: A measure transportation approach
references.venue: The Annals of Statistics,
references.authors: Bruce E. Hansen
references.title: Uniform convergence rates for kernel estimation with dependent data
references.venue: Econometric Theory,
references.authors: Wolfgang Härdle Gerard Kerkyacharian Dominique Picard Alexander Tsybakov
references.title: Wavelets, approximation, and statistical applications, volume 129 of Lecture Notes in Statistics
references.authors: Ruth Heller Yair Heller Malka Gorfine
references.title: A consistent multivariate test of association based on ranks of distances
references.authors: Jan-Christian Hütter Philippe Rigollet
references.title: Minimax estimation of smooth optimal transport maps
references.venue: The Annals of Statistics,
references.authors: Roy Jonker Anton Volgenant
references.title: A shortest augmenting path algorithm for dense and sparse linear assignment problems
references.authors: Julie Josse Susan Holmes
references.title: Measuring multivariate association and beyond
references.venue: Stat. Surv.,
references.authors: Kirthevasan Kandasamy Willie Neiswanger Jeff Schneider Barnabás Póczos Eric P. Xing
references.title: Neural architecture search with bayesian optimisation and optimal transport
references.venue: In Proceedings of the 32nd International Conference on Neural Information Processing Systems,
references.authors: L.V. Kantorovich
references.title: On a problem of Monge
references.venue: Zap. Nauchn. Sem. S.-Peterburg. Otdel. Mat. Inst. Steklov. (POMI), 312(Teor. Predst. Din. Sist. Komb. i Algoritm. Metody
references.authors: L. Kantorovitch
references.title: On the translocation of masses
references.venue: C. R. (Doklady) Acad. Sci. URSS (N.S.),
references.authors: G. Kerkyacharian D. Picard
references.title: Density estimation in Besov spaces
references.venue: Statist. Probab. Lett.,
references.authors: Ilmun Kim Sivaraman Balakrishnan Larry Wasserman
references.title: Minimax optimality of permutation tests
references.venue: arXiv preprint arXiv:2003.13208,
references.authors: Sanggyun Kim Rui Ma Diego Mesa Todd P Coleman
references.title: Efficient bayesian inference methods via convex optimization and optimal transport
references.venue: IEEE International Symposium on Information Theory,
references.authors: Diederik P. Kingma Max Welling
references.title: Auto-Encoding Variational Bayes
references.venue: In 2nd International Conference on Learning Representations,
references.authors: J. Komlós P. Major G. Tusnády
references.title: An approximation of partial sums of independent RV’s, and the sample DF
references.venue: II. Z. Wahrscheinlichkeitstheorie und Verw. Gebiete,
references.authors: Gil Kur Fuchang Gao Adityanand Guntuboyina Bodhisattva Sen
references.title: Convex regression in multidimensions: Suboptimality of least squares estimators
references.venue: arXiv preprint arXiv:2006.02044,
references.authors: Thibaut Le Gouic Jean-Michel Loubes
references.title: Existence and consistency of wasserstein barycenters
references.venue: Probability Theory and Related Fields,
references.authors: Tong Li Ming Yuan
references.title: On the optimality of gaussian kernel based nonparametric tests against smooth alternatives
references.venue: arXiv preprint arXiv:1909.03302,
references.authors: Wenbo Li Ricardo H Nochetto
references.title: Quantitative stability and error estimates for optimal transport plans
references.venue: IMA Journal of Numerical Analysis,
references.authors: Yujia Li Kevin Swersky Rich Zemel
references.title: Generative moment matching networks
references.venue: In International Conference on Machine Learning,
references.authors: Giulia Luise Alessandro Rudi Massimiliano Pontil Carlo Ciliberto
references.title: Differential properties of sinkhorn approximation for learning with Wasserstein distance
references.venue: In Proceedings of the 32nd International Conference on Neural Information Processing Systems,
references.authors: Valentina Masarotto Victor M Panaretos Yoav Zemel
references.title: Procrustes metrics on covariance operators and optimal transportation of gaussian processes
references.venue: Sankhyā A,
references.authors: Robert J. McCann
references.title: Existence and uniqueness of monotone measure-preserving maps
references.venue: Duke Math. J.,
references.authors: Gonzalo Mena Jonathan Niles-Weed
references.title: Statistical bounds for entropic optimal transport: Sample complexity and the central limit theorem
references.venue: Advances in Neural Information Processing Systems,
references.authors: Quentin Mérigot
references.title: A multiscale approach to optimal transport
references.venue: In Computer Graphics Forum,
references.authors: Quentin Merigot Boris Thibert
references.title: Optimal transport: discretization and algorithms. Handbook of Numerical Analysis 22 – Geometric PDES (arXiv
references.authors: Yves Meyer
references.title: Ondelettes et opérateurs
references.venue: I. Actualités Mathématiques. [Current Mathematical Topics]. Hermann, Paris,
references.authors: Shakir Mohamed Balaji Lakshminarayanan
references.title: Learning in implicit generative models
references.venue: In Proceedings of the International Conference in Learning Representations,
references.authors: Alfred Müller
references.title: Integral probability metrics and their generating classes of functions
references.venue: Adv. in Appl. Probab.,
references.authors: È.A. Nadaraja
references.title: On non-parametric estimates of density functions and regression
references.venue: Teor. Verojatnost. i Primenen.,
references.authors: Jonathan Niles-Weed Philippe Rigollet
references.title: Estimation of Wasserstein distances in the spiked transport model
references.venue: arXiv preprint arXiv:1909.07513,
references.authors: Nadia Oudjane Christian Musso
references.title: L2-density estimation with negative kernels
references.venue: ISPA
references.authors: Nicolas Papadakis Gabriel Peyré Edouard Oudet
references.title: Optimal transport with proximal splitting
references.venue: SIAM Journal on Imaging Sciences,
references.authors: Emanuel Parzen
references.title: On estimation of a probability density function and mode
references.venue: Ann. Math. Statist.,
references.authors: Gabriel Peyré Marco Cuturi
references.title: Computational optimal transport: With applications to data science
references.venue: Foundations and Trends® in Machine Learning,
references.authors: Svetlozar T Rachev
references.title: The monge–kantorovich mass transference problem and its stochastic applications
references.venue: Theory of Probability & Its Applications,
references.authors: Aaditya Ramdas Nicolás García Trillos Marco Cuturi
references.title: On Wasserstein two-sample testing and related families of nonparametric tests
references.venue: Entropy,
references.authors: Sebastian Reich
references.title: A dynamical systems framework for intermittent data assimilation
references.venue: BIT Numerical Mathematics,
references.authors: Sebastian Reich
references.title: A nonparametric ensemble transform method for bayesian inference
references.venue: SIAM Journal on Scientific Computing,
references.authors: Thomas Rippl Axel Munk Anja Sturm
references.title: Limit laws of the empirical Wasserstein distance: Gaussian distributions
references.venue: J. Multivariate Anal.,
references.authors: Tim Salimans Han Zhang Alec Radford Dimitris Metaxas
references.title: Improving GANs using optimal transport
references.venue: arXiv preprint arXiv:1803.05573,
references.authors: Filippo Santambrogio
references.title: Optimal transport for applied mathematicians
references.venue: Birkäuser, NY,
references.authors: Vivien Seguy Bharath Bhushan Damodaran Rémi Flamary Nicolas Courty Antoine Rolet Mathieu Blondel
references.title: Large-scale optimal transport and mapping estimation
references.venue: In Proceedings of the International Conference in Learning Representations,
references.authors: Hongjian Shi Mathias Drton Fang Han
references.title: Distribution-free consistent independence tests via center-outward ranks and signs
references.venue: Journal of the American Statistical Association,
references.authors: Hongjian Shi Marc Hallin Mathias Drton Fang Han
references.title: Rate-optimality of consistent distribution-free tests of independence based on center-outward ranks and signs
references.venue: arXiv preprint arXiv:2007.02186,
references.authors: Bernard W. Silverman
references.title: Weak and strong uniform consistency of the kernel estimate of a density and its derivatives
references.venue: Ann. Statist.,
references.authors: Max Sommerfeld Axel Munk
references.title: Inference for empirical Wasserstein distances on finite spaces
references.venue: J. R. Stat. Soc. Ser. B. Stat. Methodol.,
references.authors: Bharath K. Sriperumbudur Kenji Fukumizu Arthur Gretton Bernhard Schölkopf Gert R.G. Lanckriet
references.title: On the empirical estimation of integral probability metrics
references.venue: Electron. J. Stat.,
references.authors: Bharath K. Sriperumbudur Arthur Gretton Kenji Fukumizu Bernhard Schölkopf Gert R.G. Lanckriet
references.title: Hilbert space embeddings and metrics on probability measures
references.venue: J. Mach. Learn. Res.,
references.authors: Sanvesh Srivastava Volkan Cevher Quoc Dinh David Dunson
references.title: Wasp: Scalable bayes via barycenters of subset posteriors
references.venue: In Artificial Intelligence and Statistics,
references.authors: Sanvesh Srivastava Cheng Li David B Dunson
references.title: Scalable bayes via barycenter in wasserstein space
references.venue: The Journal of Machine Learning Research,
references.authors: Zhengyu Su Yalin Wang Rui Shi Wei Zeng Jian Sun Feng Luo Xianfeng Gu
references.title: Optimal mass transport for shape matching and comparison
references.venue: IEEE transactions on pattern analysis and machine intelligence,
references.authors: G.J. Székely M.L. Rizzo N.K. Bakirov
references.title: Measuring and testing dependence by correlation of distances
references.venue: Ann. Statist.,
references.authors: M. Talagrand
references.title: The transportation cost from the uniform measure to the empirical measure in dimension ≥ 3
references.venue: Ann. Probab.,
references.authors: Hans Triebel
references.title: Theory of function spaces. III, volume 100 of Monographs in Mathematics
references.venue: Birkhäuser Verlag, Basel,
references.authors: Alexandre B. Tsybakov
references.title: Introduction to nonparametric estimation
references.authors: A.W. van der Vaart
references.title: Asymptotic statistics, volume 3 of Cambridge Series in Statistical and Probabilistic Mathematics
references.authors: Aad W. van der Vaart Jon A. Wellner
references.title: Weak convergence and empirical processes
references.authors: Cédric Villani
references.title: Topics in optimal transportation, volume 58 of Graduate Studies in Mathematics
references.venue: American Mathematical Society,
references.authors: Cédric Villani
references.title: Optimal transport, volume 338 of Grundlehren der Mathematischen Wissenschaften [Fundamental Principles of Mathematical Sciences
references.venue: Old and new
references.authors: Gilbert G. Walter
references.title: Approximation of the delta function by wavelets
references.venue: J. Approx. Theory,
references.authors: Jonathan Weed Quentin Berthet
references.title: Estimation of smooth densities in Wasserstein distance
references.venue: In Conference on Learning Theory,
sections.heading: 1 Introduction
sections.text: Given two random variables X ∼ µ and Y ∼ ν, where µ, ν are probability measures on Rd, d ≥ 1, the problem of finding a “nice" map T0(·) such that T0(X) ∼ ν has numerous applications in machine learning such as domain adaptation and data integration [34, 35, 38, 48, 61, 112], dimension reduction [12, 66, 90], generative models [60, 81, 88, 110], to name a few. Of particular interest is the case when T0(·) is obtained by minimizing a cost function, a line of work initiated by Gaspard Monge [97] in 1781 (see (1.1) below), in which case T0(·) is termed an optimal transport (OT) map and has applications in shape matching/transfer problems [29, 47, 107, 121], Bayesian statistics [46, 75, 80, 108], econometrics [15, 28, 45, 50, 54], nonparametric statistical inference [39– 41, 113, 114]; also see [111, 128, 129] for book-length treatments on the subject. In this paper, we will focus on the OT map obtained using the standard squared Euclidean cost function, i.e.,
  T0 := argmin T :T#µ=ν
  E‖X − T (X)‖2, (1.1)
  where T#µ = ν means T (X) ∼ ν for X ∼ µ. The estimation of T0 has attracted a lot of interest in recent years due to its myriad applications (as stated above) and interesting geometrical properties (see [19, 56, 91] and Definition 1.1 below). In practice, the main hurdle in constructing estimators for T0 is that the explicit forms of the measures µ, ν are unknown; instead only random samples
  X1, . . . , Xm ∼ µ and Y1, . . . , Yn ∼ ν
  35th Conference on Neural Information Processing Systems (NeurIPS 2021).
  are available. A natural strategy in this scenario is to estimate T0 using T̃m,n, where T̃m,n is computed as in (1.1) with µ and ν replaced by µ̃m and ν̃n which are empirical approximations of µ and ν based on X1, . . . , Xm and Y1, . . . , Yn respectively (see Definition 1.2). Such estimators are often called plug-in estimators and have been used extensively; see [7, 30, 67, 93, 94, 102, 116].
  The main goal of this paper is to study the rates of convergence of general plug-in estimators of T0 under a unified framework. We show that when µ̃m and ν̃n are chosen as µ̂m and ν̂n respectively, where µ̂m and ν̂n are the standard empirical distributions supported on m and n atoms, i.e.,
  µ̂m := 1
  m m∑ i=1 δXi and ν̂n := 1 n n∑ j=1 δYj , (1.2)
  T̃m,n (appropriately defined using Definition 1.2) converges at a rate of m−2/d + n−2/d for d ≥ 4 in the sense of (1.8). This rate happens to be minimax optimal under minimal smoothness assumptions (see [72, Theorem 6]) but suffers from the curse of dimensionality. We next show that, if µ and ν are known to admit sufficiently smooth densities, it is possible to apply kernel or wavelet based smoothing techniques on µ̂m and ν̂n to obtain plug-in estimators that mitigate the aforementioned curse of dimensionality.
  Our next contribution pertains to the estimation of W 22 (µ, ν) (the squared Wasserstein distance), see (1.3) below, a quantity of independent interest in statistics and machine learning with applications in structured prediction [51, 89], image analysis [18, 59], nonparametric testing [16, 106], generative modeling [10, 96], etc. In this paper, we also obtain rates of convergence for plug-in estimators W 22 (µ̃m, ν̃n) of W 2 2 (µ, ν). We show that kernel smoothing µ̂m and ν̂n can be used to obtain plug-in estimators of W 22 (µ, ν) that mitigate the curse of dimensionality as opposed to a direct plug-in approach using µ̂m and ν̂n (as used in [30, Theorem 2]). This provides an answer to the open question of estimating W 22 (µ, ν) when µ, ν admit smooth densities laid out in [30].
sections.heading: 1.1 Background on optimal transport
sections.text: In this section, we present some basic concepts and results associated with the OT problem that will play a crucial role in the sequel. Let Pac(Rd) denote the set of all Lebesgue absolutely continuous probability measures on Rd andP2(Rd) be the set of probability measures with finite second moments. Then the 2-Wasserstein distance (squared) between µ, ν ∈ P2(Rd) is defined as:
  W 22 (µ, ν) := min π∈Π(µ,ν)
  ∫ ‖x− y‖2 dπ(x, y), (1.3)
  where Π(µ, ν) is the set of probability measures on Rd×Rd with marginals µ and ν. The optimization problem in (1.3) is often called the Kantorovich relaxation (see [76, 77]) of the optimization problem in (1.1). The existence of a minimizer in (1.3) follows from [129, Theorem 4.1].
  Proposition 1.1 (Brenier-McCann polar factorization theorem, see [91, 128]). Suppose µ ∈ Pac(Rd). Then there exists a µ-a.e. (almost everywhere) unique function T0(·) : Rd → Rd, which is the gradient of a real-valued d-variate convex function, say ϕ0(·) : Rd → R, such that T0#µ = ν. Further, the distribution defined as π(A×B) = µ(A ∩ (T0)−1(B)) for all Borel sets A,B ⊆ Rd is the unique minimizer in (1.3) provided µ, ν ∈ P2(Rd).
  Definition 1.1 (OT map and potential function). The function T0 : Rd → Rd in Proposition 1.1 which satisfies T0#µ = ν will be called the OT map from µ to ν. A convex function ϕ0(·) in Proposition 1.1 satisfying∇ϕ0 = T0 will be termed an OT potential.
  The next and final important ingredient is the alternate dual representation of (1.3) which gives: 1
  2 W 22 (µ, ν) =
  1
  2
  ∫ ‖x‖2 dµ(x) + 1
  2
  ∫ ‖y‖2 dν(y)−min
  f∈F Sµ,ν(f), where (1.4) Sµ,ν(f) = ∫ f dµ+ ∫ f∗ dν. (1.5)
  Here F denotes the space of convex functions on Rd which are also elements of L1(µ) and f∗(·) is the standard Legendre-Fenchel dual defined as:
  f∗(x) := sup y∈Rd [y>x− f(y)], for x ∈ dom(f). (1.6)
sections.heading: 1.2 Estimating OT map via barycentric projection
sections.text: Recall the setting from the Introduction. Let µ̃m, ν̃n ∈ P2(Rd). Here µ̃m, ν̃n need not be absolutely continuous and can be very general. Intuitively, µ̃m and ν̃n can be viewed as some empirical approximation of µ and ν respectively.
  Example 1.2 (Simple choices of µ̃m and ν̃n). Let X1, . . . , Xm i.i.d.∼ µ and Y1, . . . , Yn i.i.d.∼ ν; in which case a natural choice would be to set µ̃m = µ̂m and ν̃n = ν̂n where µ̂m and ν̂n are the empirical distributions of X1, . . . , Xm and Y1, . . . , Yn respectively, as defined in (1.2). This is the standard choice adopted in the discrete-discrete Kantorovich relaxation; see [104, Section 2.3]. Another popular choice is µ̃m = µ̂m, ν̃n = ν or µ̃m = µ, ν̃n = ν̂n. This is the semi-discrete Kantorovich problem and is popular when one of the measures is fully specified; see [26, 55].
  A natural way to estimate T0(·), as defined in (1.1), would be to approximate it using the OT map from µ̃m to ν̃n. However as µ̃m and ν̃n may not be elements of Pac(Rd), Proposition 1.1 does not apply and an OT map may not exist from µ̃m to ν̃n. Such is the case in Example 1.2 in the discrete-discrete case when m 6= n. To circumvent this issue, we leverage the notion of barycentric projections (see [3, Definition 5.4.2]) defined below:
  Definition 1.2 (Barycentric projection). Define the set
  Γ̃min := argmin π∈Π(µ̃m,ν̃n)
  ∫ ‖x− y‖2 dπ(x, y).
  The optimization problem above is the plug-in analog of the optimization problem on the right hand side of (1.3). Given any γ ∈ Γ̃min, define the barycentric projection of γ as the conditional mean of y given x under γ, i.e.,
  T̃m,n(x) ≡ T̃ γm,n(x) := ∫ y y dγ(x, y)∫ y dγ(x, y) , for x ∈ supp (µ̃m) . (1.7)
  In general, Γ̃min need not be a singleton which is why we index the barycentric projection T̃ γm,n(·) by γ ∈ Γ̃min. Note that T̃ γm,n(·) need not be a transport map; however, if an OT map exists then it must be equal to T̃ γm,n(·) (µ̃m-a.e.). Our goal is to obtain stochastic upper bounds for
  sup γ∈Γ̃min ∫ ∥∥T̃ γm,n(x)− T0(x)∥∥2 dµ̃m(x). (1.8) In addition, our proof techniques also yield rates of convergence for∣∣W 22 (µ̃m, ν̃n)−W 22 (µ, ν)∣∣. (1.9) In this paper, we will focus on d ≥ 2. Due to the canonical ordering of R, the case d = 1 can be handled easily using the classical Hungarian embedding theorem [82].
sections.heading: 1.3 Contributions
sections.text: 1. We provide a new and flexible stability estimate Theorem 2.1 which yields a unified approach
  to obtaining rates of convergence for general plug-in estimators of the OT map T0(·). Unlike existing stability estimates, Theorem 2.1 holds for the barycentric projection (which is the same as the OT map when it exists) and does not require any smoothness assumptions on µ̃m, ν̃n or T̃ γm,n(·); also see Remark 2.1 for a comparison with the existing literature.
  2. In Sections 2.1 and 2.2, we use Theorem 2.1 to bound (1.8) and (1.9):
  • In Section 2.1, we show that in both the discrete-discrete and semi-discrete Kantorovich relaxation problems (see Example 1.2), the rate of convergence of (1.8) is m−2/d + n−2/d for d ≥ 4 when T0 is assumed to be Lipschitz (see Theorem 2.2), which is the minimax rate (see [72, Theorem 6]). To the best of our knowledge, rates of convergence for these natural estimators weren’t previously established in the literature.
  • In Section 2.2 and Appendix A, we show that the curse of dimensionality in the above rates can be mitigated provided µ and ν admit (uniform) Sobolev smooth densities (see Section 2.2) or Besov smooth densities (see Appendix A). In Section 2.2, our plug-in estimator is obtained by choosing µ̃m (and ν̃n) as the convolution of µ̂m (and ν̂n) and a smooth kernel with an appropriate bandwidth. Under this choice, the rate of convergence in (1.8) is m−( s+2 d ∧ 1 2 ) +
  n−( s+2 d ∧ 1 2 ), where s denotes the degree of Sobolev smoothness (see Theorem 2.5). Clearly, if 2(s + 2) ≥ d, the rate of convergence becomes dimension-free and mitigates the curse of dimensionality. We also show the same rates of convergence mentioned above hold for (1.9) (see e.g., Proposition 2.6) which makes a strong case in favor of incorporating smoothness in the construction of plug-in estimators as was conjectured in [30]. In Appendix A, our plug-in estimator is obtained using natural wavelet based density estimators. The rate of convergence in (1.8) turns out to be n− 1+s d+2s where s denotes the degree of Besov smoothness (see Theorem A.1). Note that by choosing s large enough, the exponent in the rate can be made arbitrarily close to 1/2, thereby reducing the curse of dimensionality.
  3. In Section 2.3, we use a discretization technique from [131] to construct discrete approximations to the smoothed µ̃m and ν̃n from the previous paragraph that in turn yield computable plug-in estimators for T0 (provided one can sample from µ̃m and ν̃n) that also achieve the same statistical guarantees as the smoothed plug-in estimator from Section 2.2 (see Theorem 2.7). However the number of atoms required in the discretizations and correspondingly the computational complexity increases with the degree of smoothness; this highlights a statistical and computational trade-off.
  4. We provide implications of our results in popular applications of OT such as estimating the barycenter of two multivariate probability distributions (see Theorem B.1 in Appendix B.1) and in nonparametric independence testing (see Theorem B.3 in Appendix B.2).
sections.heading: 1.4 Related work
sections.text: Many recent works have focused on obtaining consistent estimators of T0 using the plug-in principle, see [26, 55] (in the semi-discrete problem) and [41, 68, 132] (in the discrete-discrete problem). In [55], the authors studied the rate of convergence of the semi-discrete optimal transport map from ν (absolutely continuous) to µ̂m. This paper complements the aforementioned papers by studying the rates of convergence for general plug-in estimators in a unified fashion. In two other papers [9, Theorem 1.1] and [87, Section 4], the authors use a “Voronoi tessellation" approach to estimate T0, however the rates obtained in this paper, even in the absence of smoothness, are strictly better than those in [9, 87]. Perhaps the most closely related paper to ours would be [67]. In [67], the author uses variational techniques to arrive at stability estimates while we exploit the Lipschitz nature of the OT map (see Definition 1.1). Further the rates in this paper have exponents s+2d ∧ 1 2 which are strictly better than the exponents s+22(s+2)+d obtained in [67, Proposition 1] under the same smoothness assumptions (Sobolev type of order s, see Definition 2.4). In another line of work [72], the authors use theoretical wavelet based estimators (not of the plug-in type) of T0 to obtain nearly minimax optimal rates of convergence. However these estimators, by themselves, are not transport maps between two probability measures, which makes them harder to interpret. In contrast, our focus is on obtaining rates of convergence for plug-in estimators, which are transport maps between natural approximations of µ and ν. Such plug-in type strategies are a lot more popular in computational OT [7, 30, 67, 93, 94, 102, 116].
  In terms of obtaining rates of convergence for (1.9), some attempts include [109, 116] where parametric rates are obtained when µ, ν are known to be finitely supported or are both Gaussian. In a related problem, bounds for W 22 (µ̂m, µ) were obtained in [6, 42, 49, 100, 123, 131]. Using these bounds, for m = n, it is easy to get a n−1/d rate of convergence for (1.9). This rate was recently improved to n−2/d in [30] under no smoothness assumptions. Our rates coincide with the n−2/d rate from [30] under no smoothness assumptions. But further, we show in this paper that the curse of dimensionality in the above rate can be mitigated by incorporating smoothness into the plug-in procedure.
sections.heading: 2 Main results
sections.text: Recall ϕ0(·) from Definition 1.1. The following is our main result.
  Theorem 2.1 (Stability estimate). Suppose that µ, ν ∈ Pac(Rd) ∩ P2(Rd) and µ̃m, ν̃n ∈ P2(Rd). Assume that T0(·) (as defined in (1.1)) is L-Lipschitz (L > 0). Then,
  sup γ∈Γ̃min
  ∫ ‖T̃ γm,n(x)− T0(x)‖2 dµ̃m(x) ≤ Lmax {∣∣∣∣ ∫ Ψ∗µ̃m,ν̃n d(ν̃n − ν†m)∣∣∣∣, ∣∣∣∣ ∫ Ψ∗µ̃m,ν†m d(ν̃n − ν†m) ∣∣∣∣}
  + 2L ∫ ϕ∗0(y) d(ν̃n − ν†m)(y), (2.1)
  where ν†m := T0#µ̃m, ϕ ∗ 0(·) is defined as in (1.6), and with S·,·(·) defined as in (1.5), Ψµ̃m,ν̃n(·) := argminf∈F Sµ̃m,ν̃n(f), Ψµ̃m,ν†m(·) := argminf∈F Sµ̃m,ν†m(f), and D denotes the space of realvalued convex functions on Rd.
  The proof of Theorem 2.1 (see Appendix C.1) starts along the same lines as the proof of the curvature estimate in [56, Proposition 3.3]. This is followed by some careful manipulations of W 22 (·, ·) (as in (1.3)) and an application of the conditional version of Jensen’s inequality, see (C.3). The final step of the proof uses the dual representation in (1.4) with techniques similar to some intermediate steps in the proof of [92, Proposition 2] and [30, Lemma 3].
  Remark 2.1 (Comparison with other stability estimates). Theorem 2.1 provides some important advantages to existing stability estimates in the literature. One of the earliest results in this direction can be found in [56, Proposition 3.3] but their bound involves a push-forward constraint which makes it hard to use for rate of convergence analysis. A bound similar to Theorem 2.1 is presented in [55, Lemma 5.1] but there the authors assume the existence of an OT map from µ̃m to ν̃n. Therefore, it does not apply to the discrete-discrete problem where µ̃m = µ̂m and ν̃n = ν̂n with m 6= n. Overcoming all these limitations is an important contribution of Theorem 2.1 and allows us to deal with popular plug-in estimators all in one go. The stability estimate in [72, Proposition 10] on the other hand requires µ̃m, ν̃n to be sufficiently smooth and hence it does not hold for discrete-discrete or semi-discrete plug-in estimators (see Example 1.2). Further their result requires all the measures involved to be compactly supported unlike the much milder requirements of Theorem 2.1. However, a shortcoming of Theorem 2.1 is that it is hard to obtain rates faster than n−1/2 using it directly, whereas [72] can obtain rates arbitrarily close to n−1. This is a price we pay for analyzing natural and popular plug-in estimators as opposed to the (more intractable) wavelet based estimators in [72].
  Remark 2.2 (How to use Theorem 2.1 to obtain rates of convergence?). Note that the second term on the right hand side of (2.1), under appropriate moment assumptions, is Op(m−1/2 +n−1/2) (free of dimension) by a direct application of Markov’s inequality. We therefore focus on the first term. By (1.5), Ψ∗µ̃m,ν̃n(·), Ψ ∗ µ̃m,ν † m
  (·) ∈ F . Further, by Caffarelli’s regularity theory [20–22], depending on the “smoothness" of µ̃m, ν̃n, it can be shown that there exists a further class of functions Fs (see Remarks 2.3 and 2.6) such that Ψ∗µ̃m,ν̃n(·), Ψ ∗ µ̃m,ν † m
  (·) ∈ F ∩ Fs. Thus, we can bound the first term on the right hand side of (2.1) as:
  max {∣∣∣∣ ∫ Ψ∗µ̃m,ν̃n d(ν̃n − ν†m)∣∣∣∣, ∣∣∣∣ ∫ Ψ∗µ̃m,ν†m d(ν̃n − ν†m) ∣∣∣∣} ≤ sup
  f∈F∩Fs ∣∣∣∣ ∫ f d(ν̃n − ν†m)∣∣∣∣. (2.2) The right hand side of (2.2) can now be bounded using the corresponding Dudley’s entropy integral bounds using empirical process techniques; see [126, Lemmas 19.35-19.37].
  To conclude, the two main steps in our strategy are identifying the family of functions Fs and computing Dudley’s entropy integral. Further, the more the smoothness of µ̃m, ν̃n, the smaller is the class of functions Fs and smaller the supremum on the right hand side of (2.2). This shows why better rates can be expected under smoothness assumptions.
sections.heading: 2.1 Natural non-smooth plug-in estimators
sections.text: In this case, we discuss the rates of convergence for the discrete-discrete problem and the semi-discrete problem, where no smoothness is available on µ̃m and ν̃n.
  Theorem 2.2. Suppose that T0(·) is L-Lipschitz, ν is compactly supported and E exp(t‖X1‖α) <∞ for some t > 0, α > 0.
  (Discrete-discrete): Set µ̃m = µ̂m and ν̃n = ν̂n. Then the following holds:
  sup γ∈Γ̃min
  ∫ ‖T̃ γm,n(x)− T0(x)‖2 dµ̃m(x) = Op ( r (m,n) d × (log (1 + max{m,n})) td,α ) , (2.3)
  where r(m,n)d :=  m−1/2 + n−1/2 for d = 2, 3, m−1/2 log (1 +m) + n−1/2 log (1 + n) for d = 4, m−2/d + n−2/d for d ≥ 5,
  (2.4)
  and
  td,α :=  (4α)−1(4 + ((2α+ 2dα− d) ∨ 0)) for d < 4, (α−1 ∨ 7/2)− 1 for d = 4, 2(1 + d−1) for d > 4.
  The same bound holds for |W 22 (µ̃m, ν̃n)−W 22 (µ, ν)| without assuming T0(·) is Lipschitz. (Semi-discrete): Set µ̃m = µ, ν̃n = ν̂n or µ̃m = µ̂m, ν̃n = ν. Then the left hand side of (2.3) is Op(r (n,n) d × (log (1 + n))td,α) or Op(r (m,m) d × (log (1 +m))td,α) respectively.
  A stronger result can be proved if both µ and ν are compactly supported.
  Corollary 2.3. Consider the setting from Theorem 2.2 and assume further that µ is compactly supported. Then, with r(m,n)d defined as in (2.4), we have:
  E [ sup
  γ∈Γ̃min
  ∫ ‖T̃ γm,n(x)− T0(x)‖2 dµ̃m(x) ] ≤ Cr(m,n)d ,
  for some constant C > 0, in both the discrete-discrete and semi-discrete settings from Theorem 2.2.
  A brief description of the proof technique of Theorem 2.2 using Theorem 2.1 is provided in Remark 2.3 below, and the actual proof is presented in Appendix C.1.
  Remark 2.3 (Proof technique). The proof of Theorem 2.2 proceeds via the strategy outlined in Remark 2.2. We first show that Fs (see Remark 2.2) can be chosen as a certain sub-class of convex functions which are in L2(ν). We then use Dudley’s entropy integral type bounds which in turn requires the bracketing entropy [126, Page 270] of Fs, recently proved in [83, Equation 26]. This strategy is slightly different from that used in the proof of [30, Theorem 2], where the authors assume that µ is compactly supported whereas we only assume the finiteness of E exp(t‖X1‖α) for some t > 0, α > 0. The compactness assumption on µ allows one to further restrict Fs to the class of Lipschitz functions. This additional restriction does not seem to be immediate without the compactness assumption.
  As discussed in Section 1.3, the exponents obtained in Theorem 2.2 are minimax optimal, up to multiplicative logarithmic factors, under bare minimal smoothness assumptions (see [72, Theorem 6]). To the best of our knowledge, rates for the discrete-discrete case for m 6= n and those for the semi-discrete case were not known previously in the literature. Our rates are also strictly better than those (for different estimators, based on space tessellations) obtained in [9, 87] and require less stringent assumptions than those in [30]. In the next section, we show how smoothness assumptions can be leveraged to mitigate the curse of dimensionality in Theorem 2.2.
sections.heading: 2.2 Smooth kernel based plug-in estimator: mitigating the curse of dimensionality
sections.text: In this section, we focus on kernel based density estimators for the probability densities associated with µ and ν (see [57, 58, 99, 103, 115]). We will show, using Theorem 2.1, that the corresponding estimators of T0(·) achieve (near) dimension-free rates under sufficient smoothness assumptions. We first introduce the Sobolev class of functions which we will exploit in this subsection to construct estimators that achieve rates of convergence which mitigate the curse of dimensionality under sufficient smoothness.
  Definition 2.4 (Uniform Sobolev class of functions). Let Ω ⊆ Rd and f(·) be uniformly continuous on Ω and admits uniformly continuous derivatives up to order s on Ω for some s ∈ N. For any m := (m1, . . . ,md) ∈ Nd, let
  ∂mf := ∂
  ∂m1x1 . . .
  ∂
  ∂mdxd f, |m| := d∑ i=1 mi.
  For any k ≤ s, we further define, ‖f‖Ck(Ω) := ∑ |m|≤k ‖∂mf‖L∞(Ω).
  The space Cs(Ω) is defined as the set of functions f(·) for which ‖f‖Ck(Ω) <∞ for all k ≤ s.
  For this subsection, assume that µ and ν admit Sobolev smooth densities fµ(·) and fν(·) in the uniform norm (see Definition 2.4 above). Given Ω ⊆ Rd and s ∈ N, let Cs(Ω) denote the set of Sobolev smooth functions on Ω of order s.
  Assumption (A1) (Regularity of the densities). Suppose that
  1. fµ and fν are supported on compact and convex subsets of Rd, say X and Y respectively. 2. There exists s,M > 0 such that fµ(·) ∈ Cs(X ;M) and fν(·) ∈ Cs(Y;M) where Cs(X ;M) is
  the space of real valued functions supported on X such that for all f(·) ∈ Cs(X ;M), we have M−1 ≤ f(x) ≤ M for all x ∈ X and ‖f‖Cs(X ) ≤ M . Here ‖·‖Cs(X ) is the standard uniform Sobolev norm as defined in Definition 2.4. The space Cs(Y;M) is defined analogously.
  We now define our estimators for fµ(·) and fν(·) using the standard kernel density estimation technique (see [125, Section 1.2]). Set
  f̂µ(x) := 1
  mhdm m∑ i=1 Kd ( Xi − x hm ) , (2.5)
  for some bandwidth parameter hm > 0 and d-variate kernel Kd(·). We assume that Kd(·) is the dfold product of univariate kernels, i.e., there exists a kernelK(·) such that for u = (u1, . . . , ud) ∈ Rd, Kd(u) = ∏d i=1K(ui). We define f̂ν(·) similarly with the same univariate kernel and bandwidth.
  Assumption (A2) (Regularity of the kernel). Assume that K(·) is a symmetric, bounded, s+ 1 times differentiable kernel on Rd with all s+ 1 derivatives bounded and integrable. Further, suppose that K(·) is of order 2s+ 2, i.e.,∫ ujK(u) du = 1(j = 0), for j = {0, 1, 2, . . . , 2s+1}, and ∫ |u|2s+2|K(u)| du <∞.
  The above assumptions on K(·) are standard for estimating smooth densities and their derivatives of different orders in the kernel density estimation literature; see e.g. [4, 57, 58, 69, 125]. There are several natural ways to construct kernels satisfying Assumption (A2), see [125, Section 1.2.2]; an example is also provided in Example 2.4 below.
  Example 2.4 (Example of a kernel satisfying Assumption (A2)). Let ψm(·) be the m-th Hermite polynomial on R (see [84]). Then the kernel function defined as
  K(u) := 2s+2∑ m=0 ψm(0)ψm(u) exp(−u2/2)
  satisfies Assumption (A2).
  It is evident from Assumption (A2) that K(·) may take some negative values, in which case, f̂µ(·) (respectively f̂ν(·)) may not be a probability density. Consequently the barycentric projection (see Definition 1.2) between f̂µ(·) and f̂ν(·) is not well-defined. We get around this by projecting
  f̂µ(·) and f̂ν(·) on an appropriate space of “smooth" probability densities (see (2.6)), via an integral probability metric (see Definition 2.5 below; also see [98, 105, 117] for examples, computational procedures and applications of such metrics).
  Definition 2.5 (Integral probability metric). Given a classH of bounded functions on Rd and two probability densities g1(·) and g2(·) on Rd, the integral probability metric/distance between g1(·) and g2(·) with respect toH is defined as
  dIP(g1, g2;H) := sup ψ(·)∈H ∣∣∣∣ ∫ ψ(x)(g1(x)− g2(x)) dx∣∣∣∣. Sufficient conditions onH for dIP(·, ·;H) to be a metric on the space of probability measures (not on the space of probability densities as they can be altered on set of Lebesgue measure 0 without altering the underlying probability measures) on Rd have been discussed in [98]. Observe that the measure dIP(g1, g2;H) is well defined even when g1(·) and g2(·) are not probability densities.
  In Theorem 2.5 below, we useH = Cs+2(X ,M ′). Note that any function in Cs+2(X ,M ′) can be extended to a function in Cs+2(Rd;M ′) (see [72, Theorem 23] and [124, Theorem 1.105]). The fact that this choice of F results in a metric follows from the argument in [98, Page 8].
  We are now in a position to describe the projection estimators for fµ(·) and fν(·), and the rates achieved by the corresponding plug-in estimator.
  Theorem 2.5. Assume that T0(·) is L-Lipschitz and fµ, fν are Lebesgue densities satisfying Assumption (A1). Also suppose that K(·) satisfies Assumption (A2). Define hm := m− 1 d+2s logm , hn := n − 1d+2s log n and T := ∫ |Kd(u)| du+ 1. Fix any M ′ > 0. Consider any probability density f̃M ′
  µ (·) ∈ Cs(X ;TM) (where M is defined as in Assumption (A1)) which satisfies
  dIP
  ( f̃M ′
  µ , f̂µ;C s+2(X ;M ′) ) ≤ inf f(·)∈Cs(X ;TM) f≥0, ∫ f=1 dIP ( f̂µ, f ;C s+2(X ;M ′) ) + r (m,n) d,s (2.6)
  where r(m,n)d,s is defined as in (2.7) and dIP(·, ·;Cs+2(X ;M ′)) is the integral probability metric defined in Definition 2.5. We define f̃M ′
  ν (·) analogously as in (2.6) with X , f̂µ(·) replaced by Y , f̂ν(·). Then the following conclusions hold.
  1. SetM ′ := 8(1+TM). If µ̃m and ν̃n are the probability measures corresponding to the probability densities f̃M ′ µ (·) and f̃M ′ ν (·), then the following holds for some constant C > 0:
  E [ sup
  γ∈Γ̃min
  ∫ ‖T̃ γm,n(x)− T0(x)‖2 dµ̃m(x) ] ≤ Cr(m,n)d,s ,
  where r(m,n)d,s :=  m−1/2 + n−1/2 for d < 2(s+ 2), m−1/2 (log (1 +m)) d + n−1/2 (log (1 + n))
  d for d = 2(s+ 2), m− s+2 d + n− s+2 d for d ≥ 2(s+ 2). (2.7)
  The same bound also holds for E|W 22 (µ̃m, ν̃n)−W 22 (µ, ν)|.
  2. f̂µ(·) satisfies
  lim n→∞ max
  { P ( ‖f̂µ‖Cs(X̃ ) ≥ TM ) ,P ( sup x∈X̃ |f̂µ(x)− fµ(x)| ≥ ε )} = 0 (2.8)
  for any ε > 0, where X̃ is any compact subset of X o. The same conclusion holds for f̂ν(·) with X replaced by Y .
  In Theorem 2.5, we have shown that the plug-in estimator for T0(·) using f̃M ′ µ (·) and f̃M ′
  ν (·) (with M ′ = 8(1 + TM)) achieves rates that mitigate the curse of dimensionality under sufficient smoothness. In fact, f̃M ′
  µ (·) can be viewed as an approximate minimizer of dIP(f̂µ, ·;Cs+2(X ,M ′)) over an appropriate class of Sobolev smooth probability densities. This is carried out because f̂µ(·) by itself may not be a probability density.
  Further note that µ̃m, ν̃n as specified in Theorem 2.5 are both smooth, and consequently Γ̃min is a singleton and the supremum in Theorem 2.5 can be dropped. A brief description of the proof technique for Theorem 2.5 is presented in Remark 2.6 below and the actual proof is given in Appendix C.1.
  Remark 2.6 (Proof technique). The proof of Theorem 2.5 proceeds along the same lines as Remark 2.3. We first show that Fs (see Remark 2.2) can be chosen as a certain subset of Cs+2(Y◦). We then use Dudley’s entropy integral type bounds which in turn requires the bracketing entropy [126, Page 270] of the class of compactly supported Sobolev smooth functions which can be found in [127, Corollary 2.7.2].
  We now explain the implications of both the parts of Theorem 2.5 in the following two remarks.
  Remark 2.7 (Mitigating the curse of dimensionality). Theorem 2.5 shows that, under enough smoothness, i.e., when 2(s+ 2) > d, both the upper bounds for (1.8) and (1.9) are Op(n−1/2). This shows that, for large dimensions, provided µ and ν admit smooth enough densities, it is possible to construct plug-in estimators that mitigate the curse of dimensionality. Note that a similar estimator was analyzed in [67, Proposition 1] when m = n. However, the rates obtained in Theorem 2.5 are strictly better than those in [67, Proposition 1]. For m = n, when d < 2(s + 2), [67] obtained a rate of n− s+2 2(s+2)+d which is worse than n−1/2 obtained in Theorem 2.5. For the other regimes, [67] obtains rates (up to log factors) of n−1/4 and n− 1
  (s+2)(d+2(s+2)) which are both worse than the respective rates of n−1/2 and n− s+2 d in Theorem 2.5.
  Remark 2.8 (Computational aspects of Theorem 2.5). Note that f̃M ′ µ (·) (with M ′ = 8(1 + TM)) is hard to compute whereas f̂µ(·) is computable easily in linear time. Note that if f̂µ(·) itself were a probability density in Cs(X ;TM), then we would have f̂µ = f̃M ′
  µ . While Theorem 2.5 does not establish that, it does come close in part 2, from which we can easily derive the following:
  lim n→∞
  P(f̂µ(·) /∈ Cs(X̃ ;TM)) = 0.
  The above shows that f̂µ(·) is indeed bounded below by (TM)−1 on X̃ (any compact subset of the interior of X ), and additionally belongs to Cs(X̃ ;TM) with probability converging to 1. This leads us to conjecture that the natural density version of f̂µ(·), i.e.,
  max{f̂µ(·), 0}∫ max{f̂µ(x), 0} dx
  should serve as a good proxy for f̃M ′
  µ (·) and lead to rates of convergence that mitigate the curse of dimensionality. From a computational perspective, the density specified above is easy to simulate from using an accept-reject algorithm without computing the integral in the denominator (see [101, Algorithm 4.3]). However, our current proof technique does not provide rates of convergence for the above density estimator based on f̂µ(·).
  Another important implication of Theorem 2.5 is the bound obtained on |W2(µ̃m, ν̃n)−W2(µ, ν)| when µ 6= ν. We first present the result and then describe the implication.
  Proposition 2.6. Consider the setting in Theorem 2.5. Then, provided µ 6= ν, the following holds:
  |W2(µ̃m, ν̃n)−W2(µ, ν)| = Op(r(m,n)d,s ).
  Proposition 2.6 (see Appendix C.1 for a proof) shows an interesting distinction between the µ 6= ν case and the µ = ν case. For µ = ν, the best possible exponent is n− 1+s 2s+d for d ≥ 3 (see [131, Theorem 3] where the result was established under more general Besov smoothness assumptions). On the contrary, when µ 6= ν, Proposition 2.6 establishes a rate of n− s+2d for the Wasserstein distance which is strictly better than the minimax achievable rate mentioned above when µ = ν. This observation complements [30, Corollary 1] where the authors make a similar remark for the special case of s = 0.
sections.heading: 2.3 Discretized plug-in estimator under smoothness assumptions
sections.text: In Section 2.1, we discussed how smoothness can be incorporated into the plug-in procedure to get faster rates of convergence. Such plug-in estimators are popular in the computational OT literature (see [7, 8, 25, 36]). However, even after f̃µ(·) ≡ f̃M ′ µ (·), f̃ν(·) ≡ f̃M ′
  ν (·) are calculated, T̃ γm,n as in Theorem 2.5 cannot be computed explicitly from data if f̃µ(·) and f̃ν(·) are continuous densities. This is in contrast to T̃ γm,n from Theorem 2.2 in the discrete-discrete case which is explicitly computable using a standard linear program, but achieves worse rates of convergence. This is not unexpected. Thanks to the no free lunch principle, better statistical accuracy is naturally accompanied by heavier computational challenges. Therefore, our goal here is to construct estimators, under smoothness assumptions as in Section 2.2, which are computable in polynomial time (with complexity increasing with smoothness) provided f̃µ(·) and f̃ν(·) can be sampled from, and also attain rates that mitigate the curse of dimensionality.
  Construction: We will illustrate the discretized estimator using the kernel based estimator from Section 2.2. Similar results also hold for the wavelet based estimator from Appendix A. Recall the kernel density estimators f̃µ(·) and f̃ν(·) (see (2.6)). Sample M ≥ 1 random points from both f̃µ(·) and f̃ν(·). Let µ̂m,M and ν̂n,M denote the standard empirical measures on the M points sampled from f̃µ(·) and f̃ν(·) respectively. Finally construct T̃m,n ≡ T̃ γm,n as in Definition 1.2 with µ̃m = µ̂m,M and ν̃n = ν̂n,M . It should be pointed out that a similar construction was also used in [131, Section 6] for estimating probability densities under the Wasserstein loss. Based on this construction, the main result of this section is as follows:
  Theorem 2.7. Consider the setting in Theorem 2.5 and the same construction of T̃ γm,n as above. For simplicity, let’s also assume m = n. Accordingly set M = n s+2 2 . Then Γ̃min is a singleton and consequently the following conclusion holds for some constant C > 0:
  E [∫ ‖T̃m,n(x)− T0(x)‖2 dµ̃m(x) ] ≤ Cr(n,n)d,s .
  The same rates also hold for E|W 22 (µ̃m, ν̃n)−W 22 (µ, ν)|.
  The proof of Theorem 2.7 is given in Appendix C.1. Once the empirical measures µ̂m,M and ν̂n,M have been obtained, an explicit computation of T̃m,n as described above requires O(M3) = O(n 3(s+2) 2 ) steps using the Hungarian algorithm, see [73]. This highlights the statistical versus computational trade-off, i.e., in order to mitigate the curse of dimensionality in convergence rates by exploiting smoothness, the computational complexity gets progressively worse by polynomial factors in n. It should be mentioned that (approximate) algorithms faster than the Hungarian algorithm stated above, can be found in [1, 36, 53] to name a few. Due to space constraints, we avoid a detailed discussion on this.
  In the above construction, sampling from the smoothed kernel densities f̃µ(·) and f̃ν(·) is crucial. If we would simply draw M bootstrap samples from the empirical distributions µ̂m and ν̂n, the rates of convergence wouldn’t improve from those observed in Theorem 2.2 no matter how large M is.
  Acknowledgments and Disclosure of Funding We would like to thank the reviewers for their constructive suggestions that greatly helped improve the quality of the paper. The third author is supported by NSF Grant DMS-2015376.
title: Rates of Estimation of Optimal Transport Maps using Plug-in Estimators via Barycentric Projections