Tung Duy Luu
email: duy-tung.luu@ensicaen.fr

Jalal Fadili
email: jalal.fadili@greyc.ensicaen.fr

Christophe Chesneau
email: christophe.chesneau@unicaen.fr

PAC-Bayesian risk bounds for group-analysis sparse regression by exponential weighting

Keywords: High-dimensional regression, exponential weighted aggregation, sparse learning, group-analysis sparsity, sparse oracle inequality, frame, forward-backward Langevin Monte-Carlo AMS subject classifications. 62G07 62G20

In this paper, we consider a high-dimensional non-parametric regression model with fixed design and i.i.d. random errors. We propose an estimator by exponential weighted aggregation (EWA) with a group-analysis sparsity promoting prior on the weights. We prove that our estimator satisfies a sharp group-analysis sparse oracle inequality with a small remainder term ensuring its good theoretical performances. We also propose a forward-backward proximal Langevin Monte-Carlo algorithm to sample from the target distribution (which is not smooth nor log-concave) and derive its convergence guarantees. In turn, this allows us to implement our estimator and validate it on some numerical experiments.

Introduction

Problem statement

Let us briefly present our statistical context. Assume that the given data (x i , Y i ), i = 1, . . . , n, is generated according to the non-parametric regression model

Y i = f (x i ) + ξ i , i ∈ {1, . . . , n}, (1.1)
where x 1 , . . . , x n are deterministic in an arbitrary set X, f : X → R is the unknown regression function and (ξ 1 , . . . , ξ n ) are random errors. (1.1) is equivalently written in vector form Y = f + ξ. Assume that there exists a dictionary H = f j : X → R, j ∈ {1, . . . , p} such that f is well approximated by a linear combination of elements in H, i.e., there exists θ = ( θ 1 , . . . , θ p ) ∈ R p such that f θ = p j=1 θ j f j is a suitable approximation of f . The f j 's are known and may be either fixed atoms in a basis or pre-estimators. In the context of high dimension, the cardinality of H is much larger than the sample size (i.e., p n). Thus, in such a setting, the classical least-squares to estimate θ is obviously not applicable.

The idea of aggregating elements in a dictionary has been introduced in machine learning to combine different techniques (see [START_REF] Littlestone | The weighted majority algorithm[END_REF][START_REF] Vovk | Aggregating strategies[END_REF]) with some procedures such as bagging [START_REF] Breiman | Bagging predictors[END_REF], boosting [START_REF] Freund | Boosting a weak learning algorithm by majority[END_REF][START_REF] Schapire | The strength of weak learnability[END_REF] and random forests [1, 4-6, 8, 32]. In the recent years, there has been a flurry of research on the use of the concept of sparsity in various areas including statistics and machine learning in high dimension. The idea is that even if the cardinality of H is very large, the number of effective elements in the dictionary is much smaller than the sample size. Namely, the number of nonzero components of θ is assumed to be much smaller than n. This makes it possible to build an estimate f θ with good provable performance guarantees under appropriate conditions on the dictionary and noise.

Overview of previous work 1.2.1. Oracle inequalities

This type of guarantees dates back, for instance, to the work [START_REF] Donoho | Ideal spatial adaptation by wavelet shrinkage[END_REF][START_REF] Donoho | Adapting to unknown smoothness via wavelet shrinkage[END_REF][START_REF] Donoho | Wavelet shrinkage: asymptopia[END_REF] on orthogonal wavelet thresholding estimators. Oracle inequalities (according to the terminology introduced in, e.g., [START_REF] Donoho | Adapting to unknown smoothness via wavelet shrinkage[END_REF]), which are at the heart of our work, quantify the quality of an estimator compared to the best possible one among a collection of estimators. Formally, let g : X → R, and denote g n = 1

n n j=1 g 2 (x j ).

The performance of an estimator f is measured by its averaged squared error, i.e.,

R( f ) = f -f 2 n .
We aim to find an estimator f that mimics as much as possible the performance of the best model of aggregation in a given class Θ (in the probabilistic sense). This idea is expressed in the following type of inequalities:

E f -f 2 n ≤ C inf θ∈Θ f θ -f 2 n + ∆ n,p (θ) , (1.2)
where C ≥ 1 and the remainder term ∆ n,p (θ) depends on the performance of the estimator, the complexity of θ, the dimension p and the sample size n. Such type of inequality is called balanced oracle inequality. Inequalities of type (1.2) are well adapted under the sparsity scenario. Namely, the complexity of θ in the remaining term is characterized by the sparsity parameters (like the number of its non-zero components), in which case these inequalities are called sparse oracle inequalities (SOI).

An estimator with good oracle properties would correspond to C is close to 1 (ideally, C = 1, in which case the inequality is said "sharp"), ∆ n,p (θ) is small even if n p and decreases rapidly to 0 as n → +∞. Besides, the choice of Θ is crucial: on the one hand, a non suitable choice can lead a large bias term in (1.2). On the other hand, if Θ is too complex, the remainder term becomes large. Then, a suitable choice for Θ must achieve a good bias-complexity trade-off.

In the literature, there are mainly two approaches to provide aggregated estimators in high dimension under the sparsity assumption: Penalization and Exponential Weighted Aggregation (EWA). Given Y = y, The penalization approach considers the minimization problem

min θ∈Θ y -f θ 2 n + pen(θ),
where pen : R p → R + is a sparsity promoting penalty function, see, e.g., [START_REF] Bühlmann | Statistics for high-dimensional data: methods, theory and applications[END_REF]. Our work focuses on EWA approach that we briefly describe now.

Exponential Weighted Aggregation (EWA)

Let (Λ, A) be a space equipped with a σ-algebra and

F Λ = { f λ : X → R : λ ∈ Λ}
be a given collection (F Λ is called dictionary of aggregation) where λ → f λ (x) is measurable ∀x ∈ X. The functions f λ may be deterministic or random. The aggregators depend on the nature of f λ if the latter is random. Otherwise, the aggregators are defined via the probability measure

µ n (dλ) = exp -n Y -f λ 2 n /β π(dλ) Λ exp -n Y -f ω 2 n /β π(dω) , (1.3)
where β > 0 called temperature parameter and π called prior which is a probability measure on Λ. Then, we define the EWA aggregate by

f n (x) = Λ f λ (x)µ n (dλ). (1.4)
This idea was initially proposed in [START_REF] Kivinen | Exponentiated gradient versus gradient descent for linear predictors[END_REF][START_REF] Littlestone | The weighted majority algorithm[END_REF][START_REF] Vovk | Aggregating strategies[END_REF] with a uniform prior on a finite set Λ. Observe that (1.4) can be interpreted as the Bayesian posterior conditional mean in the regression model but only when the noise is Gaussian noise N(0, σ 2 I n ) and provided that β = 2σ 2 and the prior density is π.

In the literature, the works in [START_REF] Dalalyan | Pac-bayesian bounds for the expected error of aggregation by exponential weights[END_REF][START_REF] Dalalyan | Aggregation by exponential weighting, sharp pac-bayesian bounds and sparsity[END_REF][START_REF] Dalalyan | Aggregation by exponential weighting and sharp oracle inequalities[END_REF][START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF] consider deterministic dictionaries. These papers proposed several PAC-Bayesian1 type of oracle inequalities under different assumptions. Especially, the assumptions in [START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF] depend only on the noise and turns out to be fulfilled for a large class of noise. This serves to construct, for a suitable prior and dictionary, a SOI with a remainder term of order O( θ 0 ln(p)/n), which scales linearly with the sparsity level and increases in p only logarithmically.

The random dictionary case are tackled in [START_REF] Nemirovski | Topics in non-parametric statistics[END_REF]. The initial idea is to obtain two independent samples from the initial sample by randomization or sample splitting (see [START_REF] Lecué | Simultaneous adaptation to the margin and to complexity in classification[END_REF][START_REF] Rigollet | Linear and convex aggregation of density estimators[END_REF][START_REF] Yang | Aggregating regression procedures to improve performance[END_REF]). The first sample is used to construct the preestimators, and the aggregation is performed on the second sample conditionally on the first one. However this idea does not work when the observations are not i.i.d. Several authors have proposed exponentially aggregating linear pre-estimators without splitting, and with discrete priors on the weights. Typical cases of linear pre-estimators are orthogonal projectors on all possible linear subspaces that are in the model set (e.g., in the sparsity context, linear subspaces spanned by the standard basis restricted to supports of increasing size). This was introduced in [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF]. More recent work such as [START_REF] Dalalyan | Sharp oracle inequalities for aggregation of affine estimators[END_REF] generalizes the idea where the pre-estimators are affine and the priors are continuous.

A shortcoming of EWA is its suboptimality in deviation. In particular, the work in [16, Section 2] has shown that the EWA leads a suboptimal remainder term for oracle inequalities in probability.

Generalization of sparsity assumption

Analysis sparsity Let q ≥ p and D ∈ R q×p be a linear analysis operator. The analysis sparsity assumption means that D θ 0 n. A typical example is total variation [START_REF] Rudin | Nonlinear total variation based noise removal algorithms[END_REF] where the operator D corresponds to "finite differences" (i.e., ( D θ) 1 = θ 1 , ( D θ) j = θ j -θ j-1 , ∀ j ≥ 2). Another example is the fused Lasso [START_REF] Tibshirani | Sparsity and smoothness via the fused lasso[END_REF] where D is a positive combination of the identity and finite differences.

Group sparsity

Group sparsity corresponds to saying that the aggregator θ is block sparse, see Section 4.1 for formal details. Group sparsity is at the heart of the group Lasso and related methods [START_REF] Huang | Variable selection in nonparametric additive models[END_REF][START_REF] Koltchinskii | Sparse recovery in large ensembles of kernel machines[END_REF][START_REF] Meier | The group lasso for logistic regression[END_REF][START_REF] Meier | High-dimensional additive modeling[END_REF][START_REF] Ravikumar | Sparse additive models[END_REF][START_REF] Yuan | Model selection and estimation in regression with grouped variables[END_REF]]. In the EWA context, the group sparsity prior is considered in [START_REF] Rigollet | Sparse estimation by exponential weighting[END_REF] as an application of the aggregation of orthogonal projectors.

Contributions

Our main contributions are summarized as follows:

• We propose an EWA estimator, with a deterministic dictionary, under a group-analysis sparsity prior (see Section 1.2.3). More precisely, we assume that D is a frame, and thus is not necessarily invertible unlike previous work. In a finite space, that equivalents to the fact that D D is invertible. In addition, our prior class (see (4.2)) is much more general than previously proposed ones [START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF] which recovered as very special cases. It also allows more flexibility to enhance the performance of EWA. This prior class is parameterized through a function g : R + → R + that satisfies mild conditions (see Assumptions (G.1)-(G.4)).

• We establish a group-analysis SOI where the remainder term depends on the number of active groups in D θ and on the function g (see Theorem 5.1)

• For an appropriate choice of g which is well-adapted to the group-analysis sparsity scenario, we exhibit a groupanalysis SOI where, as expected, the remainder term scales as O D θ 0,G ln(L)/n , where D θ 0,G is the number of active groups in D θ, and L is the total number of groups (see Corollary 5.2). This rate coincides with the classical one O( θ 0 ln(p)/n) under the sparsity scenario, i.e., D = I p and L = p.

• We also propose a forward-backward proximal Langevin Monte-Carlo (LMC) algorithm to sample from the target distribution (6.3) (which is not smooth nor log-concave), and establish several of its properties in a general setting. In turn, this allows us to efficiently implement our EWA estimator with the proposed prior. We validate this algorithm on some numerical examples.

Paper organization

Necessary notations and some preliminaries are first introduced in Section 2. Section 3 reminds the PAC-Bayesian type oracle inequalities proposed in [START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF] which are a classical starting point in literature for EWA in the deterministic case. In Section 4, we describe our EWA procedure after specifying the aggregation dictionary and our prior family. In Section 5, we establish our main results, namely group-analysis SOI. Section 6 is devoted to the forward-backward proximal LMC algorithm that implements EWA, and the numerical experiments on several numerical settings are described in Section 7. The proofs of all results are collected in Section 9.

Notations and Preliminaries

Before proceeding, let us introduce some notations and definitions.

Vectors and matrices

For a d-dimensional Euclidean space R d , we endow it with its usual inner product •, • and associated norm • 2 . I d is the identity matrix on R d . For q > 0 and x ∈ R d , we also define x q = d j=1 x j q 1/q with the usual adaptation x ∞ = max j∈{1,...,d} x j . It is the q norm for q ≥ 1, and quasi-norm for q ∈]0, 1[. x 0 is the 0 pseudo-norm which counts the number of non-zero elements in x. Let M ∈ R d×d symmetric positive definite, we denote

•, • M = •, M• and • M its associated norm. Of course, • M and • 2 are equivalent.
For a matrix M ∈ R d×r , we set σ(M) = (σ 1 (M), . . . , σ r (M)) ∈ R r be the vector of singular values of M in nonincreasing order. Note that, when M is symmetric semi-definite positive, σ(M) is also the ordered vector or positive eigenvalues of M. We denote M the spectral norm of M. For I ⊂ {1, . . . , d}, x ∈ R d and M ∈ R d×d , x I is the subvector whose entries are those of x restricted to the indices in I, and M I the submatrix whose rows and columns are those of M indexed by I. Let us denote vec : R d×d → R d 2 the vectorization operator of a matrix. For any matrix M, M denotes its transpose. For any square matrix M, det(M) is its determinant.

Definition 2.1 (Frame).

A matrix M ∈ R d×r is a frame if there exist two constants ν and κ with ν ≥ κ > 0, called frame bounds, such that the generalized Parseval relation is satisfied, i.e., κ

x 2 2 ≤ M x 2 2 ≤ ν x 2 2 , ∀x ∈ R d .
By the Courant-Fischer theorem, Definition 2.1 is equivalent to the fact that κ (resp. ν) is a lower (resp. upper) bound of the eigenvalues of M M . Moreover, since κ > 0, we have that M M is bijective and M is injective. The frame is said tight when κ = ν. Typical examples of (tight) frames that have been used in statistics are translation invariant wavelets [START_REF] Coifman | Translation invariant de-noising[END_REF], ridgelets [START_REF] Candès | Ridgelets: Estimating with Ridge Functions[END_REF] and curvelets [START_REF] Candès | New tight frames of curvelets and optimal representations of objects with piecewise C 2 singularities[END_REF] (example of groups and what they represent for wavelets/ridgelets/curvelets in applications are discussed in [START_REF] Chesneau | Stein block thresholding for image denoising[END_REF]). Let M ∈ R d×r be the canonical dual frame associated to M, i.e., M = (M M ) -1 M. We know that

M M = I d (2.1) and 1 κ ≥ σ 1 M M ≥ . . . ≥ σ d M M ≥ 1 ν . (2.2)
Note that we focus on the canonical dual frame for the sake of simplicity. In fact, our exposition in the rest of the paper remains unchanged if any other dual frame is used instead of the canonical one.

Functions For a function

f : R d → R, we define f = ( f (x 1 ), . . . , f (x d )) its d-sample vector form, f d = 1 d d j=1 f 2 (x j ) 1/2 and f ∞ = sup x∈R d f (x) . For a function f : R d → R ∪ {+∞}, its effective domain is dom( f ) = x ∈ R d : f (x) < +∞ and f is proper if f (x) > -∞ for all x and dom( f ) ∅ as is the case when it is finite-valued. A function f : R d → R ∪ {-∞, +∞} is lower semi continuous (lsc) at x 0 if lim inf x→x 0 f (x) ≥ f (x 0 ). For a differentiable function f , ∇ f is its (Euclidean) gradient.
For a set Ω, I Ω is its characteristic function, i.e., 1 if the argument is in Ω and 0 otherwise. Define sgn : R → {-1, 1} be the sign function, i.e., sgn(x) = 1 when x ≥ 0, and sgn(x) = -1 otherwise. Let x ∈ R, define x be the stands of integer part of x. Recall the Gamma function

Γ :]0, +∞] →]0, +∞], Γ : t → ∞ 0 x t-1 exp (-x)dx.
Definition 2.2 (Proximal mapping and Moreau envelope). Let M ∈ R d×d be symmetric positive definite. For a proper lsc function f and γ > 0, the proximal mapping and Moreau envelope in the metric M are defined respectively by

prox M γ f (x) def = Argmin w∈R d 1 2γ w -x 2 M + f (w) M,γ f (x) def = inf w∈R d 1 2γ w -x 2 M + f (w)
prox M γ f here is a set-valued operator since the minimizer, if it exists, is not necessarily unique. When M = I d , we simply write prox γ f and γ f .

Useful integration formulas

The following lemmas contain useful formula used throughout the paper. Lemma 2.1 ([33, 3.251.11]

). Let p, γ, ν, η > 0. If γ/ν < η + 1 we have ∞ 0 x γ-1 (p + x ν ) η+1 dx = 1 νp η+1-γ/ν Γ(γ/ν)Γ(1 + η -γ/ν) Γ(1 + η) ,
otherwise this integral is not definite.

Lemma 2.2 (Cartesian to spherical coordinates [START_REF] Fang | Symmetric multivariate and related distributions[END_REF]). Let d ≥ 1 and a mapping h : R + → R such that u → h( u 2 ) is measurable in R d . We then have

R d h( u 2 )du = C d ∞ 0 x d-1 h(x)dx,
where

C d = 2π d/2 /Γ(d/2
) is the surface area of the d-dimensional Euclidean ball of radius 1.

The following lemma provides an efficient change of variables formula, which will be a key tool in the proof of our general group-analysis SOI (see Theorem 5.1). Lemma 2.3. Let Θ ⊆ R d be a measurable set. For any M ∈ R r×d corresponding to the analysis operator of a frame of R d , let u : R r → R such that the mapping x → u(Mx) is measurable on Θ. We have

Θ u(Mx)dx = 1 det(M M ) MΘ u(v)dv (2.3)
provided either u is non-negative valued or the integral on the left converges.

Though quite natural, proving Lemma 2.3 rigorously requires nontrivial arguments from geometric measure theory; see Section 9.1.

PAC-Bayesian type oracle inequalities

This section recalls a PAC-Bayesian type oracle inequality which holds for the EWA procedure of type (1.3)-(1.4) with any deterministic aggregation dictionary, any prior and a large class of noises. Such type of oracle inequalities was introduced in [START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF] for i.i.d. noise. In the present paper, we adapt it to the non i.i.d. case. Indeed, let us start with the two following assumptions.

(P.1) The noise vector ξ = (ξ 1 , . . . , ξ n ) has zero mean.

(P.2) For any γ > 0 small enough, there exist a probability space and two random variables ξ and ζ defined on this probability space satisfying the three following points:

(a) ξ has the same distribution as ξ.

(b) ξ + ζ has the same distribution as (1 + γ)ξ and the conditional expectation satisfies E ζ|ξ = 0.

(c) There exist t 0 ∈]0, ∞] and a bounded Borel function v : R p → R + such that

lim sup γ→0 sup (t,a)∈R p ×R p ( t 2 ,a)∈[0,t 0 ]×supp(ξ ) ln E exp t ζ |ξ = a t 2 2 γv(a) ≤ 1.
where supp(ξ ) is the support of the distribution of ξ .

Assumption (P.2) is based on [23, Assumption N], and can be shown to be fulfilled for a large class of noise.

Proposition 3.1. Assume that ξ has zero mean. Assumption (P.2) is fulfilled with t 0 = ∞ when

• ξ is a Gaussian random vector with covariance matrix Σ, with t 0 = ∞ and v(a) ≡ Σ ;

• ξ is a Laplace random vector with covariance Σ, with t 0 < 2/ Σ and v(a) ≡ Σ / 1t 2 0 Σ 2 /2 ;

• ξ is a bounded symmetric random vector, i.e., Pr

{|ξ i | ≤ B i } for some B ∈ R n , with t 0 = ∞ and v(a) = a 2 ≤ B 2 .
Besides, let H ∈]0, +∞] such that sup

(λ,λ )∈Λ 2 f λ -f λ 2 ≤ H. (3.1)
Note that (3.1) is always satisfied since H is allowed to be infinite. However, for the sake of sharpness in our theoretical results, we wish to choose H as small as possible. We are now ready to state the PAC-Bayesian type oracle inequalities.

Theorem 3.1. Let Assumptions (P.1) and (P.2) be satisfied with some function v and let (3.1) holds. Then for any prior π, any probability measure p on Λ and any

β ≥ max(4 v ∞ , 2H/t 0 ) or β ≥ 4 v ∞ when H = ∞, t 0 = ∞, we have E f n -f 2 n ≤ Λ f -f λ 2 n p(dλ) + β KL(p, π) n , (3.2)
where f n is the aggregate defined in (1.4) and KL(p, π) = Λ ln (p(dλ)/π(dλ))p(dλ) is the Kullback-Leibler divergence.

The proof of Theorem 3.1 is a mild adaptation of the original one in [23, Section 2], where we used directly Assumption (P.2)-3 in the vector ζ instead of splitting it into ζ i,i∈{1,...,n} (that are no longer i.i.d.).

Related work

The work of [START_REF] Dalalyan | Aggregation by exponential weighting, sharp pac-bayesian bounds and sparsity[END_REF] proposed three types of oracle inequalities which are similar to (3.2) under different assumptions. The first type (see [START_REF] Dalalyan | Aggregation by exponential weighting, sharp pac-bayesian bounds and sparsity[END_REF]Theorem 1]) holds under a restrictive condition on the noise. The second (see [START_REF] Dalalyan | Aggregation by exponential weighting, sharp pac-bayesian bounds and sparsity[END_REF]Theorem 2]) involves conditions depending on the noise and also on the dictionary. The last (see [START_REF] Dalalyan | Aggregation by exponential weighting, sharp pac-bayesian bounds and sparsity[END_REF]Theorem 4]) works for all symmetric noises without conditions on the dictionary. However, an additional term appears in the remainder term which has a low rate for some types of noise. Therefore, Theorem 3.1 (with Assumption (P.2)) is a good trade-off between these types of oracle inequalities.

Moreover, there exist some related forms of (3.2) in different frameworks. For example, when ξ i ∼ N(0, σ 2 i ), i = 1, . . . , n, the following aggregate was proposed in [START_REF] Dalalyan | Sharp oracle inequalities for aggregation of affine estimators[END_REF]:

f = Λ f λ p(dλ), p(dλ) = exp -n β r λ π(dλ) Λ exp -n β r ω π(dω)
, where f λ , λ ∈ Λ are affine estimators satisfying some conditions imposed in [20, Theorem 1] which yield the definition of r λ , λ ∈ Λ. This aggregate satisfies oracle inequalities defined therein which are the counterparts of (3.2) for the aggregation of estimators. In addition, in the case of random design (i.e., x 1 , . . . , x n are random and i.i.d.), the works in [START_REF] Dalalyan | Mirror averaging with sparsity priors[END_REF] constructed a mirror averaging aggregate to obtain a generalized type of oracle inequalities where the performance is measured by any loss instead of the averaged square loss.

EWA estimator

Group-analysis sparsity

We now describe formally what is intended by group-analysis sparsity, which measures group sparsity of the image of a vector with an analysis linear sparsifying transform. Let q ≥ p. We partition the index set {1, . . . , q} into L nonoverlapping groups/blocks of indices {G } 1≤ ≤L such that

G = L =1 G = {1, . . . , q} and G ∩ G k = ∅, ∀l k.
For the sake of simplicity, and without loss of generality, the groups are assumed to have the same size Card G = G ≥ 1 and the total number of blocks L is supposed to be an integer. A vector v ∈ R q can be divided into L vectors v G ∈ R G which are the restrictions of v to the coordinates indexed by G . We define

v s,G = L =1 v G s 2 1/s , s > 0, which is a norm for s ≥ 1, with v ∞,G = max ∈{1,...,L} v G 2 .
It is a quasi-norm for s ∈]0, 1[. v 0,G counts the number of active (i.e., non-zero) groups in v. With these notations, the group-analysis sparsity assumption is formalized as follows.

(H.1) Let D ∈ R p×q , there exists θ ∈ R p such that f θ = f and D θ 0,G n.

In plain words, Assumption (H.1) says that the number of active groups of D θ is much smaller than the sample size. Note that this is a strict notion of group-analysis sparsity, and a weaker one could be also considered where most D θ G are nearly zero. We also impose the following assumption on D.

(H.2) D is a frame with frame bounds ν ≥ κ > 0.
Let us now introduce some applications in literature in which our sparsity context is mentioned.

Example 4.1 (2-D piecewise constant image). Let θ 0 ∈ R √ p×
√ p be a 2-D piecewise constant image. In this framework, a popular analysis operator is the isotropic total variation called D TV (see [START_REF] Rudin | Nonlinear total variation based noise removal algorithms[END_REF]). Namely, let D c : R

√ p× √ p → R √ p× √ p and D r : R √ p× √ p → R √ p×
√ p the finite difference operators along, respectively, the columns and rows of an image, with Neumann boundary conditions. We define D TV as The isotropic total variation prior on θ promotes sparsity of vec [D TV (θ 0 )] i, j 2 1≤i, j≤ √ p . By defining the set of groups by G = (i, j)∈{1,...,s p 0 } 2 {(i, j, 1), (i, j, 2)}, one immediately realizes that measuring sparsity of the above vectorized form is equivalent to group sparsity of D TV (θ 0 ) with groups of size 2 along the third dimension.

D TV : θ ∈ R √ p× √ p → vec (vec(D r (θ)), vec(D c (θ))) ∈ R 2p .
Example 4.2 (Signal with overlapping groups). Consider θ 0 ∈ R p generated from L groups which overlap. The analysis operator acts as a group extractor (see [START_REF] Chen | An efficient proximal-gradient method for general structured sparse learning[END_REF][START_REF] Peyré | Group sparsity with overlapping partition functions[END_REF]). In this framework,

D D = diag(B 1 B 1 , . . . , B L B L ) where B ∈ R q ×p ,
for ∈ {1, . . . , L}, is a countable collection of localization operators, and then q = L =1 q . Since the localization operators are injective, D D is bijective. Thus D is a frame in view of Courant-Fisher theorem.

To design an aggregation by exponential weighting, two ingredients are essential: the aggregation dictionary and the prior which promotes group-analysis sparsity. We specify them below.

Choice of dictionary

Let X ∈ R n×p where X i, j = f j (x i ), where we recall that f j ∈ H is a known function (atom) in the deterministic dictionary H.

(H.3) X is normalized such that all the diagonal entries of X X/n are 1.

Now, let us introduce our dictionary of aggregation:

F Θ = f θ = L p j=1 θ j f j : θ ∈ Θ = θ ∈ R p : D θ a a,G ≤ R , (4.1)
where a ∈]0, 1], R ∈]0, +∞] and L : R → R is twice continuously differentiable and known function that depends on the regression problem: e.g., L (x) = e x for the exponential regression, L (x) = e x /(e x + 1) for the logistic regression and L (x) = x for the linear regression. This dictionary of aggregation is similar to the one proposed in [START_REF] Dalalyan | Aggregation by exponential weighting and sharp oracle inequalities[END_REF][START_REF] Dalalyan | Mirror averaging with sparsity priors[END_REF][START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF].

However, the set of indices is modified to adapt to the group-analysis sparsity and the exponent a is varied in ]0, 1] instead of a fixed a = 1. The bound H in (3.1) for F Θ in (4.1) is established in the following result.

Proposition 4.1. Let F Θ = { f θ : θ ∈ Θ} defined in (4.1) with some R > 0, a ∈]0, 1] and L : R → R twice continuously differentiable. Let Assumption (H.2) holds for some κ > 0. We get that

sup (θ,θ )∈Θ 2 f θ -f θ 2 ≤ 2 max x∈B (x) 2 , where B = x ∈ R n : x 2 ≤ X R 1/a κ -1/2 and : x ∈ R n → (L (x 1 ), . . . , L (x n )). From Proposition 4.1, one can choose H = 2 max x∈B (x) 2 . Remark 4.1. By choosing H = 2 max x∈B (x) 2 ,
H depends on X and then on n under Assumption (H.3). So β ≥ max(4 v ∞ , 2H/t 0 ) also depends on n. In this case, ξ must satisfy Assumption (P.2) with t 0 = ∞. In view of Proposition 3.1, we can consider ξ as a Gaussian or a bounded symmetric noise.

Choice of prior 4.3.1. Main assumptions

Recall that the goal is to find a prior leading an oracle inequality with a small remainder term while promoting group-analysis sparsity. In order to promote sparsity, it is well-known that the prior is usually expected to be symmetric and sharply peaked around its mode (the origin) while having tails heavier than merely exponential [START_REF] Müller | Bayesian Inference in Wavelet-Based Models[END_REF]. This is the rationale behind our general prior which takes the form

π(dθ) = 1 C α,g,R L =1 exp -α a [D θ] G a 2 g [ D θ] G 2 I Θ (θ)dθ, (4.2)
where α > 0 and g satisfies the following requirements:

(G.1) Boundedness: g : R + → R + is a bounded function such that g 0, and θ → g

[ D θ] G 2 is measurable on R p , ∀ ∈ {1, . . . , L}.
(G.2) Integrability:

R p L =1 g [ D u] G 2 du < ∞. (G.3) Moment condition: R p [ D u] G 2 2 L k=1 g [ D u] G k 2 du < ∞, ∀ ∈ {1, . . . , L}.
(G.4) Growth condition: there exist λ ≥ 0 and h : R

+ → R + such that g( t -t * 2 )/g( t 2 ) ≤ h( t * 2 ) λ , ∀(t, t * ) ∈ R G × R G .
The exponential part can be viewed as a generalized Gaussian on the group analysis coefficients with shape parameter a and scale parameter α -1 . The choice of a ∈]0, 1] favors sparsity. The choice of α will be made clear when discussing the remainder terms in our group-analysis SOI (see for instance Remark 5.1). The role of the function g is (at least) twofold. The first is precisely to capture heavier tails than exponential. Second, it allows more flexibility to adjust to the group sparsity scenario and optimize the performance of EWA (see remarks hereafter for further discussion and examples).

Discussion of the assumptions

Assumption (G.4) controls the growth of the function g. The growth function h will impact the remainder term in the main group-analysis SOI stated in Theorem 5.1, and more precisely the term Ω D µ n ,n,L,λ (θ) therein. Assumptions (G.1)-(G.3) play a prominent role in controlling the key constant K D a,g > 0 that is involved in the construction of our general group-analysis SOI in Theorem 5.1. The following remark formalizes the existence of this constant and relates it to the integrability and moment conditions on g.

∈]0, 1], there exists K D a,g ∈]0, ∞[ such that, ∀ ∈ {1, . . . , L}, R p [ D u] G 2a 2 L k=1 g [ D u] G k 2 du R p L k=1 g [ D v] G k 2 dv ≤ K D a,g . (4.3)
Proof. From Assumption (G.3) and the fact that g 0, one can show that (4.3) holds for a = 1. Moreover, since g 0 and g satisfies Assumption (G.2), it holds that

u → L =1 g [ D u] G 2 du R p L k=1 g [D v] G k 2
dv is a probability measure. Therefore, (4.3) holds for any a in ]0, 1] by Hölder's inequality.

It is then legitimate to wonder whether there exists a simpler condition on g which implies Assumptions (G.2) and (G. 3). The answer is affirmative as we now state. When the operator D is invertible, K D a,g takes an even simpler and explicit form, and moreover, (4.4) is necessary for g to obey Assumption (G.3). This is the subject of Remark 4.

K D a,g = ∞ 0 x G-1+2a g(x)dx ∞ 0 z G-1 g(z)dz . (4.5)
Proof. The proof follows by combining Lemmas 2.3 and 2.2, i.e.,

R p [ D u] G 2a 2 L k=1 g [ D u] G k 2 du R p L k=1 g [ D v] G k 2 dv = R G u 2a 2 g( u 2 )du R G g( v 2 )dv L-1 R G g( w 2 )dw L = R G u 2a 2 g( u 2 )du R G g( v 2 )dv = ∞ 0 x G-1+2a g(x)dx ∞ 0 z G-
R p [ D u] G 2 2 L k=1 g( [ D u] G k 2 )du = C L G ∞ 0 z G+1 g(z)dz det(D D ) ∞ 0 w G-1 g(w)dw L-1
.

Examples of g

Let us now discuss some choices of g.

Example 4.3. Consider g : R + → R + defined by

g(x) = 1 (τ 2 + x 2 ) 2 , τ > 0.
This choice of g yields a prior that specializes to the one in [START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF] for the individual sparsity scenario, i.e., with D = I p , G = 1 and a = 1.

Example 4.4. Consider g : R + → R + defined by

g(x) = 1 (τ b + x b ) c , where τ > 0, b ∈]0, 1] and c > (2 + G)/b.

The choice of c guarantees the validity of Assumptions (G.2) and (G.3).

Thanks to the parameters b and c, this choice of g offers more flexibility than the one in the previous example. This allows for example to optimize the performance of EWA by tuning these parameters for the particular dataset at hand.

Group-analysis sparse oracle inequality

Once a suitable dictionary and prior are chosen according to the above, the EWA is performed via (1.3)-(1.4). Our goal now is to provide a theoretical guarantee for the aggregates by constructing a group-analysis SOI. First of all, based on PAC-Bayesian type oracle inequalities in Section 3, we establish our first main result: a group-analysis SOI for the dictionary (4.1) and the prior (4.2) with a function g obeying Assumptions (G.1)-(G.4). Theorem 5.1 (General group-analysis sparse oracle inequality). Let G ≥ 1, X satisfying Assumption (H.3), and D satisfying Assumption (H.2) with κ > 0. Let Assumptions (P.1) and (P.2) be satisfied with some function v, (3.1) holds and β ≥ max (4 v ∞ , 2H/t 0 ). For some a ∈]0, 1], take the dictionary (4.1) and the prior (4.2) with g satisfying Assumptions (G.1)-(G.4). Let K D a,g , as defined in (4.3), and assume that R > 3 K D a,g L. Then the following groupanalysis SOI holds,

E f n -f 2 n ≤ inf Θ D µn ,L,R f θ -f 2 n + Φ D µ n ,n,L (θ) + Ω D µ n ,n,L,λ (θ) + Ψ D µ n ,L,p , (5.1
)

with                            Θ D µ n ,L,R = θ ∈ R p : D θ a a,G ≤ R -3 K D a,g L Φ D µ n ,n,L (θ) = β n 1 + 3 K D a,g Lα a + α a D θ a a,G Ω D µ n ,n,L,λ (θ) = λβ n L =1 ln h D θ G 2 Ψ D µ n ,L,p = 2κ -1 K D 1,g exp 3 K D a,g Lα a pC f,L ,
and

C f,L = L 2 ∞ + L ∞ L ∞ + f ∞ .
Before proceeding, we pause to make a few important remarks.

Remark 5.1. The group-analysis SOI (5.1) is sharp. It depends on several parameters as discussed below.

(i) The parameter R appears in the dictionary. Namely, the EWA estimator f n mimics the best aggregate f θ for all possible weights belonging to

θ ∈ R p : D θ a a,G ≤ R -3 K D a,g L . Then R must be sufficiently large to cover the "good" model f θ in Assumption (H.1). Moreover, since R > 3 K D a,g L, R ∼ K D a,g L is the smallest choice we can make to reduce the rate of Φ D µ n ,n,L (θ) as D θ a a,G ≤ R.
(ii) The parameter α is used to cancel the effect of L K D a,g in the remainder terms. By choosing α ≤ (3K

K D a,g ) -1/a , Φ D µ n ,n,L (θ) ≤ βn -1 (1 + 3 K D a,g α a + α a R) ∼ n -1 and Ψ D µ n ,L,p ≤ 2eC f,L κ K D 1,g p.
(iii) The parameter K D a,g and the function h depend on the choice of g. They respectively control the rate of

Ψ D µ n ,L,p
and Ω D µ n ,n,L,λ (θ). In what follows, let us state the consequences of Theorem 5.1 with the choices of g in Example 4.3 and 4.4. Especially, we will discuss the rate of Ω D µ n ,n,L,λ (θ) and Ψ D µ n ,L,p .

We first consider the prior (4.2) in Example 4.3, under the individual sparsity scenario ( D = I p , G = 1) and the choice a = 1 (i.e., Θ = θ ∈ R p : θ 1 ≤ R ). This is the setting considered in [START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF]. We obtain the following SOI as a corollary of our main result.

Corollary 5.1. Let X satisfying Assumption (H.3), D = I p and fix G = 1. Suppose that Assumptions (P.1) and (P.2) hold with some function v, (3.1) holds and β ≥ max (4 v ∞ , 2H/t 0 ). Fix a = 1, take the dictionary (4.1) and the prior (4.2) with g defined in Example 4.3 and α ≤ 1/(3pτ). Assume that R > 3pτ. Choosing τ 2 ∼ 1/(pn) and R ∼ pτ, SOI (5.1) holds with

Θ D µ n ,L,R = θ ∈ R p : θ 1 ≤ R -3pτ , Ω D µ n ,n,L,λ (θ) ∼ θ 0 ln(p) n , and
Ψ D µ n ,L,p ∼ 1 n .
The order of Ω D µ n ,n,L,λ (θ) is the classical rate under the sparsity scenario. This scaling is similar as in [START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF] with the same prior. However, the following remark shows that this prior is not adapted in the group-analysis case for any group size strictly larger than 1. To get an explicit control of the remainder term, it is instructive to have a closed-form of K D a,g . This can be done for instance when D is invertible, see (4.5). The obtained group-analysis SOI is stated as follows. For a ∈]0, 1], let

K D a,g = Γ((2a+G)/b)Γ(c-(2a+G)/b) Γ(G/b)Γ(c-G/b)
, and set α ≤ 1/ 3τ a K D a,g L

1/a
. Choosing τ 2 ∼ 1/(pn) and R ∼ Lτ a , the group-analysis

SOI (5.1) holds with Θ D µ n ,L,R = θ ∈ R p : D θ a a,G ≤ R -3τ a K D a,g L , Ω D µ n ,n,L,λ (θ) ∼ D θ 0,G ln(L) n , and
Ψ D µ n ,L,p ∼ 1 n .
By Assumption (H.1), D θ 0,G is small when θ = θ (with R must be sufficiently large to cover θ). Thus, D θ 0,G ln(L) is small compared to n. Under the sparsity scenario, the order of Ω D µ n ,n,L,λ (θ) becomes O( θ 0 ln(p)/n) which is the same rate as the aggregate with g in Example 4.3.

Forward-Backward proximal LMC algorithm

The goal of this section is to implement our EWA estimator with the probability measure (1.3) via a novel forwardbackward proximal Monte-Carlo algorithm based on the Langevin diffusion (coined FB-LMC).

Let us consider a linear regression problem where L (x) = x, and thus 2

f n = X θ n ,
where X ∈ R n×p is the design matrix (see Section 4.2). Recall the EWA estimator

θ n = R p θµ n (θ)dθ,
with the measure

µ n (θ) ∝ exp          -y -Xθ 2 2 β          π(θ). (6.1)
2 Generalization to non-linear link functions is possible from a practical point of view. However, in this case, the convergence guarantees of the discrete Langevin-based MCMC sampling scheme would be much more intricate even for a simple exponential link function. One of the main difficulties lies in that the data loss θ → y -L (Xθ) 2 2 is not necessarily convex which prevents us from showing that the Langevin diffusion is geometrically ergodic. This is an open question that we leave to a future work.

Computing θ n corresponds to an integration problem which becomes very involved to solve analytically or even numerically in high-dimension. A classical alternative is to approximate it via a Markov chain Monte-Carlo (MCMC) method which consists in sampling from µ n by constructing an appropriate Markov chain whose stationary distribution is µ n , and to compute sample path averages based on the output of the Markov chain. The theory of MCMC methods is based on that of Markov chains on continuous state space. As in [START_REF] Dalalyan | Sparse regression learning by aggregation and langevin monte-carlo[END_REF], we here use the Langevin diffusion process; see [START_REF] Roberts | Exponential Convergence of Langevin Distributions and Their Discrete Approximations[END_REF].

The Langevin diffusion

Continuous dynamics A Langevin diffusion L in R p , p ≥ 1 is a homogeneous Markov process defined by the stochastic differential equation (SDE)

dL(t) = 1 2 ψ(L(t))dt + dW(t), t > 0, L(0) = l 0 , (6.2)
where ψ = ∇ ln(µ), µ is everywhere non-zero and a suitably smooth target density function on R p , W is a p-dimensional Brownian process and l 0 ∈ R p is the initial value. Under mild assumptions, the SDE (6.2) has a unique strong solution and, moreover, L(t) has a stationary distribution with density precisely µ [51, Theorem 2.1]. L(t) is therefore interesting for sampling from µ. In particular, this opens the door to approximating integrals R p f (θ)µ(θ)dθ by the average value of a Langevin diffusion, i.e., 1 T T 0 f (L(t))dt for a large enough T . Under additional assumptions on µ and f in a proper functional class, the expected squared error of the approximation can be controlled [START_REF] Xuerong | Stochastic differential equations and applications[END_REF].

Forward Euler discretization

In practice, in simulating the diffusion sample path, we cannot follow exactly the dynamic defined by the SDE (6.2). Instead, we must discretize it. A popular discretization is given by the forward (Euler) scheme, which reads

L k+1 = L k + δ 2 ψ(L k ) + √ δZ k , t > 0, L 0 = l 0 ,
where δ > 0 is a sufficiently small constant discretization step-size and {Z k } k are i.i.d. ∼ N(0, I p ). The average value L k .

It is then tempting to approximate θ n by applying this discretization strategy to the Langevin diffusion with µ n in (6.1) as the target density. However, quantitative consistency guarantees of this discretization require µ (hence ψ) to be sufficiently smooth, which limits their applicability in our context. To cope with this difficulty, a few works have recently proposed to replace ln(µ) with a smoothed version (typically involving the Moreau-Yosida regularization/envelope, see Definition 2.2) [START_REF] Durmus | Non-asymptotic convergence analysis for the Unadjusted Langevin Algorithm[END_REF][START_REF] Durmus | Sampling from convex non continuously differentiable functions, when Moreau meets Langevin[END_REF][START_REF] Pereyra | Proximal markov chain monte carlo algorithms[END_REF]. In [START_REF] Durmus | Sampling from convex non continuously differentiable functions, when Moreau meets Langevin[END_REF][START_REF] Pereyra | Proximal markov chain monte carlo algorithms[END_REF] for instance, the authors proposed proximal-type algorithms to sample from possibly non-smooth log-concave densities µ using the forward Euler discretization and the Moreau-Yosida regularization. In [START_REF] Pereyra | Proximal markov chain monte carlo algorithms[END_REF] 3 ,ln(µ) is replaced with its Moreau envelope, while in [START_REF] Durmus | Sampling from convex non continuously differentiable functions, when Moreau meets Langevin[END_REF], it is assumed thatln(µ) = F + G, F is convex Lipschitz continuously differentiable, and G is a proper closed convex function replaced by its Moreau envelope. In both these works, convexity plays a crucial role to get quantitative convergence guarantees and thus cannot be applied to our prior. Proximal steps within MCMC methods have been recently proposed for some simple (convex) signal processing problems [START_REF] Chaari | A hamiltonian monte carlo method for non-smooth energy sampling[END_REF], though without any guarantees.

Forward-Backward proximal Langevin MC algorithms

Consider the prior π in (4.2) with g is given in Example 4.4 and H = +∞. Then, µ n is neither differentiable nor log-concave. To overcome these difficulties, we will exploit the structure of µ n and some arguments from variational analysis [START_REF] Rockafellar | Variational analysis[END_REF]. For ease of notation, in the following, we denote with the same symbol the measure and its density with respect to the Lebesgue measure. Thus µ n reads

µ n (θ) ∝ exp (-(V(θ)) (6.3)
where

V def = F β + W λ • D , F β (θ) = y -Xθ 2 2 /β and W λ (u) = L =1 w λ u G 2 with w λ : x ∈ [0, +∞[ → α a x a + c ln(τ b + x b ), (6.4)
and

w λ is parameterized by λ = (a, b, c, α, τ) ∈]0, 1]×]0, 1]×]2 + G/b, +∞[×R + × R +, * .
We start by collecting some important properties on the function w λ and its proximal operator, as well as their implications. We denote

V γ = F β + ( γ W λ ) • D .
prox γw λ (x) =        0 if x ≤ γw λ (0 + ), t -γw λ prox γw λ (x) if x > γw λ (0 + ). (ii) For any x ≥ 0, 0 ≤ prox γw λ (x) ≤ x. (iii) prox γW λ (u) = prox γw λ u G 1 2 u G 1 u G 1 2 , . . . , prox γw λ u G L 2 u G L u G L 2 . (iv) prox γW λ •D is Lipschitz continuous. (v) ∇( γ W λ • D ) = γ -1 D • (I q -prox γW λ ) • D .

It is a uniformly bounded and Lipschitz continuous operator and

∇V γ is a Lipschitz continuous mapping.

(vi) Assume that γ ∈ 0, min

τ 2 /c, β/ 2 X 2 . Let M γ def = I p -(2γ/β)X X which is a symmetric positive definite matrix. Then, ∇ M γ ,γ V = γ -1 M γ I p -prox γW λ •D (I p -γ∇F β
) which is a Lipschitz continuous mapping.

We now describe two Langevin MC (LMC) sampling algorithms, originally proposed in [START_REF] Luu | Sampling from non-smooth distribution through Langevin diffusion[END_REF], that are based on forward-backward proximal splitting and establish their guarantees for the penalty w λ . In the rest of the section, we will fix a = b = 1.

Semi-Forward-Backward LMC (Semi-FBLMC)

Assume that γ ∈]0, τ 2 /c[. Define the following SDE with the Moreau-Yosida regularized version of W λ dL(t) = - 1 2 ∇V γ (L(t))dt + dW(t), t > 0. (6.5)
Inserting Lemma 6.1(v) into (6.5), the Euler discretization of (6.5) reads

L k+1 = L k - δ 2 ∇F β (L k ) - δ 2γ D D L k -prox γW λ ( D L k ) + √ δZ k , t > 0, L 0 = l 0 . (6.6)

Forward-Backward LMC (FBLMC)

Assume that γ ∈ 0, min τ 2 /c, β/ 2 X

2

. One can consider an alternative version of the SDE (6.5) with the Moreau-Yosida regularized version of V in the metric M γ (see Lemma 6.1(vi)), i.e.,

dL(t) = - 1 2 ∇ M γ ,γ V • M -1/2 γ (L(t))dt + M 1/2 γ dW(t), t > 0. (6.7)
By the change of variable U(t) = M -1/2 γ L(t) we get the following SDE

dU(t) = - 1 2 M -1 γ ∇ M γ ,γ V(U(t))dt + dW(t), t > 0. (6.8)
In view of Lemma 6.1(vi), the Euler discretization of (6.8) is given by

U k+1 = (1 -δ 2γ )U k + δ 2γ prox γW λ • D (U k -γ∇F β (U k )) + √ δZ k , t > 0, U 0 = l 0 . ( 6
µ γ (θ) ∝ exp -V γ (θ) .
The same also holds for (6.7) with the corresponding invariant measure. The following claim is a consequence of [40, Proposition 3.1] and Lemma 6.1. In the sequel, ν TV stands for the total variation norm of a signed measure ν. Proposition 6.1. Let µ γ be the invariant measure of either (6.5) or (6.7). Then, µ γ -µ n TV → 0 as γ → 0.

We consider the Semi-FBLMC discretization (6.6) of the Langevin diffusion (6.5). Let Lδ be the continuous-time extension of the scheme (6.6)

Lδ (t) def = L 0 - 1 2 t 0 ∇V γ ( L(t))dt + t 0 dW(s),
where L(t) = L k for t ∈ [kδ, (k + 1)δ[. We write P t,γ L (l 0 , Ω) = Pr {L(t) ∈ Ω|L(0) = l 0 } for all Borel sets Ω and initial condition l 0 . Similarly, we denote P t,δ,γ Lδ (l 0 , Ω) = Pr Lδ (t) ∈ Ω| Lδ (0) = l 0 . The superscripts stress the dependence on the parameters.

We also define the sample average Lδ,T,γ = δ T T/δ -1 k=0

L k , and the EWA estimate

θ γ = R p θµ γ (θ)dθ.
Theorem 6.1. The following holds:

(i) lim γ→0 lim T →+∞ lim δ→0 P T,δ,γ Lδ (l 0 , •) -µ n TV = 0.
(ii) Suppose that l 0 = 0, then lim T →+∞ lim δ→0 E Lδ,T,γθ γ 2 = 0.

To prove these claims, we first show that the Langevin diffusion in (6.5) is uniformly geometrically ergodic. Then we follow standard arguments, invoking Lemma 6.1, Proposition 6.1, the Girsanov formula and Pinsker inequality 4 . Similarly to Theorem 6.1, convergence guarantees of the FBLMC discretization scheme (6.9) can be established. We omit the details here for the sake of brevity.

Numerical experiments

In this section, some numerical experiments are conducted to illustrate and validate the numerical performance of the proposed EWA estimator. We consider a linear regression problem

Y = Xθ 0 + ξ,
where ξ is the noise, X is the design matrix. θ 0 ∈ R p is the unknown regression vector of interest assumed to obey Assumption (H.1). Since H = +∞, we then have to choose a distribution on the noise ξ such that β is independent of H. Such type of distributions is specified in [23, Section 2]. For our implementation, we assume ξ i.i.d. ∼ N(0, σ 2 ). The noise level σ is chosen according to the simulated θ 0 .

The parameters of EWA were chosen as prescribed in our theoretical analysis. For instance, the temperature parameter is set to β = 4σ 2 , the parameters of the prior:

a = b = 1, c > (2 + G)/b, τ ∼ 1/(pn) and α ≤ 1/ 3 K D a,g L 1/a
, where K D a,g is given in Corollary 5.3. The number of iterations N and the step-size δ are chosen respectively large and small enough to guarantee convergence and discretization consistency of the algorithm.

Following the philosophy of reproducible research, all the code implementing our EWA algorithm and reproducing the experiments of this paper are made publicly available for download at https://github.com/luuduytung/ GroupAnalyseEWAToolbox.

1-D signal recovery under group sparsity

In this example, we set D = I p , which corresponds to the classical group sparsity. The design matrix is drawn uniformly at random from the Rademacher ensemble, i.e., its entries are i.i.d. variates valued in {-1, 1} with equal probabilities. The non-zero entries of θ 0 are equal to 1 and we denote S = θ 0 0 the sparsity level of θ 0 . Two types of sparsity behavior are considered: individual sparsity where G θ 0 = 1; group structured sparsity with G θ 0 = 4. Besides, the positions of the non-zero/active entries (for G θ 0 = 1) or groups (for G θ 0 = 4) are chosen randomly uniformly on {1, . . . , p}.

The experiments are performed by fixing p = 128, and taking S ∈ 2 {2,...,7} , n ∈ 2 {3,...,7} , step-size δ = 4σ 2 /(np) and integration time T = 3500. The parameters in the prior are chosen to minimize the remainder term in the oracle inequality (5.1). For each (S , n), and each value of G θ 0 , N rep = 20 instances of the problem suite (X, θ 0 , Y) are generated, and EWA is applied with a chosen G and the other parameters as detailed above. The estimation quality/success is then assessed by

π S ,n = 1 N rep N rep j=1 I θ ( j,S ,n) n -θ ( j,S ,n) 0 n ≤ , (7.1)
where > 0 (we choose = 0.4) and θ ( j,S ,n) n (resp. θ ( j,S ,n)

0

) corresponds to θ n (resp. θ 0 ) in the j-th replication of (S , n). S /p and n/p are respectively normalized measures of sparsity and problem indeterminacy. We get a twodimensional phase space (S /p, n/p) ∈ [0, 1] 2 describing the difficulty of a problem instance, i.e., problems are easier as one moves up (more measurements) and to the left (sparser θ 0 ). Phase diagrams plotting π S ,n in (7.1) as a function (S /p, n/p) were widely advocated by Donoho and co-authors for 1 minimization [START_REF] Donoho | Observed universality of phase transitions in high-dimensional geometry, with implications for modern data analysis and signal processing[END_REF]. Such diagrams often have an interesting two-phase structure (as displayed in Figures 1(a)-(d), brighter color indicate better success), with phases separated by a specific curve, called phase transition curve. Thus, a good estimator is intended to have a large bright area which indicates its good performance at a wider range of (S , n).

Figure 1(a) (resp. (b)) shows the phase diagrams when G θ 0 = 1 and G = 1 (resp. G = 4) in EWA. In this case, the phase transition curve for G = 1, the correct group size, is slightly better that with G = 4. The situation reverses for Figures 1(c)-(d) where G θ 0 = 4, and one observes that the success area is significantly better using G = 4 than G = 1. This is expected as it reveals better performance of EWA when used with the choice G = G θ 0 . This is also confirmed by visual inspection of Figures 1(c')-(d'), where we plotted instances of recovered vectors θ ( j,S ,n) n when (S , G θ 0 ) ∈ {4, 8} × {1, 4} and n/p = 1/2. EWA was again applied with G = 1 and G = 4 in each case. Large spurious entries appear outside the true support when the group size is not correctly chosen, though the impact is less important for G = 1.

It is worth observing that S /p = θ 0 0,G G θ 0 /p. As far as the expected phase transition curve is concerned, one has from Corollary 5.3 that it is expected to occur for

n/p = C θ 0 0,G G θ 0 /p ln(p/G θ 0 )/G θ 0 = C S /p ln(p/G θ 0 )/G θ 0
for some constant C > 0 depending on . That is, the phase transition curve is linear (p and G θ 0 are fixed for each diagram), which is confirmed by visual inspection of Figures 1(a)-(d), where the overlaid blue line is the fitted linear phase transition curve.

2-D image recovery under analysis group-sparsity

In the second numerical experiment, θ 0 is a 2-D image which is a matrix in R 160×160 (a close-up of the known Shepp-Logan phantom, see Figure 2(a)). Thus vec(θ 0 ) is vector in R p with p = 160 2 , and our goal is to recover

θ 0 with values in [0, 1] from Y = X vec(θ 0 ) + ξ,
where ξ ∼ N(0, σ 2 I n ) and X ∈ R n×p is again random whose entries are i.i.d. from the Rademacher distribution. Since the targeted image is piecewise-constant, a popular prior is the so-called isotropic total variation [START_REF] Rudin | Nonlinear total variation based noise removal algorithms[END_REF] which is described in Example 4.1. It turns out that this can be cast in our analysis group-sparsity framework as a special case.

In this experiment, we use Semi FBLMC Algorithm to compute the EWA estimator. The results are depicted in Figure 2. σ = 0.525 in this experiment, the number of observations is n = 9p/16 = 14400, and we have D TV (θ 0 ) 0,G = 1376 n. A notable property of the EWA estimate is that it does not suffer from the stair-casing effect, unlike total variation minimization.

Conclusion

In this paper, we proposed a class of EWA estimators constructed from a novel and versatile family of priors which promotes analysis group-sparsity, where the analysis operator corresponds to a frame. Its quality is guaranteed by establishing a sharp SOI with a small remainder term in high-dimension. We also described a forward-backward proximal LMC algorithm, which is an implementation of EWA and can be viewed as a Forward Euler discretization of a Langevin diffusion involving the Moreau-Envelope of the potential in a proper metric. We derived convergence guarantees of this discretization. The performance of the estimator was illustrated on some numerical experiments which support our theoretical findings. There are still open problems that we leave to a future work. More precisely, one direction is to investigate how to remove the frame assumption. Another one would be to derive further/better quantitative convergence bounds of the proposed discretization.

Proofs

Proofs of Section 2 Proof of Lemma 2.3 Consider the linear mapping

M : x ∈ R d → Mx ∈ R r , r ≥ d.
The Jacobian matrix of this mapping is obviously M for any x ∈ R d . Since M is a frame, it is injective, hence so-called d-regular (see [START_REF] Maly | Lectures on change of variables in integral[END_REF]Section 1.5]). In particular, det(M M ) > 0. Thus combining [41, Theorems 1.12 and 3.4] and the Cauchy-Binet formula [41, Theorem 3.3]), we have the change of variables formula

Θ u(Mx)dx = R r x∈Θ∩ ω : Mω=v u(Mx)dv det(M M ) = Im(M) x∈Θ∩ ω : Mω=v u(Mx)dv det(M M ) .
Using once again that M is a frame, i.e., it is bijective on its image Im(M), the result follows. This concludes the proof.

Proofs of Section 3 Proof of Proposition 3.1

• Gaussian noise: Let ξ ∼ N(0, Σ). We set ζ ∼ N(0, (2γ+γ 2 )Σ). Thus conditions (a) and (b) in Assumption (P.2) are verified. We check now condition (c).

Let t ∈ x ∈ R n : x 2 ∈ [-∞, +∞] , u = t 2γ + γ 2 Σ 1/2 and ∼ N(0, I n ), we get that E e t ζ = E e u = n j=1 E {e u j j } = e u 2 2 2 ≤ e 1 2 √ 2γ+γ 2 Σ 1/2 t 2 2 ≤ e γ+ γ 2 2 Σ t 2 2 .
Thus, let a ∈ R n and v(a) ≡ Σ , we then get

ln E e t ζ |ξ = a t 2 2 γv(a) ≤ 1 + γ 2 → γ→0 1 ≤ 1.
• Laplace noise: Let ξ ∼ L(0, Σ), i.e., its associated characteristic function is ϕ ξ (t) = 1 + t Σt/2 -1 , we choose ζ according to the distribution associated to the characteristic function

ϕ ζ (t) = (1 + γ) -2 1 + 2γ + γ 2 1 + (1 + γ) 2 t Σt/2 .
For any t ∈ R n , we get that

ϕ ξ+ζ (t) = ϕ ξ (t)ϕ ζ (t) = 1 1 + (1 + γ) 2 t Σt/2 = ϕ (1+γ)ξ (t). (9.1)
Thus, ζ + ξ has the same distribution as (1 + γ)ξ. We also obtain E {ζ|ξ} = E {ζ} = (-i)∇ϕ ζ (0) = 0. It suffices to check condition (c) of Assumption (P.2). We know that

E e t ζ |ξ = E e t ζ = ϕ ζ (-it) = 1 (1+γ) 2 1 + 2γ + γ 2 1 -(1 + γ) 2 t Σt/2 .
Using Taylor's formula, we have

ln E e tζ |ξ = γt Σt 1 -t Σt/2 + O(γ 2 ). Thus, let t ∈ x ∈ R n : x 2 ≤ t 0 , a ∈ R n and v(a) ≡ Σ 1-t 2 0 Σ 2 /2
, we get

ln E e t ζ |ξ = a t 2 γv(a) → γ→0 t Σt/(1 -t Σt/2) t 2 2 Σ /(1 -t 2 0 Σ /2) ≤ 1 -t 2 0 Σ /2 1 -t Σt/2 ≤ 1 -t 2 0 Σ /2 1 -t 2 2 Σ /2 ≤ 1.
We get two last inequalities under the condition 1t 2 0 Σ /2 > 0 equivalent t 0 < 2/ Σ .

• Bounded symmetric noise: Let ξ are symmetric and Pr

ξ i ≤ B i = 1 for some B ∈ R n , we set ζ = (ζ 1 , . . . , ζ n ) such that ζ i = (1 + γ) ξ i sgn(sgn(ξ i ) -(1 + γ)U i ) -ξ i , U i ∼ U([-1, 1]) for any i ∈ {1, . . . , n}. Using [23, Equation (22)
], for any t ∈ R n and a ∈ x ∈ R n :

x i ∈ [-B i , B i ], ∀i ∈ {1, . . . , n} , we get that E e t ζ |ξ = a = n j=1 E e t i ζ i |ξ i = a i = e -t a e (1+γ)t a 2 + γ 2 + 2γ + e -(1+γ)t a γ 2 + 2γ . (9.2)
From (9.2) and the symmetry of ξ, we obtain

E e t (ζ+ξ) = E E e t (ζ+ξ) |ξ = e (1+γ)t ξ 2 + γ 2 + 2γ + e -(1+γ)t ξ γ 2 + 2γ = E e (1+γ)t ξ .
Thus, ζ + ξ has the same distribution as (1 + γ)ξ. Since E {ζ|ξ = a} equals to the gradient of E e t ζ |ξ = a at t = 0, from (9.2) we have then

E {ζ|ξ = a} = 0, ∀a ∈ [-B, B].
It suffices to check the condition (c) of Assumption (P.2). Owing to [19, Lemma 3] and [23, Equation ( 22)], we get that ln

E e t i ζ i |ξ i =a i ≤ (t i a i ) 2 γ(1+γ). Thus, let t ∈ x ∈ R n : x 2 ∈ [-∞, ∞] and v(a) = a 2 2 , we get that ln E e t ζ |ξ = a t 2 2 γv(a) = n i=1 ln E e t i ζ i |ξ i = a i t 2 2 γ a 2 2 ≤ n i=1 t 2 i a 2 i γ(1 + γ) t 2 2 γ a 2 2 ≤ 1 + γ → γ→0 1 ≤ 1.

Proofs of Section 4

Proof of Proposition 4.1 Let θ ∈ Θ and i ∈ {1, . . . , n}. Setting u θ i = p j=1 θ j f j (x i ) and u θ = (u θ 1 , . . . , u θ n ) , and by virtue of (2.1), (2.2), (4.1) and the fact that a ∈]0, 1], we have

u θ 2 = Xθ 2 ≤ X D D θ 2 ≤ X D θ a,G √ κ ≤ X R 1/a √ κ ,
which in turn implies u θ ∈ B. Therefore, for any θ, θ

∈ Θ 2 , f θ -f θ 2 = u θ - u θ 2 ≤ 2 max x∈B (x) 2 .
Proof of Lemma 4.1 Let us first check the integrability condition (G.2). By Lemmas 2.3 and 2.2, we obtain

R p L =1 g( [ D u] G 2 )du = Im( D ) L =1 g( v G 2 )dv det(D D ) ≤ R q L =1 g( v G 2 )dv det(D D ) = R G g( u 2 )du L det(D D ) = C L G ∞ 0 z G-1 g(z)dz L det( D D )
.

Since G ≥ 1 and g : R + → R + , by (4.4), we get

R p L =1 g( [ D u] G 2 )du ≤ C L G 1 0 z G-1 g(z)dz + ∞ 1 w G-1 g(w)dw L det(D D ) ≤ C L G sup z∈[0,1] g(z) + ∞ 1 w G+1 g(w)dw L det(D D ) < ∞. (9.3)
Therefore, g satisfies Assumption (G.2). Now, we check the moment condition (G.3). Using similar arguments to the bound (9.3), we have

R p [D u] G 2 2 L k=1 g( [ D u] G k 2 )du ≤ R q v G 2 2 L k=1 g( v G k 2 )dv det(D D ) = R G u 2 2 g( u 2 )du R G g( v 2 )dv L-1 det(D D ) = C L G ∞ 0 z G+1 g(z)dz ∞ 0 w G-1 g(w)dw L-1 det(D D ) < ∞, (9.4)
whence we conclude that g satisfies Assumption (G.3). Remind the prior π(dθ) from (4.2), where

Proofs of

Θ = θ ∈ R p : D θ a a,G ≤ R . Let r L = 3 K D a,g L, Θ p 0 D = θ ∈ R p : D θ -D θ * a a,G ≤ r L and θ * ∈ θ ∈ R p : D θ a a,G ≤ R -3 K D a,g L = R -r L . (9.5)
We define the probability measure

p 0 D (dθ) = 1 C L dπ dθ (θ -θ * ) I Θ p 0 D (θ)dθ,
where C L > 0 is the normalization factor for p 0 D . Since r L < R, θ ∈ Θ p 0 D implies that θθ * ∈ Θ. Therefore,

p 0 D (dθ) = 1 C L L =1 exp -α a [ D θ -D θ * ] G a 2 g [ D θ -D θ * ] G 2 I Θ (θ -θ * )I Θ p 0 D (θ)dθ = 1 C L L =1 exp -α a [ D θ -D θ * ] G a 2 g [ D θ -D θ * ] G 2 I Θ p 0 D (θ)dθ.
For any i ∈ {1, . . . , n}, with X i = f 1 (x i ), . . . , f p (x i ) , one can write f θ (x i ) = L p j=1 θ j f j (x i ) = L X i θ . Taylor-Lagrange formula then gives us

( f θ (x i ) -f (x i )) 2 ≤ ( f θ * (x i ) -f (x i )) 2 + C f,L X i (θ -θ * ) 2 + 2( f θ * (x i ) -f (x i ))L X i θ * X i (θ -θ * ), (9.6)
where

C f,L = L 2 ∞ + L ∞ L ∞ + f ∞ .
By summing over i from 1 to n, normalizing by 1/n, taking the integral in Θ w.r.t. p 0 D , inequality (9.6) becomes

Θ f θ -f 2 n p 0 D (dθ) ≤ f θ * -f 2 n + C f,L R p 1 n n i=1 X i (θ -θ * ) 2 p 0 D (dθ) + 2 n n i=1 ( f θ * (x i ) -f (x i ))L X i θ * X i Θ (θ -θ * )p 0 D (dθ). (9.7)
Note that, the right term of inequality (9.7) corresponds to a sum of three components. In the following, we keep the first component and treat the other two.

Let us first show that the last component vanishes. Indeed, let θ ∈ Θ p 0 D , from (9.5) and the fact that a ∈]0, 1], we have

D θ a a,G = L =1 [ D θ] G a 2 ≤ L =1 [ D θ -D θ * ] G 2 + [D θ * ] G 2 a ≤ D θ -D θ * a a,G + D θ * a a,G ≤ r L + D θ * a a,G ≤ R. Then θ ∈ θ ∈ R p : D θ a a,G ≤ R = Θ.
Therefore, we have the embedding

Θ p 0 D ⊆ Θ. (9.8)
In what follows, we denote B a a,G (x) = z ∈ R q : z a a,G ≤ x , ∀x > 0 for brevity. By (9.8), property (2.1), Lemma 2.3 and symmetry of B a a,G (r L ) ∩ Im( D ), we obtain

Θ (θ -θ * )p 0 D (dθ) ∝ Θ∩Θ p 0 D (θ -θ * ) L =1 exp -α a [ D θ -D θ * ] G a 2 g [D θ -D θ * ] G 2 dθ = Θ p 0 D (θ -θ * ) L =1 exp -α a [ D θ -D θ * ] G a 2 g [D θ -D θ * ] G 2 dθ = D det(D D ) B a a,G (r L )∩Im(D ) z L =1 exp -α a z G a 2 g z G 2 d z = 0, (9.9)
which is the desired claim. We now bound the second term in the right hand side of (9.7). Define

p 0 (du) = 1 C L det(D D ) L =1 exp -α a u G a 2 g u G 2 I Im( D )∩B a a,G (r L ) (u)du. (9.10)
One can see that p 0 coincides with the probability measure p 0 D on R p via a change of variables of type (2.3). So, p 0 is a probability measure on R q . For any i, j ∈ {1, . . . , L}, i j, by a change of variables, we get

R q u G i u G j p 0 (du) = - R q u G i u G j p 0 (du), so R q u G i u G j p 0 (du) = 0. (9.11)
For any j ∈ {1, . . . , L}, as all groups have the same size, we have

R q u G j u G j p 0 (du) = R q u G 1 u G 1 p 0 (du). (9.12)
We obtain

R p 1 n n i=1 X i (θ -θ * ) 2 p 0 D (dθ) ((2.1) and Lemma 2.3) = 1 n R q X Du X Dup 0 (du) = 1 n R q tr u D X X Du p 0 (du) = 1 n tr X D X D R q
uu p 0 (du) ((9.11) and (9.12)

) = 1 n L =1 tr X D X D G R q u G 1 u G 1 p 0 (du) (Von Neumann's trace inequality) ≤ 1 n L =1 G j=1 σ j X D X D G σ j R q u G 1 u G 1 p 0 (du) ≤ 1 n R q u G 1 2 2 p 0 (du) L =1 tr X D X D G = 1 n tr X D X D R q u G 1 2 2 p 0 (du). ( 9
tr X D X D n ≤ p j=1 σ j X X n σ j D D ≤ σ 1 D D p j=1 σ j X X n ≤ p κ . (9.14)
Putting together (9.13) and (9.14), we get the bound

C f,L R p 1 n n i=1 X i (θ -θ * ) 2 p 0 D (dθ) ≤ C f,L p κ R q u G 1 2 2 p 0 (du). (9.15)
Thanks to (9.9) and (9.15), inequality (9.7) becomes

Θ f θ -f 2 n p 0 D (dθ) ≤ f θ * -f 2 n + C f,L p κ R q u G 1
E f n -f 2 n ≤ f θ * -f 2 n + C f,L p κ R q u G 1 2 2 p 0 (du) + β KL(p 0 D , π) n . (9.17)
To complete the proof, it remains to bound the last two terms in the right hand side of (9.17). This is the goal of the following lemma. Lemma 9.1. Consider the same framework as the one in Theorem 5.1, we have

R q u G 1 2 2 p 0 (du) ≤ 2K D 1,g e r L α a , (9.18)
and

KL(p 0 D , π) ≤ 1 + r L α a + λ L =1 ln h D θ * G 2 + α a D θ * a a,G . (9.19)
With r L = 3 K D a,g L, it follows from (9.17) and Lemma 9.1 that

E f n -f 2 n ≤ f θ * -f 2 n + β n 1 + 3 K D a,g Lα a + α a D θ * a a,G + λβ n L =1 ln h D θ * G 2 + 2K D 1,g e 3 √ K D a,g Lα a pC f,L κ .
According to (9.5), this completes the proof of Theorem 5.1.

Proof of Lemma 9.1

To prove Lemma 9.1, we need an intermediate result.

Lemma 9.2. Let s > L K D a,g . The following inequality holds

1 T u∈R p : D u a a,G >s L =1 g( [ D u] G 2 )du ≤ L 2 K D a,g s -L K D a,g 2 ,
where

T = R p L =1 g( [ D u] G 2 )du.
Proof of Lemma 9.2. Let U be a random vector in R p with density u → 1

T L =1 g( [ D u] G 2 )
, where T < ∞ by Assumption (G.2). By Chebyshev inequality, we have

1 T u∈R p : D u a a,G >s L =1 g( [ D u] G 2 )du = Pr L =1 [ D U] G a 2 > s = Pr L =1 [ D U] G a 2 -E [ D U] G a 2 > s -L =1 E [ D U] G a 2 ≤ E L =1 [ D U] G a 2 -E [ D U] G a 2 2 /(s -L =1 E [ D U] G a 2 ) 2 = var L =1 [ D U] G a 2 /(s -L =1 E [ D U] G a 2 ) 2 ≤ E L =1 [ D U] G a 2 2 /(s -L =1 E [ D U] G a 2 ) 2 . (9.20)
Next, by Cauchy-Schwartz inequality and Remark 4.2, we obtain

E L =1 [ D U] G a 2 2 ≤ E L L =1 [ D U] G 2a 2 ≤ L 2 K D a,g (9.21)
and by Jensen inequality

s -L =1 E [ D U] G a 2 ≥ s -L =1 E [ D U] G 2a 2 ≥ s -L K D a,g > 0. (9.22)
Thus, combining (9.20), (9.21) and (9.22), we get

1 T u∈R p : D u a a,G >s L =1 g( [ D u] G 2 )du ≤ L 2 K D a,g s -L K D a,g 2 .
We now turn to the proof of Lemma 9.1

Proof of Lemma 9.1. Let us begin by the proof of inequality (9.18). We have

R q u G 1 2 2 p 0 (du) = 1 C L det(D D ) B a a,G (r L )∩Im( D ) u G 1 2 2 L =1 e -α a u G a 2 g u G 2 du ≤ 1 C L det(D D ) B a a,G (r L )∩Im( D ) u G 1 2 2 L =1 g u G 2 du. (9.23)
In the following, we show inequality (9.18) by bounding the right term of inequality (9.23). By Lemma 2.3 and Remark 4.2, we get

B a a,G (r L )∩Im( D ) u G 1 2 2 L =1 g( u G 2 )du det(D D ) R p L =1 g( [ D u] G 2 )du ≤ Im( D ) u G 1 2 2 L =1 g( u G 2 )du det(D D ) R p L =1 g( [ D u] G 2 )du = R p [ D u] G 1 2 2 L =1 g( [D u] G 2 )du R p L =1 g( [D u] G 2 )du ≤ K D 1,g . Then 1 det(D D ) B a a,G (r L )∩Im(D ) u G 1 2 2 L =1 g( u G 2 )du ≤ K D 1,g T. (9.24)
We now bound C L -1 . By a change of variables, we obtain

C L -1 =          Θ p D 0 L =1 e -α a [ D θ-D θ * ] G a 2 g [ D θ -D θ * ] G 2 dθ          -1 =         u∈R p : D u a a,G ≤r L e -α a D u a a,G L =1 g( [ D u] G 2 )du         -1 ≤ e r L α a         u∈R p : D u a a,G ≤r L L =1 g( [D u] G 2 )du         -1 . Since r L = 3 K D a,g L > K D a,g L, Lemma 9.2 gives us C L -1 ≤ e r L α a         T         1 - 1 T u∈R p : D u a a,G >r L L =1 g( [ D u] G 2 )du                 -1 ≤ e r L α a T -1             1 - L 2 K D a,g (r L -L K D a,g ) 2             -1 = e r L α a T -1 1 - 1 4 -1 ≤ 2e r L α a T -1 . ( 9
, π) = R p ln p 0 D (dθ) π(dθ) p 0 D (dθ) = R p ln             C α,g,R C L L =1 e -α a [ D θ-D θ * ] G a 2 g [ D θ -D θ * ] G 2 L =1 e -α a [ D θ] G a 2 g [ D θ] G 2             p 0 D (dθ) = R q ln             C α,g,R C L L =1 e α a t G a 2 g t G -t * G 2 e α a t G -t * G a 2 g t G 2             p 0 (d t) = ln C α,g,R C L + α a L =1 R q t G a 2 -t G -t * G a 2 p 0 (d t) + L =1 R q ln         g t G -t * G 2 g t G 2         p 0 (d t),
where p 0 is a probability measure in R q defined in (9.10). We know that t * = D θ * , according to the fact that t

G a 2 -t G -t * G a 2 ≤ t * G a 2
and Assumption (G.4), we get

KL(p 0 D , π) ≤ ln C α,g,R C L + α a D θ * a a,G + λ L =1 ln h D θ * G 2 . (9.26)
Now, it remains to bound ln(C α,g,R /C L ). Remind that C α,g,R is the normalization factor of π, and thus

C α,g,R = Θ L =1 exp -α a [ D θ] G a 2 g [ D θ] G 2 dθ ≤ R p L =1 g [ D θ] G 2 dθ = T.
Combining this with the bound of C L -1 in (9.25), we obtain

ln C α,g,R C L ≤ r L α a + ln(2) ≤ 1 + r L α a . (9.27)
Inserting (9.27) into (9.26), we get inequality (9.19). This completes the proof.

Proofs of corollaries

Proof of Corollary 5.1 Let γ = 3, ν = 2 and η = 1. We have γ/ν < η + 1 so that Lemma 2.1 applies. We thus obtain

∞ 0 z G+1 g(z)dz = ∞ 0 z 2 (z 2 + τ 2 ) 2 dz < ∞.
From Lemma 4.1, g satisfies Assumptions (G.2) and (G.3). Moreover, taking h(t) = 1 + t/τ and λ = 4, for all (t, t * ) ∈ R 2 , we have by Young's inequality

g(|t -t * |) g(|t|) = τ 2 + t 2 τ 2 + (t -t * ) 2 2 = 1 + 2τ(t -t * )t * /τ + t * 2 τ 2 + (t -t * ) 2 2 ≤ 1 + |t * | τ + t * 2 τ 2 2 ≤ h(|t * |) λ .
Therefore, g satisfies Assumptions (G.1)-(G.4) for G = 1. Owing to Remark 4.3 and Lemma 2.1, we obtain

K D 1,g = ∞ 0 x 2 (τ 2 +x 2 ) 2 dx ∞ 0 1 (τ 2 +y 2 ) 2 dy = τ 2 .
We are now in position to apply Theorem 5.1 with D = I p (then q = p), G = 1 (then L = q), a = 1 and α ≤ 1/(3pτ) to conclude. Namely, since τ 2 ∼ (pn) -1 and R ∼ pτ, we get that Ψ D µ n ,L,p ≤ 2eC f,L τ 2 p ∼ n -1 , and

Ω D µ n ,n,L,λ (θ) = 4β n p j=1 ln(1 + |θ j |/τ) ≤ 4β n θ 0 ln(1 + R/τ) ∼ θ 0 ln(p) n .
This completes the proof.

Proof of Corollary 5.2 Let γ = 2 + G, ν = b and η = c -1. We have γ/ν < η + 1 and thus Lemma 2.1 applies, whence we obtain

∞ 0 x G+1 g(x)dx = ∞ 0 x G+1 (τ b + x b ) c dx < ∞. From Lemma 4.1, g satisfies Assumptions (G.2) and (G.3). Recall that b ∈]0, 1]. Taking h(x) = 1 + (x/τ) b and λ = c, for all (t, t * ) ∈ R G × R G , we have g( t-t * 2 ) g( t 2 ) = τ b + t b 2 τ b + t-t * b 2 c ≤ τ b + t-t * b 2 + t * b 2 τ b + t-t * b 2 c ≤ 1 + t * b 2 τ b + t-t * b 2 c ≤ h( t * 2 ) λ .
Therefore, g satisfies Assumptions (G.

K D a,g = ∞ 0 r G-1+2a (τ b +r b ) c dr ∞ 0 q G-1 (τ b +q b ) c dq = τ 2a Γ 2a+G b Γ c -2a+G b Γ G b Γ c -G b = K D a,g τ 2a .
Since τ 2 ∼ (pn) -1 and R ∼ Lτ a , we get that

Ψ D µ n ,L,p ≤ 2C f,L K D 1,g eκ -1 pτ 2 ∼ n -1 ,
and

Ω D µ n ,n,L,λ (θ) = cβ n L =1 ln         1 +       [ D θ] G 2 τ       b         ≤ cβ n D θ 0,G ln        1 + R 1/a τ b        ∼ D θ 0,G ln(L) n .
This ends the proof. (ii) For γ ≤ τ 2 /c, it is immediate to see that the function u defined in (i) is nondecreasing on [0, +∞[, and thus so is prox γw λ . It then follows that prox γw λ (x) ≥ prox γw λ (0) = 0 for any x ∈ [0, +∞[. Since w λ is increasing, the second inequality follows.

(iii) See [START_REF] Luu | Sampling from non-smooth distribution through Langevin diffusion[END_REF]Lemma 7.1].

(iv) Denote the function x ∈]0, +∞[ → v(x) = x r + w λ (x) for r > 0. For x, z ≥ 0, we have

(v(z) -v(x))(z -x) = (z -x) 2 /(r) -c (z -x) 2 (τ + z)(τ + x) ≥ (z -x) 2 (r -1 -cτ -2 ) ≥ 0, ∀r ∈]0, τ 2 /c].
This shows that v is nondecreasing, or equivalently, that w λ (x) + 1 2r x 2 is convex (i.e., w λ is semi-convex). Thus

u 2 2 2r + W λ (u) = L =1          u G 2 2 2r + w λ u G 2          = L =1 (•) 2 2r + w λ • • 2 (u G ) .
Since x → x 2 2r + w λ (x) is convex and increasing, and the norm is also convex, we deduce that

• 2 2
2r + W λ is convex, and so is (v) [START_REF] Luu | Sampling from non-smooth distribution through Langevin diffusion[END_REF]Lemma 3.2], which applies thanks to (i), yields the expression of prox γW λ . ∇( γ W λ • D ) is Lipschitz continuous thanks to (iii), which, together with Lipschitz continuity of ∇F β yields that of ∇V γ . Let us turn to uniform boundedness of ∇( γ W λ • D ). For any θ ∈ R p , denote u = D θ. Thus, using (i) and (iii) and that w λ is a decreasing function on ]0, +∞[, we have ∀θ ∈ R p , 2 R 2 + 2β -1 X y 2 R + η 2 /16 and R = max(r, 1 + 2p/η). In turn, L(t) in (6.5) satisfies (9.28).

∇( γ W λ • D )(θ) 2 2 ≤ γ -2 ν L =1 u G 2 -prox γw λ u G 2 2 ≤ γ -2 ν L =1        γw λ prox γw λ u G 2 if u G 2 > γw λ (0 + ) γw λ (0 + ) otherwise ≤ γ -1 νLw λ (0 + ) = γ -

Proof of Lemma 9.3. Elementary derivations give

A H(θ) = -η/4 H(θ)

1 + θ 2 2             ∇F β (θ), θ + ∇( γ W λ )(D θ), D θ -η/4 θ 2 2 1 + θ 2 2 -p + θ 2 2 1 + θ 2 2            
.

Denote u = D θ. In view of Lemma 6.1(ii), (iii) and (v), we have

∇( γ W λ )(u), u = γ -1 u 2 2 -prox γW λ (u), u = γ -1 L =1 u G 2 u G 2 -prox γw λ u G 2 ≥ 0. (9.31)
This, together with (9.30) yields the bound

A H(θ) ≤ -η/4H(θ) ηR -p 1 + R -η/4 I B c R (θ) + η/4H(R) η/4 + d + 2β -1 X 2 R 2 + 2β -1 X y 2 R I B R (θ),
whence we deduce the desired claim.

Let us now turn to the proof of the theorem.

(i) By the triangle inequality, we have P T,δ,γ Lδ (l 0 , •) -µ n TV ≤ P T,δ,γ Lδ (l 0 , •) -P T,γ L (l 0 , •)

TV + P T,γ L (l 0 , •) -µ γ TV + µ γ -µ n TV .

In view of Lemma 6.1(v) and [START_REF] Dalalyan | Theoretical guarantees for approximate sampling from smooth and log-concave densities[END_REF]Lemma 2] and the Pinkser inequality, the first term in the right hand side goes to 0 as δ → 0. The second term converges to 0 as T → +∞ thanks to uniform geometric ergodicity (see Lemma 9.3), where we apply (9.28) with f ∞ ≤ 1. The last term vanishes as γ → 0 by virtue of Proposition 6.1.

(ii) In the following, C is any positive constant that does not depend on δ and T . Let δ ∈ 0, 2β/ X 2 . Since E {Z k } = 0 and ∇( γ W λ ) • D is uniformly bounded, we have

E {L k+1 } 2 ≤ E L k -δ/2∇F β (L k ) 2 + Cδ ≤ E (I p -δ/2∇F β )(L k ) -(I p -δ/2∇F β )(0) 2 + C + ∇F β (0) 2 /2 δ ≤ E {L k } + (C + X y 2 /2)δ ≤ (C + X y 2 /2)T, ∀k[0, T/δ ],
where we used the fact that F β is a differentiable convex function whose gradient is 2 X 2 /β-Lipschitz, and thus, I p -δ∇F β is non-expansive for the prescribed choice of δ [START_REF] Bauschke | Convex analysis and monotone operator theory in Hilbert spaces[END_REF]. In addition, by independence of Z k from L k and Y, and in view of (9.31) and Lipschitz continuity of ∇F β , we have

E L k+1 2 2 ≤ E L k -δ/2∇V γ (L k ) 2 2 + δp ≤ E L k -δ/2∇F β (L k ) + δ/2∇F β (0) 2 2 -δ L k -δ/2(∇F β (L k ) -∇F β (0)), ∇( γ W λ • D )(L k ) + Cδ 2 + δp ≤ E L k
E Lδ,T,γ -θ γ 2 ≤ E Lδ,T,γ - LT δ 2 + E LT,R δ -L T,R 2 + E L T,R -θ R γ 2 + 1 T T 0 E Lδ (t) 2 I B c R ( Lδ (t)) dt + B c R θ 2 µ γ (θ)dθ.
For the last two terms, we have the bounds This completes the proof.

1 T T 0 E Lδ (t) 2 I B c R ( Lδ (t)) dt ≤ R -1 1 T T 0 E Lδ (t)

By vectorizing θ

  in Assumption (H.1), D can be considered as the matrix version of the linear operator D TV , called D TV . Here, D TV = [ D r D c ] where D r ∈ R p×p (resp. D c ∈ R p×p ) is the matrix counterpart of D r (resp. D c ). With Neumann boundary conditions, D r and D c are bijective implying injectivity of D TV . Thus, D is a frame in view of Courant-Fisher theorem.

Remark 4 . 2 .

 42 Let G ≥ 1 and D ∈ R p×q satisfying Assumption (H.2). For any function g satisfying Assumptions (G.1)-(G.3), and any a

Lemma 4 . 1 .

 41 Let G ≥ 1 and D ∈ R p×q satisfying Assumption (H.2). Suppose that g satisfies Assumption (G.1) and∞ 0 z G+1 g(z)dz < ∞.(4.4)Then Assumptions (G.2) and (G.3) are in force.

  3 and Remark 4.4. Remark 4.3. Let G ≥ 1 and D ∈ R p×p be invertible. For any function g satisfying Assumptions (G.1)-(G.3), and for any a ∈]0, 1], one can choose K D a,g in (4.3) as

Remark 5 . 2 .Corollary 5 . 2 .

 5252 Suppose that G ≥ 2, and let γ = G + 2, ν = 2 and η = 1. We have γ/ν ≥ η + 1, and thus Lemma 2.1 yields ∞ 0 x G+1 (τ 2 + x 2 ) -2 dx is not definite. Consequently, condition (4.4) is not fulfilled with g defined in Example 4.3 when G ≥ 2.According to Remark 4.4, Remark 5.2 implies that Assumption (G.3) is not fulfilled for g in Example 4.3 when the group size G ≥ 2 and D invertible. Thus one cannot construct a group-analysis SOI from Theorem 5.1 to guarantee the quality of the corresponding estimator. Overcoming this limitation was yet another motivation behind the choice of g in Example 4.4, which turns out to work well under the group-analysis sparsity scenario. In a nutshell, an aggregate with g in Example 4.4 exhibits the group-analysis SOI defined in the following corollary with any G ≥ 1, any D ∈ R p×q satisfying Assumption (H.2) and any a ∈]0, 1]. Let X satisfying Assumption (H.3), G ≥ 1 and D satisfying Assumption (H.2) with κ > 0. Let Assumptions (P.1) and (P.2) be satisfied with some function v, (3.1) holds and β ≥ max (4 v ∞ , 2H/t 0 ). Take the dictionary (4.1) and the prior (4.2) with a ∈]0, 1], α ≥ 0 and g defined in Example 4.4. We get that g satisfies Assumptions (G.1)-(G.4). Then, let K D a,g as defined in (4.3), and assume that R > 3 K D a,g L. Then the group-analysis SOI (5.1) holds, with λ = c and h(x) = 1 + (x/τ) b .

Corollary 5 . 3 .

 53 Consider the same framework as Corollary 5.2 with D invertible.

  t)dt can then be naturally approximated via the Riemann sum δ T T/δ -1 k=0

Lemma 6 . 1 .

 61 The function w λ in (6.4) is bounded from below, increasing and continuously differentiable on ]0, +∞[. Fix a = b = 1. Then, for any γ ∈]0, τ 2 /c[, the following holds, (i) prox γw λ is single-valued on [0, +∞[ and is given by

1 S = 4 , 1 S = 4 , 1 S = 8 , 1 S = 8 , 4 Figure 1 :

 1414181841 Figure 1: (a)-(d): Phase diagrams of EWA for D = I p , the color bar ranges from dark (π S ,n = 0) to bright (π S ,n = 1). The blue line is the fitted phase transition curve. (a')-(d'): Examples of vectors θ ( j,S ,n) n recovered by EWA with n/p = 1/2, two sparsity levels S = 4 and S = 8 and two group sizes G θ 0 = 1 and G θ 0 = 4.

Figure 2 :

 2 Figure 2: (a): Original close-up of Shepp-Logan phantom image. (b): Image recovered by EWA with δ = 2 • 10 -8 and T = 10 4 . (c) Profiles of a row extracted from each image.

Section 5 9. 4 . 1 .

 541 Proof of Theorem 5.1

2 2

 2 p 0 (du). (9.16) Now, inserting (9.16) into Theorem 3.1 (with p = p 0 D ), we arrive at

9. 5 .

 5 Proofs of Section 6 9.5.1. Proof of Lemma 6.1 w λ is clearly increasing, bounded from below by w λ (0), and continuously differentiable (in fact even C ∞ ) on ]0, +∞[. (i) The expression of the proximal mapping follows from [40, Lemma 7.2] provided that it is single-valued. To show the latter statement, it is sufficient to prove that the function u : x ∈ [0, +∞[ → x + γw λ (x) = x + γc/(τ + x) + γα has a unique minimizer occurring at 0. One can see that u admits a local maximum at x = -√ γcτ [0, +∞[ and a local minimum at x = √ γcτ. Thus, the problem min x∈[0,+∞[ u(x) has a unique solution at 0 when x ≤ 0 equivalent to γ ≤ τ 2 /c.

• 2 2 2r

 2 + W λ • D for r = r/ν (recall ν from (H.2)). It then follows from [40, Lemma 5.3] that prox γW λ •D is Lipschitz continuous for any γ ∈]0, r].

E. 30 ) 9 . 3 .

 3093 { f (L(t))|L 0 = l 0 } -R p f (θ)µ γ (θ)dθ ≤ KH(l 0 )ρ t (9.28)for some K < +∞ and ρ ∈ [0, 1[ and any functionf satisfying | f | H ≤ 1. Denote A def = -∇V γ , ∇ + ∆ is the µ γ symmetric natural operator. By [51, Theorem 2.2], L(t) is H-uniformly geometrically ergodic if H is a Lyapunov function for the generator A , i.e. H is a C 2 function and there exist δ > 0, ϑ ≥ 0 and some R > 0 such that for all θ ∈ R p A H(θ) ≤ -δH(θ) + ϑI B R (θ)(9.29)whereB R def = θ ∈ R p : θ 2 ≤ R . We denote B c R def = R p \ B R .To prove geometric ergodicity, we use the following fact. Since F β is differentiable, convex and e -F β is integrable, it follows from [2, Lemma 2.2] that there exist η > 0 and r > 0 such that for all θ ∈ B c r ,∇F β (θ), θ ≥ η θ 2 . (9Lemma The drift condition (9.29) is satisfied with the Lyapunov function H(θ) = exp η/4 1 + θ 2 2 , δ = η 2 /16, ϑ = η/4H(R) η/4 + d + 2β -1 X

2 2 - 2 2+ 2 2+ 2 +

 2222 δ L k -δ/2(∇F β (L k ) -∇F β (0)), ∇( γ W λ • D )(L k ) + Cδ 2 + δp ≤ E L k Cδ 2 E {L k } 2 + Cδ 2 + δp ≤ E L k Cδ 2 E {L k } 2 + Cδ 2 + δp ≤ E L k 2 Cδ 2 T + Cδ 2 + δp ≤ CδT 2 + CδT + pT = O(T 2 ), ∀k[0, T/δ ].(9.32)For T > 0 and R > 0, denote LT t)I B R ( Lδ (t))dt, L T,R def = 1/T T 0 L(t)I B R (L(t))dt and θ R γ def = B Rθµ γ (θ)dθ. The triangle and Jensen inequalities yield

γ 2 ≤

 2 .g., R = 1/T 3 , these terms converge to 0 as T → ∞. Using (9.32) and arguing as in[23, Step 1, Proposition 2], we haveE Lδ,T,γ -LT δ 2 ≤ Cδ(1 + δT 2 ) → 0 as δ → 0.29Using Girsanov formula and Pinsker inequality, as in[23, Step 2, Proposition 2], the distribution of Lδ (t)t∈[0,T ]converges to that of {L(t)} t∈[0,T ] in total variation as δ → 0. In turn,E LT,R δ -L T,R2→ 0 as δ → 0. Using Lemma 9.3 with f (θ) = θ i I B R (θ), and arguing as in[23, Step 4, Proposition 2], we have E L T,Rθ R CT -1/2 → 0 as T → ∞.

  Observe that prox γW λ • D in (6.9) is no separable in general. Owing to Assumption (H.2), one can show quite immediately thatprox γW λ • D = D • prox D D γW λ +ι Im(D ) • D , where ι Im( D ) (u) = 0 if u ∈ Im(D )and +∞ otherwise, and D D is indeed definite positive on Im(D ). This means that, unless D is orthogonal, W λ • D does not have an easy-to-compute expression, but rather necessitates to solve an optimization subproblem. Thus, from a computational perspective, for a general frame D, Semi-FBLMC is more

	.9)
	Remark 6.1.

efficient. For the case where D is invertible, one can operate a simple change of variable u = D θ and replace F β with F β • D -1 , and W λ • D with W λ . In this case, FBLMC can be applied efficiently. 6.3. Convergence guarantees By virtue of Lemma 6.1(v)-(vi), standard results (see, e.g., [51, Theorem 2.1]) show that for any initial point L(0) such that E L(0) 2 2 < ∞, (6.5) has a unique solution which is strongly Markovian, E L(t) 2 2 < ∞ for all t > 0, and L admits an (unique) invariant measure having the density µ γ

  1)-(G.4) with any G ≥ 1. Applying Theorem 5.1, we conclude the proof.

	Proof of Corollary 5.3 Since g satisfies Assumptions (G.2), (G.3) and D is invertible, by Remark 4.3 and Lemma 2.1,
	we get

  1 νL(α + c/τ). (vi) Combine (i), [40, Lemma 3.2 and Lemma 6.1], (iii) and that M γ is positive definite. 9.5.2. Proof of Theorem 6.1 Uniform geometric ergodicity We start by showing that L(t) in (6.5) is H-uniformly geometrically ergodic, where H : R p → [1, +∞[ is a measurable function. Let | f | H

def = sup θ∈R p | f (θ)| H(θ) .

Recall that H-uniform geometric ergodicity requires that for all l 0

PAC stands for Probably Approximately Correct. A PAC-bound is probably correct as it not a deterministic guarantee allowing a small probability that the estimator does not behave well. It is also approximately correct as it is tolerant to an inexact performance of the estimator.

The author however applied it to problems whereln(µ) = F + G. But the gradient of the Moreau envelope of a sum, which amounts to computing the proximity operator ofln(µ) does not have an easily implementable expression even if those of F and G do.

We have made no effort to sharpen the constants and the rates, and for instance their dependence on the dimension. This is beyond the scope of this paper.

Acknowledgements. This work was supported by Conseil Régional de Basse-Normandie and partly by Institut Universitaire de France.