Ai Security Projects

This project focuses on ai security projects using modern AI and machine learning techniques. The content below is adapted from research literature and practical implementation notes.

Minimax convergence rate for estimating the Wasserstein

barycenter of random measures on the real line

Jérémie Bigot1 , Raúl Gouet2 , Thierry Klein3 & Alfredo López4

Institut de Mathématiques de Bordeaux et CNRS (UMR 5251)1

Université de Bordeaux

Depto. de Ingenierı́a Matemática and CMM (CNRS, UMI 2807)2

Universidad de Chile

ENAC- Ecole nationale de l’aviation civile

et Institut de Mathématiques de Toulouse et CNRS (UMR 5219)3

Université de Toulouse

CSIRO Chile International Centre of Excellence4

November 28, 2021

Abstract

This paper is focused on the statistical analysis of probability measures ν 1 , . . . , ν n on

R that can be viewed as independent realizations of an underlying stochastic process. We

consider the situation of practical importance where the random measures ν i are absolutely continuous with densities f i that are not directly observable. In this case, instead of the densities, we have access to datasets of real random variables (Xi,j )1≤i≤n; 1≤j≤pi organized in the form of n experimental units, such that Xi,1 , . . . , Xi,pi are iid observations sampled from a random measure ν i for each 1 ≤ i ≤ n. In this setting, we focus on first-order statistics methods for estimating, from such data, a meaningful structural mean measure.

For the purpose of taking into account phase and amplitude variations in the observations, we argue that the notion of Wasserstein barycenter is a relevant tool. The main contribution of this paper is to characterize the rate of convergence of a (possibly smoothed) empirical

Wasserstein barycenter towards its population counterpart in the asymptotic setting where

both n and min1≤i≤n pi may go to infinity. The optimality of this procedure is discussed from the minimax point of view with respect to the Wasserstein metric. We also highlight the connection between our approach and the curve registration problem in statistics. Some numerical experiments are used to illustrate the results of the paper on the convergence rate of empirical Wasserstein barycenters.

Keywords: Wasserstein space; Fréchet mean; Barycenter of probability measures; Functional

data analysis; Phase and amplitude variability; Smoothing; Minimax optimality.

AMS classifications: Primary 62G08; secondary 62G20.

1 Introduction

In this paper, we are concerned with the statistical analysis of a set of absolutely continuous measures ν 1 , . . . , ν n on the real line R, with supports included in (a possibly unbounded) interval Ω ⊂ R, that can be viewed as independent copies of an underlying random measure ν. In this setting, it is of interest to define and estimate a mean measure ν0 of the random probability measure ν. The notion of mean or averaging depends on the metric that is chosen to compare elements in a given data set. In this work, we consider the Wasserstein metric dW associated to

the quadratic cost for the comparison of probability measures and we define ν0 as the population

Wasserstein barycenter of ν, given by

ν0 = arg min E d2W (ν, µ) , µ∈W2 (Ω)

where the above expectation is taken with respect to the distribution of ν, and W2 (Ω) denotes the space of probability measures with support included in Ω and with finite second moment. A Wasserstein barycenter corresponds to the Fréchet mean [Fré48] that is an extension of the usual Euclidean mean to non-linear metric spaces. Throughout the paper, the population mean measure ν0 is also referred to as the structural mean of ν, which is a terminology borrowed from curve registration (see [ZM11] and references therein).

Data sets leading to the analysis of absolutely continuous measures appear in various re-

search fields. Examples can be found in neuroscience [WS11], demographic and genomics studies [Del11, ZM11], economics [KU01], as well as in biomedical imaging [PM15]. Nevertheless, in such applications, one does not directly observe raw data in the form of absolutely continuous measures. Indeed, we generally only have access to random observations sampled from different distributions that represent independent subjects or experimental units. Thus, we propose to study the estimation of the structural mean measure ν0 (the popula-

tion Wasserstein barycenter) from a data set consisting of independent real random variables (Xi,j )1≤i≤n; 1≤j≤pi organized in the form of n experimental units, such that (conditionally on ν i ) the random variables Xi,1 , . . . , Xi,pi are iid observations sampled from the measure ν i with density f i , where pi denotes the number of observations for the i-th subject or experimental unit. The main purpose of this paper is to propose nonparametric estimators of the structural mean measure ν0 and to characterize their rates of convergence with respect to the Wasserstein

metric in the asymptotic setting, where both n and min1≤i≤n pi may go to infinity.

1.1 Main contributions

Two types of nonparametric estimators are considered in this paper. The first one isPgiven by the empirical Wasserstein barycenter of the set of measures ν̃ 1 , . . . , ν̃ n , with ν̃ i = p1i pj=1 i δXi,j for 1 ≤ i ≤ n. This estimator will be referred to as the non-smoothed empirical Wasserstein barycenter. Alternatively, since the unknown probability measures ν i are supposed to be abso- lutely continuous, a second estimator is based on a preliminary smoothing step which consists in

using standard kernel smoothing to construct estimators f̂ i of the unknown densities f i for each 1 ≤ i ≤ n. Then, an estimator of ν0 is obtained by taking the empirical Wasserstein barycenter ˆ R of the measures ν̂i , . . . , ν̂n , with ν̂i (A) := A fi (x)dx, A ⊂ R measurable. We refer to this class

of estimators as smoothed empirical Wasserstein barycenters whose smoothness depend on the choice of the bandwidths in the preliminary kernel smoothing step. The rates of convergence of both types of estimators are derived for their (squared) Wasser- stein risks, defined as their expected (squared) Wassertein distances from ν0 , and their optimality is discussed from the minimax point of view. Finally, some numerical experiments with simulated data are used to illustrate these results.

1.2 Related work in the literature

The notion of barycenter in the Wasserstein space, for a finite set of n probability measures supported on Rd (for any d ≥ 1), has been recently introduced in [AC11] where a detailed characterization of such barycenters in terms of existence, uniqueness and regularity is given using arguments from duality and convex analysis. However, the convergence (as n → +∞) of such Wasserstein barycenters is not considered in that work. In the one dimensional case (d = 1), computing the Wasserstein barycenter of a finite set of

probability measures simply amounts to averaging (in the usual way) their quantile functions.

In statistics, this approach has been referred to as quantile synchronization [ZM11]. In the

presence of phase variability in the data, quantile synchronization is known to be an appropriate alternative to the usual Euclidean mean of densities to compute a structural mean density that is more consistent with the data. Various asymptotic properties of quantile synchronization are studied in [ZM11] in a statistical model and asymptotic setting similar to that of this paper with min1≤i≤n pi ≥ n. However, other measures of risk than the one in this paper are considered in

[ZM11], but the optimality of the resulting convergence rates of quantile synchronization is not discussed. The results of this paper are very much connected with those in [PZ16] where a new frame- work is developed for the registration of multiple point processes on the real line for the purpose of separating amplitude and phase variation in such data. In [PZ16], consistent estimators of the structural mean of multiple point processes are obtained by the use of smoothed Wasser-

stein barycenters with an appropriate choice of kernel smoothing. Also, rates of convergence of such estimators are derived for the Wasserstein metric. The statistical analysis of multiple point processes is very much connected to the study of repeated observations organized in sam- ples from independent subjects or experimental units. Therefore, some of our results in this paper on smoothed empirical Wasserstein barycenters are built upon the work in [PZ16]. Nev- ertheless, novel contributions include the derivation of an exact formula to compute the risk

of non-smoothed Wasserstein barycenters in the case of samples of equal size, and new upper bounds on the rate of convergence of the Wasserstein risk of non-smoothed and smoothed em- pirical Wasserstein barycenters, together with a discussion of their optimality from the minimax point of view.

The construction of consistent estimators of a population Wasserstein barycenter for semi-

parametric models of random measures can also be found in [BK16] and [BLGL15], together with a discussion on their connection to the well known curve registration problem in statistics

[RL01, WG97].

1.3 Organization of the paper

In Section 2, we first briefly explain why using statistics based on the Wasserstein metric is a relevant approach for the analysis of a set of random measures in the presence of phase and amplitude variations in their densities. Then, we introduce a deformable model for the registration of probability measures that is appropriate to study the statistical properties of empirical Wasserstein barycenters. The two types of nonparametric estimators described above are finally introduced at the end of Section 2. The convergence rates and the optimality of

these estimators are studied in Section 3. Some numerical experiments with simulated data are proposed in Section 4 to highlight the finite sample performances of these estimators. Section

5 contains a discussion on the main contributions of this work and their potential extensions.

The proofs of the main results are gathered in a technical Appendix. Finally, note that we use bold symbols f , ν, . . . to denote random objects (except real random variables).

2 Wasserstein barycenters for the estimation of the structural

mean in a deformable model of probability measures

2.1 The need to account for phase and amplitude variations

To estimate a mean measure from the data (Xi,j )1≤i≤n; 1≤j≤pi , a natural approach is the fol- lowing one. In a first step, one uses the Xi,j ’s to compute estimators f̂ 1 , . . . , f̂ n (e.g. via kernel smoothing) of the unobserved density functions f 1 , . . . , f n of the measures ν 1 , . . . , ν nP . Then, an estimator of a mean density might be defined as the usual Euclidean mean f̄ n = n ni=1 f̂ i , 1

which is also classically referred to as the cross-sectional mean in curve registration.PAt the level of measures, it corresponds to computing the arithmetical mean measure ν̄ n = n1 ni=1 ν̂ i . The Euclidean mean f̄ n is to the Fréchet mean of the f̂ i ’s with respect to the usual squared distance in the Hilbert space L2 (Ω) of square integrable functions on Ω. Therefore, it only accounts for linear variations in amplitude in the data. However, as remarked in [ZM11], in many applica-

tions, it is often of interest to also incorporate an analysis of phase variability (i.e. time warping) in such functional objects, since it may lead to a better understanding of the structure of the data. In such settings, the use of the standard squared distance in L2 (Ω) to compare density functions ignores a possible significant source of phase variability in the data. To better account for phase variability in the data, it has been proposed in [ZM11] to introduce the so-called method of quantile synchronization as an alternative to the cross sectional

⊕ mean f̄ n . It amounts to computing the mean measure ν ⊕ n (and, if it exists, its density f n ) whose quantile function is n

− 1X −

F̄ n = Fi , (2.1)

where F − i denotes the quantile function of the measure ν i with density f i .

The statistical analysis of quantile synchronization, as studied in [ZM11], complements the

quantile normalization method originally proposed in [BIAS03] to align density curves in mi- croarray data analysis. This method is therefore appropriate for the registration of density functions and the estimation of phase and amplitude variations as explained in details in [PZ16].

Let us now assume that ν 1 , . . . , ν n are random elements taking values in the set of absolutely continuous measures contained in W2 (Ω). In this setting, it can be checked (see e.g. Proposition

2.1 below) that quantile synchronization corresponds to computing the empirical Wasserstein

barycenter of the random measures ν 1 , . . . , ν n , namely n

1X 2

ν⊕ n = arg min dW (ν i , µ). µ∈W2 (Ω) n i=1

Therefore, the notion of averaging by quantile synchronization corresponds to using the Wasser-

stein distance dW to compare probability measures, which leads to a notion of measure averaging that may better reflect the structure of the data than the arithmetical mean in the presence of phase and amplitude variability.

To illustrate the differences between using Euclidean and Wasserstein distances to account

for phase and amplitude variation, let us assume that the measures ν 1 , . . . , ν n have densities f 1 , . . . , f n obtained from the following location-scale model: we let f0 be a density on R having a finite second moment and, for (ai , bi ) ∈ (0, ∞) × R, i = 1, . . . , n a given sequence of independent random variables, we define

f i (x) := a−1 −1 i f0 ai (x − bi ) , x ∈ R, 1 ≤ i ≤ n. (2.2)

The sources of variability of the densities from model (2.2) are the variation in location along the x-axis, and the scaling variation. In Figure 1(a), we plot a sample of n = 100 densities from model (2.2) with f0 being the standard Gaussian density, ai ∼ U ([0.8, 1.2]) and bi ∼ U ([−2, 2]), where U ([x, y]) denotes the uniform distribution on the interval [x, y]. In this numerical experiment, there is more variability in phase (i.e. location) than in amplitude (i.e. scaling), which can also

be observed at the level of quantile functions as shown by Figure 1(b). In the location-scale model (2.2), it can be checked, e.g. using the quantile averaging formula (2.1), that the empirical Wasserstein barycenter ν ⊕ n is the probability measure with density

f⊕ −1 n (x) = ān f0 ān (x − b̄n ) ,

where ān = n1 ni=1 ai and b̄n = n1 ni=1 bi . Hence, if we assume that E(a1 ) = 1 and E(b1 ) = 0,

P P

it follows that d2W (ν ⊕ ⊕ n , ν0 ) converges almost surely to 0 as n → ∞, meaning that ν n is a consis- tent estimator of ν0 as shown by Figure 1(f). On the contrary, the arithmetical mean measure ν̄ n is clearly not a consistent estimator of ν0 , as it can be observed in Figure 1(d).

Remark 2.1. It is clear that, in the above location-scale model, one may easily prove that f ⊕ n converges almost surely to f0 as n → ∞ for various distances between density functions as illustrated by Figure 1(e). However, in this paper, we restrict our attention to the problem of how the structural mean measure ν0 can be estimated from empirical Wasserstein barycenters with respect to the Wasserstein distance between probability measures. Showing that the density (if it it exists) of such estimators converges to the density f0 of ν0 is not considered in this work.

−6 −6 −4 −2 0 2 4 6 0.0 0.2 0.4 0.6 0.8 1.0

(a) Densities f 1 , . . . , f n sampled from a location- (b) Quantile functions F − − 1 , . . . , F n of f 1 , . . . , f n scale model

−6 −4 −2 0 2 4 6 0.0 0.2 0.4 0.6 0.8 1.0

(e) Density f ⊕ n by quantile synchronization (f) Quantile function of the Wasserstein barycen- ter ν ⊕ ⊕ n with density f n

sampled from the location-scale model (2.2) with f0 the standard Gaussian density, ai ∼ U([0.8, 1.2]) and bi ∼ U([−2, 2]). (c,d) The solid-black curves are the Euclidean mean f̄ n and its quantile function. (e,f) The solid-red curves are the structural mean f⊕n given by quantile synchronization and the quantile function of the empirical Wasser- stein barycenter ν ⊕n . In all the figures, the dashed-blue curves are either the density f0 or its quantile function in the location-scale model (2.2).

2.2 Barycenters in the Wasserstein space

Let Ω be an interval of R, that is possibly unbounded. We let W2 (Ω) be the set of probability measures over (Ω, B(Ω)), with finite second moment, where B(Ω) is the σ-algebra of Borel subsets of Ω. We also denote by W2ac (Ω) the set of measures ν ∈ W2 (Ω) that are absolutely continuous with respect to the Lebesgue measure dx on R. The cumulative distribution function (cdf) and the quantile function of ν are denoted respectively by Fν and Fν− .

Definition 2.1. The quadratic Wasserstein distance dW in W2 (Ω) is defined by

Z 1

d2W (µ, ν) := (Fµ− (α) − Fν− (α))2 dα, for any µ, ν ∈ W2 (Ω). (2.3)

It can be shown that W2 (Ω) endowed with dW is a metric space, usually called Wasserstein space. For a detailed analysis of W2 (Ω) and its connection with optimal transport theory, we refer to [Vil03]. A probability measure ν in W2 (Ω) is said to be random if it is sampled from a distribution P on (W2 (Ω), B (W2 (Ω)), where B (W2 (Ω)) is the Borel σ-algebra generated by the topology induced by the distance dW .

Definition 2.2 (Square-integrability). The random measure ν is said to be square-integrable if

E(d2W (µ, ν)) = d2W (µ, ν)dP(ν) < +∞

W2 (Ω)

for some (thus for every) µ ∈ W2 (Ω).

Definition 2.3 (Population and empirical Wasserstein barycenters). Let ν be a W2 (Ω)-valued

square integrable random measure with distribution P. The population Wasserstein barycenter of ν is defined as the minimizer of Z µ 7→ d2W (µ, ν)dP(ν) over µ ∈ W2 (Ω).

W2 (Ω)

The empirical Wasserstein barycenter of ν1 , . . . , νn ∈ W2 (Ω) is defined as the minimizer of n

1X 2

µ 7→ dW (µ, νi ) over µ ∈ W2 (Ω). n i=1

Remark 2.2. In the whole paper, we assume that the model of random probability measure is well defined in the sense that all applications from an abstract probability space to W2 (Ω) are measurable and hence we can apply Fubini’s theorem. For an example of a rigorous model satisfying this kind of assumptions we refer to [BK16].

Proposition 2.1. Let ν ∈ W2 (Ω) be a square-integrable random measure then

(i) The exists a unique barycenter ν0 of ν.

R1

(iii) Var(ν) := E d2W (ν, ν0 ) = 0 Var (Fν− (α)) dα.

Proof. Points (i) and (ii) are consequences of Proposition 4.1 in [BGKL15]. Let us prove (iii).

From (2.3) and Fubini’s theorem, we have

Z 1 Z 1 h Z 1 − − 2 − − 2 i Var Fν− (α) dα. Var (ν) = E Fν (α) − Fν0 (α) dα = E Fν (α) − Fν0 (α) dα = 0 0 0

2.3 A deformable model of probability measures

Let ν be a W2 (Ω)-valued square integrable random probability measure. We use the notation F and F − to denote the cumulative distribution function (cdf) and the quantile function of the random measure ν. Let us also denote by ν0 the barycenter of ν (the existence and unicity of ν0 is ensured thanks to Proposition 2.1) and by ν 1 , . . . , ν n independent copies of ν. In this paper, we consider a deformable model of random probability measures satisfying the following assumptions:

Assumption 2.1. ν ∈ W2ac (Ω), a.s. and is a square integrable random probability measure in the sense of Definition 2.2.

Assumption 2.3. For each 1 ≤ i ≤ n, conditionally on ν i , the observations Xi,1 , . . . , Xi,pi are iid random variables sampled from ν i , where pi ≥ 1 is a known integer.

Remark 2.3. Since ν is square integrable, it follows from Proposition 2.1 that

Z 1

− − Var F − (α) dα < +∞. E F (α) = F0 (α), for all α ∈]0, 1[, and 0 ≤ (2.4)

It should be remarked that similar assumptions are considered in [PZ16] to characterize a pop- ulation barycenter in W2 (Ω) for the purpose of estimating phase and amplitude variations from the observations of multiple point processes. For examples of parametric models satisfying the

Assumptions 2.1-2.3, we refer to [BK16] and [BLGL15]. The main restriction of this deformable

model is that ν0 is assumed to be absolutely continuous.

2.4 Non-smoothed empirical barycenter

To estimate the structural mean measure ν0 from the data (Xi,j )1≤i≤n; 1≤j≤pi , a first approach consists in computing straightaway P i the barycenter of the empirical measure ν̃ 1 , . . . , ν̃ n where p = (p1 , . . . , pn ), and ν̃ i = p1i pj=1 δXi,j (δa denotes the Dirac mass at point x ∈ Ω). The ¯non-smoothed empirical barycenter is thus defined as

1X 2

ν̂ n,p = arg min dW (ν̃ i , µ). (2.5) ¯ µ∈W2 (Ω) n i=1

In the case where p1 = p2 = . . . = pn = p, we have the following procedure for computing the non-smoothed empirical barycenter. For each 1 ≤ i ≤ n, we denote by Xi,1 ∗ ≤ X∗ ≤ . . . ≤ X∗ i,2 i,p the order statistics corresponding to the i-th sample of observations (Xi,j )1≤j≤p , and we define n

1X ∗

X̄j∗ = Xi,j , for all 1 ≤ j ≤ p. n i=1

Thanks to Proposition 2.1, the quantile function of the empirical Wasserstein barycenter is the

average of the quantile functions of ν̃ 1 , . . . , ν̃ n , and thus we obtain the formula p 1X ν̂ n,p = δX̄j∗ . (2.6) p j=1

Note that we use the notation ν̂ n,p instead of ν̂ n,p to denote the non-smoothed empirical barycen- ¯ ter in the case p1 = p2 = . . . = pn = p.

2.5 Smoothed empirical barycenter

An alternative approach is to use a smoothing step to obtain estimated densities and then h1 hn compute the barycenter. In a first step, to obtain estimators f̂ 1 , . . . , f̂ n of f 1 , . . . , f n , we use kernel smoothing, where h1 , . . . , hn are positive bandwidth parameters that may be different for each subject or experimental unit. In this paper, to analyze the convergence of smoothed empirical barycenter in W2 (Ω), we shall investigate a non-standard choice for the kernel function

that has been proposed in [PZ16]. In Section 3, we give a precise definition of the resulting estimators based on the work in [PZ16]. However, at this point, it is not necessary to go into such details. Then, in a second step, an estimator of ν0 is given by ν̂ hn,p , with p = (p1 , . . . , pn ), defined as the measure whose quantile function is given by ¯ ¯

− 1X −

F̂ h (α) = F̂ i (α), α ∈ [0, 1], (2.7)

− hi where F̂ i denotes the quantile function of the density f̂ i for each 1 ≤ i ≤ n. If for each hi 1 ≤ i ≤ n, we denote by ν̂ hi i the measure with density f̂ i , then by Proposition 2.1, one has that ν̂ hn,p is also characterized as the following smoothed empirical Wasserstein barycenter ¯ n

1X 2

ν̂ hn,p = arg min dW (ν̂ hi i , µ). (2.8) ¯ µ∈W2 (Ω) n i=1

3 Convergence rate for estimators of the population Wasserstein

barycenter In this section, we discuss the rates of convergence of the estimators ν̂ n,p and ν̂ hn,p , that are ¯ ¯ respectively characterized by equations (2.5) and (2.8). Some of the results presented below

are using the work in [BL14] on a detailed study of the variety of rates of convergence of an empirical measure on the real line toward its population counterpart in the Wasserstein metric. Then, we discuss the optimality of these estimators from the minimax point of view following the guidelines in nonparametric statistics to derive optimal rates of convergence (see e.g. [Tsy09] for an introduction to this topic).

3.1 Non-smoothed empirical barycenter in the case of samples of equal size

Let us first characterize the rate of convergence of ν̂ n,p , in the specific case where samples of observations per unit are of equal size, namely when p1 = p2 = . . . = pn = p. In what follows, we let Y1 , . . . , Yp be iid random variables sampled from the population

P mean measure ν0 (inde-

pendently of the data (Xi,j )1≤i≤n; 1≤j≤p ), and we denote by µp = 1p pk=1 δYk the corresponding empirical measure.

Theorem 3.1. If Assumptions 2.1, 2.2 and 2.3 are satisfied and if p1 = p2 = . . . = pn = p, then the estimator ν̂ n,p satisfies

Z 1 p p Z j/p

1 1 X X 2

E d2W (ν̂ n,p , ν0 ) = Var F − (α) dα + Var Yj∗ + E Yj∗ − F0− (α) dα, n 0 pn (j−1)/p j=1 j=1

Z 1 p

1 1−nX Var F − (α) dα + Var Yj∗ + E d2W (µp , ν0 ) , = (3.1) n 0 pn j=1

where Y1∗ ≤ Y2∗ ≤ . . . ≤ Yp∗ denote the order statistics of the sample Y1 , . . . , Yp .

Theorem 3.1 provides exact formulas to compute the rate of convergence (for the expected

squared Wasserstein distance) of ν̂ n,p . Formula (3.1) relies on the computation of the variances of the order statistics of iid variables Y1 , . . . , Yp sampled from the population mean measure ν0 , and on the computation of the rate of convergence of E d2W (µp , ν0 ) . We discuss below some

examples where equality (3.1) may be used to derive a sharp rate of convergence for ν̂ n,p .

The case where ν0 is the uniform distribution on [0, 1]. In this setting, it is known (see e.g. Section 4.2 in [BL14]) that p j(p − j + 1) X p

Var Yj∗ = Var Yj∗ =

and thus . (p + 1) (p + 2) 6(p + 1) j=1

Moreover, from Theorem 4.7 in [BL14], it follows that E d2W (µp , ν0 ) = 6p . Therefore, thanks to equality (3.1), we obtain that

1 1 1−n 1 Z 2

Var F − (α) dα +

E dW (ν̂ n,p , ν0 ) = +

n 0 6n(p + 1) 6p

Z 1

1 − 1 1 1 = Var F (α) dα + + . (3.2) n 0 6 n(p + 1) p(p + 1)

Equality (3.2) thus shows that, when ν0 is the uniform distribution on [0, 1], the rate of conver- gence of ν̂ n,p is of the order 1 1 1

E d2W (ν̂ n,p , ν0 ) ≍ +

+ 2, (3.3) n np p and that this rate is sharp.

The case where ν0 is the one-sided exponential distribution. From Theorem 4.3 in

[BL14], one has that p p

1 X 2X

Var Yj∗ ≤ E d2W (µp , ν0 ) ≤ Var Yj∗ , (3.4) 2p p j=1 j=1

for any distribution ν0 ∈ W2 (Ω). Therefore, combining the above inequalities with (3.1), it follows that p 1 1 1+n X Z 2 Var F − (α) dα + Var Yj∗ .

E dW (ν̂ n,p , ν0 ) ≤ (3.5)

Now (using e.g. Remark 6.13 in [BL14]) one has that if ν0 is the one-sided exponential distribution (with density e−x for x ≥ 0) then p p

X X 1

Yj∗ Var = ∼ log(p) as p → +∞. j j=1 j=1

Therefore, there exist a constant c > 0 such that

1 log(p)

Z E d2W (ν̂ n,p , ν0 ) ≤ Var F − (α) dα + c 1 + (3.6) n 0 n p for all sufficiently large p. Hence, when ν0 is the exponential distribution the above inequalities 1 1 1 show that the rate of convergence of ν̂ n,p is of the order O n + log(p) np + p .

The case where ν0 is the standard Gaussian distribution. Deriving a sharp rate of con- vergence for ν̂ n,p using inequalities (3.1) combined with (3.4) requires computing the variances of the order statistics of iid random variables. To the best of our knowledge, obtaining a sharp estimate for Var Yj∗ for any 1 ≤ j ≤ p remains a difficult task except for specific distributions. Nevertheless, if ν0 is assumed to be a log-concave measure, then it is possible to use the results

in Section 6 of [BL14] which provide sharp bounds on the variances of order statistics for such probability measures. For example, if ν0 is the standard Gaussian distribution, then by Theorem 4.3 and Corollary

6.14 in [BL14] we obtain that there exist two constants c1 , c2 > 0 such that

p log(log(p)) 1X log(log(p))

Var Yj∗ ≤ c2

c1 ≤ . p p p j=1

Therefore, combining the above upper bound with (3.5), one finally has that

1 log(log(p))

Z 2 − E dW (ν̂ n,p , ν0 ) ≤ Var F (α) dα + c2 +1 . (3.7) n 0 n p when ν0 is the standard Gaussian distribution. In this setting, the rate of convergence is thus 1 1 1 of the order O n + log(log(p)) np + p .

Upper 2 bounds in more general cases. If one is interested in deriving an upper bound on E dW (ν̂ n,p , ν0 ) for a larger class of measures ν0 ∈ W2 (Ω) (e.g. Pp beyond ∗the log-concave case), another approach is as follows. Noting that the term 1−n pn j=1 Var(Yj ) in equality (3.1) is negative, a straightforward consequence of Theorem 3.1 is the following upper bound 1 1 Z E d2W (ν̂ n,p , ν0 ) ≤ Var F − (α) dα + E d2W (µp , ν0 ) .

(3.8) n 0 Then, thanks to inequality (3.8), to derive the rate of convergence of ν̂ n,p , it remains to control the rate of convergence of the empirical measure µp to ν0 for the expected squared Wasserstein distance. This issue is discussed in detail in [BL14]. In particular, the work in [BL14] describes a variety of rates for the expected distance E d2W (µp , ν0 ) , from the standard one O p1 to slower rates. For example, by Theorem 5.1 in [BL14], the following upper bound

E d2W (µp , ν0 ) ≤

J2 (ν0 ), (3.9)

p+1 where the so-called J2 -functional is defined as

F0 (x)(1 − F0 (x))

Z J2 (ν0 ) = dx,

Ω f0 (x)

where f0 is the density of ν0 , and F0 denotes its cdf. Therefore, provided that J2 (ν0 ) is finite, the empirical measure µp converges to ν0 at the rate O p . Hence, using inequality (3.9), we have:

Corollary 3.1. Suppose that Assumptions 2.1, 2.2 and 2.3 are satisfied. Then, the estimator

ν̂ n,p satisfies 1 1 2 Z 2

Var F − (α) dα +

E dW (ν̂ n,p , ν0 ) ≤ J2 (ν0 ). (3.10)

By Corollary 3.1, if J2 (ν0 ) < +∞, then it follows that ν̂ n,p converges to ν0 at the rate O n1 + p1 . Hence, in the setting where p ≥ n, ν̂ n,p converges at the classical parametric rate O n1 , provided that J2 (ν0 ) < +∞. The case p ≥ n is usually refereed to as the dense case in the

literature on functional data analysis (see e.g. [LH10] and references therein) which corresponds to the situation where the number of observations per unit/subject is larger than the sample size n of functional objects. In thesparse case (when p < n), the non-smoothed Wasserstein barycenter converges at the rate O p , provided that J2 (ν0 ) < +∞.

Remark 3.1. Whenν0 is the uniform distribution on [0, 1] one has that J2 (ν0 ) < +∞, but we have shown that E d2W (ν̂ n,p , ν0 ) ≍ n1 + np + p12 . Hence, in this setting, ν̂ n,p converges at the √ parametric rate O n1 provided that p ≥ n, which is a dense regime condition weaker than

To conclude this discussion on the rate of convergence of the non-smoothed Wasserstein

barycenter in the case of samples of equal size, we study in more detail the control of the rate of convergence of the term E d2W (µp , ν0 ) in inequality (3.8). As pointed out in many works (see for

example [dBGU05, BL14] and the references therein) the fact that the functional J2 (ν0 ) is finite or not is the key point to control the convergence of the empirical measure µp to the population measure ν0 in the Wasserstein space. Some known facts concerning J2 are the following.

1. If J2 (ν0 ) < +∞ then ν0 is supported on an interval of R and its density is a.e. strictly

2. If ν0 is compactly supported with a density bounded away from zero or with a log-concave

3. If the density of ν0 is of the form Cα e−|x| then J2 (ν0 ) is finite if and only if α > 2. In

particular, J2 (ν0 ) = +∞ for the Gaussian distribution.

Some further comments can be made in the case where ν0 is a Gaussian distribution. In this setting, one has that J2 (ν0 ) = +∞ and the rate of convergence of E d2W (µp , ν0 ) to zero is slower than O 1p . Indeed, from Corollary 6.14 in [BL14], if ν0 is the standard Gaussian distribution, then there exist two constants c1 , c2 > 0 such that

log(log(p)) log(log(p)) ≤ E d2W (µp , ν0 ) ≤ c2 c1 . (3.11) p p

Hence, using again inequality (3.8) combined with the above upper bound, we have:

Corollary 3.2. Suppose that Assumptions 2.1, 2.2 and 2.3 are satisfied. If ν0 is the standard

Gaussian distribution, then the estimator ν̂ n,p satisfies

Z 1

1 log(log(p)) E d2W (ν̂ n,p , ν0 ) ≤ Var F − (α) dα + c , (3.12) n 0 p

for some numerical constant c > 0.

Hence by Corollary 3.2, if p is sufficiently large with respect to n (namely when p ≥ n log(log(p))), then ν̂ n,p also converges at the classical parametric rate O n when ν0 is the standard Gaussian distribution. Remark 3.2. Following the work of [BL14], if ν0 is a log-concave distribution, then one may obtain rates of convergence for E d2W (ν̂ n,p , ν0 ) that are slower than the standard O p1 rate

(e.g. for beta or exponential distributions). Moreover, it is also possible to considerer for any

q ≥ 1 and for any probability measure µ on the real line (with density f and distribution function

F ) the functional

Jq (µ) = dx

R f (x)q−1

in order to control the rate of convergence of the empirical measure to µ for the q-Wasserstein distance.

3.2 Non-smoothed empirical barycenter in the general case

Let us now consider the general situation where the pi ’s are possibly different. The result below gives an upper bound on the rate of convergence of ν̂ n,p where p = (p1 , . . . , pn ). ¯ ¯ Theorem 3.2. Suppose that Assumptions 2.1, 2.2 and 2.3 are satisfied. Then, s

Z 1 n q

−1/2 1X

Var F − (α) dα + E d2W (ν̃ i , ν i ) ,

E dW (ν̂ n,p , ν0 ) ≤ n

¯ 0 n i=1

1 Pp i

where ν̃ i = pi j=1 δXi,j for each 1 ≤ i ≤ n

For the random measure ν, we define the random variable

F (x)(1 − F (x))

Z J2 (ν) = dx.

Ω f (x)

Since the ν i ’s are independent copies of ν by applying inequality (3.9), it follows that q p −1/2 E d2W (ν̃ i , ν i ) ≤ 2E [J2 (ν)]pi .

Hence, from Theorem 3.2, we finally obtain the following upper bound on the rate of convergence for the non-smoothed empirical barycenter Corollary 3.3. Suppose that Assumptions 2.1, 2.2 and 2.3 are satisfied. If J2 (ν) has a finite expectation, then s n

Z 1 !

h i 1 X −1/2

E dW (ν̂ n,p , ν0 ) ≤ n−1/2

Var F − (α) dα + 2E [J2 (ν)]

pi . ¯ 0 n i=1

−1/2 From Corollary 3.3, one has that if min1≤i≤n pi ≥ n (dense case), then n1 ni=1 pi P ≤ n−1/2 , and thus, the non-smoothed empirical barycenter converges of the parametric rate n−1/2 (provided that E [J2 (ν)] < +∞), namely s  h i Z 1 Var F − (α) dα + 2E [J2 (ν)) n−1/2 . p

E dW (ν̂ n,p , ν0 ) ≤  (3.13)

Remark 3.3. Knowing if J2 (ν) has a finite expectation is in general a difficult task. But, if we assume that the density f of ν is bounded below by a non-random positive constant then (obviously) E [J2 (ν)] < +∞.

3.3 The case of smoothed empirical barycenters

In this section, we assume that Ω = [0, 1] and we discuss the rate of convergence of smoothed empirical barycenters ν̂ hn,p (note that the following results hold if Ω is any compact interval). ¯ To choose an appropriate kernel function to study the convergence rate of the estimator ν̂ hn,p , ¯ we follow the proposal madeR in2 [PZ16]. We let ψ be a positive, smooth and symmetric density on the real line, such that R x ψ(x)dx = 1. We also denote by Ψ the cdf of the density ψ and,

for a bandwidth parameter h > 0, we let ψh (x) = h1 ψ hx . Then, for any y ∈ [0, 1] and h > 0, we denote by µyh the measure supported on [0, 1] whose density fµy is defined as h

fµy (x) = ψh (x − y) + 2b2 ψh (x − y)11{x−y>0} + 2b1 ψh (x − y)11{x−y<0} + 4b1 b2 , x ∈ [0, 1], (3.14) h

where b1 = 1 − Ψ ((1 − y)/h) and b2 = Ψ (−y/h). Then, for each 1 ≤ i ≤ n, we construct a hi kernel density estimator of f i by defining f̂ i as the density associated to the measure pi

1 X Xi,j

ν̂ hi i = µhi , (3.15) pi j=1

Frequently Asked Questions

What is this project about?

This project covers practical implementation and research aspects of the topic using AI/ML techniques.

This section provides additional detailed analysis and supporting information derived from the research paper content to ensure comprehensive coverage of the topic with expanded discussion on key concepts, methods, and findings presented in the source material. The practical implications and implementation considerations are further elaborated to support complete understanding of the subject matter across different application domains and use cases.

Speaker Verification Xvector

Ai Security Projects

Minimax convergence rate for estimating the Wasserstein

Jérémie Bigot1 , Raúl Gouet2 , Thierry Klein3 & Alfredo López4

Institut de Mathématiques de Bordeaux et CNRS (UMR 5251)1

Université de Bordeaux

Depto. de Ingenierı́a Matemática and CMM (CNRS, UMI 2807)2

Universidad de Chile

ENAC- Ecole nationale de l’aviation civile

Université de Toulouse

CSIRO Chile International Centre of Excellence4

November 28, 2021

Abstract

R that can be viewed as independent realizations of an underlying stochastic process. We

Wasserstein barycenter towards its population counterpart in the asymptotic setting where

Keywords: Wasserstein space; Fréchet mean; Barycenter of probability measures; Functional

1 Introduction

Wasserstein barycenter of ν, given by

Data sets leading to the analysis of absolutely continuous measures appear in various re-

1.1 Main contributions

1.2 Related work in the literature

In statistics, this approach has been referred to as quantile synchronization [ZM11]. In the

The construction of consistent estimators of a population Wasserstein barycenter for semi-

[RL01, WG97].

1.3 Organization of the paper

5 contains a discussion on the main contributions of this work and their potential extensions.

2 Wasserstein barycenters for the estimation of the structural

2.1 The need to account for phase and amplitude variations

− 1X −

F̄ n = Fi , (2.1)

The statistical analysis of quantile synchronization, as studied in [ZM11], complements the

2.1 below) that quantile synchronization corresponds to computing the empirical Wasserstein

1X 2

Therefore, the notion of averaging by quantile synchronization corresponds to using the Wasser-

To illustrate the differences between using Euclidean and Wasserstein distances to account

P P

2.2 Barycenters in the Wasserstein space

Definition 2.1. The quadratic Wasserstein distance dW in W2 (Ω) is defined by

Z 1

Definition 2.2 (Square-integrability). The random measure ν is said to be square-integrable if

E(d2W (µ, ν)) = d2W (µ, ν)dP(ν) < +∞

W2 (Ω)

Definition 2.3 (Population and empirical Wasserstein barycenters). Let ν be a W2 (Ω)-valued

W2 (Ω)

1X 2

Proposition 2.1. Let ν ∈ W2 (Ω) be a square-integrable random measure then

R1

From (2.3) and Fubini’s theorem, we have

2.3 A deformable model of probability measures

Remark 2.3. Since ν is square integrable, it follows from Proposition 2.1 that

Z 1

Assumptions 2.1-2.3, we refer to [BK16] and [BLGL15]. The main restriction of this deformable

2.4 Non-smoothed empirical barycenter

1X 2

1X ∗

Thanks to Proposition 2.1, the quantile function of the empirical Wasserstein barycenter is the

2.5 Smoothed empirical barycenter

− 1X −

F̂ h (α) = F̂ i (α), α ∈ [0, 1], (2.7)

1X 2

3 Convergence rate for estimators of the population Wasserstein

3.1 Non-smoothed empirical barycenter in the case of samples of equal size

P mean measure ν0 (inde-

Z 1 p p Z j/p

1 1 X X 2

Z 1 p

Theorem 3.1 provides exact formulas to compute the rate of convergence (for the expected

Var Yj∗ = Var Yj∗ =

Var F − (α) dα +

E dW (ν̂ n,p , ν0 ) = +

Z 1  

E d2W (ν̂ n,p , ν0 ) ≍ +

The case where ν0 is the one-sided exponential distribution. From Theorem 4.3 in

1 X 2X

E dW (ν̂ n,p , ν0 ) ≤ (3.5)

X X 1

Therefore, there exist a constant c > 0 such that

1 log(p)

6.14 in [BL14] we obtain that there exist two constants c1 , c2 > 0 such that

Var Yj∗ ≤ c2

Z 1