Yu Golubev

D Ostrovski

Concentration Inequalities for the Exponential Weighting Method *

. The estimators withing this family are aggregated based on the exponential weighting method and the performance of the aggregated estimate is measured by the excess risk controlling deviation of the square losses from the oracle risk. Based on natural statistical properties of ordered smoothers, we propose a novel method for obtaining concentration inequalities for the exponential weighting method.

Introduction and main results

This paper deals with recovering an unknown vector µ ∈ R n from the noisy observations Y i = µ i + σξ i , i = 1, 2, . . . , n,

where ξ is a standard white Gaussian noise, i.e., ξ i are i.i.d. Gaussian random variables with Eξ i = 0 and Eξ 2 i = 1. For the sake of simplicity it is assumed that the noise level σ > 0 is known. It is also assumed that n is large. Let us emphasize that all results in this paper can be extended to the case n = ∞ provided that µ ∈ 2 .

Notice also that in spite of the obvious simplicity of this statistical model, it can cover a very wide class of statistical problems ranging from regression estimation to inverse problems (see, e.g., [START_REF] Kneip | Ordered linear smoothers[END_REF], [START_REF] Golubev | On universal oracle inequalities related to high dimensional linear models[END_REF]).

In this paper, the quality of an estimate μ(Y ) is measured by the mean square risk ρ(μ, µ) = E µ μ(Y ) -µ 2 ,

where E µ stands for the expectation with respect to the measure P µ generated by the observations (1.1) under given µ, and

• is the standard Euclidean norm in R n .
In what follows it is assumed that we have at our disposal a family of linear estimates μh

k (Y ) = h k Y k , h ∈ H, (1.2)
where H is a given set of vectors in R n . This set consists of vectors with specific properties that will be described later on.

Our goal is to construct a good estimate of µ with the help of the family of linear estimates μh (Y ), h ∈ H. To explain what estimates might be viewed as good, let us notice that for given h, the risk of μh (Y ) is computed as follows :

R(µ, h) def = E µ μh (Y ) -µ 2 = (1 -h) • µ 2 + σ 2 h 2 , (1.3)
where here and in what follows • means the component-wise multiplication of vectors in R n , i.e, z = x • y means that z k = x k y k , k = 1, . . . , n.

Of coarse, we want to find a method of combining estimates {μ h (Y ), h ∈ H} that provides an estimate with the minimal risk. If the family of estimates is reach enough, then we would be satisfied with estimators whose risks are close to the value

r(µ, H) = min h∈H (1 -h) • µ 2 + σ 2 h 2 .
This risk is often called oracle risk since it coincides with the risk of the following psudo-estimate

μ * (Y ) = h * (µ) • Y, h * (µ) = arg min h∈H (1 -h) • µ 2 + σ 2 h 2 . (1.4)
Obviously, this object is not a real statistical estimate since it depends on the unknown vector µ.

A rather natural approach to constructing new estimates with the help of the available ones is to compute their convex combination μw (Y ) = h∈H w h (Y )μ h (Y ), (1.5) where w h (Y ) are some nonnegative numbers called weights such that h∈H w h (Y ) = 1.

To simplify unessential technical details, here and below H is assumed to be finite.

The main issue in this method is evidently related to the selection of weights w h (Y ). In this paper, we focus on the so-called exponential weights defined as follows :

w h (Y ) = π h Z(Y ) exp - R(Y, h) 2βσ 2 , (1.6)
where

Z(Y ) def = h∈H π h exp - R(Y, h) 2βσ 2 , (1.7) R(Y, h) def = μh (Y ) 2 -2 Y, μh (Y ) + 2σ 2 h∈H h k , (1.8)
is the empirical counterpart of R(µ, h) and β, π h are positive numbers. It is clear that without loss of generality one may assume that h∈H π h = 1.

To explain heuristically how these weights can be obtained, let us simplify the initial statistical problem assuming that we have at our disposal the auxiliary observations

Y k = µ k + σξ k , k = 1, 2 . . . , n, (1.9)
where ξ is independent of ξ in (1.1) a standard white Gaussian noise. These observations will be used to select weights

w h (Y ) in the estimate μw (Y, Y ) = h∈H w h (Y )μ h (Y ), (1.10)
where w h (Y ) ≥ 0 and h∈H w h (Y ) = 1.

Apparently the approach to constructing new estimates with the help of the auxiliary sample were firstly studied and developed by A. Nemirovsky and, independently, by A. Catoni (see [START_REF] Nemirovski | Topics in non-parametric statistics[END_REF][START_REF] Catoni | Statistical learning theory and stochastic optimization[END_REF]). Subsequently, these methods have been adapted for a wide class of statistical models ( see, for example, [START_REF] Yang | Combining different procedures for adaptive regression[END_REF], [START_REF] Rigollet | Linear and convex aggregation of density estimators[END_REF], [START_REF] Lecué | Simultaneous adaptation to the margin and to complexity in classification[END_REF], [START_REF] Rigolet | Sparse estimation by exponential weighting[END_REF]).

Let us assume that Y is frozen. The idea in choosing w h (•) in (1.10) is related to the so-called roughness penalty approach yielding the weights wh (Y ) = arg min

w h 1 2 2 Y - h∈H w h μh (Y ) 2 + βK(w h , π h ) , (1.11)
where

K(w h , π h ) = h∈H w h log w h π h
is the Kulback-Leibler divergence between probability distributions w h and the a priori weights π h . Since there is no explicit formula for wh (Y ), to simplify their computing, let us replace in (1.11) the distance between the data and the estimate

Y - h∈H w h μh (Y )
2 by the following upper bound obtained by the Jensen inequality

Y - h∈H w h μh (Y ) 2 ≤ h∈H w h Y -μh (Y ) 2 .
Thus we arrive at

w h (Y ) = arg min w h 1 2σ 2 h∈H w h Y -μh (Y ) 2 + βK(w h , π h ) .
Now, these weights can be easily computed with a simple algebra

w h (Y ) = π h Z(Y ) exp - μh (Y ) 2 -2 Y , μh (Y ) 2σ 2 β , (1.12)
where

Z(Y ) = h∈H π h exp - μh (Y ) 2 -2 Y , μh (Y ) 2σ 2 β .
Obviously, the described above two samples estimation procedure results in loss of a significant statistical information contained in the data. Therefore computing estimates and their aggregation should be done, of course, with the help of the single sample. The main difficulty in implementing this idea is related to the inner product Y , μh (Y ) in (1.12). It is clear that we need to replace this inner product by something that should to be close to it, but depending only on Y . In order to understand how to implement this idea, let us look at the probabilistic structure of Y , μh (Y ) . We have

Y , μh (Y ) = n k=1 h k (µ k + σξ k )(µ k + σξ k ) = n k=1 h k µ 2 k +σ n k=1 h k µ k ξ k + σ n k=1 h k µ k ξ k + σ 2 n k=1 h k ξ k ξ k .
(1.13) Note that the second line in this equation contains random variables with zero mean for given h. If we want to make use of the single sample, we need first to look at

Y, μh (Y ) = n k=1 h k (µ k + σξ k ) 2 = n k=1 h k µ 2 k + 2σ 2 n k=1 h k +2σ n k=1 h k µ k ξ k + σ 2 n k=1 h k (ξ 2 k -1). ( 1
E Y , μh (Y ) = E Y, μh (Y ) -2σ 2 n k=1 h k .
These intuitive considerations combined with (1.12) result in the exponential weights defined by (1.6)- (1.8).

Notice that in recent years the exponential weighting method has been extensively studied and enough good upper bounds for its performance have been obtained for several statistical models [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF][START_REF] Dalayan | Sharp oracle inequalities for aggregation of affine estimators[END_REF][START_REF] Alquier | Pac-bayesian theorems for sparse regression estimation with exponential weights[END_REF][START_REF] Arias-Castro | K Variable Selection with Exponential and 0 -Penalization[END_REF].

In this paper, we study this method in the case when H is a set of ordered multipliers defined as follows : Definition 1.1. H is a set of ordered multipliers if the following properties hold:

1. h i ∈ [0, 1], i = 1, . . . , n for all h ∈ H, 2. h i+1 ≤ h i , i = 1, . . . , n for all h ∈ H, 3. if for some integer k and some h, g ∈ H, h k < g k , then h i ≤ g i for all i = 1, . . . , n.

Property 3 means that vectors in H may be naturally ordered, since for any h, g ∈ H there are only two possibilities h i ≤ g i or h i ≥ g i for all i = 1, . . . , n. Therefore the estimators defined by (1.2), where H is a set of ordered multipliers, are often called ordered smoothers [START_REF] Kneip | Ordered linear smoothers[END_REF]. Let us emphasize that ordered smoothers are very common in statistics (see for instance [START_REF] Kneip | Ordered linear smoothers[END_REF], [START_REF] Chernousova | Ordered Smoothers With Exponential Weighting[END_REF]).

We start our study of the exponential weighting method with the case β → 0. It is easy to see that in this case the limiting estimator has the following form :

μ• (Y ) = h • (Y ) • Y, where h • (Y ) = arg min h∈H R(Y, h).
(1.15)

This standard method of estimate selection with the help of the unbiased risk estimation has long been well known in statistics and goes back to [START_REF] Akaike | Information theory and an extension of the maximum likelihood principle[END_REF][START_REF] Mallows | Some comments on C p[END_REF]. A classical fact about its performance is given by the following theorem [START_REF] Kneip | Ordered linear smoothers[END_REF].

Theorem 1.1. Let H be a set of ordered multipliers. Then for any µ ∈ R n

E µ μ• (Y ) -µ 2 2 ≤ r(µ, H) + Kσ 2 1 + r(µ, H) σ 2 , (1.16)
where here and in what follows K denotes generic constants.

Equation (1.16) is often called oracle inequality since it controls the risk of μ• (Y ) via the oracle risk r(µ, H).

When β > 0, statistical analysis of the exponential weighting is more involved (see, e.g., [START_REF] Chernousova | Ordered Smoothers With Exponential Weighting[END_REF]). In order to get an oracle inequality similar to (1.16), some additional assumptions about π h and H are required.

First, we will make use of a priory weights defined by the following condition Condition 1.1.

π h def = 1 -exp - h + 1 -h 1 β , (1.17)
where h + = min{g ∈ H : g > h} and π h max = 1, h max is the maximal multiplier in H.

In this definition and below • 1 stands for 1 -norm in R n , i.e.,

h 1 = n i=1 |h i |.
Along with the above condition, we will need also the following one.

Condition 1.2. There exists a constant

K • > 0 such that h 2 -g 2 ≥ K • h 1 -g 1 (1.18)
for all h ≥ g from H.

An upper bound for mean square risk of the exponential weighting method is given by the following theorem [START_REF] Chernousova | Ordered Smoothers With Exponential Weighting[END_REF].

Theorem 1.2. Assume that H is a set of ordered multipliers, β ≥ 4, and Conditions 1.1, 1.2 hold. Then, uniformly in µ ∈ R n , E µ μβ (Y ) -µ 2 ≤ r β (µ, H) + σ 2 C(K • , β), (1.19)
where

r β (µ, H) def = r(µ, H) + 2βσ 2 log 2 + r(µ, H) σ 2 . (1.20)
Here and in what follows C(K • , β) denotes strictly positive and bounded constants depending on K • and β.

At first glance, it may seem comparing the remainder terms in Equations (1.16) and (1. [START_REF] Yang | Combining different procedures for adaptive regression[END_REF])-(1.20), that the exponential weighting should be significantly better than the classical methods of model selection. In fact, this is not exactly so, and to find out what is really going on, we need to understand first of all what methods were used in proving Theorems 1.1 and 1.2.

Theorem 1.1 results from the following concentration inequality (see Kneip [START_REF] Kneip | Ordered linear smoothers[END_REF]) :

P µ μ• (Y ) -µ ≥ r(µ, H) + Kσ + x ≤ exp - x 2 Kσ 2 , x ≥ 0. (1.21)
At the same time, Theorem 1.2 represents a fact of an entirely different class. At the very core in the proof of (1.19) in [START_REF] Chernousova | Ordered Smoothers With Exponential Weighting[END_REF] is an original method proposed in [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF] which is based on Stein's formula [START_REF] Stein | Estimation of the mean of a multivariate normal distribution[END_REF] for the unbiased risk estimate. In essence, this means that in contrast to (1.21) Equation (1.19) cannot control the deviations of the square losses μβ (Y ) -µ 2 from the oracle risk.

Note also that Theorem 1.2 holds true only for β ≥ 4. Unfortunately, it is impossible to say whether this condition is crucial for practical applications of the exponential weighting, or it is a purely mathematical constraint resulting from the method of the proof.

In order to understand more profoundly statistical properties of the exponential weighting method, we focus in this paper on the so-called excess risk ∆ β (Y, H)

def = μβ (Y ) -µ 2 -r β (µ, H) + ,
where [x] + = max(0, x) and r β (µ, H) is defined by (1.20). The next theorem, controlling moments of ∆ β (Y, H), is the main result of this paper.

Theorem 1.3. Assume that H is a set of ordered multipliers and Conditions 1.1, 1.2 hold. Then, uniformly in µ ∈ R n and in m ≥ 1,

E 1/m µ ∆ β (Y, H) m + ≤ Kσ mr β (µ, H) + Kσ 2 m + C(K • , β)σ 2 .
(

Notice that compared to majority of results related to the exponential weighting, see e.g. [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF][START_REF] Rigolet | Sparse estimation by exponential weighting[END_REF][START_REF] Dalayan | Sharp oracle inequalities for aggregation of affine estimators[END_REF][START_REF] Chernousova | Ordered Smoothers With Exponential Weighting[END_REF], the main advantage of this theorem is that it holds for any β > 0 and provides an upper bound for the excess risk that does depend neither n no the cardinality of H. It justifies, in particular, the use of the exponential weighting with β = 1 that demonstrates a good performance in practice as shown in the next section.

Combining Theorem 1.3 with the Markov inequality, one obtains the following fact similar to (1.21).

Theorem 1.4. Assume that H is a set of ordered multipliers and Conditions 1.1, 1.2 hold. Then, uniformly in µ ∈ R n and x ≥ σ,

P µ μβ (Y ) -µ ≥ r β (µ, H) + x ≤ exp - x 2 Kσ 2 + C(K • , β) .

Simulation study

To illustrate numerically Theorem 1.3, the following experiment has been carried out. Its goal was to find out how the exponential weighting with β = {0, 1, 2, 4} works in regression estimation with the help of the cubic smoothing splines. Recall that these splines are usually used in recovering an unknown smooth function f (x), x ∈ [0, 1] from the noisy observations

Y i = f (X i ) + ξ i , i = 1, . . . , n, (2.1)
where it is assumed that the design points X i belong to (0, 1) and ξ is a standard white Gaussian noise. Smoothing spline of order 2m -1 is defined by

fα (x, X, Y ) = arg min f 1 n n i=1 [Y i -f (X i )] 2 + α 1 0 [f (m) (u)] 2 du , (2.2)
where f (m) (•) denotes the derivative of order m and α ∈ A ⊂ R + is a smoothing parameter. In what follows the set of possible smoothing parameters A is assumed to be finite.

In order to show that the regression estimation with the help of the smoothing splines is equivalent to the sequence space model described by (1.1) and (1.2), consider the Demmler-Reinsch [START_REF] Demmler | Oscillation matrices with spline smoothing[END_REF] basis

ψ k (x), x ∈ [0, 1], k = 1, . . . , n having double orthogonality ψ k , ψ l n = δ kl , 1 0 ψ (m) k (x)ψ (m) l (x) dx = δ kl λ k , k, l = 1, . . . , n, (2.3)
where δ kl = 1 if k = l, and δ kl = 0 otherwise. Here and below u, v n stands for the inner product

u, v n = 1 n n i=1 u(X i )v(X i ).
Let us assume for definiteness that the eigenvalues λ k are sorted in ascending order

λ 1 ≤ • • • ≤ λ n .
With this basis, representing the underlying regression function as follows:

f (x) = n k=1 ψ k (x)µ k , (2.4)
we get from (2.1) and (2.3)

Y k = 1 n n i=1 Y i ψ k (X i ) = µ k + ε √ n ξ k , (2.5)
where µ k = f, ψ k n and ξ is a standard white Gaussian noise. So, substituting (2.4) in (2.2) and using (2.3), we arrive at

fα (x, X, Y ) = n k=1 μk ψ k (x),
where μ = arg min

µ n k=1 [Y k -µ k ] 2 + α n k=1 λ k µ 2 k .
It can be seen easily that μk = Y k 1 + αλ k and thus the spline regression model (2.1)-(2.2) is equivalent to the sequence space model defined by (1.1) and (1.2) with σ = ε/ √ n and

H = h : h k = 1 1 + αλ k , k = 1, . . . , n, α ∈ A . (2.6)
In order to simulate cubic splines, the following family of ordered multipliers was used

H = h : h k = 1 1 + [απ(k -1)] 4 , k = 1, . . . , n, α ∈ A n ,
where

A n = α : α = (1 + ε) s , s = 0, 1, . . . , log(n)/ log(1 + ε)
with ε = 0.01. The motivation of this family is due to the well-known asymptotic formula λ k (πk) 4 , k 1 (see [START_REF] Demmler | Oscillation matrices with spline smoothing[END_REF]). The simulations were organized as follows. For given A ∈ [0, 50], 40000 replications of the observations

Y k = µ k (A) + ξ k , k = 1, . . . , 300
were generated. Here µ(A) ∈ R 300 is a Gaussian vector with independent components and

Eµ k (A) = 0, Eµ 2 k (A) = A exp - k 2 2Σ 2 .
Next, the mean oracle risk

r(A, H) = E min h∈H (1 -h) • µ(A) 2 + h 2
and the absolute excess risk

∆β 1 (A) = E µ(A) -μβ (Y ) 2 -r(A, H) + ,
were computed with the help of the Monte-Carlo method. Finally, the data r(A, H), ∆β 1 (A), A ∈ [0, 50], β = 0, 1, 2, 4 are plotted on Figure 1 to illustrate graphically statistical properties of the exponential weighting method.

Looking at this picture we see that there is no universal β minimizing the excess risk uniformly in µ, but it seems that a reasonable choice would be β ≈ 1. Note also that the exponential weighting with β ∈ [START_REF] Alquier | Pac-bayesian theorems for sparse regression estimation with exponential weights[END_REF][START_REF] Catoni | Statistical learning theory and stochastic optimization[END_REF] works better compared to the classical unbiased risk estimation (β = 0) when r(µ, H)/σ 2 is not large. 3 Proofs

Auxiliary facts

In what follows it is assumed that H be a set of ordered multipliers. Basic probabilistic facts related to the ordered multipliers are described is the following lemma.

Lemma 3.1. Let ξ i be i.i.d. N (0, 1). Then for any α ∈ (0, 1/4)

P max h∈H ± n i=1 h i (ξ 2 i -1) -α n i=1 h 2 i ≥ x α ≤ exp - x K , (3.1)
P max h∈H n i=1 (1 -h i )ξ i µ i -α n i=1 (1 -h i ) 2 µ 2 i ≥ x α ≤ exp - x K , (3.2)
where K is a generic constant.

The proof of this lemma follows immediately from Lemma 2 in [START_REF] Golubev | On universal oracle inequalities related to high dimensional linear models[END_REF].

Lemma 3.2. I) If for some R, r ≥ 0 R ≤ r + min α∈(0, ] x 2 α + αR , then √ R ≤ √ r + |x| max 2, 2 .
II) If for some R, r ≥ 0 R ≤ r + min α∈(0, ) x 2 α + αr , then √ R ≤ √ r + |x| max 1, 2 .
Proof. It is easy to check with a simple algebra that min α∈(0, ]

x 2 α + αR = 2|x| √ R1 x ≤ √ R + ( -1 x 2 + R)1 x > √ R 2|x| √ R1 x ≤ √ R + 2 -1 x 2 1 x > √ R . Therefore if x > √ R, then R ≤ r + 2 -1 x 2 ,
and so,

√ R ≤ √ r + √ 2|x|/ √ . Next, when x ≤ √ R we have R ≤ r + 2|x| √ R, or, equivalently, ( √ R -|x|) 2 ≤ r + x 2 . Hence √ R ≤ √ r + 2|x|.
The proof of the second part of the lemma is quite similar.

Lemma 3.3. Let ζ be a nonnegative random variable with Eζ m < ∞. Then for any A > 0

E 1/m log m (A + ζ) ≤ log(A + E 1/m ζ m ).
Proof. Consider the following function

F (z) = max x>-A log(A + x) -zx , z ∈ R + .
It can be checked with a simple algebra that

F (z) = -log(z) + Az -1.
and therefore for any z ≥ 0 log(A + x) ≤ zx + F (z).

So, we have

E 1/m log m (A + ζ) ≤ zE 1/m ζ m + F (z)
and minimizing the right-hand side at this equation in z ≥ 0, we complete the proof.

Lemma 3.4. Let ĥ = max h ∈ H : R(Y, h) -R(Y, h • ) ≤ 2β σ 2 h 2 -h • 2 + -1 , (3.3)
where ∈ (0, 1/(2β)) and h • (Y ) is defined by (1.15). Then for any integer m ≥ 1

1 -2β E 1/(2m) µ ĥ 2m ≤ r(µ, H) σ 2 + K (m + β). (3.4)
Proof. By the definition of R(Y, h) (see (1.8))) and (3.3), we get ĥ = max h :

(1 -h) • µ 2 + σ 2 (1 -2β ) h 2 + 2σ n i=1 (1 -h i ) 2 µ i ξ i + σ 2 n i=1 (h 2 i -2h i )(ξ 2 i -1) ≤ (1 -ĥ) • µ 2 + σ 2 (1 -2β ) ĥ 2 + 2σ n i=1 (1 -ĥi ) 2 µ i ξ i + σ 2 n i=1 ( ĥ2 i -2 ĥi )(ξ 2 i -1) + 2βσ 2 .
Let us fix some γ, γ ∈ (0, 1/4). Then we can rewrite the above equation as follows:

ĥ = max h ∈ H : σ 2 (1 -2β -γ) h 2 + 2σ n i=1 (1 -h i ) 2 µ i ξ i + (1 -h) • µ 2 + σ 2 n i=1 (h 2 i -2h i )(ξ 2 i -1) + γσ 2 h 2 ≤ (1 + γ ) (1 -ĥ) • µ 2 + σ 2 (1 -2β + γ ) ĥ 2 + 2σ n i=1 (1 -ĥi ) 2 µ i ξ i -γ (1 -ĥ) • µ 2 + σ 2 n i=1 ( ĥ2 i -2 ĥi )(ξ 2 i -1) -γ σ 2 ĥ 2 + 2βσ 2 . Hence ĥ ≤ h def = max h ∈ H : σ 2 (1 -2β -γ) h 2 + min g∈H 2σ n i=1 (1 -g i ) 2 µ i ξ i + (1 -g) • µ 2 + min g∈H σ 2 n i=1 (g 2 i -2g i )(ξ 2 i -1) + γσ 2 g 2 ≤ (1 + γ ) (1 -h • ) • µ 2 + (1 -2β + γ )σ 2 h • 2 + max g∈H 2σ n i=1 (1 -g i ) 2 µ i ξ i -γ (1 -g) • µ 2 + max g∈H σ 2 n i=1 (g 2 i -2g i )(ξ 2 i -1) -γ σ 2 g 2 + 2βσ 2 .
(3.5)

Next, bounding max and min in (3.5) with the help of Lemma 3.1, we obtain that for any integer m > 0 and any γ, γ ∈ (0, 1/4) the following inequality holds

(1 -2β -γ)σ 2 E 1/m µ h 2m ≤ (1 + γ )E 1/m µ R m (µ, h • ) + K(m!) 1/m σ 2 γ + K[(m!) 1/m + β]σ 2 γ ,
where R(•, •) is defined by (1.3).

Therefore, applying Lemma 3.2, we obtain from the above equation

1 -2β σE 1/(2m) µ h 2m ≤ E 1/(2m) µ R m (µ, h • ) + K m + βσ. (3.6)
In order to control the expectation at the right-hand side in (3.6), note that the following inequality

n i=1 [1 -h • i ] 2 Y 2 i + 2σ 2 n i=1 h • i ≤ n i=1 [1 -g i ] 2 Y 2 i + 2σ 2 n i=1 g i
holds for any given g ∈ H. Hence, we have

R(µ, h • )+2σ n i=1 (1 -h • i ) 2 µ i ξ i + σ 2 n i=1 (h • i ) 2 -2h • i )(ξ 2 i -1) ≤R(µ, g) + 2σ n i=1 (1 -g i ) 2 µ i ξ i + σ 2 n i=1 (g 2 i -2g i )(ξ 2 i -1).
So, for any γ, γ ∈ (0, 1/4], we get with this equation and Lemma 3.1

E 1/m µ R m (µ, h • ) ≤ R(µ, g) + γE 1/m µ R m (µ, h • ) + γ R(µ, g) + 2σE 1/m max h∈H - n i=1 (1 -h i ) 2 µ i ξ i - γ 2σ n i=1 (1 -h i ) 2 µ 2 i m + σ 2 E 1/m max h∈H n i=1 (2h i -h 2 i )(ξ 2 i -1) -γ n i=1 h 2 i m + 2σE 1/m max h∈H - n i=1 (1 -h i ) 2 µ i ξ i - γ 2σ n i=1 (1 -h i ) 2 µ 2 i m + σ 2 E 1/m max h∈H n i=1 (2h i -h 2 i )(ξ 2 i -1) -γ n i=1 h 2 i m ≤ R(µ, g) + Kσ 2 m γ + γE 1/m µ R m (µ, h • ) + Kσ 2 m γ + γ R(µ, g).
Next, minimizing the right-hand side in g ∈ H and applying Lemma 3.2, we obtain

E 1/(2m) µ R m (µ, h • ) ≤ r(µ, H) + Kσ √ m
and substituting this inequality in (3.6), we complete the proof.

The next lemma collects some useful facts about a priory weights defined by (1.17). Let

D h = g ∈ H : h 1 ≤ g 1 ≤ h 1 + 1 .
(3.7)

Lemma 3.5. Under Condition 1.1, for any h ∈ H, the following assertions hold:

g≥h π g exp - g 1 β = exp - h 1 β , (3.8
)

g≤h π g ≤ 1 + h 2 K • β , (3.9
)

g∈D h π g ≤ 1 + 1 β , (3.10
)

g∈D h π g ≥ 1 2β exp - 1 β . (3.11)
The proof of this lemma can be found in [START_REF] Chernousova | Ordered Smoothers With Exponential Weighting[END_REF].

The following technical lemma, whose proof is given in [START_REF] Chernousova | Ordered Smoothers With Exponential Weighting[END_REF], is a cornerstone in the proof of Theorem 1.3. Lemma 3.6. Suppose {q h ≤ 1, h ∈ H} is a nonnegative sequence such that for all h ≥ h

q h ≤ exp -γ h 1 -h 1 -1 , γ > 0. Let W h = π h q h g∈H π g q g -1
and G be a subset in H.

Then h∈H W h log π h W h ≤ log h≤ h π h + exp R(G, γ) ,
where

E(G, γ) = log 2 γβe + h∈G π h + h∈G π h q h -1 8 γβe + h∈G π h . (3.12)
3.2 Proof of Theorem 1.3

The proof is based essentially on the approach proposed in [START_REF] Golubev | Exponential weighting and oracle inequalities for projection estimates[END_REF][START_REF] Chernousova | Ordered Smoothers With Exponential Weighting[END_REF]. We begin with a basic inequality for the excess risk. By Jensen's inequality we have

μβ (Y ) -µ 2 = h∈H w h (Y ) h • Y -µ 2 ≤ h∈H w h (Y ) h • Y -µ 2 .
(

This is why in what follows we will deal with

h∈H w h (Y ) h • Y -µ 2 .
In order to bound this value from above, let us express h • Y -µ 2 in terms of the unbiased risk estimate of h • Y . With a simple algebra we get

h • Y -µ 2 = h • Y -Y + σξ 2 = (1 -h) • Y 2 + σ 2 ξ 2 -2σ ξ, (1 -h) • (µ + σξ) = (1 -h) • Y 2 + 2σ 2 n i=1 h i -σ 2 ξ 2 -2σ ξ, (1 -h) • µ + 2σ 2 n i=1 h i (ξ 2 i -1). (3.14)
Denote for brevity

R(Y, h) def = (1 -h) • Y 2 + 2σ 2 n i=1 h i -σ 2 ξ 2 , r(Y, H) def = min h∈H R(Y, h), h * = arg min h∈H E R(Y, h). (3.15)
With these notations and (3.14), we obtain

h∈H w h (Y ) h • Y -µ 2 = r(Y, H) + h∈H w h (Y ) R(Y, h) -r(Y, H) + h∈H w h (Y ) 2σ 2 n i=1 h i (ξ 2 i -1) -2σ n i=1 (1 -h i )ξ i µ i . (3.16)
It is also obvious that

r(Y, H) ≤ R(Y, h * ) = (1 -h * ) • µ 2 + σ 2 h * 2 + σ 2 n i=1 (h * i 2 -2h * i )(ξ 2 i -1) + 2σ n i=1 (1 -h * i ) 2 µ i ξ i = r(µ, H) + σ 2 n i=1 (h * i 2 -2h * i )(ξ 2 i -1) + 2σ n i=1 (1 -h * i ) 2 µ i ξ i . (3.17)
Recalling the definition of the weights w h (Y ), we rewrite them as follows:

w h (Y ) = π h exp - R(Y, h) -r(Y, H) 2βσ 2 g∈H π g exp - R(Y, g) -r(Y, H) 2βσ 2
and we obtain

R(Y, h)-r(Y, H) = 2βσ 2 h∈H w h (Y ) log π h w h (Y ) -2βσ 2 log h∈H π h exp - R(Y, h) -r(Y, H) 2βσ 2 . (3.18)
Therefore, substituting Equations (3.17) and (3.18) in (3.16), we arrive at the basic inequality

h∈H w h (Y ) h • Y -µ 2 ≤ r(µ, H) + δ 1 (µ) + δ 2 (µ) + 2βσ 2 δ 3 (µ) + 2βσ 2 δ 4 (µ), (3.19)
where

δ 1 (µ) =σ 2 n i=1 (h * i 2 -2h * i )(ξ 2 i -1) + 2σ n i=1 (1 -h * i ) 2 µ i ξ i , δ 2 (µ) = h∈H w h (Y ) 2σ 2 n i=1 h i (ξ 2 i -1) -2σ n i=1 (1 -h i )ξ i µ i , δ 3 (µ) = -log h∈H π h exp - R(Y, h) -r(Y, H) 2βσ 2 , δ 4 (µ) = h∈H w h (Y ) log π h w h (Y )
.

The first term at the right-hand side of (3.19) can be easily controlled. Indeed, denote for brevity one can check easily that

Eδ 2 1 (µ) = 2σ 4 n i=1 (h * i 2 -2h * i ) 2 + 4σ 2 n i=1 (1 -h * i ) 4 µ 2 i ≤ σ 2 8σ 2 n i=1 h * i 2 + 4 n i=1 (1 -h * i ) 2 µ 2 i ≤ 8σ 2 r(µ, H).
(3.20)

Next, since by the definition of ordered multipliers (h * i ) 2 -2h * i ≤ 0, we get by a simple algebra that for any λ > 0 and so, for any x, λ > 0

log[exp(λδ 1 (µ))] = n i=1 λσ 2 (2h * i -h * 2 i ) - 1 2 log 1 + 2λσ 2 (2h * 2 i -h * i ) 2λ 2 σ 2 (1 -h * i ) 4 µ 2 i 1 + 2λσ 2 (2h * 2 i -h * i ) ≤ n i=1 λ 2 σ 4 (2h * 2 i -h * i ) 2 + 2λ 2 σ 2 (1 -h * i ) 4 µ 2 i ≤ 4λ 2 σ 2 r(µ, H).
x m ≤ exp(λmx -m) λ m .
Thus, by (3.22)

E δ 1 (µ) σ 8r(µ, H) m + ≤ exp λ 2 m 2 2 -m log(λ) -m .
Minimizing the right-hand side at this equation in λ > 0, we arrive at

E δ 1 (µ) σ 8r(µ, H) m + ≤ exp - m 2 + m 2 log(m)
and thus

E 1/m [δ 1 (µ)] m + ≤ σ Kmr(µ, H). (3.23)
The last three terms at the right-hand side in (3.19) are bounded from above with the help of more delicate arguments. Since h∈H w h (Y ) = 1, we obtain that for any γ ∈ (0, 1/4)

δ 2 (µ) ≤2σ h∈H w h (Y ) σ n i=1 h i (ξ 2 i -1) -σγ n i=1 h 2 i + + 2σ h∈H w h (Y ) n i=1 (1 -h i )ξ i µ i - γ σ n i=1 (1 -h i ) 2 µ 2 i + + 2γ h∈H w h (Y )R(µ, h) ≤2σ 2 max h∈H n i=1 h i (ξ 2 i -1) -γ n i=1 h 2 i + + 2σ max h∈H n i=1 (1 -h i )ξ i µ i - γ σ n i=1 (1 -h i ) 2 µ 2 i + + 2γ h∈H w h (Y )R(µ, h).
Hence, using Lemma 3.1, we obtain from this equation

E 1/m µ [δ 2 (µ)] m + ≤ 2γE 1/m µ h∈H w h (Y )R(µ, h) m + Kmσ 2 γ (3.24) Our next step is to relate h∈H w h (Y )R(µ, h) and h∈H w h (Y ) µ -h • Y 2 . Since h∈H w h (Y ) µ -h • Y 2 = h∈H w h (Y )R(µ, h) + 2σ h∈H w h (Y ) (1 -h i ) -(1 -h i ) 2 µ i ξ i + σ 2 h∈H w h (Y )h 2 i (ξ 2 i -1), we get h∈H w h (Y )R(µ, h) ≤ h∈H w h (Y ) µ -h • Y 2 + 2σ h∈H w h (Y ) (1 -h i ) 2 -(1 -h i ) µ i ξ i + σ 2 h∈H w h (Y )h 2 i (1 -ξ 2 i ).
With this equation, using the same arguments as in proving (3.24), we arrive at

(1 + 2γ)E 1/m h∈H w h (Y )R(µ, h) m ≤E 1/m h∈H w h (Y ) µ -h • Y 2 m + Kmσ 2 γ .
Therefore by Lemma 3.2 we have

E 1/m µ h∈H w h (Y )R(µ, h) m 1/2 ≤ E 1/m µ h∈H w h (Y ) µ -h • Y 2 m 1/2 + σ √ Km.
Next, substituting this equation in (3.24) and using that (x+y) 2 ≤ 2x 2 +2y 2 we arrive at the following upper bound

E 1/m µ [δ 2 (µ)] m + ≤ 4γE 1/m µ h∈H w h (Y ) µ -h • Y 2 m + Kmσ 2 γ . ( 3

.25)

Let us now consider δ 3 (µ). We have by (3.8)

log h∈H π h exp - R(Y, h) -r(Y, H) 2βσ 2 ≥ log h≥h • (Y ) π h exp - R(Y, h) -r(Y, H) 2βσ 2 = log h≥h • (Y ) π h exp - (1 -h) • Y 2 -[1 -h • (Y )] • Y 2 2βσ 2 - 1 β n i=1 [h i -h • i (Y )] ≥ log h≥h • (Y ) π h exp - 1 β n i=1 [h i -h • i (Y )] ≥ 0 (3.26)
and hence

δ 3 (µ) ≤ 0. (3.27)
The final step in the proof is to bound from above the last term at the right-hand side of Equation (3.19), namely δ 4 (µ). We will do this the help of Lemma 3.6. First of all note that for all h > ĥ , where ĥ is defined by

(3.3), we have R(Y, h) -r(Y, H) ≥ 2β σ 2 h 2 -h • (Y ) 2 + 2βσ 2 .
In order to make use of Lemma 3.6, let us set

q h = exp - R(Y, h) -r(Y, H) 2βσ 2 = exp - R(Y, h) -R(Y, h • ) 2βσ 2 , G = h ∈ H : h • 1 ≤ h 1 < h • 1 + 1 , h = ĥ .
By Condition (1.18) we have that for all h ≥ ĥ

q h ≤ exp - 2σ 2 βK • ( h 1 -h • 1 ) 2βσ 2 -1 = exp -K • ( h 1 -h • 1 ) -1 .
(3.28)

Note also that similar to (3.26) it can be checked easily that

q h ≥ exp - 1 β , h ∈ G
and thus it follows immediately from the definition of G and (3.11) that

h∈G π h q h ≥ exp - 1 β h∈G π h ≥ 1 2β exp - 2 β
and hence

E(G, K • ) ≤ C(K • , β) 2 + 1 2 log C(K • , β) . (3.29) Denote for brevity ρ(µ, H) = 2 + r(µ, H σ 2 .
With Equations (3.29), (3.9), and Lemmas 3.6, 3.3 we obtain

E 1/m µ δ 4 (µ) -log[ρ(µ, H)] m + ≤ 2E 1/m µ log m ĥ + exp[E(G, K • )/2] βK • ρ(µ, H) ≤ 2 log 1 √ 1 -2β r(µ, H) σ 2 + K √ m + β √ 1 -2β + C(K • , β) exp C(K • , β) 2 -log K • βρ(µ, H) .
(3.30) Our next step is to minimize the right-hand side at this equation in ∈ (0, 1/(2β)). Note that that for any ≤ 1/(3β) we have

1 √ 1 -2β r(µ, H) σ 2 + K √ m + β √ 1 -2β ≤ r(µ, H) σ 2 + K m + β + K r(µ, H) σ 2 + K m + β . (3.31)
Consider the following function

Ψ(x) = min ∈[0,1/(3β)] x + C(K • , β) exp C(K • , β) 2 .
It is clear that Ψ(0) is bounded from above. It is also easy to check with = 4C(K • , β)/log(x) that for any x ≥ C(K • , β)

Ψ(x) ≤ 4C(K • , β)x log(x) + x log(x) 2 ≤ C(K • , β)x log(x) .
Therefore, combining this equation with (3.30), (3.31), we arrive at

Figure 1 :

 1 Figure 1: Absolute excess risks ∆β 1 (•) (left panel Σ = 5; right panel Σ = 50).

( 3 .

 3 21) Therefore, combining (3.20) and (3.21), we obtain that for any λ > 0 from convexity of the function exp(•) that z ≤ exp(z -1),

E 1/m µ δ 4

 4 (µ) -log[ρ(µ, H)] m + ≤ log C(K • , β) m + β . (3.32) Next, substituting (3.23), (3.25), (3.27), and (3.32) in (3.19) we arrive atE 1/m µ h∈H w h (Y ) µ -h • Y 2 m ≤ r(µ, H) + Kσ mr(µ, H) + 4γE 1/m µ h∈H w h (Y ) µ -h • Y 2 m + Kmσ 2 γ + 2βσ 2 log C(K • , β)(m + β) + 2βσ 2 log[ρ(µ, H)].

( 3 . 2 E 1 /≤+≤+≤ 1 / 2 ≤ 1 / 2 +

 3211212 33)Using this equation and minimizing the right-hand at(3.33) in γ ∈ (0, 1/4), we obtain with the help of Lemma 3.(2m)µ h∈H w h (Y ) µ -h • Y 2 m ≤ K √ mσ + r β (µ, H) + 2βσ 2 log C(K • , β)(m + β) 3.33), we have that for any γ ∈ (0, 1/4)E 1/m µ h∈H w h (Y ) µ -h • Y 2 -r β (µ, H) m + Kσ mr β (µ, H) + 2βσ 2 log C(K • , β)(m + β) + 4γE 1/m µ h∈H w h (Y ) µ -h • Y 2 m + Kmσ 2 γ .and substituting in this equation (3.34), we arrive atE 1/m µ h∈H w h (Y ) µ -h • Y 2 -r β (µ, H) m Kσ mr β (µ, H) + 2βσ 2 log C(K • , β)(m + β) + γ 8r β (µ, H) + 2βσ 2 log C(K • , β)(m + β) + 16Kσ2 m is to minimize the right-hand side at this equation in γ ∈ (0, 1/4). Notice thatγ • = σ √ Km 8r β (µ, H) + 2βσ 2 log C(K • , β)(m + β) + 16Kσ 2 m -1/2is the unconstrained minimizer of the right-hand side. However it is clear that γ • ≤ 1/4 and hence we have by(3.35)E 1/m µ h∈H w h (Y ) µ -h • Y 2 -r β (µ, H) m Kσ mr β (µ, H) + 2βσ 2 log C(K • , β)(m + β) + Kσ √ m r β (µ, H) + 2βσ 2 log C(K • , β)(m + β) + Kσ 2 m Kσ √ m r β (µ, H) + Kσ 2 m Kβσ 2 log C(K • , β)(m + β) ≤ Kσ mr β (µ, H) + Kσ 2 m + Kβσ2 log C(K • , β) + Kβσ 2 log(m + β) ≤ Kσ mr β (µ, H) + Kσ 2 m + Kβ[1 + log(β)]σ 2 log C(K • , β) .

( 3 .

 3 36)In deriving this inequality it was used that2 log(x) -x ≤ 2 [1 + log( )], x > 0, > 0.(3.37)Finally, Equation (3.36) with (3.13), we finish the proof.

3. 3 1 E 1 /x 2 K 2 2 ≤ exp - x 2

 3112222 Proof of Theorem 1.4 It follows from (1.22) that for any m≥ m μβ (Y ) -µ -r β (µ, H) Kσ 2 [C(K • , β) + m] r β (µ, H) .Hence, by the Markov inequality we obtainP µ μβ (Y ) -µ ≥ r β (µ, H) + x ≤ exp -m log(x) + m log E 1/m μβ (Y ) -µ -r β (µ, H) m + ≤ exp -m log(x) + m log K √ mσ + Kσ 2 [C(K • , β) + m] r β (µ, 2 σ 2 exp(1) ,we get with (3.37) that for any x ≥ σP µ μβ (Y ) -µ ≥ r β (µ, H) + x ≤ exp -x Kσ 2 + log C(K • , β) + Kσ r β (µ, H) Kσ 2 + C(K • , β) .

* This work is partially supported by Laboratory for Structural Methods of Data Analysis in Predictive Modeling, MIPT, RF Government grant, ag. 11.G34.31.0073; and RFBR research projects 13-01-12447 and 13-07-12111.


E Chernousova

Yu Golubev
email: golubev.yuri@gmail.com.

Spectral Cut-off Regularizations for Ill-posed Linear Models

Keywords: ill-posed linear model, spectral cut-off regularization, datadriven cut-off frequency, oracle inequality, minimax risk. 2010 Mathematics Subject Classification: Primary 62C99; secondary 62C10, 62C20, 62J05

This paper deals with recovering an unknown vector β from the noisy data Y = Xβ + σξ, where X is a known n × p -matrix with n ≥ p and ξ is a standard white Gaussian noise. In order to estimate β, a spectral cut-off estimate β m(Y ) with a data-driven cut-off frequency m(Y ) is used. The cut-off frequency is selected as a minimizer of the unbiased risk estimate of the mean square prediction error, i.e. m = arg min m Y -X βm (Y ) 2 + 2σ 2 m . Assuming that β belongs to an ellipsoid W, we derive uppers bounds for the maximal risk sup β∈W E β m(Y ) -β 2 and show that β m(Y ) is a rate optimal minimax estimator over W.

Introduction and main result

This paper deals with recovering an unknown vector β ∈ R p from the noisy observations

Y = Xβ + σξ, (1.1)
where X is a known n × p -matrix with n ≥ p and ξ is a standard white Gaussian noise. Let us emphasize that all results below can be extended to the case p = ∞ provided that Xβ ∈ ℓ 2 . For the sake of simplicity it is assumed also that the noise level σ is known.

The standard way to estimate β based on the observations (1.1) is to make use of the maximum likelihood estimate

β(Y ) = arg min β Y -Xβ 2 = (X ⊤ X) -1 X ⊤ Y,
where here and in what follows • stands for the standard Euclidean norm.

It is well known that the mean square risk of this method is given by

E β(Y ) -β 2 = σ 2 p i=1 1 λ(k) ,
where λ(k) are the eigenvalues of X ⊤ X, i.e.

X ⊤ Xe k = λ(k)e k , k = 1, . . . , p, (1.2)
and so, it is clear that it may be very large when X ⊤ X is ill-posed. In this case a regularization of β(Y ) is required. Nowadays, statisticians have at their disposal a very vast family of regularization methods (see e.g. [START_REF] Engl | Regularization of Inverse Problems[END_REF]). In this paper, we focus on the so-called spectral cut-off regularizations which are computed with the help of the Singular Value Decomposition.

Let e k ∈ R p be the eigenvectors of X ⊤ X (see (1.2)). Using this basis, we represent the vector of interest as follows :

where integer m is often called cut-off frequency.

In view of (1.4), statistical analysis of this method for a fixed cut-off frequency is rather simple. If the performance of βm (Y ) is measured by the prediction mean square error then one can check easily that

r(β, m) def = E X[ βm (Y ) -β] 2 = p k=m+1 λ(k)β 2 (k) + σ 2 m. (1.6)
On the other hand, for the standard mean square risk one obtains

R(β, m) def = E βm (Y ) -β 2 = p k=m+1 β 2 (k) + σ 2 m k=1 1 λ(k) . (1.7)
Thus, we see that in the both cases the risks depend on m and to get a good estimate of β we have to select properly the cut-off frequency. Since β is unknown, this selection must be data-driven. The standard approaches to the data-driven choice of m are often based on the principle of the unbiased risk estimation and go back to [START_REF] Akaike | Information theory and an extension of the maximum likelihood principle[END_REF][START_REF] Mallows | Some comments on C p[END_REF]. There are two basic methods.

The first one is related to the unbiased estimate of E X[ βm (Y ) -β] 2 and yields the following cut-off frequency :

m = arg min m Y -X βm (Y ) 2 + 2σ 2 m . (1.8)
Notice that there is a vast literature devoted to statistical analysis of this method. The most precise fact about its performance was firstly obtained in [START_REF] Kneip | Ordered linear smoothers[END_REF].

Theorem 1 Uniformly in β ∈ R p E X( β m -β) 2 ≤ r • (β) + Kσ 2 r • (β) σ 2 , (1.9)
where K is a generic constant and r • (β) is the so-called oracle risk defined by r

• (β) def = min m r(β, m).
Unfortunately, with the help of this theorem we cannot obtain good upper bounds for the risk of β m(Y ) measured by E β m(Y )-β 2 . Of course, one can bound this risk as follows

E β m(Y ) -β 2 ≤ λ -1 (p)E X( β m(Y ) -β) 2 ≈ λ -1 (p)r • (β), (1.10)
but it can be seen easily that the right-hand side in this equation may be very far from the oracle risk (see (1.7)) defined by

R • (β) def = min m R(β, m) (1.11)
when X ⊤ X is ill-posed. Therefore, in order to get good upper bounds for E β m(Y ) -β 2 , the cut-off frequency m is selected as a minimizer of the unbiased risk estimate

of E βm (Y ) -β 2 , i.e., m = arg min m β(Y ) -βm (Y ) 2 + 2σ 2 m k=1 1 λ(k)
.

(1.12)

The risk of this method is controlled by following oracle inequality (see, [START_REF] Cavalier | Risk hull method and regularization by projections of ill-posed inverse problems[END_REF]).

Theorem 2 Suppose λ(k) ≥ λ(1)k -α for some α > 0. Then uniformly in β ∈ R p E β m(Y ) -β 2 ≤ R • (β) + C(α)σ 2 R • (β) σ 2 (2α+1)/(2α+2) , (1.13)
where C(α) ≥ Kα is a constant depending on α and the oracle risk R • (β) is defined by (1.11).

Comparing (1.9) and (1.13) one can see the principal difference between these oracle inequalities: the remainder term in (1.9) doesn't depend on α, whereas the one in (1.13) goes to infinity as α → ∞. This observation confirms a suspicion that β m(Y ) may perform poorly when X ⊤ X is ill-posed. This effect can be seen in simulation, see e.g. [START_REF] Cavalier | Risk hull method and regularization by projections of ill-posed inverse problems[END_REF] or Section 2 below.

In order to improve the performance of the unbiased risk estimation method in the ill-posed case, penalties heavier than 2σ 2 m k=1 λ -1 (k) should be used, see for details [START_REF] Cavalier | Risk hull method and regularization by projections of ill-posed inverse problems[END_REF][START_REF] Golubev | On universal oracle inequalities related to high dimensional linear models[END_REF]. It is shown in these papers that the estimate β m+ (Y ) with

m+ = arg min m β(Y ) -βm (Y ) 2 + 2σ 2 m k=1 1 λ(k) + P en λ (m) (1.14)
works better than β m(Y ) provided that the additional penalty P en λ (m) is properly chosen. Unfortunately, computing good penalties P en λ (m) in (1.14) is a time-consuming numerical problem (see, e.g. [START_REF] Cavalier | Risk hull method and regularization by projections of ill-posed inverse problems[END_REF][START_REF] Golubev | On universal oracle inequalities related to high dimensional linear models[END_REF]). This is why the main goal in this paper is to improve significantly Inequality (1.10) and to show that the standard spectral cut-off estimator

β m(Y ) = m k=1 Y, e * k λ(k) e k , m = arg min m Y -X βm (Y ) 2 + 2σ
W = β : p k=1 w(k)β 2 (k) ≤ 1 .
Here w(k), k = 1, . . . , p is a positive increasing sequence.

In order to obtain such upper bounds, we will need some basic notions and facts related to the minimax estimation theory over ellipsoids. Nowadays, the minimax approach to ill-posed linear models is well-developed and literature on this topic is very vast, see e.g. [START_REF] Donoho | Nonlinear solutions of linear inverse problems by wavelet-vaguelette decomposition[END_REF][START_REF] Donoho | Neoclassical minimax problems, thresholding and adaptive function estimation[END_REF][START_REF] Ermakov | Minimax estimation of the solution of an illposed convolution type problem[END_REF][START_REF] Mair | Statistical estimation in Hilbert scale[END_REF][START_REF] Koo | Optimal rates of convergence for nonparametric statistical inverse problems[END_REF], where additional references can be found.

With W we associate the following risks :

r(W, m) def = max β∈W r(β, m), R(W, m) def = max β∈W R(β, m), and r • (W) def = min m r(W, m), R • (W) def = min m R(W, m). Proposition 1 For any integer m r(W, m) = λ(m + 1) w(m + 1) + σ 2 m, R(W, m) = 1 w(m + 1) + σ 2 m k=1 1 λ(k) .
Proof. It follows immediately from (1.6) and (1.7).

Let us emphasize that in general case the spectral cut-off estimates are not sharp minimax over ellipsoids, they are only rate optimal. Minimax estimates optimal up to a constant have been obtained in the seminal paper [START_REF] Pinsker | Optimal filtering of square integrable signals in Gaussian white noise[END_REF].

To simplify some technical details, it is assumed in what follows that

λ(k) = Lk -α and w(k) = k γ /P,
where α, γ, L, P are some positive constants. In order to make computations more transparent, let us assume also that the noise level is small, i.e. σ → 0.

Then we obtain easily

r(W, m) = LP (m + 1) -γ-α + σ 2 m
and denoting for brevity

q def = 1 γ + α + 1 ,
we get a simple algebra the minimax cut-off frequency

m r (P, σ 2 ) = arg min m P L(m + 1) -γ-α + σ 2 m = (1 + o(1)) (γ + α)P L σ 2 q .
With this cut-off frequency we obtain the following formula for the minimax risk

r • (W) = (1 + o(1)) q 1 -q 1-q σ 2 P L σ 2 q . (1.16)
The same technique is used in computing R • (W). We have

R(W, m) = P (m + 1) -γ + σ 2 L m k=1 k α
and therefore we get the following minimax cut-off frequency

m R (P, σ 2 ) = arg min m P (m + 1) -γ + σ 2 L m k=1 k α = (1 + o(1)) γP L σ 2 q
and the minimax risk

R • (W) =(1 + o(1)) γ -qγ q(α + 1) P σ 2 P L qγ =(1 + o(1)) γ -qγ q(α + 1) σ 2 L σ 2 P L -q(1+α)
.

(1.17)

Comparing (1.17) and (1.16), we see that the minimax risks r • (W) and R • (W) are related as follows :

R • (W) = (1 + o(1)) γ -qγ (γ + α) (1+α)(q-1) q(α + 1) σ 2 L r • (W) σ 2 α+1
.

(1.18)

The second important remark is that the optimal minimax cut-off frequencies have the same order, i.e., lim

σ→0 m R (P, σ 2 ) m r (P, σ 2 ) = γ γ + α 1/(1+α+γ)
.

(1.19)

In practice, we cannot make use of the minimax cut-off frequencies because they strongly depend on the ellipsoid parameters that are hardly known in practice. However, Equations (1.18) and (1.19) may be viewed as a heuristic motivation of β m(Y ). These equations show that there are strong links between the minimax risks r(W) and R(W) as well between the cut-off frequencies m R (P, σ 2 ) and m r (P, σ 2 ). So, there is a hope the cut-off frequency m(Y ), which is nearly optimal (see Theorem 1) when the risk is measured by the prediction error

E X[ β m(Y ) -β] 2 , is also good for the risk E β m(Y ) -β 2 .
The next theorem provides a mathematical justification of this conjecture.

Theorem 3 For any β ∈ W E β m(Y ) -β 2 ≤ 2σ 2 L P L σ 2 α/(α+γ) r • (β) σ 2 + K 2γ/(α+γ) + σ 2 (α + 1)L r • (β) σ 2 + Kα 2α+2 + Kσ 2 L √ 2α + 1 r • (β) σ 2 + Kα 2α+1 + Kσ 2 , (1.20)
where K is a generic constant.

With the help of this theorem one can easily check that β m(Y ) is a rate optimal asymptotically minimax estimator over W.

Theorem 4 As σ → 0 sup β∈W E β m(Y ) -β 2 ≤ (1 + o(1))C(α, γ)R • (W), (1.21)
where

C(α, γ) = γ γ/(1+γ+α) 1 + α + γ 2(α + 1) + (γ + α) (1+α)(α+γ)/(1+α+γ) .
Proof. It is clear that r • (β) ≤ r • (W) for any β ∈ W, and thus by (1.16) and (1.17) we obtain

2σ 2 L P L σ 2 α/(α+γ) r • (β) σ 2 + K 2γ/(α+γ) = (1 + o(1)) 2σ 2 L P L σ 2 α/(α+γ) r • (β) σ 2 γ/(α+γ) = 2(1 + o(1))P σ 2 LP = 2q(α + 1)γ qγ (1 + o(1))R • (W).
Similarly, by (1.18) we have

σ 2 (α + 1)L r • (β) σ 2 + Kα 2α+2 = 1 + o(1) 1 + α σ 2 L r • (β) σ 2 1+α ≤ qγ qγ (γ + α) (1+α)(1-q) (1 + o(1))R • (W)
and

σ 2 (α + 1)L r • (β) σ 2 + Kα 2α+1 ≤ o(1)R • (W).
Substituting the above equations in (1.20), we complete the proof.

Remark. The spectral cut-off method requires the SVD which may be time-consuming from a numerical viewpoint when p is large. This is why for very large linear models simpler regularization techniques are usually used. A classic example of such a regularization is the famous Phillips-Tikhonov method [START_REF] Tikhonov | Solution of Ill-posed Problems[END_REF] known in statistics as ridge regression. In order to select the regularization parameter in this method, the Generalized Cross Validation is usually used (see, e.g. [START_REF] Golub | Generalized crossvalidation as a method to choosing a good ridge parameter[END_REF][START_REF] Golub | Generalized crossvalidation for large scale problems[END_REF]). Unfortunately, in spite of a very long history of this idea, all known facts about its performance are related to upper bounds for the prediction error (see e.g. Theorem 1). In view of Theorem 3, there is a hope that bounds similar to (1.20) and (1.21) may be obtained for this method as well, and so, the use of the GCV will be justified.

Simulation study

To illustrate numerically how does the spectral cut-off method with different data-driven cut-off frequencies work, the following simulations have been carried out. For given A ∈ [0, 20], 20000 replications of the observations

Z(k) = λ(k)µ(k) + ξ ′ (k), k = 1, . . . , 300
were generated. Here ξ ′ is a standard white Gaussian noise, µ ∈ R 300 is a Gaussian vector (independent from ξ ′ ) with independent components and

Eµ(k) = 0, Eµ 2 (k) = A exp - k 2 2Σ 2 .
To simulate the spectral cut-off estimator, the unknown vector µ is estimated as follows:

μm (k) = Z(k) λ(k) 1{k ≤ m}.
With the help of the Monte-Carlo method the following risks were computed

• the oracle risk r or (A) = min m E µ -μm 2 ,
• the risk related to the cut-off frequency minimizing the unbiased estimate of E µ -μm 2 , i.e.,

r inv (A) = E µ -μ m 2 , where m = arg min m λ -1/2 • Z -μm 2 + 2 m k=1 1 λ(k) , (2.1)
• the risk related to the cut-off frequency minimizing the unbiased estimate of E √ λ(µ -μm ) 2 , i.e.,

r dir (A) = E µ -μ m 2 , where m = arg min m Z - √ λ • μm 2 + 2m , (2.2)
• the risk related to the cut-off frequency minimizing the over-penalized unbiased risk estimator, i.e. ,

r ovp (A) = E µ -μ m+ 2 ,
where

m+ = arg min m Z - √ λ • μm 2 + 2m + 2m log(m) . (2.3)
A motivation of this method and, in particular, the use of the additional penalty 2m log(m) can be found in [START_REF] Cavalier | Risk hull method and regularization by projections of ill-posed inverse problems[END_REF].

Finally, the following oracle efficiencies were computed

θ inv (A) = r or (A) r inv (A) , θ dir (A) = r or (A) r dir (A) , θ ovp (A) = r or (A) r ovp (A) .
In order to compare graphically the above methods of selection of the cut-off frequency, the data A, θ inv (A), θ dir (A), θ ovp (A) were plotted on Figures 1, 2, 3, and4.

Looking at these pictures we see that when X ⊤ X becomes ill-posed, the spectral cut-off method with the cut-off frequency from (2.2) works significantly better than the cut-off frequency minimizing the unbiased estimate of E µ -μm 2 . Let us emphasize also that in the case of ill-posed X ⊤ X the data-driven cut-off frequency from (2.3) improves the risk of the spectral cut-off estimate.

Proofs

Auxiliary facts

The proof of Theorem 3 is based on two cornerstone facts. The first one is used to control the bias of β m(Y ) and has the following form: Lemma 1 For any β ∈ W and any integer m = 1, . . . , p Proof. Let us prove (3.1) for p = ∞. Since β ∈ W we have that for any

p k=m β 2 (k) ≤ 2L -γ/(α+γ) P α/(α+γ) p k=m β 2 (k)λ(k) γ/(α+γ) . ( 3
x ≥ 1 ∞ k=m β 2 (k) ≤ mx k=m β 2 (k) + k>mx β 2 (k) ≤ 1 λ(mx) mx k=m λ(k)β 2 (k) + 1 w(xm) ≤ x α λ(m) ∞ k=m λ(k)β 2 (k) + 1 w(m)x γ = 1 w(m) x α w(m) λ(m) ∞ k=m λ(k)β 2 (k) + x -γ . (3.2)
Our next step is to minimize the right-hand side in this equation in x ≥ 1. Denote for brevity

S(m) = w(m) λ(m) ∞ k=m λ(k)β 2 (k)
and note that S(m) ≤ 1 since β ∈ W. So, we can choose x as a root of the equation x α S(m) = x -γ , or, equivalently, x = S(m) -1/(α+γ) . Hence, substituting this in (3.2), we obtain

∞ k=m β 2 (k) ≤ 2 w(m) w(m) λ(m) ∞ k=m λ(k)β 2 (k) γ/(α+γ) =2[w(m)] -α/(α+γ) [λ(m)] -γ/(α+γ) ∞ k=m λ(k)β 2 (k) γ/(α+γ) =2P α/(α+γ) L -γ/(α+γ) ∞ k=m λ(k)β 2 (k) γ/(α+γ)
.

Equation (3.1) follows now from this inequality if we assume that

β(k) = 0, k > p.
The following lemma plays an essential role in controlling the variance of β m(Y ).

Lemma 2 For any β and any µ > 0

E exp µ r(β, m) ≤ exp µ r • (β) + Kσµ + Kσ 2 µ 2 , (3.3)
where K is a generic constant.

Proof. It follows immediately from the definition of m (see (1.4), (1.5), and (1.8)) that for any given positive integer m the following inequality holds

- m i=1 Z 2 (i) + 2σ 2 m ≤ - m i=1 Z 2 (i) + 2σ 2 m
or, equivalently,

r(β, m) ≤ r(β, m) + σ 2 m i=1 [ξ ′2 (i) -1] + 2σ p i= m+1 λ(i)β(i)ξ ′ (i) -σ 2 m i=1 [ξ ′2 (i) -1] -2σ p i=m+1 λ(i)β(i)ξ ′ (i).
We can easily derive from this inequality the following one

(1 -ǫ)r(β, m) ≤ r(β, m) + σ 2 sup k≥1 k i=1 [ξ ′2 (i) -1] -ǫk + 2σ sup k≥1 p i=k+1 λ(i)β(i)ξ ′ (i) - ǫ 2σ p i=k+1 λ(i)β 2 (i) + σ 2 m i=1 [1 -ξ ′2 (i)] -2σ p i=m+1 λ(i)β(i)ξ ′ (i), (3.4)
where ǫ ∈ (0, 1). The last line in this equation can be easily controlled since for any µ > 0

E exp µ σ 2 m i=1 1 -ξ ′2 (i) -2σ p i=m+1 λ(i)β(i)ξ ′ (i) ≤ exp σ 2 mµ - m 2 log 1 + 2µσ 2 + 2σ 2 µ 2 p i=m+1 λ(i)β 2 (i) ≤ exp σ 4 mµ 2 + 2σ 2 µ 2 p i=m+1 λ(i)β 2 (i) ≤ exp 2σ 2 µ 2 r(β, m) . (3.5)
To bound from above sup k≥1 at the right-hand side in (3.4), we make use of the following inequalities (see Lemma 2 in [START_REF] Golubev | Exponential weighting and oracle inequalities for projection estimates[END_REF])

P max k≥1 k i=1 ξ ′2 (i) -1 -U (ǫ)k ≥ x ≤ exp(-ǫx), P max k≥1 p i=k ξ ′ (i)β(i) λ(i) - ǫ 2 p i=k β 2 (i)λ(i) ≥ x ≤ exp(-ǫx),
where ǫ ∈ R + and

U (ǫ) = - ǫ + log[1 -2ǫ] + /2 ǫ .
With the help of these inequalities and a simple algebra one can check easily that for any ǫ ∈ (0, 1/4)

E exp λσ 2 sup k≥1 k i=1 [ξ ′2 (i) -1] -ǫk ≤ 1 [1 -Kλσ 2 /ǫ] + (3.6)
and

E exp 2λσ sup k≥1 p i=k+1 λ(i)β(i)ξ ′ (i) - ǫ 2σ p i=k+1 λ(i)β 2 (i) ≤ 1 [1 -Kλσ 2 /ǫ] + . (3.7)
It follows from the Hölder inequality that for any random variables ζ 1 , ζ 2 , ζ 3 with bounded exponential moments the following inequality holds

E exp λζ 1 + λζ 2 + λζ 3 ≤ 3 i=1 E exp 3λζ i 1/3 .
With this equation and (3.4)-(3.7) we obtain

E exp (1 -ǫ)λr(β, m) ≤ exp λr(β, m) + Kσ 2 λ 2 r(β, m) - 2 3 log 1 - Kλσ 2 ǫ + or equivalently, E exp λr(β, m) ≤ exp λr(β, m) + ǫλr(β, m) 1 -ǫ + Kσ 2 λ 2 r(β, m) (1 -ǫ) 2 - 2 3 log 1 - Kλσ 2 (1 -ǫ)ǫ + .
Minimizing the right-hand side at this equation in m and recalling that

r • (β) def = min m r(β, m),
we arrive at the following inequality

E exp λ[r(β, m) -r • (β)] ≤ exp Kσ 2 λ 2 r • (β) + Kǫλr • (β) - 2 3 log 1 - Kλσ 2 ǫ + , (3.8)
that holds for any ǫ ∈ (0, 1/2). The next step is to minimize the right-hand side at (3.8) in ǫ ∈ (0, 1/2). Suppose

λ ≤ 1 4Kσ 2 .
(3.9)

Since r • (β)/σ 2 ≥ 1, therefore

ǫ = Kλσ 2 + σ 4 r • (β) ≤ 1 2 ,
and substituting this ǫ in (3.8), we obtain

E exp λ[r(β, m) -r • (β)] ≤ exp Kσ 2 λ 2 r • (β) + Kσλ r • (β) + 2 3 log 1 + 4Kλσ r • (β) ≤ exp Kσ 2 λ 2 r • (β) + Kσλ r • (β) .
(3.10)

Let ζ be a nonnegative random variable. Then for any µ, λ > 0 we have

E exp µ ζ = ∞ 0 exp(µ √ x -λx + λx) dP{ζ ≤ x} ≤ exp max x≥0 µ √ x -λx ∞ 0 exp(λx) dP{ζ ≤ x} ≤ exp µ 2 4λ + log E exp(λζ) .
and combining this equation with (3.9) and (3.10), we have

E exp µ r(β, m) ≤ exp min λ≤1/(4Kσ 2 ) µ 2 4λ + λr + Kσ r • (β)λ + Kσ 2 r • (β)λ 2 .
(3.11)

Note that choosing λ = 1 4Kσ 2 , we obtain from (3.11) that for any

µ ∈ R + E exp µ r(β, m) ≤ exp Kµ 2 σ 2 + K r • (β) σ 2 + K r • (β) σ 2 . (3.12)
On the other hand, if µ is small, namely,

µ ≤ 2 r • (β) 4Kσ 2 , then with λ = µ 2 r • (β)
we get from (3.11)

E exp µ r(β, m) ≤ exp µ r • (β) + Kµ 2 σ 2 + Kµσ . (3.13)
Finally notice that if

µ ≥ 2 r • (β) 4Kσ 2 , then µ r • (β) + Kµ 2 σ 2 + Kµσ ≥ Kµ 2 σ 2 + K r • (β) σ 2 + K r • (β) σ 2
and thus in view of (3.12), Equation (3.13) holds for any µ ∈ R + .

The following simple technical lemma is used together with Lemma 2 to control the variance of β m(Y ).

Lemma 3 Let η ≥ 1 be a random variable such that

E exp(λ √ η) ≤ exp(Aλ + Bλ 2 )
for any λ > 0. Then for any p ≥ 1

Eη p ≤ [A + (2p -1)B] 2p . (3.14)
Proof. Let us consider the following function f (x) = log 2p (x), x ≥ 0 and compute its second derivative

f ′′ (x) = - 2p log 2p-1 (x) x 2 + 2p(2p -1) log 2p-2 (x) x 2 = 2p log 2p-2 (x)[(2p -1) -log(x)] x 2 . So, f (x) is concave on [exp(2p -1), ∞). Therefore noticing that η p = f [exp(λ √ η)] λ 2p
and using Jensen's inequality, we obtain that for any λ ≥ 2p -1

Eη p ≤ f [E exp(λ √ η)] λ 2p = [Aλ + Bλ 2 ] 2p λ 2p = [A + Bλ] 2p .
Minimizing the right-hand side at this equation in λ ≥ 2p -1 we finish the proof.

In order to control the variance of β m(Y ), we will need to bound from above Eζ( m), where

ζ(m) def = m i=1 ξ ′2 (i) -1 λ(i) .
We will do this with the help of following fact. Let η(t) be a separable zero mean random process on R + . Denote for brevity

∆ η (u, v) = η(u) -η(v).
Lemma 4 Let σ 2 (u), u ∈ R + , be a continuous strictly increasing function with σ 2 (0) = 0. Then for any µ > 0,

log E exp µ max 0<u≤t ∆ η (u, t) σ(t) ≤ log(2) √ 2 √ 2 -1 + + max 0<u<v≤t max |z|≤ √ 2/( √ 2-1) log E exp zµ ∆ η (u, v) ∆σ (v, u) , (3.15)
where ∆σ

(v, u) = |σ 2 (v) -σ 2 (u)|.
Proof. It is similar to the one of Dudley's entropy bound (see, e.g., [START_REF] Van Der Vaart | Weak convergence and empirical processes[END_REF]) and can be found in [START_REF] Golubev | On universal oracle inequalities related to high dimensional linear models[END_REF].

Lemma 5 Assume that λ(k)/λ(1) ≥ 2 -k . Then for any ǫ ∈ (0, 1/4)

E max m≥1 ζ(m) -ǫEζ 2 (m) + ≤ K ǫ .
Proof. Let us apply Lemma 4 for the random process η

(k) = ζ(k + k 0 ) - ζ(k 0 ), k = 0, . . . , t and σ 2 (k) = E[ζ(k + k 0 ) -ζ(k 0 )] 2 .
Using a Taylor expansion, we obtain from (3.15)

E exp µσ -1 (t) max k=1,...,t η(k) ≤ K exp(Kµ 2 ), µ ∈ (0, µ • ), (3.16)
where

µ • = √ 2 -1 4 √ 2 . Let m ǫ k = max m : Eζ 2 (m) ≤ 2 k ǫ 2 .
Then we have

E max m≥1 ζ(m) -ǫEζ 2 (m) + ≤ E p k=0 max m∈[m ǫ k ,m ǫ k+1 ) ζ(m) -ǫEζ 2 (m) + ≤ p k=0 E ζ(m ǫ k ) + max m∈[m ǫ k ,m ǫ k+1 ) [ζ(m) -ζ(m ǫ k )] -ǫEζ 2 (m ǫ k ) + .
(3.17)

Denote for brevity

Σ ǫ k = 2 m ǫ k+1 k=m ǫ k +1 1 λ 2 (k) 1/2 . and notice that 2 k/2 4ǫ ≤ Σ ǫ k ≤ 2 k/2 ǫ . ( 3

.18)

To control the right-hand side at (3.17), we make use of the following inequality :

E ζ -x} + ≤ µ -1 exp(-µx)E exp(µζ), x > 0,
+ = Σ ǫ k E ζ(m ǫ k ) Σ ǫ k + 1 Σ ǫ k max m∈[m ǫ k ,m ǫ k+1 ) [ζ(m) -ζ(m ǫ k )] -ǫ Eζ 2 (m ǫ k ) Σ ǫ k + ≤ KΣ ǫ k exp -Kǫ Eζ 2 (m ǫ k ) Σ ǫ k ≤ Kǫ -1 2 k/2 exp -K2 k/2 .
Substituting this equation in (3.17), we finish the proof.

Proof of Theorem 3

From the definition of the spectral cut-off estimator (see (1.5)) and (1.4) it follows that

E β m(Y ) -β 2 = E p k= m+1 β 2 (k) + σ 2 E m k=1 1 λ(k) +σ 2 E m i=1 ξ ′2 (i) -1 λ(i) . (3.19)
To bound from above the first term at the right-hand side of this equation, we combine Lemma 1, Jensen's inequality, and Lemma 2. So, we have Similar arguments may be are used in controlling the last term at the right-hand side of (3.19). By Lemma 5 we get that for any ǫ ∈ (0, 1/4)

E m i=1 ξ 2 (i) -1 λ(i) ≤ ǫE m i=1 1 λ 2 (i) + K ǫ .
Next, minimizing the right-hand side at this equation in ǫ ∈ (0, 1/4), we obtain with the help of (3.21), Lemma 3, and a simple algebra

E m i=1 ξ 2 (i) -1 λ(i) ≤ K E m k=1 1 λ 2 (k) 1/2 + K ≤ K L √ 2α + 1 E m2α+1 1/2 + K ≤ K L √ 2α + 1 r • (β) σ 2 + Kα 2α+1 + K.
Finally, substituting this equation, (3.20), and (3.22) in (3.19), we finish the proof.

= 1 ,

 1 )e k , where •, • is the ordinary inner product in R p , and thus we obtain Y = p k=1 β(k)Xe k + σξ. . . . , p is an orthonormal system in R n and projecting (1.3) onto these basis vectors, we arrive at the following equivalent representation of Y Z(k) def = Y, e * k = λ(k)β(k) + σξ ′ (k), k = 1, . . . , p,(1.4)where ξ ′ is a standard white Gaussian noise. Spectral cut-off estimators of β are defined by βm (Y ) e k , m = 1, . . . , p,

3 )

 3 Figure 1: Oracle efficiencies of data-driven cut-off frequencies from (2.1), (2.2), and (2.3) for λ(k) = k -1 (left panel Σ = 10; right panel Σ = 100).

Figure 2 :

 2 Figure 2: Oracle efficiencies of data-driven cut-off frequencies from (2.1), (2.2), and (2.3) for λ(k) = k -2 (left panel Σ = 10; right panel Σ = 100).

Figure 3 :

 3 Figure 3: Oracle efficiencies of data-driven cut-off frequencies from (2.1), (2.2), and (2.3) for λ(k) = k -3 (left panel Σ = 10; right panel Σ = 100).

Figure 4 :

 4 Figure 4: Oracle efficiencies of data-driven cut-off frequencies from (2.1), (2.2), and (2.3) for λ(k) = k -4 (left panel Σ = 10; right panel Σ = 100).

  that holds true for any random variable ζ and any µ > 0. With this inequality,(3.16), and (3.18) we arrive atE ζ(m ǫ k ) + max m∈[m ǫ k ,m ǫ k+1 ) [ζ(m) -ζ(m ǫ k )] -ǫEζ 2 (m ǫ k )

β 2

 2 (k) ≤2L -γ/(α+γ) P α/(α+γ) E p k= m+1β 2 (k)λ(k) γ/(α+γ)≤2L -γ/(α+γ) P α/(α+γ) Er(β, m) γ/(α+γ) ≤2L -γ/(α+γ) P α/(α+γ) r • (β) + Kσ 2γ/(α+γ) at the right-hand side in(3.19) may be bounded with the help of Lemmas 2 and 3. Since m ≤ r(W, m)/σ 2 , we obtain from (3.3) that for any µ > 0E exp µ √ m ≤ exp µ r • (β)/σ 2 + Kµ + Kµ 2 . (3.21) With the help of this equation and (3.14), we get E

  we will see below, it is possible to bound from above the maximal risk sup β∈W E β m(Y ) -β 2 with the help of sup β∈W R • (β), where W is an ellipsoid in R p defined by

		2 m	(1.15)
	has a reasonable risk.
	Unfortunately, in general case, point-wise oracle inequalities having the
	form	E β m(Y ) -β 2 ≤ R • (β) + remainder terms
	cannot be proved.
	However, as

Acknowledgments

This work is partially supported by Laboratory for Structural Methods of Data Analysis in Predictive Modeling, MIPT, RF Government grant, ag. 11.G34.31.0073; and RFBR research projects 13-01-12447 and 13-07-12111.


E Chernousova

Yu Golubev
email: golubev.yuri@gmail.com

E Krymova

Ordered Smoothers With Exponential Weighting *

. The estimators withing this family are aggregated using the exponential weighting method and the aim is to control the risk of the aggregated estimate. Based on natural probabilistic properties of the unbiased risk estimate, we derive new oracle inequalities for mean square risk and show that the exponential weighting permits to improve Kneip's oracle inequality.

Introduction and main results

This paper deals with the simple linear model

Y i = µ i + σξ i , i = 1, 2, . . . , n, (1.1)
where ξ is a standard white Gaussian noise, i.e. ξ i are Gaussian i.i.d. random variables with Eξ i = 0 and Eξ 2 i = 1. For the sake of simplicity it is assumed that the noise level σ > 0 is known.

The goal is to estimate an unknown vector µ ∈ R n based on the data Y = (Y 1 , . . . , Y n ) . In this paper, µ is recovered with the help of the following family of linear estimates

μh i (Y ) = h i Y i , h ∈ H, (1.2)
where H is a finite set of vectors in R n which will be described later on.

In what follows, the risk of an estimate μ(Y ) = (μ 1 (Y ), . . . , μn (Y )) is measured by

R(μ, µ) = E µ μ(Y ) -µ 2 ,
where E µ is the expectation with respect to the measure P µ generated by the observations from (1.1) and • , •, • stand for the standard norm and inner product in R n

x 2 = n i=1

x 2 i , x, y = n i=1

x i y i .

One can check very easily that the mean square risk of μh (Y ) is given by

R(μ h , µ) = (1 -h) • µ 2 + σ 2 h 2 ,
where x • y denotes the coordinate-wise product of vectors x, y ∈ R n , i.e. z = x • y, means that z i = x i y i , i = 1, . . . , n. So, R(μ h , µ) depends on h ∈ H and one can minimize it choosing properly h ∈ H. Very often the minimal risk r H (µ) = min cannot be used since it depends on the unknown vector µ. However, one could try to construct an estimate μH (Y ) based on the family of linear estimates μh (Y ), h ∈ H, the risk of which is close to the oracle risk. This means that the risk of μH (Y ) might be bounded from above by the so-called oracle inequality R(μ H , µ) ≤ r H (µ) + ∆H (µ), which holds uniformly in µ ∈ R n . Heuristically, this inequality assumes that the remainder term ∆H (µ) is negligible with respect to the oracle risk for all µ ∈ R n . In general, such an estimator doesn't exist, but for certain statistical models one can construct an estimator μH (Y ) (see, e.g., Theorem 1.1 below) such that:

• ∆H (µ) ≤ Cr H (µ) for all µ ∈ R n , where C > 1 is a constant.

• ∆H (µ) r H (µ) for all µ such that r H (µ) σ 2 .

It is also well-known that one can find an estimator with the above properties provided that H is not very rich. In particular, as shown in [START_REF] Kneip | Ordered linear smoothers[END_REF], this can be done for the so-called ordered smoothers. This is why this paper deals with H containing solely ordered multipliers defined as follows:

Definition 1.1. H is a set of ordered multipliers if the following properties hold

• contracting property: h i ∈ [0, 1], i = 1, . . . , n for all h ∈ H,

• decreasing property: h i+1 ≤ h i , i = 1, . . . , n for all h ∈ H,

• totally ordered property: if for some integer k and some h, g ∈ H, h k < g k , then h i ≤ g i for all i = 1, . . . , n.

The totally ordered property means that vectors in H may be naturally ordered, since for any h, g ∈ H there are only two possibilities h i ≤ g i or h i ≥ g i for all i = 1, . . . , n. Therefore the estimators defined by (1.2), where H is a set of ordered multipliers, are often called ordered smoothers [START_REF] Kneip | Ordered linear smoothers[END_REF].

Note that ordered smoothers are common in statistics. Below we give two basic examples, where these smoothers appear naturally.

Smoothing splines. They are usually used in recovering smooth regression functions f (x), x ∈ [0, 1], given the noisy observations Z = (Z 1 , . . . , Z n )

Z i = f (X i ) + εξ i , i = 1, . . . , n, (1.3)
where the design points X i belong to (0, 1) and ξ is a standard white Gaussian noise. It is well known that smoothing splines are defined by fα (x, Z) = arg min

f 1 n n i=1 [Z i -f (X i )] 2 + α 1 0 [f (m) (x)] 2 , (1.4)
where f (m) (•) denotes the derivative of order m and α > 0 is a smoothing parameter which is usually chosen with the help of the Generalized Cross Validation (see, e.g., [START_REF] Wahba | Spline Models for Observational Data[END_REF]).

In order to show that the regression estimation with the help of the smoothing splines is equivalent to the sequence space model (1.1), consider the Demmler-Reinsch [START_REF] Demmler | Oscillation matrices with spline smoothing[END_REF] basis ψ k (x), x ∈ [0, 1], k = 1, . . . , n having double orthogonality

ψ k , ψ l n = δ kl , 1 0 ψ (m) k (x)ψ (m) l (x) dx = δ kl λ k , k, l = 1, . . . , n, (1.5)
where δ kl = 1 if k = l, and δ kl = 0 otherwise. Here and below u, v n stands for the inner product

u, v n = 1 n n i=1 u(X i )v(X i ).
Let us assume for definiteness that the eigenvalues λ k are sorted in ascending order λ 1 ≤ . . . ≤ λ n .

With this basis, representing the underlying regression function as follows:

f (x) = n k=1 ψ k (x)µ k , (1.6)
we get from (1.3) and (1.5)

Y k = 1 n n i=1 Z i ψ k (X i ) = µ k + ε √ n ξ k , (1.7)
where µ k = f, ψ k n and ξ is a standard white Gaussian noise. So, substituting (1.6) in (1.4) and using (1.5), we arrive at

fα (x, Z) = n k=1 μk ψ k (x), where μ = arg min µ n k=1 [Y k -µ k ] 2 + α n k=1 λ k µ 2 k .
It can be seen easily that

μk = Y k 1 + αλ k and thus the model (1.3)-(1.4) is equivalent to (1.1)-(1.2) with σ = ε/ √ n and H = h : h k = 1 1 + αλ k , α ∈ R + . (1.8)
If we are interested in minimax regression estimates on Sobolev's classes, they can be easily constructed using the statistics from (1.7) and the following set of ordered multipliers

H = h : h k = max(1 -αλ k , 0), α ∈ R + .
See [START_REF] Nussbaum | Spline smoothing in regression models and asymptotic efficiency in L 2[END_REF] and [START_REF] Speckman | Spline smoothing and optimal rates of convergence in nonparametric regression[END_REF] for details.

Note that the Demmler-Reinsch basis is a very useful tool for statistical analysis of spline methods. However, in constructing statistical estimates, this basis is rarely used since there are very fast algorithms for computing smoothing splines (see, e.g., [START_REF] Green | Nonparametric Regression and Generalized Linear Models. A roughness penalty approach[END_REF] and [START_REF] Wahba | Spline Models for Observational Data[END_REF]).

Spectral regularizations of large linear models. Very often in linear models, we are interested in estimating Xθ ∈ R n based on the observations

Z = Xθ + σξ, (1.9)
where X is a known n × p -matrix, θ ∈ R p is an unknown vector, and ξ is a standard white Gaussian noise. It is well known that if X X has a large condition number or p is large, then the standard maximum likelihood estimate X θ0 (Z), where

θ0 (Z) = arg min θ Z -Xθ 2 = (X X) -1 X Z
may result in a large risk. In particular, if X X > 0, then

E Xθ -X θ0 2 = σ 2 p.
When p is large, this risk may be improved with the help of a regularization term. For instance, one can use the Phillps-Tikhonov regularization [START_REF] Tikhonov | Solution of Ill-posed Problems[END_REF] (often called ridge regression in statistics)

θα (Z) = arg min θ Z -Xθ 2 + α θ 2 ,
where α > 0 is a smoothing parameter. It can be seen easily that

θα (Z) = [I + α(X X) -1 ] -1 θ0 (Z).
This formula is a particular case of the so-called spectral regularizations defined as follows (see, e.g., [START_REF] Engl | Regularization of Inverse Problems[END_REF]):

θα (Z) = H α (X X) θ0 (Z),
where H α (X X) is a matrix depending on a smoothing parameter α ∈ R + and X X and admitting the following representation

H α (X X) = p k=1 h α (λ k )e k e k
where e k , k = 1, . . . , p and λ 1 ≤ λ 2 ≤ . . . ≤ λ p are eigenvectors and eigenvalues of X X, and

h α (•) is a function R + → [0, 1].
The standard way to construct an equivalent model for the spectral regularization method is to make use of the SVD. Note that

e * k = Xe k √ λ k , k = 1, . . . , p
is an orthonormal system in R n . Therefore the observations Z (see (1.9)) admit the following equivalent representation

Y k = e * k , Z = e * k , Xθ + σξ k , k = 1, . . . , p, (1.10)
where ξ is a standard white Gaussian noise. Noticing that

X θα (Z) = XH α (X X)(X X) -1 X Z,
we have

X θα (Z), e * k = p s=1 Y s XH α (X X)(X X) -1 X e * s , e * k = p s=1 Y s λ k H α (X X)(X X) -1 e s , e k = h α (λ k )Y k .
Therefore from this equation and (1.10) and we see that the spectral regularization method is equivalent to the statistical model defined by (1.1) and (1.2) with H = {h :

h k = h α (λ k ), α ∈ R + }.
Note that for the Phillps-Tikhonov method we have

h α (λ) = 1 1 + α/λ
and it is clear that the corresponding H is a set of ordered multipliers. Along with this regularization method, the spectral cut-off method and Landweber's iterations (see, e.g., [START_REF] Engl | Regularization of Inverse Problems[END_REF] for details) are typical examples of ordered smoothers.

Nowadays, there are a lot of approaches aimed to construct estimates mimicking the oracle risk. At the best of our knowledge, the principal idea in obtaining such estimates goes back to [START_REF] Akaike | Information theory and an extension of the maximum likelihood principle[END_REF] and [START_REF] Mallows | Some comments on C p[END_REF] and related to the method of the unbiased risk estimation [START_REF] Stein | Estimation of the mean of a multivariate normal distribution[END_REF]. The literature on this approach is so vast that it would be impractical to cite it here. We mention solely the following result by Kneip [12] since it plays an important role in our presentation. Denote by r(Y, μh )

def = Y -μh (Y ) 2 + 2σ 2 n i=1 h i -σ 2 n, (1.11)
the unbiased risk estimate of μh (Y ).

Theorem 1.1. Let ĥ = arg min h∈H r(Y, μh ) be a minimizer of the unbiased risk estimate. Then uniformly in µ ∈ R n ,

E µ ĥ • Y -µ 2 ≤ r H (µ) + Kσ 2 1 + r H (µ) σ 2 , (1.12)
where K is a generic constant.

Another well-known idea to construct a good estimator based on the family μh , h ∈ H is to aggregate the estimates within this family using a held-out sample. Apparently, this approach was firstly developed by Nemirovsky in [START_REF] Nemirovski | Topics in non-parametric statistics[END_REF] (see also [START_REF] Juditsky | Functional aggregation for nonparametric regression[END_REF]) and independently by Catoni (see [START_REF] Catoni | Statistical learning theory and stochastic optimization[END_REF] for a summary). Later, the method was extended to several statistical models (see, e.g., [START_REF] Yang | Combining different procedures for adaptive regression[END_REF], [START_REF] Rigollet | Linear and convex aggregation of density estimators[END_REF], [START_REF] Lecué | Simultaneous adaptation to the margin and to complexity in classification[END_REF], [START_REF] Rigolet | Sparse estimation by exponential weighting[END_REF]).

To overcome the well-know drawbacks of sample splitting, one would like to aggregate estimators using the same observations for constructing estimators and performing the aggregation. This can be done, for instance, with the help of the exponential weighting. The motivation of this method is related to the problem of functional aggregation, see [START_REF] Rigolet | Sparse estimation by exponential weighting[END_REF]. It has been shown that this method yields rather good oracle inequalities for certain statistical models [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF], [START_REF] Dalayan | Sharp oracle inequalities for aggregation of affine estimators[END_REF], [START_REF] Rigolet | Sparse estimation by exponential weighting[END_REF], [START_REF] Alquier | Pac-bayesian theorems for sparse regression estimation with exponential weights[END_REF], [START_REF] Arias-Castro | K Variable Selection with Exponential Weights and l 0 -Penalization[END_REF].

In the considered statistical model, the exponential weighting estimate is defined as follows:

μ(Y ) = h∈H w h (Y )μ h (Y ), (1.13)
where

w h (Y ) = π h exp - r(Y, μh ) 2βσ 2 g∈H π g exp - r(Y, μg ) 2βσ 2 , β > 0,
r(Y, μh ) is the unbiased risk estimate of μh (Y ) defined by (1.11), and a priori weights {π h , h ∈ H} are non-negative and such that h∈H π h > 0. Recall that for simplicity, it is assumed here and in what follows that H is discrete and finite.

It has been shown in [START_REF] Dalayan | Sharp oracle inequalities for aggregation of affine estimators[END_REF] that for this method the following oracle inequalities hold.

Theorem 1.2. Assume that h∈H π h = 1. If β ≥ 4, then uniformly in µ ∈ R n R(μ, µ) ≤ min λ h ≥0: λ 1 =1 h∈H λ h R(μ h , µ) + 2σ 2 βK(λ, π) , R(μ, µ) ≤ min h∈H R(μ h , µ) + 2σ 2 β log 1 π h , (1.14)
where K(•, •) is the Kullback-Leibler divergence and • 1 stands for 1 -norm, i.e.,

K(λ, π) = h∈H λ h log λ h π h , λ 1 = h∈H |λ h |.
Note that for projection methods (h k ∈ {0, 1}) this theorem holds for β ≥ 2, see [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF].

It is clear that if we want to derive from (1.14) an oracle inequality similar to (1.12), then we have to chose π h = 1/(#H), where #H denotes the cardinality of H, and thus we arrive at

R(μ, µ) ≤ r H (µ) + 2σ 2 β log(#H).
This oracle inequality is good only when the cardinality of H is not very large. If we deal with H having a very large cardinality like those related to smoothing splines, this inequality is not good. To some extent, this situation may be improved, see Proposition 2 in [START_REF] Dalayan | Sharp oracle inequalities for aggregation of affine estimators[END_REF]. However, looking at the oracle inequality in this proposition, one cannot say, unfortunately, that it is always better than (1.12).

The main goal in this paper is to show that for the exponential weighting method one can obtain an oracle inequality with a smaller remainder term than the one in Theorem 1.1, Equation (1.12).

In order to attain this goal and to cover H with both small and large cardinality, we make use of the special prior weights defined as follows:

π h def = 1 -exp - h + 1 -h 1 β . (1.15)
Here h + = min{g ∈ H : g > h} and π h max = 1, where h max is the maximal multiplier in H.

Along with these weights we will need also the following condition:

Condition 1.1. There exists a constant K • ∈ (0, ∞) such that h 2 -g 2 ≥ K • h 1 -g 1 (1.16)
for all h ≥ g from H.

The next theorem, yielding an upper bound for the mean square risk of μ(Y ) defined by (1.13), is the main result of this paper.

Theorem 1.3. Assume that H is a set of ordered multipliers, β ≥ 4, and Condition 1.1 holds. Then, uniformly in µ ∈ R n ,

E µ μ(Y ) -µ 2 ≤ r H (µ) + 2βσ 2 log C 1 + r H (µ) σ 2 .
(1.17)

Here and in what follows C = C(K • , β) denotes strictly positive and bounded constants depending on K • , β.

We finish this section with some remarks regarding this theorem.

Remark 1. The condition β ≥ 4 may be improved when the ordered multipliers h ∈ H take only two values 0 and 1. In this case it is sufficient to assume that β ≥ 2 (see [START_REF] Golubev | Exponential weighting and oracle inequalities for projection estimates[END_REF]).

Remark 2. Usually Condition 1.1 may be checked rather easily. For instance, for smoothing splines and equidistant design, the set of ordered multipliers is given by (1.8) and this condition follows from the well-known asymptotic formula λ k (πk) 2m for large k (see [START_REF] Demmler | Oscillation matrices with spline smoothing[END_REF] for details). Heuristically, for small α and large n we have

h α 2 = n k=1 1 (1 + αλ k ) 2 ≈ n k=1 1 [1 + α(πk) 2m ] 2 ≈ 1 πα 1/(2m) ∞ 0 1 [1 + x 2m ] 2 dx and h α 1 = n k=1 1 1 + αλ k ≈ n k=1 1 1 + α(πk) 2m ≈ 1 πα 1/(2m) ∞ 0 1 1 + x 2m dx. (1.18)
With these equations Condition 1.1 becomes obvious. A rigorous proof of (1. [START_REF] Nemirovski | Topics in non-parametric statistics[END_REF]) is based on a non-asymptotic version of these arguments. It is technical but unfortunately cumbersome and therefore, in order not to overload the paper, we omit it.

For spectral regularizations Condition 1.1 is obvious for the spectral cutoff method. For the Phillps-Tikhonov method the proof is more involved but similar to the one for splines. Note, however, that in this case the following condition λ s ≥ Ks q , for some q > 1 is required.

Remark 3. In practice, the multipliers in H are often chosen so that

h + 1 = (1 + ) h 1
with some small > 0 and the initial condition h min 1 = 1, where h min is the minimal multiplier in H. In this case one may chose the a priori weights π h = 1 since π h from (1.15) are strictly bounded from below and it can be checked very easily that Lemma 2.1 holds true for π h = 1. This means in particular that if the smoothing parameter α in smoothing splines (1.4) takes values {(1 + ) -2mk , k = 0, 1, . . .}, then π h = 1 may be used in the exponential weighting (see (1.18)).

Remark 4. In contrast to Proposition 2 in [START_REF] Dalayan | Sharp oracle inequalities for aggregation of affine estimators[END_REF], the remainder term in (1.17) does depend neither on the cardinality of H nor n. It has the same structure as Kneip's oracle inequality in Theorem 1.1.

Remark 5. Comparing (1.17) with (1.12), we see that when

r H (µ) σ 2 ≈ 1,
then the remainder terms in (1.12) and (1.17) have the same order, namely, σ 2 . However, when r H (µ)

σ 2 1,
we get

2βσ 2 log C 1 + r H (µ) σ 2 Kσ 2 1 + r H (µ) σ 2 ,
thus showing that the upper bound for the remainder term in the oracle inequality related to the exponential weighting is better than the one in Theorem 1.1.

Remark 6. To compare the actual remainder terms in (1.17) and (1.12) and to find out what β is good from a practical viewpoint, a numerical experiment has been carried out. The goal in this experiment is to compare the exponential weighting methods for β = {0, 1, 2, 4} combined with the cubic smoothing splines for the equidistant design. In order to simulate these splines, the following family of ordered multipliers was used

H = h : h k = 1 1 + [α(k -1)] 4 , α > 0 .
The motivation of this family is due to the well-known asymptotic formula λ k (πk) 4 , k → ∞.

The simulations are organized as follows. For given A ∈ [0, 300], 100000 replications of the observations

Y k = µ k (A) + ξ k , k = 1, . . . , 400
are generated. Here µ(A) ∈ R 400 is a Gaussian vector with independent components and

Eµ k (A) = 0, Eµ 2 k (A) = A exp - k 2 2Σ 2 .
Next, the mean oracle risk

rH (A) = E min h∈H (1 -h) • µ(A) 2 + h 2
and the mean excess risk

∆β (A) = E µ(A) -μ(Y ) 2 -rH (A),
were computed with the help of the Monte-Carlo method. Finally, the data {r H (A), ∆β (A), A ∈ [0, 300]} are plotted on Figure 1 to illustrate graphically the remainder term ∆ β (r

H ) = E µ μ -µ 2 -r H (µ).
Looking at this picture we see that there is no universal β minimizing the excess risk uniformly in µ. However, it seems to us that a reasonable choice is β ≈ 1, but unfortunately, good oracle inequalities are not available for this case. Note also that the exponential weighting can provide only a moderate improvement of the risk compared to the classical unbiased risk estimation (β = 0). All methods demonstrate almost similar statistical performance. However, when r H (µ)/σ 2 is not large, the exponential weighting works usually better.

Proofs

The proof of Theorem 1.3 is based on a combination of methods for deriving oracle inequalities proposed in [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF] and [START_REF] Golubev | Exponential weighting and oracle inequalities for projection estimates[END_REF]. We indicate here solely the main steps in the proof, all details are given below. With the help of Stein's formula for the unbiased risk estimate it can be shown similar to [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF] that for β ≥ 4

E µ μ -µ 2 ≤E µ h∈H w h (Y )r(Y, μh ) ≤r H (µ) + 2βσ 2 E µ h∈H w h (Y ) log π h w h (Y ) -2βσ 2 E µ log h∈H π h exp - r(Y, μh ) -r(Y, μĥ ) 2βσ 2 , (2.1)
where ĥ = arg min h∈H r(Y, μh ).

To control the right-hand side at this equation, we make use of the ordering property of estimates μh , h ∈ H. First, we check using (??) that if π h is defined by (1.15), then

h∈H π h exp - r(Y, μh ) -r(Y, μĥ ) 2βσ 2 ≥ h≥ ĥ π h exp - r(Y, μh ) -r(Y, μĥ ) 2βσ 2 ≥ 1,
and so, the last term in Equation (2.1) is always negative.

The most difficult and delicate part of the proof is related to the average Kullback-Leibler divergence

E µ h∈H w h (Y ) log w h (Y ) π h .
To compute a good lower bound for this value, we follow the approach proposed in [START_REF] Golubev | Exponential weighting and oracle inequalities for projection estimates[END_REF]. The main idea here is to make use of the following property of the unbiased risk estimate: for any sufficiently small < 1, there exists ĥ depending on Y such that with probability 1, for all h ≥ ĥ r(Y, μh ) -r(Y, μĥ ) ≥ 2βσ 2 h 2 -ĥ 2 + 2βσ

π h ≤ 1 + ĥ 2 K • β .
More cumbersome probabilistic technique is required to prove the following upper bound (see Lemma 2.5):

E µ ĥ 2 ≤ r H (µ) (1 -2β )σ 2 + √ 1 + 2β 1 -2β √ K.
Finally, combining the above equations, we arrive at (1.17).

Auxiliary facts

The next lemma collects some useful facts about a priori weights defined by (1.15). Let

D h = g ∈ H : h 1 ≤ g 1 ≤ h 1 + 1 .
Lemma 2.1. Under Condition 1.1, for any h ∈ H, the following assertions hold:

g≥h π g exp - g 1 β = exp - h 1 β , (2.2)
g≤h π g ≤ 1 + h 2 K • β , (2.3
)

g∈D h π g ≤ 1 + 1 β , (2.4
)

g∈D h π g ≥ 1 2β exp - 1 β . (2.5)
Proof. Denote for brevity

S(h) = g≥h π g exp - g 1 -h 1 β .

Then we have

S(h) -S(h

+ ) = π h + exp - h + 1 -h 1 β × g≥h + π g exp - g 1 -h + 1 β - g≥h + π g exp - g 1 -h + 1 β = π h -1 -exp - h + 1 -h 1 β S(h + ).
Therefore in view of the definition of π h , it is clear that if S(h max ) = 1, then S(h) = S(h + ), thus proving (2.2).

To prove (2.3), note that

π g ≤ g + 1 -g 1 β . (2.6)
This inequality follows from (1.15) and from the inequality inequality 1exp(-x) ≤ x. Hence, by Condition (1.16) we obtain

g≤h π g ≤ 1 + g≤h - π g ≤ 1 + 1 β g≤h - g + 1 -g 1 = 1 + h 1 -h min 1 β ≤ 1 + h 2 -h min 2 K • β ,
where h min is the minimal element in H.

The same arguments can be used in proving (2.4).

In order to check (2.5), denote by g h be the maximal element in D h . Then there are two possibilities

• g h 1 ≤ h 1 + 1/2, • g h 1 > h 1 + 1/2.
In the first case, g + h 1 -g h 1 ≥ 1/2, and thus by

g∈D h π g ≥ π g h ≥ 1 -exp - g + h 1 -g h 1 β ≥ 1 -exp - 1 2β . (2.7)
In the second case, where h 1 + 1/2 < g h 1 ≤ h 1 + 1, we obtain by a Taylor expansion that for any g < g h

π g ≥ g + 1 -g 1 β exp - g h 1 -h 1 β ≥ g + 1 -g 1 β exp - 1 β ,
and thus

g∈D h π g ≥ g h 1 -h 1 β exp - 1 β ≥ 1 2β exp - 1 β .
This equation together with (2.7) ensures (2.5) since it can be checked with a simple algebra that

1 2β exp - 1 β ≤ 1 -exp - 1 2β , β > 0.
The following lemma is a cornerstone in the proof of Theorem 1.3.

Lemma 2.2. For β ≥ 4 the risk of μ(Y ) is bounded from above as follows:

E µ μ(Y ) -µ 2 ≤ E µ h∈H w h (Y )r(Y, μh ).
Proof. It is based essentially on the method proposed in [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF]. Unfortunately, we cannot use directly Corollary 2 in [START_REF] Leung | Information theory and mixing least-squares regressions[END_REF] because it holds only for h k , k = 1, . . . , n taking values 0 and 1. In the case of ordered smoothers h k belongs to the interval [0, 1] and we will see below that this fact results in the condition β ≥ 4.

Recall that the unbiased risk estimates for μi (Y ) and μh i (Y ) are computed as follows (see, e.g. [START_REF] Stein | Estimation of the mean of a multivariate normal distribution[END_REF])

r(Y i , μi ) = [μ i (Y ) -Y i ] 2 + 2σ 2 ∂ μi (Y ) ∂Y i -σ 2 , r(Y i , μh i ) = [μ h i (Y ) -Y i ] 2 + 2σ 2 h i -σ 2 .
(2.8)

Since h∈H w h = 1, we have [μ i (Y ) -Y i ] 2 = h∈H w h (Y )[μ i (Y ) -Y i ] 2 = h∈H w h (Y )[μ i (Y ) -μh i (Y ) + μh i (Y ) -Y i ] 2 = h∈H w h (Y )[μ i (Y ) -μh i (Y )] 2 + h∈H w h (Y )[μ h i (Y ) -Y i ] 2 + 2 h∈H w h (Y )[μ i (Y ) -μh i (Y )][μ h i (Y ) -Y i ] = h∈H w h (Y )[μ i (Y ) -μh i (Y )] 2 + h∈H w h (Y )[μ h i (Y ) -Y i ] 2 + 2 h∈H w h (Y )[μ i (Y ) -μh i (Y )][μ h i (Y ) -μi (Y ) + μi (Y ) -Y i ] = - h∈H w h (Y )[μ i (Y ) -μh i (Y )] 2 + h∈H w h (Y )[μ h i (Y ) -Y i ] 2 .
(2.9)

From the definition of μ(Y ) we obviously get

∂ μi (Y ) ∂Y i = h∈H w h (Y ) ∂ μh i (Y ) ∂Y i + h∈H ∂w h (Y ) ∂Y i μh i (Y )
and combining this equation with (2.9) (see also (2.8)), we arrive at

r(Y i , μi ) = [μ i (Y ) -Y i ] 2 + 2σ 2 ∂ μi (Y ) ∂Y i -σ 2 = h∈H w h (Y ) [μ h i (Y ) -Y i ] 2 + 2σ 2 ∂ μh i (Y ) ∂Y i -σ 2 -[μ i (Y ) -μh i (Y )] 2 + 2σ 2 ∂ log(w h (Y )) ∂Y i μh i (Y ) = h∈H w h (Y )r(Y i , μh i )+ + h∈H w h (Y ) -[μ i (Y ) -μh i (Y )] 2 + 2σ 2 ∂ log[w h (Y )] ∂Y i μh i (Y ) = h∈H w h (Y )r(Y i , μh i ) + h∈H w h (Y ) -[μ i (Y ) -μh i (Y )] 2 + 2σ 2 ∂ log[w h (Y )] ∂Y i [μ h i (Y ) -μi (Y )] .
(2.10)

In deriving the above equation it was used that h∈H w h (Y ) = 1 and hence

h∈H ∂w h (Y ) ∂Y i = h∈H w h (Y ) ∂ log w h (Y ) ∂Y i = 0.
To control the second sum at the right hand side of (2.10) we make use of the following equation

log[w h (Y )] = - r(Y, μh ) 2βσ 2 + log(π h ) -log g∈H π g exp - r(Y, μg ) 2βσ 2 . (2.11) Therefore h∈H w h (Y ) ∂ log w h (Y ) ∂Y i [μ h i (Y ) -μi (Y )] = - 1 2βσ 2 h∈H w h (Y ) ∂ r(Y, μh ) ∂Y i [μ h i (Y ) -μi (Y )].
Substituting in the above equation (see (1.11)

) ∂ r(Y i , μh i ) ∂Y i = 2(1 -h i ) 2 Y i , we obtain h∈H w h (Y ) ∂ log w h (Y ) ∂Y i [μ h i (Y ) -μi (Y )] = - 1 βσ 2 Y 2 i h∈H w h (Y )[h i -1] 2 [h i -hi ], (2.12)
where hi =

h∈H w h (Y )h i . Next noticing that (1 -h i ) 2 = (1 -hi ) 2 + ( hi -h i ) 2 + 2(1 -hi )( hi -h i ), we have -Y 2 i h∈H w h (Y )(h i -1) 2 (h i -hi ) = Y 2 i (1 -hi ) 2 h∈H w h (Y )( hi -h i ) +Y 2 i h∈H w h (Y )( hi -h i ) 2 ( hi -h i + 2 -2 hi ) = 2Y 2 i h∈H w h (Y )( hi -h i ) 2 1 - h i + hi 2 ≤ 2 h∈H w h (Y )[μ i (Y ) -μh i (Y )] 2 .
Combining this equation with (2.9)-(2.12), we finish the proof.

Lemma 2.3. Suppose {q h ≤ 1, h ∈ H} is a nonnegative sequence such that for all h ≥ h

q h ≤ exp -γ h 1 -h 1 -1 , γ > 0. Let W h = π h q h g∈H π g q g -1
and G be a subset in H. Then

H(W, π) def = h∈H W h log π h W h ≤ log h≤ h π h + exp[R(γ)] ,
where

R(γ) = log 2 γβe + h∈G π h + h∈G π h q h -1 8 γβe + h∈G π h . (2.13)
Proof. Decompose H onto two subsets

Q = {h ≥ h} ∪ G, P = H \ Q
and denote for brevity

P = h∈P π h q h , Q = h∈Q π h q h .
By convexity of log(x) we obtain

H(W, π) = P P + Q h∈P π h q h P log (P + Q)/P q h /P + Q P + Q h∈Q π h q h Q log (P + Q)/Q q h /Q ≤ P P + Q log P + Q P + Q P + Q log P + Q Q + P P + Q log h∈P π h + 1 P + Q h∈Q π h q h log 1 q h + Q log(Q) . (2.14)
Next, note that f (x) = x log(1/x) is an increasing function, when x ∈ (0, e -1 ], and max x∈[0,1] f (x) = e -1 . Therefore, if

q h ≤ exp -γ( h 1 -h 1 ) -1 , then q h log 1 q h ≤ exp -γ( h 1 -h 1 ) -1 γ( h 1 -h 1 ) + 1 ,
and we get

h∈Q π h q h log 1 q h ≤ 1 e h∈G π h + 1 e h≥ h π h exp -γ( h 1 -h 1 ) × γ( h 1 -h 1 ) + 1 .
We continue this equation with (2.6) as follows:

h∈Q π h q h log 1 q h ≤ 1 e h∈G π h + 1 βe h≥ h exp -γ( h 1 -h 1 ) × γ h 1 -h 1 + 1 h + 1 -h 1 . (2.15)
In order to bound from above the right-hand side at this equation, let us index the elements in {h ∈ H : h ≥ h} denoting them as {h k , k ≥ 0}, so that h k+1 1 ≥ h k 1 and h 0 = h. Then, denoting for brevity

S i = h i 1 -h 1 ,
we can rewrite the sum at the right hand side in (2.15) as follows:

h≥ h exp -γ( h 1 -h 1 ) γ h 1 -h 1 + 1 h + 1 -h 1 = i≥0 exp -γS i γS i + 1 (S i+1 -S i ).
(2.16) To bound from above the right hand side at this equation, let us check that max

S k ,k≥1 i≥0 exp -γS i (S i+1 -S i ) ≤ 2 γ , (2.17)
where max is computed over all nondecreasing sequences. Solving the equation

∂ ∂S k i≥0 exp -γS i [S i+1 -S i ] + = 0,
we obtain with a simple algebra

S k+1 -S k = exp γ(S k -S k-1 ) -1 γ . Hence exp(-γS k )(S k+1 -S k ) = exp[-γS k-1 ] -exp[-γS k ]
γ and summing up these equations, we arrive at (2.17). Next notice that for any z > 0 we have (1 + x) ≤ z exp(x/z). With this inequality and (2.17) we obtain for any z >

1 i≥0 exp -γS i γS i + 1 (S i+1 -S i ) ≤ z i≥0 exp[-γ(1 -1/z)S i ](S i+1 -S i ) ≤ 2z 2 (z -1)γ .
It can be seen easily that the minimun of the right hand side at this equation is attained at z = 2. So, we obtain

i≥0 exp -γS i γS i + 1 (S i+1 -S i ) ≤ 8 γ . ( 2
H(W, π) ≤ max x∈[0,1] -x log(x) -(1 -x) log(1 -x) +(1 -x) log h∈P π h + xρ , (2.22)
where

ρ def = log(Q) + 1 Q h∈Q π h q h log 1 q h ≤ log 2 γβe + h∈G π h + h∈G π h q h -1 8 γβe + h∈G π h = R(γ).
It is seen easily that the minimizer x * of the right-hand side at (2.22) is a solution to the following equation .

Therefore from (2.22) we get

H(W, π) ≤ log h∈P π h -log(1 -x * ) -x * log x * 1 -x * + log h∈P π h -ρ = log h∈P π h -log(1 -x * ) = log h∈P π h + exp(ρ) ≤ log h< h π h + exp(ρ) .
Lemma 2.4. Let ξ i be i.i.d. N (0, 1) and H be a set of ordered multipliers. Then for any α ∈ (0, 1/4)

E max h∈H ± n i=1 (h 2 i -2h i )(ξ 2 i -1) -α n i=1 h 2 i ≤ K α , E max h∈H n i=1 (1 -h i ) 2 ξ i µ i -α n i=1 (1 -h i ) 2 µ 2 i ≤ K α ,
where K is a generic constant.

Proof. It follows from Lemma 2 in [START_REF] Golubev | On universal oracle inequalities related to high dimensional linear models[END_REF].

Lemma 2.5. Let

ĥ = max h : r(Y, μh ) -rH (Y ) ≤ 2β σ 2 h 2 -ĥ 2 + 2βσ 2 , (2.23)
where ∈ (0, 1/(2β)) and rH (Y ) = min h∈H r(Y, μh ). Then

E µ ĥ 2 ≤ r H (µ) (1 -2β )σ 2 + √ 1 + 2β 1 -2β √ K. (2.24)
Proof. By the definition of r(Y, μh ), see also (1.11) and (2.23), we get ĥ = max h :

(1 -h) • µ 2 + σ 2 (1 -2β ) h 2 + 2σ n i=1 (1 -h i ) 2 µ i ξ i + σ 2 n i=1 (h 2 i -2h i )(ξ 2 i -1) ≤ (1 -ĥ) • µ 2 + σ 2 (1 -2β ) ĥ 2 + 2σ n i=1 (1 -ĥi ) 2 µ i ξ i + σ 2 n i=1 ( ĥ2 i -2 ĥi )(ξ 2 i -1) + 2βσ 2 .
Let us fix some γ, γ > 0. Then we can rewrite the above equation as follows:

ĥ = max h : σ 2 (1 -2β -γ) h 2 + 2σ n i=1 (1 -h i ) 2 µ i ξ i + (1 -h) • µ 2 +σ 2 n i=1 (h 2 i -2h i )(ξ 2 i -1) + γσ 2 h 2 ≤ (1 + γ ) (1 -ĥ) • µ 2 + σ 2 (1 -2β + γ ) ĥ 2 +2σ n i=1 (1 -ĥi ) 2 µ i ξ i -γ (1 -ĥ) • µ 2 +σ 2 n i=1 ( ĥ2 i -2 ĥi )(ξ 2 i -1) -γ σ 2 ĥ 2 + 2βσ 2 . Therefore ĥ ≤ h def = max h : σ 2 (1 -2β -γ) h 2 + min g∈H 2σ n i=1 (1 -g i ) 2 µ i ξ i + (1 -g) • µ 2 + min g∈H σ 2 n i=1 (g 2 i -2g i )(ξ 2 i -1) + γσ 2 g 2 ≤(1 + γ ) (1 -ĥ)µ 2 + (1 -2β + γ )σ 2 ĥ 2 + max g∈H 2σ n i=1 (1 -g i ) 2 µ i ξ i -γ (1 -g) • µ 2 + max g∈H σ 2 n i=1 (g 2 i -2g i )(ξ 2 i -1) -γ σ 2 g 2 + 2βσ 2 .
Next, bounding max and min in this equation with the help of Lemma 2.4, we arrive at

(1 -2β )σ 2 E µ h 2 - Kσ 2 γ -γσ 2 E µ h 2 -Kσ 2 ≤ E µ (1 -ĥ) • µ 2 + (1 -2β )σ 2 E µ ĥ 2 + Kσ 2 γ + γ E µ (1 -ĥ) • µ 2 + Kσ 2 γ + γ σ 2 E µ ĥ 2 + 2βσ 2 .
Maximizing the left-hand side in γ ∈ (0, 1/4) and minimizing the right-hand side in γ ∈ (0, 1/4), we obtain with a simple algebra

σ 2 (1 -2β )E µ h 2 - K 1 -2β 1/2 2 ≤ E µ (1 -ĥ) • µ 2 + σ 2 E µ ĥ 2 + σ √ K 2 + 2β Kσ 2 1 -2β + (2β + K)σ 2 .
Combining this equation with ĥ 2 ≤ h 2 we arrive at

(1 -2β )E µ ĥ 2 ≤ σ -1 E µ (1 -ĥ) • µ 2 + σ 2 E µ ĥ 2 + √ 1 + 2β √ 1 -2β √ K + 2β + K.
(2.25)

To control the expectation at the right-hand side in (2.25), we note that for any given g ∈ H the following inequality

n i=1 [1 -ĥi ] 2 Y 2 i + 2σ 2 n i=1 ĥi ≤ n i=1 [1 -g i ] 2 Y 2 i + 2σ 2 n i=1 g i holds. Denoting for brevity ρ(h) def = (1 -h) • µ 2 + σ 2 h 2 ,
we rewrite this equation as follows:

ρ( ĥ)+2σ n i=1 (1 -ĥi ) 2 µ i ξ i + σ 2 n i=1 ( ĥ2 i -2 ĥi )(ξ 2 i -1) ≤ρ(g) + 2σ n i=1 (1 -g i ) 2 µ i ξ i + σ 2 n i=1 (g 2 i -2g i )(ξ 2 i -1).
So, for any γ > 0, we get with this equation and Lemma 2.4

E µ ρ( ĥ) ≤ρ(g) + γE µ ρ( ĥ)

+ 2σE µ max h∈H - n i=1 (1 -h i ) 2 µ i ξ i - γ 2σ n i=1 (1 -h i ) 2 µ 2 i + σ 2 E µ max h∈H n i=1 [2h i -h 2 i ](ξ 2 i -1) -γ n i=1 h 2 i ≤ρ(g) + Kσ 2 γ + γE µ ρ( ĥ).
Next, minimizing the right-hand side in g ∈ H and γ ∈ (0, 1/4), we obtain

E µ ρ( ĥ) ≤r H (µ) + 2σ KE µ ρ( ĥ) 1/2 + Kσ 2
or, equivalently,

E µ ρ( ĥ) 1/2 - √ Kσ 2 ≤ r H (µ) + Kσ 2 .
This yields obviously

E µ (1 -ĥ) • µ 2 + σ ĥ 2 1/2 ≤ r H (µ) + 2σ √ K.
Finally, substituting this inequality in (2.25), we get (2.24). , where ĥ = arg min h∈H r(Y, μh ). Therefore = log g≥ ĥ π g exp -

Proof of

(1 -g) • Y 2 -(1 -ĥ) • Y 2 2βσ 2 - 1 β n i=1
[g i -ĥi ] ≥ log g≥ ĥ π g exp -

1 β n i=1 [g i -ĥi ] ≥ 0.
(2.27) Our next step is to bound from above the second term at the right-hand side of Equation (2.26). We note that for all h > ĥ , where ĥ is defined by (2.23), r(Y, μh ) -rH (Y ) ≥ 2β σ 2 h 2 -ĥ 2 + 2βσ 2 .

Let q h = exp -r(Y, μh ) -rH (Y ) 2βσ 2 and G = h ∈ H : ĥ 1 ≤ h 1 < ĥ 1 + 1 .

By Condition (1.16) we have that for all h ≥ ĥ q h ≤ exp -2σ 2 βK • ( h 1 -ĥ 1 ) 2βσ 2 -1 = exp -K • ( h 1 -ĥ 1 ) -1 .

So, we obtain with Lemma 2.3, Jensen's inequality, and (2.3)

E µ h∈H w h (Y ) log π h w h (Y ) ≤ E µ log h≤ ĥ π h + exp[R(K • )] ≤ E µ log 1 + ĥ 2 βK • + exp[R(K • )] ≤ log 1 + E µ ĥ 2 βK • + E µ exp[R(K • )] ≤ 2 log 1 + E µ ĥ 2 βK • 1/2 + E µ exp[R(K • )]

h∈HR

  (μ h , µ)is called the oracle risk.Obviously, the oracle estimateµ * (Y ) = h * • Y, where h * = arg min h∈H R(μ h , µ),

Figure 1 :

 1 Figure 1: Mean excess risks ∆β (r H ) (left panel Σ = 5; right panel Σ = 50).

2 ≥

 2 control the right-hand side at (2.26) with the last term. Ordering the elements in H, we obtain by (2.2) log g∈H π g exp -r(Y, μg ) -r(Y, μĥ ) 2βσ log g≥ ĥ π g exp -r(Y, μg ) -r(Y, μĥ ) 2βσ 2

  Theorem 1.3

	By (2.11) we have by
	log[w h (Y )] =	1 2σ 2 β	r(Y,	μĥ ) -	1 2σ 2 β	r(Y, μh ) + log π h
		-log	g∈H	π g exp -	r(Y, μg ) -r(Y, 2βσ 2	μĥ )

  With this equation, bounding from above E µ ĥ 2 with the help of (2.24), we get from (2.26)-(2.28)To finish the proof of the theorem, it remains to minimize the right hand side at this equation in . Assuming ≤ 1/(3β), we obtain So, combining (2.29) and (2.30) with Lemma 2.2, we complete the proof of (1.17) since E µ rH (Y ) ≤ r H (µ).

	where R(•) is defined by (2.13). It is clear that Ψ(0) is bounded from above. It is also easy to check with
	Note that similar to (2.27) it can be checked easily that = 4C/log(x) that for any x ≥ 2
					Ψ(x) ≤	q h ≥ exp -4Cx log(x) +	1 β x log(x) , h ∈ G 2 ≤	Cx log(x)	.
	and thus it follows immediately from the definition of G and (2.4)-(2.5) that
				h∈G	π h q h ≥ exp -	1 β		h∈G	π h ≥	1 2β	exp -	2 β
	and hence		R(K • ) ≤	C 2	+	1 2	log	C .
	E µ	h∈H	w h (Y )r(Y, μh ) ≤ E µ + rH (Y ) + 4βσ 2 log 1 + 2β 1 -2β √ K + exp √ √	1 1 -2β C 2 + 1 2	r H (µ) σ 2 βK • log C .	(2.29)
		√	1 1 -2β	r H (µ) σ 2 βK •	+	√ 1 -2β 1 + 2β	√	K + exp	C 2	+	1 2	log	C
				≤	r H (µ) σ 2 βK •	+	2β 1 -2β √		r H (µ) σ 2 βK •	+	√ 1 + 2β 1 -2β	√ K
													+ exp	C 2	+	1 2	log	C	(2.30)
					≤	r H (µ) σ 2 βK •	+ 3 (1 + 2β)K + 4β	r H (µ) σ 2 βK •
													+ exp	C 2	+	1 2	log	C .	(2.28)
	Let
				Ψ(x) = min ∈[0,1/(3β)]		x +	C exp	C 2	.	1/2 ,

* This work is partially supported by Laboratory for Structural Methods of Data Analysis in Predictive Modeling, MIPT, RF Government grant, ag. 11.G34.31.0073; and RFBR research projects 13-01-12447 and 13-07-12111.