Manuel Loth
email: manuel.loth@inria.fr

Philippe Preux
email: philippe.preux@inria.fr

Samuel Delepoulle
email: delepoulle@lil.univ-littoral.fr

Christophe Renaud
email: renaud@lil.univ-littoral.fr

ECON: a Kernel Basis Pursuit Algorithm with Automatic Feature Parameter Tuning, and its Application to Photometric Solids Approximation

This paper introduces a new algorithm, namely the Equi-Correlation Network (ECON), to perform supervised classification, and regression. ECON is a kernelized LARS-like algorithm, by which we mean that ECON uses an l 1 regularization to produce sparse estimators, ECON efficiently rides the regularization path to obtain the estimator associated to any regularization constant values, and ECON represents the data by way of features induced by a feature function. The originality of ECON is that it automatically tunes the parameters of the features while riding the regularization path. So, ECON has the unique ability to produce optimally tuned features for each value of the constant of regularization. We illustrate the remarkable experimental performance of ECON on standard benchmark datasets; we also present a novel application of machine learning in the field of computer graphics, namely the approximation of photometric solids.

Introduction

We consider supervised learning, classification as well as regression problems: given a set of N examples (x i , y i ) i∈{1,...,N } , x i ∈ D ⊂ R P , we wish to predict the label y(x) of any data x ∈ D as accurately as possible. If y i ∈ R, this is a regression problem, whereas if y i ∈ C with C a discrete set of values, this is a supervised classification task. y i is a noisy measurement of some function y at point x i ; the set of examples is the only knowledge we have about this function y that we try to estimate at any point x ∈ D.

Each data x is supposed to be a real-valued vector, thus x ∈ R P . Then, each data x may be represented in a feature space as (φ 1 (x), ...φ M (x)), where each feature function φ m associate a real number to a data (φ m : D → R, ∀m ∈ {1, ..M }). Each φ m is generally a non-linear combination of the original attributes of the data. Typically, the number of features M is very big 1 . We denote Φ the set of feature functions, that is Φ ≡ {φ 1 , ...φ m , ..., φ M }, with M ≫ P . From a geometrical point of view, this can be seen as representing data in an other, generally huge, space, named the "feature space". Our goal is to build an accurate predictor ŷ, linear in the feature space, of the form w 0 + k=K k=1 w k φ k , which is as sparse as possible, that is, K is as small as possible (where φ k ∈ Φ, and w 0 , w k ∈ R). To this end, we define a loss function L(y(x), ŷ(x)) which quantifies the cost of mis-predicting the label y(x) of a data x; in this paper, we use the quadratic loss function L(y(x), ŷ(x)) ≡ ||y(x)ŷ(x)|| 2 . Accordingly, we measure the error E(ŷ) as the average over the whole domain of definition of the data of this loss: E(ŷ) ≡ D (y(x)ŷ(x)) 2 dx. We use the empirical error to estimate this loss: E emp (ŷ) ≡ i (y iŷ(x i )) 2 . To meet the sparsity expectation, we regularize the problem by minimizing the loss function with the l 1 -norm of ŷ, that is l 1 (ŷ) ≡ k≥1 |w k |. We end-up with the well-known LASSO problem of finding y * :

y * ≡ arg min ŷ∈Y i (y i -ŷ(x i )) 2 + λ k |w k | (1)
where λ ∈ R + is the regularization constant, and Y the set of all biased, linear combinations of the elements of Φ.

The optimal solution to the LASSO problem depends on λ, and we note ŷ(λ) this optimal estimator for a given value of λ. In practice, the value of λ is difficult to set. A common practice to tune it is to try different values, and pick the best. However, efficient algorithms have been proposed to overcome this difficulty by computing all the solutions of the LASSO problem for potentially all values of λ [START_REF] Efron | Leastangle regression[END_REF][START_REF] Guigue | Kernel basis pursuit[END_REF]. These algorithms, LARS and kernel basis pursuit, provide a (noncountable) collection of solutions, for each single value of λ, thus ŷ(λ). The set of all couples (λ, ŷ(λ)) is called the "regularization path".

The features φ m often have parameters themselves; for instance, Gaussian multivariates are often used, for which one has to set the center and the covariance matrix. Setting these parameters is yet an other pratical problem which has also been studied [START_REF] Raykar | Very fast optimal bandwidth selection for univariate kernel density estimation[END_REF][START_REF] Wang | The kernel path in kernelized lasso[END_REF][START_REF] Han | Sequential kernel density approximation and its application to real-time visual tracking[END_REF].

The resolution of both problems at once has not yet been much studied. Rosset et al. [START_REF] Rosset | l1 regularization in infinite dimensional feature spaces[END_REF] have studied this issue from a fundamental point of view, but the practical issues are very loosely tackled. This is precisely the aim of this paper. Hence, we propose an algorithm that builds the regularization path and tunes automatically the feature parameters at the same time. In the sequel, we will first introduce the necessary material on LARS and kernel basis pursuit algorithms in Sec. 2. Then, we present our algorithm named the "equi-correlated network" (ECON) in Sec. 3. Then, we provide an experimental study of ECON based on standard datasets in Sec. 4: there, we show that we obtain "state-ofthe-art" performance on supervised classification, and regression problems, in particular meeting the performance of SVM and boosting. In Sec. 5, we present an original application of this work to the field of computer graphics, a field which is mostly a virgin field of applications for machine learning. Finally, we conclude and discuss future work.

The LARS and Kernel Basis Pursuit

We wish to solve the LASSO problem (cf. eq. ( 1)). The Least-Angle Regression (LARS) algorithm [START_REF] Efron | Leastangle regression[END_REF] solves the LASSO in the case where the data are represented by their original attributes, that is, are not shattered into a feature space. The case where the data are shattered into a feature space is dealt by the kernel basis pursuit (KBP) [START_REF] Guigue | Kernel basis pursuit[END_REF]. Clearly, both LARS and KBP may be seen as essentially the same algorithm 2 , and we use the term LARS as a generic term to designate an algorithm to build the regularization path of a LASSO problem.

Let us denote ζ the function we wish to minimize in the LASSO:

ζ(ŷ) = E emp (ŷ) + λ l 1 (ŷ).
To build the regularization path, the LARS draws from the homotopy idea: we wish to minimize E emp (ŷ) which 2 this being said without wishing to diminish the merits of the authors of the Kernel Basis Pursuit algorithm. is a non trivial task. To the opposite, we know how to minimize l 1 (ŷ): simply take K = 0; then, l 1 (ŷ) = 0, and ŷ is simply ŷ = w 0 . This solution corresponds to λ = +∞; so, ŷ(+∞) = w 0 = 1 N i=N i=1 y i . Starting from this solution, the LARS looks for the largest value of λ for which K is 1, that is the value of λ at which it pays to use one of the features to obtain a smaller error than merely predicting w 0 for any data. At this point, one feature enters the expansion ŷ; we denote this value of λ by λ 1 , and for this value, the estimator has the form: ŷ(λ 1 ) = w 0 + w 1 λ 1 φ i1 , with φ i1 , a certain feature drawn among the set of M (potential) features, which is the one that precisely minimizes ζ; please note that for all values λ > λ i1 , ŷ(λ) = w 0 . This φ i1 enters the set of "active features" A, and leaves the set of potential features P. Then, the second iteration of the LARS starts which will provide the value λ 2 < λ 1 at which a second feature enters the expansion, thus yields:

ŷ(λ 2 ) = w 0 + w 1 λ 2 φ i1 + w 2 λ 2 φ i2 . Please, note that ∀λ 2 < λ ≤ λ 1 , ŷ = w 0 + w(λ)φ i1 , where w(λ) ≡ w i2 + wi 1 -wi 2 λ1-λ2 λ.
As demonstrated in the seminal paper on the LARS [START_REF] Efron | Leastangle regression[END_REF], the search for φ i2 may be done greedily; furthermore, once λ i has been computed for the i th iteration, we can compute the next λ i+1 < λ i ; at the i + 1 th iteration, either a potential feature enters the active set and is added to the expansion, or an active feature leaves the expansion and gets back to the set potential features, which means that in this case, the expansion shrinks. Furthermore, between two values of λ corresponding to two subsequent iterations i and i + 1, the value of the weights of the active features vary linearly; so, it is only necessary to keep track of the values of these weights at points where a feature becomes active, or when an active feature is inactivated. The fact that the weights of active features are not constant led us to use this cumbersome notation w i λ to denote the dependence of its value on the value of λ.

A crucial point here is that each time a feature becomes active, the associated weight is not set so as to correct the error as much as this feature could do it, but so that the next feature to enter the expansion, at the next iteration, will be as much correlated with the residual error as all other active features (thus the "equi-correlation" in ECON name)). This is an important property of the LARS, but tricky, and we refer, again, the interested reader to the seminal paper [START_REF] Efron | Leastangle regression[END_REF].

The LARS may iterate either until λ gets to 0, but this would yield non sparse solutions since at this point, the complexity of the estimator would not be penalized any longer, and all features would be used. In practice, one uses a stopping criterion, such as until the error measured on a test set reaches a minimum.

The computational complexity of the LARS is linear in the number of potential features and quadratic in the number of active features which tends to remain small as will be seen in the experimental section. The LARS is sketched in Algorithm 1 Sketch of the LARS algorithm.

Require: A regression problem defined by N examples (x i , y i ), x i ∈ D ⊂ R P , y i ∈ R. Require: A set of M potential features Φ. λ ← +∞.
Set of active units: A ← ∅. Set of potential units: P ← Φ\A.

Compute the bias:

w 0 ← 1 N i y i Number of active units: K ≡ |A|. Let ŷ ≡ w 0 + k=K k=1 w k φ k . while stopping criterion not fullfilled do
Compute the residual r on the training set (i.e., r i ← y iŷ(x i )).

if the sign of the weight of an active feature has changed at the previous iteration "update active features" step then remove it from A, and put it back into P, and K --; else Best feature selection: select φ K+1 ≡ φ * the feature among P which is most correlated with r.

Update the active features: update the weights of all active features, Include the new feature to the active set: compute its weight w K+1 , add φ K+1 to A, remove it from P. Udpate: Update λ, and set K ++. end if end while Algorithm 13 .

Description of ECON

Since all potential features are tested to enter the expansion at each iteration of the LARS, needless to point out that if the number of potential features is large, each iteration will be costly (even if the complexity is linear in the number of potential units), and in any case, should remain finite. However, the features generally have some hyperparameters, and, as done in the kernel basis pursuit, one has to select a finite set of hyper-parametrizations a priori to apply the LARS. In practice, the choice of this set is not easy. We want to go beyond this limitation, and be able to consider all possible hyper-parametrizations, which typically define a continuous space, for instance of dimension P + P 2 for multivariate Gaussian features. To that end, we propose the equi-correlated network algorithm (ECON) that precisely aims at selecting the best feature to enter the expansion, along with its best hyper-parametrization.

The principle of ECON

The key difference between ECON and the LARS lies in the function that is minimized in order to select the best feature to enter the ŷ at each iteration, and how the minimization is performed. Otherwise, ECON retains the same principle of iterating over λ, starting from +∞, and down to a certain value which corresponds to the stopping criterion used in the algorithm. In the LARS, the set of potential features Φ is finite. Let r ∈ R N the vector of residual of an estimator ŷ, that is the i th component of r is r i = ŷ(x i )-y i . At each iteration, the LARS selects the feature φ * ∈ P ⊂ Φ which maximizes the correlation between φ and r. Since Φ is finite, this is merely done by enumerating P. Now, let us turn to features that have p hyper-parameters. We denote them as θ ∈ Θ ⊂ R p , so that features are really φ θ , where the θ is a free parameter which value has to be set. So, now, we really look for the optimal combination (φ * , θ * ), and there is thus a non-countable number of potential features to consider.

Technically, at the step which selects φ * in Algorithm 1, the LARS solves the following minimization problem:

min{ min i∈1..N λ -x i , r 1 -x i , ∆r + , min i∈1..N λ + x i , r 1 + x i , ∆r + }
where < ., . > denotes the inner product, ∆r is the latest correction to the residual, and:

(x) + ≡ x if x ≥ 0 +∞ if x < 0
Instead of that, ECON solves the following minimization problem:

ξ ≡ min{min θ∈Θ ξ + (θ), min θ∈Θ ξ -(θ)}
where:

ξ + (θ) ≡ λ + φ θ , r 1 + φ θ , ∆r + , ξ -(θ) ≡ λ -φ θ , r 1 -φ θ , ∆r + where φ θ = φ θ (x 1 ), . . . , φ θ (x N )
T . If φ is continuous and differentiable, ξ + and ξ -are continuous and differentiable everywhere except at rare unfeasible points, and ξ is also continuous, and looses differentiability only at the frontiers where arg min(ξ + (θ), ξ -(θ)) changes, that is precisely at each iteration of ECON (which exactly corresponds to where λ is computed in the LARS). Thus, the main task becomes to minimize, at each step, a relatively smooth function of R p (smooth but not convex).

Assuming the features are multivariate Gaussians, we know that ECON is optimal since the results of [START_REF] Rosset | l1 regularization in infinite dimensional feature spaces[END_REF] apply straightforwardly.

Solving the minimization problem

To exemplify things, let us assume without loss of generality, that the features are multivariate Gaussians, that is φ θ=(µ,C) ≡ e (x-µ) T C(x-µ) in which µ ∈ R P , and the covariance matrix C ∈ R P ×P . We have to find µ * and C * that minimizes < φ µ,C , r >. There is no method to solve this minimization problem analytically, so we have to resort to a numerical minimization. Using a numerical minimization instead of an analytical solution of the minimization problem leads to sub-optimal solutions. This seems a strong weakness, so we would like to discuss this point.

We use DiRect [START_REF] Jones | Lipschitzian optimization without the lipschitz constant[END_REF], a global optimization algorithm which is guaranteed to converge asymptotically towards the global optimum of the function 4 .

Regarding the issue of sub-optimal solutions, two cases may happen: DiRect returns a non optimal point which is close, to the optimum; or, DiRect returns a point which is not even close to the global optimum, but from a local one. In the first case, though we do not do it (yet), postprocessing may be done by way of a gradient descent to get to the global optimum. In the second case, this means that ECON makes active a certain feature instead of an other one which is better. Then, the missed feature may be very well found by DiRect in a subsequent iteration of ECON; actually, the odds for this do increase as a feature is missed because the missed optimum will get more and more attractive (in terms of the energy landscape), so that it will quickly be caught; furthermore, the sub-optimal features that have been included instead of the best feature may well be removed in subsequent ECON iterations, since ECON removes active features that are detected as no longer useful. 4 We use the implementation of DiRect available in the NL-OPT library [START_REF] Johnson | The NLopt nonlinear-optimization package[END_REF].

Experimental study

In this section we provide experimental data regarding the performance of ECON. First, we compare ECON with published results of state of the art algorithms on standard supervised classification, and regression datasets. Then, we present an original application to image synthesis.

Methodology

We compare the performance of ECON with published results on support vector machine (SVM) [START_REF] Rätsch | Soft margins for adaboost[END_REF], Relevance Vector Machine (RVM) [START_REF] Tipping | Sparse bayesian learning and the relevance vector machine[END_REF], various boosting algorithms investigated by Rätsch [START_REF] Rätsch | Soft margins for adaboost[END_REF] and the Fisher Kernel Discrimant (KFD) [START_REF] Mika | Fisher discriminant analysis with kernels[END_REF] as reported on the website [START_REF] Rätsch | Benchmark repository[END_REF], the kernel basis pursuit [START_REF] Guigue | Kernel basis pursuit[END_REF], and the LOGREG-LASSO algorithm [START_REF] Roth | Generalizd LASSO[END_REF]. These published results provide an interesting basis of comparison; SVM are well-known to perform extremely well in supervised classification tasks in terms of accuracy, though not being sparse; RVM is meant to be sparser than SVM; boosting algorithms are known to perform very well too for supervised classification; KFD was reported as also very good for classification. In the LARS family which is meant to provide sparse estimators, kernel basis pursuit and LOGREG-LASSO both perform very well on classification, and regression tasks. This all means that an algorithm that reaches the same kind of performance as the best results available in these publications may be considered as performing very well, as long as it uses the same amount of computational resources.

Basically, we split each dataset into a training set, and a test-set. We use the same amount of data in each set as in the cited publications. We perform 100 such splits on each dataset, and thus perform 100 runs on each dataset. This provide performance on which statistics may be performed.

To cope with the published results, we compare the mean-square error; its standard deviation provides a measure of the variability of the precision of predictions; importantly, the number of terms in the estimator is also measured, to be related to the support vectors of SVM, or the number of weak learners in boosting.

In the performance reported below on supervised classification, and regression tasks, we perform a fixed amount of iterations (500), and then retain the best MSE measured on the test set during these iterations. The figures below are averaged over the 100 training-testing splits. On the photometric solids representation task presented in Sec. 5, the stopping criterion has been investigated further, and is thus discussed in this section. 500 iterations is in general far too much iterations to obtain the best test error. Actually, the test error first decreases along the iterations, and then begins to increase. The iteration at which the test error is minimal is the iteration at which to stop, and this can effectively be used as the stopping criterion.

Interestingly, while ECON may in principle add a new feature at each iteration, we observe that the number of features increases sub-linearly, and finally saturates after a certain number of iterations. After that, features keep on being removed and added, keeping the same number of features K in ŷ. This is a very good property of ECON that helps prevent overfitting. The value at which K saturates varies from a dataset to an other; it may remain very small (10 for the Titanic dataset for instance). We think that this value provides some interesting information about the complexity of the dataset; this issue has yet to be investigated more thoroughly.

Supervised classification

ECON is essentially a regression algorithm. To deal with supervised binary classification, we merely encode the class as ±1; to predict the label of a data, we consider the sign of the prediction for this data: sign(ŷ(x)).

Table 1 presents the results: they are striking: on 5 datasets out of 10, ECON obtains the best results in accuracy; ECON improves the best accuracy by 3.5 % on these 5 datasets.

On the other 5 datasets, one (Ringnorm) seems to exhibit an anomaly since there is an order of magnitude between ECON performance and other algorithms; this dataset does have a training set which is much smaller that the test set (15 times smaller), but the situation is the same for the Titanic dataset where ECON performs the best, and Twonorm where the difference in performance is not so large (by far); furthermore, the dimention P of the data is the same for Ringnorm and Twonorm, so that we conclude that this is also not the reason for this order of magnitude. So, let alone this anomalous Ringnorm, on the other 4 datasets on which ECON is not the best, ECON is only beaten by a small factor; on Thyroid, and Twonorm, the performance of ECON is better than that of SVM, and of the same order as boosting (please, refer to [START_REF] Rätsch | Benchmark repository[END_REF] for these figures not reported here).

So despite being an algorithm for regression, and definitely not relying on any principle of margin maximization, ECON obtains "state of the art" performances on supervised classification problems.

Furthermore, ECON provides sparse estimators, as expected. They are much sparser than those provided by SVM, and a bit less sparse than those provided by LOGREG-LASSO, and RVM, but more accurate.

In more details, table 1 provides the following information about classification results on UCI benchmark datasets. For each dataset, we provide the size of the training set N train and the size of the test set N test , as well as the number of attributes of data P , and the results of different algo-rithms. The column entitled "best from [START_REF] Rätsch | Benchmark repository[END_REF]" provides the best results available from Rätsch's data available on his website [START_REF] Rätsch | Benchmark repository[END_REF] and discussed in [START_REF] Rätsch | Soft margins for adaboost[END_REF][START_REF] Mika | Fisher discriminant analysis with kernels[END_REF]; they compared 7 algorithms: RBF-networks with a fixed number of Gaussian units, Adaboost, and variants (LP Reg -AdaBoost, QP Reg -AdaBoost, AdaBoost Reg ), SVM, and Kernel Fisher Discriminant; we only provide the best results among these (highly performing) algorithms. SVM uses Gaussian kernels. For each algorithm, we provide a pattern of the form x(y)/z where x denotes the average MSE, y its standard deviation, and z is the mean number of terms in the best ŷ, that is the best value of K (this figure is provided when it is relevant and available). The best performance is highlighted with boldface font. (There is an annoyance here that the results provided in [START_REF] Rätsch | Benchmark repository[END_REF] are not always compatible with those provided in [START_REF] Roth | Generalizd LASSO[END_REF]; actually, the results provided here show a better accuracy than the latter reference; so, in short, they are more difficult to beat.)

Regression problems

For regression problems, we illustrate ECON performance on a few standard benchmark datasets in this section. Then, we delve into the real application presented in this paper, the representation of photometric solids, into the next section.

Again, we use the same datasets as Rosset [START_REF] Roth | Generalizd LASSO[END_REF], using the same training, and testing sets. Some results are presented in table 2. We compare the performance of ECON with those published in [START_REF] Roth | Generalizd LASSO[END_REF], concerning the Support Vector Machine, the Relevance Vector Machine, and the LOGREG-LASSO. For each algorithm, we provide the accuracy measured on a test set. In the table, boldface font highlights the best performance for each dataset.

We note the very good performance of ECON on Friedmann's F1 and F3 functions. On Abalone, though not the best, ECON performs quite well with regards to other algorithms. On Boston housing, it is a little bit behind SVM and RVM, though still very competitive with most other algorithms and published results. Noticeably, the estimators produced by ECON are quite sparse again.

Application to Photometric Solids representation

In this section, we present and discuss the application of ECON to a real application, originating from computer graphics.

Photometric Solids

In computer graphics, one of the main challenges of the last two decades has been to develop efficient algorithms for One key to these simulations lies in an elaborate description of the light sources as they play the central role of any lighting (figure 1 illustrates the effect of different light sources onto the final scene lighting) 5 . Furthermore, the direction of the light emission has to be carefully taken into account when accuracy is required.

Numerous sorts of lighting devices exist: these range from basic incandescent light bulbs, neon tubes, ... to laser beams; we may also consider particles of dust, or droplets that reflect light, as sources of light. Goniophotometers are used to record the light intensity of the device for a large number of emitting directions [START_REF] Ashdown | Near-Field Photometry: Measuring and Modeling Complex 3-D Light Sources[END_REF]. This set of directions and intensities is called a Photometric Solid (PS). These photometric solids can then be used in place of the accurate description of the entire device (under some assumptions) to obtain the intensity emitted from the device along any direction (see figure 2). 5 Images have been computed with the Indigo photorealistic renderer (http://www.indigorenderer.com).

One of the main drawbacks of PS lies in the large number of sampling directions they require when complex luminaries have to be used: sampling the directions for each degree of angle in spherical coordinates system generates 180 × 360 = 64800 samples that have to be stored and managed during all the simulation process.

Approximating PS with ECON

The current practice is to perform linear interpolation based on measurements. Dealing with several light sources, each being finely sampled, leads to the use of a large amount of memory in any lighting simulation software. We thus studied the use of ECON in order to approximate the PS.

The first step is to provide data to ECON: a subset of the (θ, φ, L) values stored in the PS are provided to the learning algorithm (with (θ, φ) being the lighting direction from the light source expressed in spherical coordinates, and L the light intensity (luminance) along this direction). Classically 80% of the data are used during the learning stage, the remaining 20% being used for error computation. A quadratic error threshold is used for the ANN, convergence being assumed when the error threshold is less than 5%. In ECON, we investigated this issue a lot because the estimator that minimizes the MSE on a test set is not the most satisfying estimator; the reason is that the data are far from being uniformly distributed in the domain (θ, φ) because of the way the measurements are acquired, and that merely minimizing the MSE leads to overfit the data. We have tested various approaches to this model selection problem, in particular Akaike information criterion (AIC). After a lot of inquiry, we found that a good stopping heuristic is based on the second order difference of the MSE at two subsequent iterations, which is a traditional hint to detect a good model (intuitively, we detect an elbow in the learning curve).

ECON provides an estimator ŷ that is used in the lighting software: when an intensity value along a certain lighting direction is required, this direction is input to the estimator ŷ to approximate the PS in this direction.

Some results

The use of ECON solutions has been implemented into a global illumination framework. The results we obtained can be evaluated according to two criterions: on one hand, the absolute error between PS and ECON; on the other hand, we can estimate the quality of the generated images by comparing their content to images generated by more the classical approach. In the end, the second criterion is really what matters here: whether the final image is what we, humans, expect and consider as realistic.

The error should be measured by comparing the values generated by an ECON solution to physical measures. But these measures are generally not available, and only the captured values stored into the PS can be used. In table 3, we report the Root Mean Square Error (RMSE) between a spiral arc interpolation and the approximations provided by a neural network (ANN) [START_REF] Delepoulle | Light source storage and interpolation for global illumination: a neural solution[END_REF], and by ECON. In each case, the difference over all pixels is computed to obtain the RMSE. The reference image use arc spiral interpolation. All images are rendered with path-tracing with estimation of the next event [START_REF] Lafortune | Mathematical Models and Monte Carlo Algorithms for Physically Based Rendering[END_REF]. The RMSE is generally small for both ANNs and ECON solutions. ECON consistently provides a better accuracy than ANN 6 .

From a visual perception criterion, Fig. 3 displays two images that have been computed respectively with spiral arc interpolated PS, and with a solution provided by ECON. The 3D scenes that are visible into these images have been voluntarily chosen as simple as possible in order to allow any observer to compare them easily. It appears that no visible difference is visible between the two images, confirming the interest in the solutions obtained by ECON.

Conclusion and future directions of work

In this paper, we have presented the equi-correlated network (ECON) which is a featurized LARS-like algorithm. Its unique ability is that while riding the regularization path, it automatically tunes the parameters of the features used to represent the data. On standard datasets, we have shown that ECON obtains state of the art performance on supervised classification, and regression problems. Then, we have illustrated its use in a novel field of application in computer graphics, aiming at providing a compact representation of photometric solids.

Regarding ECON itself, we have the feeling that its current implementation may be still slightly improved by studying even more carefully the optimizer that tunes the parameters of the features. We have yet observed large improvements by working on it, and we feel that some more refinements are possible.

From a data mining point of view, as we have noted, the number of terms in the estimator built by ECON saturates after a certain amount of iterations. This is likely to convey interesting information about the complexity of the dataset; the location, and shape of the kernels being used by ECON also provide some information that has to be exploited: ECON does not yield a black box estimator.

From the point of view of computer graphics, the use of ECON for approximating photometric solids appears to provide results that are both numerically, and visually convincing. Furthermore, the memory requirements are considerably reduced as compared to the classical use of PS.

Figure 1 .

 1 Figure 1. Illustration of the effect of the lighting devices. See text for explanations.

Figure 2 .

 2 Figure 2. A 3D representation of two Photometric Solids: each one represents the light distribution from the location of the light device that is symbolized by the red point. Some lighting directions are represented with red arrows.

Figure 3 .

 3 Figure 3. See text for comments.

Table 1 . Supervised classification problems. See text for explanations.

 1

	Dataset Banana	best from [12] 2 LP Reg -AB: 10.73(0.43)/18 10.8/11.4 RVM N train , N test P 4900, 400	LOGREG-LASSO 10.7(0.5)/21.8	ECON 10.1(1.9)/42.5(10.6)
	Breast-cancer	200, 77	9	KFD: 24.77(4.63)	29.9/6.3	26.1(4.6)/17.3	23.3(4.3)/21.3
	Diabetes	468, 300	8	KFD: 23.21(1.63)	NA	23.5(1.9)/7.6	22.7(1.8)/28.0
	Flare-solar	666, 400	9	SVM: 32.43(1.82)	NA	33.3(1.6)/6.4	32.5(1.7)/8.4
	German	700, 300	20	SVM: 23.61(2.07)	22.2/12.5	23.63(2.3)/22.3	23.2(2.1)/72.1
	Heart	170, 100	13	SVM: 15.95(3.26)	NA	16.0(3.1)/10.7	16.2(3.5)/28.3
	Ringnorm	400, 7000	20	KFD: 1.49(0.12)	NA	1.8(0.3)/12.3	10.4(3.0)/53.0
	Thyroid	140, 75	5	KFD: 4.2(2.07)	NA	4.8(2.3)/6.3	4.6(2.3)/21.6
	Titanic	150, 2051	3	SVM: 22.42(1.02)	23.0/93.7	22.9(1.2)/5.0	22.0(0.8)/5.5
	Twonorm	400, 7000	20	KFD: 2.61(0.15)	NA	2.6(0.2)/8.7	2.9(0.4)/56

Table 2 . Regression problems. See text for explanations.

 2

	DATASET	N TRAIN , N TEST	P	SVM	RVM	LOGREG-LASSO	ECON
				AVE. MSE	AVE. MSE	AVE. MSE	AVE. MSE
	FRIEDMAN'S F1	240, 1000	10	2.92/116.2	2.80/59.4	2.84/73.5	2.16/18
	FRIEDMAN'S F2	240, 1000	4	4140/110.3	3505/6.9	3808/14.2	4062.46/10
	FRIEDMAN'S F3	240, 1000	4	0.0202/106.5 0.0164/11.5	0.0192/16.4	0.0139/24
	ABALONE	1000, 2000	8	4.37/972	4.35/12	4.35/19	4.59/63 ± 0.21/22
	BOSTON HOUSING	481, 25	13	8.04	7.46	NA	11.51/87
	simulating light propagation and light/matter interactions.
	Accurate approaches have been studied that now allow both
	computer scientists, and computer artists to generate accu-
	rate simulations, or photorealistic, images [9][6].

Table 3 . RMSE for images rendered with ECON, and a neural network approach (a mul- tilayer perceptron). Lighting of the reference image has been computed using a spiral arc interpolation of the PS. Lighting of each im- age use path tracing algorithm.

 3

	PS dataset	ANN	ECON
	2518T1EF 1919.61 (2.93 %) 1340.89 (2.04 %)
	8013H1EN 2348.77 (3.58 %) 1176.28 (1.79 %)
	8816H1PN 2675.73 (4.08 %) 3164.57 (2.11 %)
	LTL14123 3098.26 (4.73 %) 1722.58 (2.63 %)

the kernel method informed reader should not think that we assume here the existence of a kernel function, and so on: a feature function is merely a non linear combination of the attributes of the data.

We use the following notations: i indexes examples, k indexes active features; a ≡ b means "a is defined by the concept b", whereas a ← b is the assignment operator of procedural programming language, that is: the value of the variable named a gets the current value of b; in this case, the value of a does not change whether the value of b changes later, or not; in the former case (≡), the value of a is always the same as the value of b. Boldfaces, as in x, indicates vectors.

The photometric data originate from the Ledalite (http://www. ledalite.com/), and Lithonia (http://www.lithonia.com/) databases.

Acknowledgment

Part of the experiments were carried out using the Grid'5000 experimental testbed, an initiative from the French Ministry of Research, INRIA, CNRS and RE-NATER and other contributing partners.

* ML acknowledges the support of an INRIA-Nord-Pas de Calais PhD grant during this work.


Olivier Wintenberger

DEVIATION INEQUALITIES FOR SUMS OF WEAKLY DEPENDENT TIME SERIES

Keywords: AMS 2000 subject classifications: Primary 60F99; Secondary 60G10 Exponential inequality, weak dependence, coupling scheme, dynamical systems, Bernoulli schifts

In this paper we give new deviation inequalities of Bernstein's type for the partial sums of weakly dependent time series. The loss from the independent case is studied carefully. We give non mixing examples such that dynamical systems and Bernoulli shifts for whom our deviation inequalities hold. The proofs are based on the blocks technique and different coupling arguments.

Introduction

The aim of this paper is to extend the deviation inequality of Bernstein's type from the independent case to some weakly dependent ones. We consider a sample (X 1 , . . . , X n ) of a stationary process (X t ) in a metric space (X , d). Considering the set F of 1-Lipschitz functions from X to [-1/2, 1/2], we are interested by the deviation of the partial sum S(f ) = n i=1 f (X i ) for any f ∈ F assuming that E(f (X i )) = 0. If the X i are independent and if σ 2 k (f ) = k -1 Var( k i=1 f (X i )), the classical Bernstein inequality gives the deviation estimate, see Bennett [START_REF] Bennett | Probability inequalities for the sum of independant random variables[END_REF]:

(1.1) P S(f ) ≥ 2nσ 2 1 (f )x + x/6 ≤ e -x for all x ≥ 0. This inequality reflects the gaussian approximation of the tail of S(f ) for small values x. And for large values of x, it reflects the exponential approximation of the tail of S(f ). This deviation inequality is very useful in statistics, see for example the monographs of Catoni [START_REF] Catoni | Statistical Learning Theory and Stochastic Optimization[END_REF] and of Massart [START_REF] Massart | Concentration Inequalities and Model Selection[END_REF].

To extend such deviation inequality to the dependent cases, a tradeoff between the sharpness of the estimates and the generality of the context has to be done. Estimates as sharp as in the independent cases (up to constants) are obtained for Markov chains in Lezaud [START_REF] Lezaud | Chernoff-type bound for finite markov chains[END_REF], Joulin and Ollivier [START_REF] Joulin | Curvature, concentration, and error estimates for markov chain monte carlo[END_REF] under granularity. Bertail and Clemenon [START_REF] Bertail | Sharp bounds for the tail of functionals of markov chains. to appear Probability Theory and its applications[END_REF] obtain a deviation inequality for recurrent Markov chains. There exists C > 0 such that for all M > 0 and all x > 0:

P(S(f ) ≥ C( nσ 2 (f ) T x + M x)) ≤ e -x + nP(T 1 ≥ M ),
where the T i are the iid regeneration times and σ 2 T (f ) = E(T ) -1 Var( T 1 i=1 f (X i )). Up to a constant, it is the limit variance in the CLT of S(f ), more natural than σ 2 1 (f ) in (1.1). This primitive estimate of the tail is natural as, through the splitting technique of Nummelin [START_REF] Nummelin | A splitting technique for harris recurrent markov chains[END_REF], the partial sums S(f ) are sums of iid sequences of blocks of size T i . If the regeneration times are bounded, then up to different constants the same estimate than in the iid case is obtained. If the regeneration times admit finite exponential moments, fixing M ≈ ln n Adamczak [START_REF] Adamczak | A tail inequality for suprema of unbounded empirical processes with applications to markov chains[END_REF] obtains estimates of the deviations with a constant C > 0:

P S(f ) ≥ C nσ 2
T (f )x + ln nx ≤ e -x for all x ≥ 0.

A loss of rate ln n, that cannot be reduced, appears in the exponential approximation compared with the iid case, see Adamczak [START_REF] Adamczak | A tail inequality for suprema of unbounded empirical processes with applications to markov chains[END_REF] for more details.

In all these works, the strong Markov property is crucial. To bypass the Markov assumption, one way is to use dependent coefficients. Ibragimov [START_REF] Ibragimov | Some limit theorems for stationary processes[END_REF] introduced the uniformly φ-mixing coefficients. In this settings, Samson [START_REF] Samson | Concentration of measure inequalities for markov chains and φ-mixing processes[END_REF] achieves the deviation inequality (1.1) with different constants as soon as √ φ r < ∞. Less accurate results have been obtained for more general mixing coefficients than the φ-mixing ones: Viennet [START_REF] Viennet | Inequalities for absolutely regular sequences: application to density estimation[END_REF] for absolutely regular mixing and Merlevede et al. [START_REF] Merlevede | Bernstein inequality and moderate deviations under strong mixing conditions[END_REF] for geometrically strongly mixing. Recently, mixing coefficients have been extended to weakly dependent ones, see Doukhan and Louhichi [START_REF] Doukhan | A new weak dependence condition and applications to moment inequalities[END_REF] and Dedecker and Prieur [START_REF] Dedecker | New dependence coefficients: Examples and applications to statistics[END_REF]. Under the exponential decrease of these coefficients, deviation inequalities for S(f ) with a loss in the exponential approximation are given in Doukhan and Neumann [START_REF] Doukhan | A bernstein type inequality for times series[END_REF]. Merlevede et al. [START_REF] Merlevede | A bernstein type inequality and moderate deviations for weakly dependent sequences[END_REF] extends these results for the partial sum S(f ) for unbounded functions f . The dependence context of this paper is the one of the so-called ϕ-weakly dependent coefficients introduced by Rio in [START_REF] Rio | The berry-esseen theorem for weakly dependent sequences[END_REF] to extend the uniformly φ-mixing coefficients. We provide new deviation inequalities for non mixing processes, such that dynamical systems called expanding maps, see Collet et al. [START_REF] Collet | Exponential inequalities for dynamical measures of expanding maps of the interval[END_REF] and continuous functions of Bernoulli shifts. The Bernstein's deviation inequality in these non mixing contexts sharpens the existing ones. The deviation inequality is obtained by dividing the sample (X 1 , . . . , X n ) in different blocks (X i , . . . , X i+k * ), where the length k * must be carefully chosen and then by approximating non consecutive blocks by independent blocks using a coupling scheme.

The coupling scheme follows from a conditional Kantorovitch-Rubinstein duality due to Dedecker et al. [START_REF] Dedecker | Parametrized kantorovich-rubinstein theorem and application to the coupling of random variables[END_REF] and detailed in Section 2. Using this coupling argument, a new deviation inequality can be stated in Section 3:

(1.2) P S(f ) ≥ 5.8 nσ 2 k * (f )x + 1.5 k * x) ≤ e -x for all x ≥ 0, with σ 2 j (f ) = sup j≤k≤n σ 2 k (f ) for all 1 ≤ j ≤ n and k * = min{k ≥ 1; kδ k ≤ σ 2 k (f )}
, where (δ k ) only depends of the ϕ-coefficients, see condition (3.1) for more details. Unlike σ 1 (f ) in (1.1), the variance term σ 2 k * (f ) is natural as it tends to the limit variance in the CLT with k * . When the TLC holds, i.e. σ 2 k (f ) → σ 2 (f ) > 0, then the classical Bernstein's inequality (1.1) holds up to constants with σ 2 1 (f ) replaced by σ 2 (f ), see Subsection 3.3 for more details. On the opposite, if the functionals f n are such that σ 2 1 (f n ) → 0, then for exponentially decreasing rate of the ϕcoefficients, k * ≈ln(σ 2 1 (f n )) and a logarithmic loss appears. As in the recursive Markov chains case, the loss in the exponential approximation depends on the size of blocks k * . We do not know if this loss in the exponential domain may be reduced for such non uniformly φ-mixing sequences.

In many practical examples such that chains with infinite memory introduced by Doukhan and Wintenberger [START_REF] Doukhan | Weakly dependent chains with infinite memory[END_REF], Bernoulli shifts and Markov kernels, an L ∞ coupling scheme is tractable, see Section 5 for a detailed definition. In these specific cases of ϕ-weakly dependent sequences, an improved version of the deviation inequality (1.2) is given in Theorem 5.1:

P S(f ) ≥ 2 nσ 2 k * ′ (f )x + 1.34 k * ′ x) ≤ e -x for all x ≥ 0, with k * ′ = min{1 ≤ k ≤ n / nδ ′ k ≤ kx},
where (δ ′ k ) only depends of the L ∞ coupling scheme, see condition (5.1) for more details. The paper finishes with the proofs collected in Section 6.

Preliminaries: coupling and weak dependence coefficients

Let (X 1 , . . . , X n ) with n ≥ 1 be a sample of random variables on some probability space (Ω, A, P) with value in a metric space (X , d). We assume in all the sequel that for any n ≥ 1 there exists a strictly stationary process (X

(n) t ) such that (X 1 , . . . , X n ) = (X (n) 1 , . . . , X (n) n ). Let us consider F the set of measurable functions f : X → R satisfying: (2.1) |f (x) -f (y)| ≤ d(x, y), ∀(x, y) ∈ X × X and sup x∈X |f (x)| ≤ 1/2.
We denote the partial sum S(f ) = n i=1 f (X i ) and M j = σ(X t ; 1 ≤ t ≤ j) for all 1 ≤ j ≤ n. 2.1. Kantorovitch-Rubinstein duality. The technique of coupling is related with the Kantorovitch-Rubinstein duality. The duality states that given two distribution P and Q on X there exists a random couple

Y = (Y 1 , Y 2 ) with Y 1 ∼ P and Y 2 ∼ Q satisfying E(d(Y 1 , Y 2 )) = sup f ∈Λ 1 E|f (dP -dQ)| = inf Y ′ E(d(Y ′ 1 , Y ′ 2 )),
where Y ′ have the same margins than Y and Λ 1 denotes the set of 1-Lipschitz functions such that |f (x)f (y)| ≤ d(x, y).

Dedecker, Prieur and Raynaud de Fitte [START_REF] Dedecker | Parametrized kantorovich-rubinstein theorem and application to the coupling of random variables[END_REF] extend the classicalKantorovitch-Rubinstein duality in the time series framework by considering it conditionally on some event M ∈ A. Assuming that the original space Ω is rich enough, i.e. it exists a random variable U uniformly distributed over [0, 1] and independent of M, for any Y 1 ∼ P with values in a Polish space it exists a random variable Y 2 ∼ P independent of M satisfying (2.2)

E(d(Y 1 , Y 2 ) | M) = sup{|E(f (Y 1 )|M) -E(f (Y 1 ))|, f ∈ Λ 1 } a.s.
2.2. ϕ-weak dependence coefficients and coupling schemes. Let us recall the weak dependence coefficient ϕ introduced in Rio [START_REF] Rio | The berry-esseen theorem for weakly dependent sequences[END_REF] Definition 2.1. For any X ∈ X , for any σ-algebra

M of A then ϕ(M, X) = sup{ E(f (X)|M) -E(f (X)) ∞ , f ∈ F}.
Another equivalent definition is

(2.3) ϕ(M, X r )) = sup{| Cov(Y, f (X r ))|, f ∈ F and Y is M-measurable and E|Y | = 1},
see [START_REF] Dedecker | New dependence coefficients: Examples and applications to statistics[END_REF].

We will denote by (A) the specific case where X is a Polish space with sup (x,y)∈X 2 d(x, y) ≤ 1. In the case (A) we have ϕ(M, X) = τ ∞ (M, X) where τ ∞ is the coupling coefficient defined in [START_REF] Dedecker | Weak Dependence, Examples and Applications[END_REF] by the relation

τ ∞ (M, X) = sup{ E(f (X)|M) -E(f (X)) ∞ , f ∈ Λ 1 }.
This last coefficients is the essential supremum of the right hand side term of the conditional Kantorovicth-Rubinstein duality (2.2). Thus in the case (A) we get a coupling scheme directly on the variable X via the Kantorovitch-Rubinstein duality (2.2): it exists a version X * ∼ X independent of M such that

E(d(X, X * ) | M) ∞ = τ ∞ (M, X) = ϕ(M, X).
When d is the Hamming distance d(x, y) = 1 1 x =y the coefficient ϕ(M, X) coincides with the uniform mixing coefficient φ(M, σ(X)) of Ibragimov defined for 2 σ-algebra M and M ′ as:

φ(M, M ′ ) = sup M ∈M,M ′ ∈M ′ |P(M ′ | M ) -P(M ′ )|.
In more general context than (A), we have ϕ(M, X) ≤ τ ∞ (M, X) and coupling scheme directly on X is not tractable. Then we do coupling scheme on the variables f (X i ) for some function f : If the sample (X 1 , . . . X n ) is such that the coefficients ϕ(M j , X i ) are finite for 1 ≤ j < i ≤ n and if f ∈ F then the coupling scheme for f (X i ) follows from the conditional Kantorovitch-Rubinstein duality (2.2) and the relation

τ ∞ (M j , f (X i )) ≤ ϕ(M j , X i ) : There exists f (X i ) * such that f (X i ) * ∼ f (X i ) is independent of M j and E(|f (X i ) * -f (X i )| | M j ) ∞ = τ ∞ (M j , f (X i )) ≤ ϕ(M j , X i ).
In the case (A) we also have another possible coupling scheme for f (X i ), see Section 5 for practical examples:

f (X * i ) ∼ f (X i ) is independent of M j and E(|f (X * i ) -f (X i )| | M j ) ∞ ≤ E(d(X * i , X i )| | M j ) ∞ = τ ∞ (M j , f (X i )) = ϕ(M j , X i ).
2.3. Extensions on the product space X q , q > 1. To consider conditional coupling schemes of length q > 1 we need to extend the notions of weak dependence coefficients on X = (X t ) r≤t<r+q ∈ X q . It depends on the metric d q chosen for X q : Definition 2.2. For any q ≥ 1, any X ∈ X q and any σ-algebra M of A let us define the coefficients

ϕ(M, X) = sup{ E(f (X)|M) -E(f (X)) ∞ , f ∈ F q },
where F q is the set of 1-Lipschitz functions with values in [-1/2, 1/2] of X d equipped with the metric d q (x, y) = q -1 q i=1 d(x i , y i ). Let us discuss the consequences of the choice of the metric d q :

• The τ ∞ coupling coefficients on X q are defined for the metric d q , see [START_REF] Dedecker | Weak Dependence, Examples and Applications[END_REF], and for all f ∈ F:

τ ∞ (M, (f (X 1 ), . . . , f (X q ))) ≤ ϕ(M, (X 1 , . . . , X q )).

Moreover, in the case

(A) it holds τ ∞ (M, X) = ϕ(M, X). • If d is the Hamming metric, as d q (x, y) ≤ 1 1 x =y then ϕ(M, X) ≤ φ(M, σ(X)).
Thus the definition of the weakly dependent coefficients ϕ differs here from the one of Rio in [START_REF] Rio | The berry-esseen theorem for weakly dependent sequences[END_REF] where X q is equipped with d ∞ (x, y) = max 1≤i≤q d(x i , y i ).

2.4. First application: deviation inequality of Hoeffding type. This application is due to Dedecker and Prieur [START_REF] Dedecker | New dependence coefficients: Examples and applications to statistics[END_REF]. Assume that X is a Polish space such that sup x,y d(x, y) ≤ 1, i.e. we are in the case (A). Assume that the coefficients ϕ(M j , (X j+1 , . . . , X n )) are finite for all 1 ≤ j ≤ n -1 and that g :

X n → R satisfies |g(x 1 , . . . , x n ) -g(y 1 , . . . , y n )| ≤ n i=1 d(x i , y i ).
As we are in the case (A) it exists a coupling scheme (X * j+1 , . . . , X * n ) of (X j+1 , . . . , X n ) for any 1 ≤ j ≤ n -1 such that, keeping the same notation than in [START_REF] Rio | Ingalités de hoeffding pour les fonctions lipschitziennes de suites dpendantes[END_REF]:

Γ(g) = E(g(X j+1 , . . . , X n ) | M j ) -E(g(X j+1 , . . . , X n ) ∞ = E(g(X j+1 , . . . , X n ) -g(X * j+1 , . . . , X * n ) | M j ) ∞ ≤ (n -j)ϕ(M j , (X j+1 , . . . , X n )).
Applying Theorem 1 of [START_REF] Rio | Ingalités de hoeffding pour les fonctions lipschitziennes de suites dpendantes[END_REF], if E(g(X 1 , . . . , X n )) = 0 then for all x ≥ 0 it holds:

P   f (X 1 , . . . , X n ) ≥ 2 -1 n j=1 (1 + 2(n -j)ϕ(M j , (X j+1 , . . . , X n ))) 2 x   ≤ e -x .
This deviation inequality of Hoeffding type only differs from the one for independence by a constant. However, such inequalities are not as satisfactory as Bernstein ones for statistical applications.

Deviation inequality around the mean inequality

Let us give an inequality for the deviation around the mean of S(f ) = n i=1 f (X i ) for f ∈ F, with (X 1 , . . . , X n ) on the metric space (X , d) and such that there exists a non increasing sequence (δ r ) that satisfies ϕ(M j , , (X r+j , . . . , X 2r+j-1 )) ≤ δ r for all r ≥ 1.

3.1.

A deviation inequality of Bernstein type. Assume with no loss of generality that E(f (X 1 )) = 0.

Theorem 3.1. For any integer n, if there exists (δ r ) as in (3.1) then

P S(f ) ≥ 5.8 nσ 2 k * (f )x + 1.5 k * x ≤ e -x ,
where

k * = min{1 ≤ k ≤ n / kδ k ≤ σ 2 k (f )} and σ 2 k * (f ) = max{σ 2 k (f ) / k * ≤ k ≤ n}.
The proof of this Theorem is given in Subsection 6.1. We adopt the convention min ∅ = +∞ and the estimate is non trivial when rδ r → 0 and nδ n ≥ σ n (f ), i.e. for not too small values of n.

Remark that the variance term σ 2 k * (f ) is more natural than σ 2 1 (f ) in (1.1) as in the central limit theorem σ 2 k * (f ) converges to the limit variance as k * goes to infinity. Before giving some remarks on this Theorem, the next proposition give estimates of the quantity

σ 2 k (f ) = k -1 Var k i=1 f (X i ) .

The variance terms σ 2

k (f ). Under suitable assumptions on (δ r ), it is always possible to obtain rough estimates of σ 2 k (f ) in function of σ 2 1 (f ) and E|f (X 1 )|: Proposition 3.2. If the condtion (3.1) is satisfied then we have for all k ≤ n the inequality:

σ 2 k (f ) ≤ σ 2 1 (f ) + 2E|f (X 1 )| k-1 r=1 δ r .
See Subsection for a straightforward proof 6.2 of this Proposition. The estimate given in Proposition 3.2 can be rough, for example in the degenerate cases when σ 2 k (f ) tends to 0 with k. 3.3. Remarks on Theorem 3.1. The gaussian behavior around the mean is, up to a universal constant, the same than in the iid case with the more natural variance σ 2 k * (f ) instead of σ 2 1 (f ) in (1.1). However in the exponential domain the estimates given in Theorem 3.1 is sometimes less sharp than the one obtained for φ-mixing in Samson [START_REF] Samson | Concentration of measure inequalities for markov chains and φ-mixing processes[END_REF].

In the non degenerate case σ 2 k (f ) → σ 2 (f ) > 0 then k * is finite as soon as rδ r ↓ 0. The deviation inequality of Theorem 3.1 becomes similar than the one in the iid case (1.1) with the variance term σ 2 (f ) instead of σ 2 1 (f ): there exists C > 0 such that for n sufficiently large we have

P(S(f ) ≥ C( nσ 2 (f )x + x)) ≤ e -x for all x ≥ 0.
However, the estimate of the exponential behavior in Theorem 3.1 may differ from the one of the iid case. For example, for statistical issues it is often assumed that f is chosen depending on n such that σ 2 1 (f n ) → 0. Assume that rδ r is summable. Using Proposition 3.2 and Jensen's inequality we have the estimate

σ 2 1 (f n ) σ 2 1 (f n ) 1/2 . If σ 2 1 (f n ) -1/2 nδ n ↓ 0 then for n sufficiently large such that k * n = min{k ≤ n / kδ k ≤ σ 2 1 (f n ) 1/2 } exists, it holds P(S(f ) ≥ C( σ 2 k * n (f n )nx + k * n x)
) ≤ e -x for all x ≥ 0, with C > 0. As k * n ↑ ∞ there is a loss compare with the iid case (1.1). We do not know if this loss may be reduced outside the cases of uniformly mixing processes where (1.1) holds, see Samson [START_REF] Samson | Concentration of measure inequalities for markov chains and φ-mixing processes[END_REF]. This loss may be reduced when the autocorrelations are controlled, choosing a smaller size of blocks k

* n . Assume that σ 2 1 (f n ) σ 2 1 (f n ) (such relation is satisfied in the uniformly φ-mixing context). If σ 2 1 (f n ) -1 nδ n ↓ 0 then for n sufficiently large such that k * n = min{k ≤ n / kδ k ≤ σ 2 1 (f n )} exists, it holds P S(f n ) ≥ C σ 2 1 (f n )nx + k * n x ≤ e -x for all x ≥ 0, with C > 0.
The loss compare with the iid case is due to k * n ↑ ∞. More precisely

• If δ r = Cδ r for C > 0 and 0 < δ < 1 then k * n ≈ -ln(σ 2 1 (f n )), • If δ r = Cr δ for C > 0 and δ > 1 then k * n ≈ σ 2 1 (f n ) 1/(1-δ) .

Examples

We focus on non φ-mixing examples as for them the inequality (1.1) holds up to constants, see Samson [START_REF] Samson | Concentration of measure inequalities for markov chains and φ-mixing processes[END_REF]. We present dynamical systems that are known to be non φ-mixing processes but they satisfy (3.1) in the case (A). Other examples in the case (A) are presented in the Section 5 as a sharpened deviation inequality holds for them, see Theorem 5.1. We also present in this Section continuous functions of Bernoulli shifts that are examples not φ-mixing and not in the case (A) and thus cannot be treated by the approach of Section 5 and of [START_REF] Samson | Concentration of measure inequalities for markov chains and φ-mixing processes[END_REF]. 4.1. Dynamical systems. Here we are in the case X = [0; 1] and d(x, y) = |x -y|, i.e. in the (A) and then ϕ = τ ∞ . Since Andrews [START_REF] Andrews | Nonstrong mixing autoregressive processes[END_REF], dynamical systems, defined as stationary solutions of X t = T (X t+1 ) for all t are classical examples of non-mixing processes. Let us consider X t the stationary solution of

X t = 1 2 (X t-1 + ξ t )
where (ξ t ) is an iid sequence distributed as a Bernoulli(1/2). Then X t = T (X t+1 where T (x) = 2x modulo 1. Even if it is not mixing, easy computation shows that (X t ) satisfies (3.1) with rδ r = (4/9)2 -r (in fact this specific case satisfies also rδ ′ r = (4/9)2 -r , see Section 5 for more details).

More general examples of dynamical systems are studied in Collet et al. [START_REF] Collet | Exponential inequalities for dynamical measures of expanding maps of the interval[END_REF]. They obtain estimates of covariances terms, multivariate versions of (2.3), for dynamical systems called expanding maps. Then it follows the existence of C > 0 and 0 < ρ < 1 such that (3.1) is satisfied with rδ r = Cρ r , see Dedecker and Prieur [START_REF] Dedecker | New dependence coefficients: Examples and applications to statistics[END_REF] for more details.

Continuous functions of Bernoulli shifts.

Let us consider a φ-mixing stationary process (ξ t ) in some measurable space Y and a sequence (U t ) in the metric space X defined as

U t = F (ξ t-j ; j ∈ N),
where F is a measurable function. Assume that the original state space is large enough such that it exists (ξ ′ t ) distributed as (ξ t ) but independent of it. As in [START_REF] Rio | The berry-esseen theorem for weakly dependent sequences[END_REF], assume that there exists a non increasing sequence (v k ) satisfying almost surely d(F (ξ j ; j ∈ N), F (ξ k j ; j ∈ N)) ≤ v k , with the sequence (ξ k t ) satisfying ξ t = ξ k t for 0 ≤ t ≤ k and for t > k, ξ k t = ξ ′ t . Finally set X t = H(U t ) for some measurable function H : X → X and t = {1, . . . , n} and denote

w H (x, η) = sup d(x,y)≤η

d(H(x), H(y)).

Proposition 4.1. The sample (X 1 , . . . , X n ) satisfies (3.1) with

δ r = inf 1≤k≤r-1 {2φ r-k + E(3w H (U 0 , 2v k )) ∧ 1}.
See the Subsection 6.3 for the proof of this Proposition. Remark that by construction the process (X t ) is non necessarily in the case (A) 5. In the case (A) with a coupling scheme in L ∞ .

In all this section we place us in the case (A) where X is a Polish metric space and d(x, y) ≤ 1 for all x, y ∈ X . For all r ≥ 1 a coupling scheme in L ∞ for (X i ) r+j≤i<2r+j-1 , j ≥ 1, exists when we can construct (X * i ) r+j≤i<2r+j-1 distributed as (X i ) r+j≤i<2r+j-1 and independent of M j such that (5.1) sup

1≤j≤n-2r+1 2r+j-1 i=r+j d(X i , X * i ) ≤ rδ ′ r a.
s. for all r ≥ 1.

5.1.

A sharper deviation inequality of Bernstein's type. Remark that condition (5.1) with (δ ′ r ) implies condition (3.1) with δ r = δ ′ r . Then we obtain a slightly sharper deviation inequality than in Theorem 3.1: Theorem 5.1. For f ∈ F such that E(f (X 1 )) = 0 then we have for all x ≥ nδ ′ k and all 1 ≤ k ≤ n:

P(S(f ) ≥ x) ≤ exp - 2nσ 2 k (f ) k 2 h k(x -nδ ′ k ) 2nσ 2 k (f ) where h(x) = (1 + x) ln(1 + x) -x for all x ≥ 0.
Then it holds for all x ≥ 0:

P S(f ) ≥ 2 nσ 2 k * ′ (f )x + 1.34 k * ′ x) ≤ exp(-x) with k * ′ = min{1 ≤ k ≤ n / nδ ′ k ≤ kx}.
The proof of this Theorem is given in Subsection 6.4.

Let us compare this deviation inequality with the result of Theorem 3.1. In Theorem 5.1 the variance term σ 2 k (f ) sharpens σ 2 k (f ) and the normal approximation is better here. For the exponential approximation, in both Theorems losses are due to the chosen blocks sizes. As

k * ′ = min{1 ≤ k ≤ n / kδ ′ k ≤ xk 2 /n}, if kδ ′ k is decreasing as kδ k then k * ′ ≤ k * as soon as nσ 2 k (f ) ≤ k 2 x or equivalently nσ 2 k (f )x ≤ kx, i.e.
as soon as x is in the domain of the exponential approximation. Thus for the normal and the exponential approximations, he deviation inequality in Theorem 5.1 improves the one of Theorem 3.1.

A tradeoff between the generality of the context and the sharpness of the deviation inequalities is done. Even if (5.1) is less general than (3.1), it is satisfied for many examples, see below. [START_REF] Dedecker | New dependence coefficients: Examples and applications to statistics[END_REF], let us consider a stationary Markov chain (X t ) with transition kernel P satisfying, for all f ∈ Λ 1 , that P (f ) = f (y)P (x, dy) is a κ-Lipschitz function with κ < 1. Then [START_REF] Dedecker | New dependence coefficients: Examples and applications to statistics[END_REF] for more details. 5.3. Bounded chains with infinite memory. Let the sequence of the innovations (ξ t ) t∈Z be an iid process on a measurable space Y. We define X = (X t ) t∈Z as the solution of the equation (5.2)

Bounded Markov Chains. Following Dedecker and Prieur

rδ ′ r = κ r (1 + • • • + κ r ), see
X t = F (X t-1 , X t-2 , . . . ; ξ t ) a.s.,
for some bounded function F :

X (N \{0}) × Y → X satisfying the condition (5.3) d(F ((x k ) k∈N \{0} ; ξ 0 ), F ((y k ) k∈N \{0} ; ξ 0 )) ≤ ∞ j=1 a j (F )d(x j , y j ), a.s.
for all (x k ) k∈N \{0} , (y k ) k∈N \{0} ∈ X N \{0} such that there exists N > 0 as x k = y k = 0 for all k > N and with a j (F ) ≥ 0 satisfying

(5.4) ∞ j=1 a j (F ) := a(F ) < 1.
Let (ξ * t ) t∈Z be a stationary sequence distributed as (ξ t ) t∈Z , independent of (ξ t ) t≤0 and such that ξ t = ξ * t for t > 0. Let (X * t ) t∈Z be the solution of the equation

X * t = F (X * t-1 , X * t-2 , . . . ; ξ * t ), a.s.
Using similar arguments than in Doukhan and Wintenberger [START_REF] Doukhan | Weakly dependent chains with infinite memory[END_REF] we have the following result, Lemma 5.2. Under condition (5.4) there exists some bounded (by 1/2) stationary process X solution of the equation (5.2). Moreover, this solution satisfies (5.1) with

rδ ′ r = 2r-1 j=r inf 0<p≤j    a(F ) r/p + ∞ j=p a j (F )    .
As the proof of this Lemma is similar than the one in [START_REF] Doukhan | Weakly dependent chains with infinite memory[END_REF], it is omitted here.

Many solutions of econometrical models may be written as chains with infinite memory. However, the assumption of boundedness is very restrictive for practical models. 5.4. Bernoulli shifts. Solutions of the recurrence equation (5.2) may always be written as X t = H((ξ j ) j≤t ) for some measurable function H : Y N → X were (ξ t ) is an iid process called the innovations. In this very general framework, a coupling version X * t is given by X * t = H((ξ * t )) where (ξ * t ) is a stationary sequence distributed as (ξ t ), independent of (ξ t ) t≤0 and such that

ξ t = ξ * t for t > 0. If there exist a i ≥ 0 such that d(H(x), H(y)) ≤ i≥1 a i d(x i , y i ) with i≥1 a i < ∞,
and if Y is a metric space such that it exists y ∈ Y with d(ξ 1 , y) bounded a.s., then (X t ) satisfies (5.1) with

rδ ′ r = C i≥r a r
for some C > 0.

Proofs

This Section contains the proofs.

6.1. Proofs of the Theorems 3.1. This section contains the proofs of the Bernstein's type estimates on the partial sums S(f ) for f ∈ F. As in the independent case, the proofs follow the Chernoff device. We will proceed using Bernstein's block technique as in [10]. Let us denote by I j the j-th block of indices of size k, i.e. {(j -1)k + 1, jk} except the last blocks and let p be an integer such that 2p -

1 ≤ k -1 n ≤ 2p.
Let us denote by S 1 and S 2 the sums of even and odd blocks defined as

S 1 = i∈I 2j , 1≤j≤p f (X i ) and S 2 = i∈I 2j-1 , 1≤j≤p f (X i ).
From Cauchy-Schwartz inequality, it holds:

ln E[exp(tS(f ))] ≤ 1 2 (ln E exp (2tS 1 ) + ln E exp (2tS 2 )) .
Now let us treat in detail the term depending on S 1 , the same argument applies identically to S 2 . We want to prove that for any 0 ≤ t ≤ 1, choosing k = [1/t] ∧ n as in [10] it holds:

(6.1) ln E(exp(tS(f ))) ≤ 4nt 2 (2(e -2)σ 2 k (f ) + ekδ k ). Denoting L m = ln E(exp(2t i∈I 2j , 1≤j≤m f (X i )))
for any 1 ≤ m ≤ p, we do a recurrence on m remarking that ln E(exp(2tS 1 )) = L p . From Holder inequality, we have for any 2 ≤ m ≤ p -1 the inequalities:

exp(L m+1 ) -exp(L m ) exp(L 1 ) ≤ exp(L m ) E exp 2t i∈I 2(m+1) f (X i ) | M 2mk -E exp 2t i∈I 2(m+1) f (X i ) ∞ ≤ exp(L m ) E exp 2t i∈I 2(m+1) f (X i ) -exp 2t i∈I 2(m+1) f (X i ) * | M 2mk ∞ , where exp 2t i∈I 2(m+1) f (X i ) *
is a coupling version of the variable exp 2t i∈I 2(m+1) f (X i ) , independent of M 2mk . From the definition of the coupling coefficients τ ∞ , we know that

E exp 2t i∈I 2(m+1) f (X i ) -exp 2t i∈I 2(m+1) f (X i ) * | M 2mk ∞ ≤ τ ∞ M 2mk , exp 2t i∈I 2(m+1) f (X i ) .
As i∈I 2(m+1) f (X i ) is bounded with k/2, then u → exp(2tu) is a Lipschitz function with constant 2kt exp(kt) with respect to d k and bounded with exp(kt) for all t ≥ 0. We then deduce that for n -1 < t ≤ 1, choosing k = [1/t] ∧ (n -1) and under condition (3.1) we have

τ ∞ M 2mk , exp 2t i∈I 2(m+1) f (X i ) ≤ 2kte kt ϕ(M 2mk , (X i ) i∈I 2(m+1) ) ≤ 2eδ k .
Collecting this inequalities, we achieve that exp(L m+1 ) ≤ exp(L m )(exp(L 1 ) + 2eδ k ).

The classical Bennett's inequality on i∈I 2 f (X i ) gives the estimates exp(L 1 ) ≤ 1+4σ 2 k (f )/k(e ktkt -1) and as kt ≤ 1 we obtain

L m+1 ≤ L m + ln 1 + 4(e -2)σ 2 k (f ) + 2ekδ k k ≤ L m + 4(e -2)σ 2 k (f ) + 2ekδ k k .
The p steps of the recurrence leads to the desired inequality

ln E(exp(2tS 1 )) ≤ 2p 2(e -2)σ 2 k (f ) + ekδ k k .
As the same inequality holds for S 2 we obtain (6.1) for n -1 < t ≤ 1 remarking that 2pk -1 ≤ 4nt 2 . For t ≤ n -1 , classical Bennett inequality on S 1 gives ln E(exp(2tS 1 )) ≤ 4σ 2 n (f )/n(e ntnt -1).

Remarking that e ntnt -1 ≤ (nt) 2 k≥0 (nt) k /(k + 2)! and (k + 2)! ≥ 23 k we derive that e ntnt -1 ≤ 2 -1 (nt) 2 k≥0 3 -k ≤ 3/4(nt) 2 for nt ≤ 1. Then collecting thes bounds, for t ≤ n -1 it holds ln E(exp(2tS 1 )) ≤ 3nσ 2 n (f )t 2 ≤ 4nt 2 (2(e -2)σ 2 n (f ) + enδ n ). The same holds for S 2 and then (6.1) follows for 0 ≤ t ≤ n -1 and then for all 0 ≤ t ≤ 1.

Note that for k ≥ k * we have σ 2 k (f ) ≤ σ 2 k * (f ) and kδ k ≤ σ 2 k (f ) by definition. From (6.1) we achieve ln E(exp(tS(f ))) ≤ Knσ 2 k * (f )t 2 , for 0 ≤ t ≤ k * -1 , with K = 4(3e-4). Follow the Chernoff's device, i.e. using ln P(S(f ) ≥ x) ≤ ln E(exp(tS(f )))-tx and optimizing in 0 ≤ t ≤ k * -1 , we obtain

P(S(f ) ≥ x) ≤ exp - x 2 2Knσ 2 k * (f ) 1 1 k * x≤2Knσ 2 k * (f ) + exp Knσ 2 k * (f ) k * 2 - x k * 1 1 k * x>2Knσ 2 k * (f ) .
Easy calculation yields for all x ≥ 0

P(S(f ) ≥ 2Knσ 2 k * (f )1 1 k * 2 x≤2Knσ 2 k * (f ) + (k * t + k * -1 Knσ 2 k * (f ))1 1 k * 2 x>2Knσ 2 k * (f ) ≤ e-x. A rough bound k * t + k * -1 Knσ 2 k * (f ) ≤ 3k * x/2 for k * 2 x > 2Knσ 2 k * (f )
k i=1 f (X i ) = k Var(f (X 1 )) + 2 k-1 r=1 (k -r) Cov(f (X 1 ), f (X r+1 )).
Now let us consider the coupling scheme f (X r+1 ) * distributed as f (X r+1 ) but independent of M 1 . Then from Holder inequality it holds

Cov(f (X 1 ), f (X r+1 )) = E(E(f (X r+1 ) -f (X r+1 ) * | M 1 )f (X 1 )).
But as f (X r+1 )f (X r+1 ) * ≤ δ r conditionally to M 0 we get the desired result. 6.3. Proof of Proposition 4.1. We adapt the proof of [START_REF] Rio | The berry-esseen theorem for weakly dependent sequences[END_REF]. We are interested in estimated the coefficients ϕ(M j , (X r+j , . . . , X 2r-1+j )) for any (j, r) satisfying 1 ≤ j ≤ j + r ≤ 2r -1 + j ≤ n. Let us fix (j, r) and denote (ξ k t ) a sequence such that ξ k t = ξ t for all t ≥ r + jk > j and

ξ k t = ξ ′ t otherwise. Denote U k t = F (ξ k t-j ; j ∈ N) and X k t = H(U k t ). For any f ∈ F, we have (6.2) f (X r+j , . . . , X 2r-1+j ) -f (X k r+j , . . . , X k 2r-1+j ) ≤   1 r 2r-1+j i=r+j d(X i , X k i )   ∧ 1.
By definition of the modulus of continuity and as d(

U k i , U i ) ≤ v k for any r + j ≤ i ≤ 2r -1 + j, we have d(X i , X k i ) = d(H(U i ), H(U k i )) ≤ w H (U k i , v k ).
Remarking that r -1 2r-1+j i=r+j w H (U k i , v k ) ∧1 is a measurable function of ((ξ ′ t ) t<r+j-k , (ξ t ) t≥r+j-k ) bounded by 1, it holds from the definition of the φ-mixing coefficients:

E   r -1 2r-1+j i=r+j w H (U k i , v k ) ∧ 1 / M j   ≤ φ r-k + E   r -1 2r-1+j i=r+j w H (U k i , v k ) ∧ 1   . Using again that d(U k i , U i ) ≤ v k , then w H (U k i , v k ) ≤ 2w H (U i , 2v k )
. By stationarity of (U t ), we obtain

E   r -1 2r-1+j i=r+j w H (U k i , v k ) ∧ 1   ≤ E(2w H (U 0 , 2v k )) ∧ 1.
So combining these inequalities we obtain for all 1 ≤ k ≤ r -1:

(6.3) E f (X r+j , . . . , X 2r-1+j ) -f (X k r+j , . . . , X k 2r-1+j ) M j ∞ ≤ φ r-k + E(2w H (U 0 , 2v k )) ∧ 1.
Using again the definition of the φ-mixing coefficients, as f is bounded by 1 it holds

(6.4) E f (X k r+j , . . . , X k 2r-1+j ) M j -E f (X k r+j , . . . , X k 2r-1+j ) ∞ ≤ φ r-k .
Finally, using again (6.2) and that d(X i , X k i ) ≤ w H (U i , v k ), by stationarity of (U t ) we obtain (6.5)

Ef (X r+j , . . . , X 2r-1+j ) -Ef (X k r+j , . . . , X k 2r-1+j ) ≤ E(w H (U 0 , v k )) ∧ 1. The result of the Proposition 4.1 follow from the definition of the ϕ-coefficients, the inequalities (6.3), (6.4) and (6.5). To prove (6.6), let us use the L ∞ -coupling scheme and (5.1) to derive for all 1 ≤ m ≤ p:

i∈I 2m f (X i ) - i∈I 2m f (X * i ) ∞ ≤ i∈I 2(m+1) d(X i , X * i ) ∞ ≤ kδ ′ k ,
where, as in Subsection 6.1, |I j | = k for all 1 ≤ j ≤ 2p with 2p -1 ≤ nk -1 ≤ 2p. Then, for all t ≥ 0 we have:

exp   2t i∈I 2m f (X i )   ≤ e 2tkδ ′ k exp   2t i∈I 2m f (X * i )   a.s.
for all 1 ≤ m ≤ p. In particular, by independence of (X * i ) i∈I 2m with M 2i-1 and by stationary we deduce that and the inequality (6.6) follows remarking that 4pk -1 ≤ 2nk -2 and 2(p -1)k ≤ n.

E   exp   2t i∈I 2m f (X i )   | M 2(m-1)   ≤ e 2tkδ ′ k E   exp   2t
For the Bernstein's type inequality, we use (6.6), the series expansion of the function exp(x)x -1 and that k! ≥ 23 k-2 for k ≥ 2 to derive: ln(E(exp(tS(f ))) ≤ nσ 2 k (f )t 2 1 -(k/3)t + nδ ′ k t for all t ≥ 0.

With the same notation than in [START_REF] Massart | Concentration Inequalities and Model Selection[END_REF], for x ≥ nδ ′ k the Chernoff device leads to:

P(S(f ) ≥ x) ≤ exp 2nσ 2 k (f ) (k/3) 2 h 1 (k/3)(x -nδ ′ k ) 2nσ 2 k (f )
, where h 1 (x) = 1 + x -√ 1 + 2x for all x ≥ 0. Then for all x ≥ 0 we have

P(S(f ) ≥ x + nδ ′ k ) ≤ exp 2nσ 2 k (f ) (k/3) 2 h 1 (k/3)x 2nσ 2
k (f ) and the desired result follows as h -1 1 (x) = √ 2x + x for all x ≥ 0.

6. 4 .

 4 Proof of Theorem 5.1. Let us keep the same notation than in the proof of Theorem 3.1. The Benett's type deviation inequality follows classically from the Chernoff device applies with the estimate: (6.6) ln(E(exp(tS(f ))) ≤ 2nσ 2 k (f ) k 2 (exp(kt)kt -1) + nδ ′ k t for all t ≥ 0.

1 ≤

 1 m ≤ p. Applying this inequality for m = p we have E exp(2tS 1 ) = E Let us do the same reasoning recursively on m = p -1, . . . , 2 to obtain finally ln E exp(2tS 1 ) ≤ 2(p -1)kδ ′ k t + p ln E kt)kt -1)

  leads to the result of the Theorem.6.2. Proof of Proposition 3.2. We have the classical decompositionVar

Acknowledgments. The author is grateful to Jérôme Dedecker for his helpful comments.


S Kawka

B A Van Tiggelen

Quantum Electrodynamics of Casimir Momentum: Momentum of the Quantum Vacuum?

The electromagnetic vacuum is known to have energy. It has been recently argued that the quantum vacuum can possess momentum, that adds up to the momentum of matter. This "Casimir momentum" is closely related to the Casimir effect, in which case energy is exchanged. In previous theory it was treated semi-classically. We present a non-relativistic quantum theory for the linear momentum of electromagnetic zero-point fluctuations, considering an harmonic oscillator subject to crossed, quasi-static magnetic and electric and coupled to the quantum vacuum. We derive a contribution of the quantum vacuum to the linear pseudo-momentum and give a new estimate for the achievable speed. Our analysis show that the effect exists and that it is finite.

Casimir energy refers to the electromagnetic (EM) energy that shows up when dielectric or metallic objects interact with the quantum vacuum. It is undoubtedly one of the most fascinating phenomena in physics, with a rich history in the 20Th century. Casimir forces become important on sub-micron scales and are thus believed to play an important role in nano-optics [1]. Casimir energy has been the subject of many speculations, such as its role in sonoluminescence [START_REF] Schwinger | Proc. Natl. Acad. Sci[END_REF] or in the cosmological constant problem [START_REF] Birell | Quantum Fields in Curved Space[END_REF].

The standard Casimir effect refers to the reduction of the EM zero-point energy when two ideal metallic plates approach [START_REF] Casimir | [END_REF]. Other well-known phenomena related to Casimir energy are Van Der Waals and Casimir-Polder forces between neutral atoms [5], the Lifshitz forces between dielectric media, and arguably the most famous among all, the Lamb shift of atomic levels. Shortly after its observation by Lamb in 1947 [6], Bethe explained the Lamb shift by the change in EM vacuum energy caused by the interaction of the atom with the quantum vacuum [7,[START_REF] Milonni | The Quantum Vacuum[END_REF]. The Lamb shifts in light atoms are now understood to be basically nonrelativistic QED phenomena, although full relativistic theory, including the contribution of several percents due to vacuum polarization, is necessary to come to the extraordinary agreement with experiment, unprecedented in physics. For the two-photon 1S-2S transition in atomic hydrogen, the shift is known up to several cycles [START_REF] Zimmermann | [END_REF].

Energy and momentum are naturally related by relativity. The search for "Casimir momentum" seems therefore obvious. In 2004 Feigel [10] proposed a quantum correction to the momentum of dielectric media exposed to static electric and magnetic fields. In this case, classical electrodynamics provides the following expression for the linear momentum of a neutral, polarizable object with mass M ,

Q class = M v -α(0)E 0 × B 0 (1)
which is conserved in time, even if the external electric field E 0 is varied slowly in time. Here α(0) is the static polarizability, with the dimension of a volume. The semiclassical theory of Ref. [10] predicts a strongly diverging contribution of the quantum vacuum to Eq. ( 1), quite similar to the one encountered for Casimir energy. Fortunately, spatial gradients of Casimir energy -observable as forces -are often found not to diverge. Momentum however is an observable parameter and the divergence does pose a problem. It has been suggested that UV divergences are not physical and should disappear into the values attributed to physical observables, such as inertial mass, electric charge or cosmological constant [11]. If this is true it is not evident that the prediction of "Casimir momentum" found in Ref. [10] will survive or be measurable. An obvious next question is what physical observable will then absorb the UV divergence of Casimir momentum. In this work we provide first answers to these questions. We use the method of mass renormalization first employed by Bethe and Kramers that results in a finite Casimir momentum of simple quantum objects.

It is now realized that Casimir momentum emerges quite generally in so-called bi-anisotropic media, in which also magnetic fields can induce an electric polarization [12,13]. Except in media exposed to external EM fields, bi-anisotropy also occurs in moving dielectric media. This follows from the relativistic transformations of EM fields, and shall be discussed elsewhere [14]. In general, like spin, bi-anisotropic behavior and Casimir momentum can be viewed as "remnants" of special relativity [15] in non-relativistic theory that often suffices to describe phenomena quantitatively.

We consider here the following system: a 3D harmonic oscillator -composed of two particles with opposite charge q 1 = +e and q 2 = -e and masses m i -exposed to crossed, homogeneous static EM fields E 0 , B 0 which constitute our bi-anisotropic object, coupled to the EM vacuum. In the Schrödinger picture the Hamiltonian is given by

H = 2 i=1 1 2m i p i -q i A t (r i ) 2 -q i E 0 .r i + V (r) + kǫ ω k a † kǫ a kǫ + 1 2 (2)
where A t = A 0 + A is the total vector potential containing a contribution from the external, static, classical magnetic field, described by the classical vector potential A 0 (r) = 1 2 B 0 × r, and the quantum operator A of the EM field. The EM bath will be treated in the Coulomb gauge. We will use R = (m 1 r 1 +m 2 r 2 )/M and r = r 1 -r 2 for the center of mass position and the inter-particle distance, with conjugate momenta P and p, respectively; M = m 1 + m 2 and µ = (1/m 1 + 1/m 2 ) -1 are the total and reduced mass. For a harmonic oscillator we can then write V (r) = 1 2 µω 2 0 r 2 . As we are looking for momentum, we notice the existence of a conserved pseudo-momentum K that commutes with H [START_REF] Cohen-Tannoudji | Photons et atomes[END_REF], even when the electric field is varied in time. The electric field is a parameter that can be varied experimentally. This momentum has contributions from both atom and radiation,

K = P + e 2 B 0 × r + kǫ k a † kǫ a kǫ + 1 2 (3)
It is a pseudo-momentum as it is a constant of the motion only if the external magnetic field B 0 is time independent and homogeneous, as we will assume here. Note that we are interested in the change of the total kinetic momentum P kin of the oscillator in the presence of vacuum, which is not equal to K, neither to P. Yet, because K is conserved, even for a slowly time-dependent electric field E 0 , it is the appropriate momentum to look at. In the presence of magnetic fields the kinetic and conjugate momenta operators are related by P kin = P-e∆A t with e∆A t = eA t (r 1 ) -eA t (r 2 ), so that we can obtain

K = P kin + eB 0 × r + e∆A + kǫ k a † kǫ a kǫ + 1 2 (4)
The operator e∆A = eA(r 1 ) -eA(r 2 ) guarantees the gauge-invariant contribution of the "longitudinal" vacuum field to the pseudo-momentum, in terms of the vector potential A(r) quantized as usual inside a quantization volume V , A(r) = kǫ A k ǫ a kǫ e ikr + a † kǫ e -ikr . We will be obliged to go beyond the electric dipole approximation -in which e∆A would be neglected-to treat the high wave numbers of zero-point fluctuations accurately. The last term in Eq. ( 4) stems from the "transverse" electromagnetic field in the vicinity of the atom [START_REF] Cohen-Tannoudji | Photons et atomes[END_REF].

We wish to express the expectation value K = Ψ 0 |K|Ψ 0 , of the pseudo-momentum in the total ground state |Ψ 0 in terms of the two relevant vectors: the kinetic momentum M v of the oscillator and the magneto-electric vector E 0 × B 0 . We will use perturbation theory in the coupling between the magnetoelectric oscillator and the EM field, where the small parameter of the expansion is the fine structure constant α. To lowest order, only the emission and subsequent re-absorption of one virtual photon will contribute at this one-loop level of the theory, which will thus be second order perturbation theory. In order to facilitate perturbation theory, we split the Hamiltonian (2) up as follow,

H = H 0 + H F + W (5)
H 0 = 2 i=1 1 2m i p i -q i A 0 (r i ) 2 -q i E 0 .r i + 1 2 µω 2 0 r 2 H F = kǫ ω k a † kǫ a kǫ + 1 2 W = 2 i=1 - q i m i (p i -q i A 0 (r i ))A(r i ) + q 2 i 2m i A(r i ) 2
The operator W represents the perturbation of the quantum vacuum on the atom. The term A 2 can be disregarded because it does not couple field and matter, and at this order its contribution to K will vanish due to basic selection rules. The photon field in free space described by H F is well known. Finally, the quantum-mechanics of a 3D harmonic oscillator exposed to crossed, homogeneous static EM fields, described by H 0 , was discussed in detail and non-perturbationally by Dippel etal [16]. H 0 and H F act in different Hilbert-spaces and the basis will be the direct product of their eigenvectors. As we are interested in the center-of-mass motion of the 3D harmonic oscillator, one key point is the distinction between the center-of-mass and the internal coordinates. Without the quantum vacuum, the pseudo-momentum K reduces to,

Q = P kin + eB 0 × r (6)
which commutes with the atomic Hamiltonian. Q is the quantum-mechanical operator corresponding to the classical pseudo-momentum Q class in Eq. ( 1). For the purpose of this work it is convenient to choose eigenfunctions that simultaneously diagonalize the atomic Hamiltonian and the pseudo-momentum Q, labeled by the eigenvalue Q 0 , as Q appears directly in K as seen by inserting Eq. (4) into Eq. ( 6). It will be sufficient to ignore all contributions other than on those linear in either E 0 and B 0 , as indicated by the sign ≍. In this approximation, the magneto-electric oscillator is unitary equivalent to an isotropic harmonic oscillator as expressed by,

|n, Q 0 ≍ exp i Q 0 • R exp - i 2 (B 0 × r) • R exp - i p 0 • r exp - i p • r 0 |φ n (7)
n = (n x , n y , n z ), n i = 0, 1, • • • denotes the quantum levels of the oscillator. The first two exponentials on the right denote translational momentum of the center of mass, and governed by the conjugate momentum P. The last pair of exponentials eliminate the static electric field from the picture, with the eigenfunctions of the oscillator shifted out of the center of mass over a distance r 0 = e -1 α(0)(E 0 + Q 0 × B 0 /M ), and the reduced momentum shifted by p 0 ≍ (2M ) -1 (m 2 -m 1 )α(0)(E 0 ×B 0 ) with α(0) = e 2 /µω 2 0 the static polarizability of the oscillator. This last feature is reminiscent of bi-anisotropic activity, and will generate a dominant contribution to K. Note that it vanishes for m 1 = m 2 . Due to the static magnetic field, the oscillator states |φ n are in principle still anisotropic and even in B 0 . The anisotropy is estimated by the small parameter eB 0 a/( /a) ≈ 10 -5 , with a the atomic size. This anisotropy constitutes corrections nonlinear in the applied fields to the final result for the total momentum. We can therefore neglect it.

Upon taking the quantum-expectation value of Eq. ( 6) for the atomic ground state, which still ignores the quantum vacuum, reveals that the eigenvalue Q 0 is just equal to the classical expression (1). The total energy of the oscillator in the ground state E 0 ≍ 3 2 ω 0 + Q 2 0 /2M is minimal when Q 0 = 0, i.e. for a finite kinetic momentum.

In the absence of the interaction with the quantum vacuum, the eigenstates are just the direct products

|n, Q 0 , n k = |n, Q 0 ⊗ |n k , with unperturbed energies E nQ0n k = E nQ 0 + k ω k (n k + 1 2 ).
Here n k is the occupation of the EM Fock states with photon momentum k. The ground state |Ψ 0 follows from second-order perturbation in the coupling W to the quantum vacuum ( ′ avoids zeros in the denominator),

|Ψ 0 = |0, Q 0 , 0 + lQn ′ W lQn,0Q00 E 0Q00 -E lQn |l, Q, n + lQn ′ sQ ′ m ′ W lQn,sQ ′ m W sQ ′ m,0Q00 (E 0Q00 -E lQn )(E 0Q00 -E sQ ′ m ) |l, Q, n + lQn ′ |W lQn,0Q00 | 2 (E 0Q00 -E lQn ) 2 |l, Q, n
Only the emission and subsequent re-absorption of one virtual photon is considered, which generates a tempo-rary recoil momentum Q 0k of the oscillator,

K = Q 0 (8) +e 2 2Re lkǫ A 2 k ǫ φ 0 | (e ik(r+r0) m 2 M -e -ik(r+r0) m 1 M ) |φ l Ω * l,0 E 0Q00 -E l(Q0-k)1 k +B 0 ×2Re lskǫ e 3 A 2 k φ 0 | r + r 0 |φ l Ω l,s Ω * s,0 (E 0Q00 -E lQ00 )(E 0Q00 -E s(Q0-k)1 k ) +B 0 × lskǫ e 3 A 2 k φ s | r + r 0 |φ l Ω * l,0 Ω 0,s (E 0Q00 -E l(Q0-k)1 k )(E 0Q00 -E s(Q0-k)1 k )
The second term on the right-hand-side is due to e∆A in Eq. ( 4), the two last stem from eB 0 × r. We introduced the matrix element Ω l,s = φ l | Ω |φ s of the operator

Ω = eǫ • [B 0 × (r + r 0 )] e ik(r+r0) m 2 M m 1 - e -ik(r+r0) m 1 M m 2 - Q 0 M • ǫ e ik(r+r0) m 2 M -e -ik(r+r0) m 1 M -(p -p 0 ) • ǫ e ik(r+r0) m 2 M m 1 + e -ik(r+r0) m 1 M m 2
and A 2 k = /2ε 0 V kc familiar from quantum optics. Ω stems directly from the development of W applied to a one photon transition.

Three kinds of contributions to K can be identified. The first class is proportional to Q 0 (see the middle term in Ω l,s ) that survives even in the absence of external fields. This term can be seen to affect the inertial mass of the atom, by typically 3 2 ω0 c 2 in accordance with the equivalence principle of energy and inertia [14]. Its UV divergence can be absorbed into the total mass M in the same way as will be discussed below for the ME divergences.

The second class, represented by the two last contributions in Eq. ( 8), are actually QED contributions to the induced electrical dipole moment Ψ 0 |er|Ψ 0 of the oscillator, that find their way to the total momentum via the classical expression (1). It is straightforward to calculate these corrections -they actually do not diverge and they are relatively small -but we note that if an experimental value for α(0) is used to evaluate the "classical" contribution, these terms are automatically included. In this sense they do not constitute a genuine "Casimir momentum".

The term e∆A in Eq. ( 4) is a genuine contribution of the vacuum radiation to the pseudo-momentum. It will be seen to generate a momentum linear in E 0 × B 0 by means of the second term in Eq. ( 8) and the third term in Ω. The following calculation will focus on this contribution. It suffers from a UV divergence, that can be eliminated by exactly the same mass regularization as applied by Bethe in his calculation of the Lamb shift [7]. In particular, this procedure establishes that the reduced mass featuring in the static polarizability α(0) = e 2 /µω 2 0 will be replaced by the "observed" reduced mass. It is straightforward to show that the so-called transverse electromagnetic momentum, represented by the last term of Eq. ( 4), does not contribute a net momentum. This follows from selection rules and spatial symmetry.

We apply the closure relation

l |φ l φ l | (E(k) + E l ) -1 = (E(k) + H ho ) -1
in the second contribution to K. This generates two terms involving exponentials with opposite phases that cancel. The contribution from the first term in Ω to K is

α(0)E 0 × B 0 4α 3π 2 m i kdk 2 k 2 2mi + ck + O ω 0 M c 2 (9)
with i = 1, 2. Here we neglect in the denominator the Doppler terms Q 0 • k and p • k, generated by the identity e ikr H ho (p)e -ikr = H ho (pk). They provide a finite correction of order ω0 Mc 2 ∼ 10 -8 . The leading contribution diverges logarithmically. In the Bethe theory for the Lamb shift [7] exactly the same kind of divergency was encountered. The two diverging mass-like terms δm i = 4α 3π 2 kdk( 2 k 2 2mi + ck) -1 -here with recoil effects included -stem from the QED coupling of the free particles 1 and 2 with the quantum vacuum and are therefore naturally interpreted to be part of their intrinsic, observable masses [START_REF] Milonni | The Quantum Vacuum[END_REF]. Two other diverging contributions are generated by the term α(0)E 0 × B 0 contained in Q 0 , see (1), in the expression for Ω and add up to -δm 1 /Mδm 2 /M = -δM /M . Adding up all diverging terms we obtain

µ -1 (δm 1 /m 1 + δm 2 /m 2 -δM /M ) = µ -2 δµ = -δ(1/µ)
Since the static polarizability is proportional to 1/µ, these UV-divergent corrections all disappear into the factor α(0) of Eq. ( 1), which thus becomes defined in terms of the observed reduced mass µ * :

Q 0 -δ( 1 µ ) e 2 ω 2 0 E 0 × B 0 = M v - e 2 µ * ω 2 0 E 0 × B 0
All other terms generated by Eq. ( 8) are finite. In particular, the term involving p 0 in Ω generates the following contribution to K:

K 1 = α(0)E 0 × B 0 m 2 -m 1 2M 4α 3π × lim δ↓0 ∞ δ dk k k 2 /2 + ckm 2 / - k k 2 /2 + ckm 1 / = -α(0)E 0 × B 0 4α 3π m 1 -m 2 M log m 1 m 2 . ( 10
)
The UV divergency cancels out, and the integral is finite. Of course part of the k-integral enters the relativistic regime k > m i c in which the present theory is not valid. However, by subtracting and adding 2/k to both integrands in Eq. ( 10) reveals two terms whose range of integration is typically m i c/ . This wavenumber was used by Bethe to cut off his nonrelativistic theory for the Lamb shift [7]. We will thus adopt the final result in Eq. ( 10) as a "reasonable" nonrelativistic estimate for the Casimir momentum. It is nevertheless clear that a relativistic theory is required to get the complete picture for Casimir momentum.

All other cross-terms in Eq. ( 8) contain oscillating exponential factors and converge rapidly for k > 1/a, i.e. stay in the nonrelativistic regime . They generate a Casimir momentum that is typically a factor ω 0 /µc 2 ∼ α smaller. We obtain,

K 2 = α(0)E 0 × B 0 α ω 0 µc 2 × - 14 15 √ π + 2 3 √ π ∆m M 2 + 8 3 √ π µ M (11)
When both masses are equal, K 2 becomes the sole contribution to Casimir momentum, with a relative correction of order α 2 . For m 1 ≫ m 2 , K 1 dominates. Since it is independent on details of the force between the two particles, it is tempting to apply Eq. ( 10) to the hydrogen atom. With ω 0 = 10 eV, m 1 = m p and m 2 = m e , E 0 = 10 5 V/m and B 0 = 17 T, we find for the velocity associated with the classical contribution (1) v cl ≈ 5 µm/s, and a QED correction of 2 % in the same direction; K 2 yields a negligible correction of 0.01 %.

In conclusion, we have presented a non-relativistic quantum electrodynamic theory for the total of an harmonic oscillator, subject to external classical fields, and coupled to the electromagnetic quantum vacuum. The most important conclusions of this work are that Casimir momentum exists and that its UV divergences are renormalizable. The theory shows it to be basically a nonrelativistic quantity, but that relativistic corrections are likely to be significant, much like in the Lamb shift problem. To our knowledge nor the classical contribution to magneto-electric momentum, neither the QED correction have ever been observed.

We are indebted to Geert Rikken for many helpful discussions. This work was supported by the contract PHO-TONIMPULS ANR-09-BLAN-0088-01.