Paul Rochet
email: rochet@math.univ-toulouse.fr

A Threshold Regularization Method for Inverse Problems

Keywords: Inverse problems, regularization, oracle inequalities, hard thresholding Subject Class. MSC-2000 : 62G05, 62G08

A number of regularization methods for discrete inverse problems consist in considering weighted versions of the usual least square solution. However, these so-called filter methods are generally restricted to monotonic transformations, e.g. the Tikhonov regularization or the spectral cut-off. In this paper, we point out that in several cases, non-monotonic sequences of filters are more efficient. We study a regularization method that naturally extends the spectral cut-off procedure to nonmonotonic sequences and provide several oracle inequalities, showing the method to be nearly optimal under mild assumptions. Then, we extend the method to inverse problems with noisy operator and provide efficiency results in a newly introduced conditional framework.

Introduction

We are interested in recovering an unobservable signal x 0 , based on noisy observations of the image of x 0 through a linear operator A. The observation y satisfies the following relation y(t) = Ax 0 (t) + ε(t), where ε(.) is a random process representing the noise. This problem is studied in [START_REF] Cavalier | Oracle inequalities for inverse problems[END_REF], [START_REF] Per | The use of the L-curve in the regularization of discrete ill-posed problems[END_REF], [START_REF] Loubes | l 1 penalty for ill-posed inverse problems[END_REF] and in many applied fields such as medical imaging in [START_REF] Natterer | The mathematics of computerized tomography[END_REF] or seismography in [START_REF] Scales | Robust methods in inverse theory[END_REF] for instance. When the measured signal is only available at a finite number of points t 1 , ..., t n , the operator A must be replaced by a discrete version A n : x → (Ax(t 1 ), ..., Ax(t n )) ′ , leading to a discrete linear model

y = A n x 0 + ε,
with y ∈ R n . Difficulties in estimating x 0 occur when the problem is ill-posed, in the sense that small perturbations in the observations induce large changes in the solution. This is caused by an ill-conditioning of the operator A n , reflected by a fast decay of its spectral values b i . In such problems, the least square solution, although having a small bias, is generally inefficient due to a too large variance. Hence, regularization of the problem is required to improve the estimation. A large number of regularization methods are based on considering weighted versions of the least square estimator. The idea is to allocate low weights λ i , or filters, to the least square coefficients that are highly contaminated with noise, thus reducing the variance, at the cost of increasing the bias at the same time.

The most famous filter-based method is arguably the one due to Tikhonov (see [START_REF] Tikhonov | Solutions of ill-posed problems[END_REF]), where a collection of filters is indirectly obtained via a minimization procedure with ℓ 2 penalization. Tikhonov filters are entirely determined by a parameter τ that controls the balance between the minimization of the ℓ 2 norm of the estimator and the residual. Another well spread filter method that will be given a particular interest in this paper, is the spectral cut-off discussed in [START_REF] Bissantz | Convergence rates of general regularization methods for statistical inverse problems and applications[END_REF], [START_REF] Engl | of Mathematics and its Applications[END_REF] and [START_REF] Per | The truncated SVD as a method for regularization[END_REF]. One simply considers a truncated version of the least square solution, where all coefficients corresponding to arbitrarily small eigenvalues are removed. Thus, spectral cut-off is associated to binary filters λ i , equal to 1 if the corresponding eigenvalue b i exceeds in absolute value a certain threshold τ , and 0 otherwise.

A common feature of spectral cut-off and Tikhonov regularization is the predetermined nature of the filters λ i , defined in each case as a fixed non-decreasing function f (τ, .) of the eigenvalues b 2 i , and where only the parameter τ is allowed to depend on the observations. However, in many situations, non-monotonic sequences of filters may provide a more efficient estimation of x 0 . Actually, optimal values for λ i generally depend on both the noise level, which is determined by the eigenvalue b i , and the component, say x i , of x 0 in the direction associated to b i . A restriction to monotonic collections of filters turns out to be inefficient in situations where the coefficients x i are uncorrelated to the spectral values b i of the operator A n .

Regularization methods involving more general classes of filters have also been treated in the literature. In [START_REF] Cavalier | Oracle inequalities for inverse problems[END_REF], the authors study a general procedure known as unbiased risk estimation, that applies to arbitrary classes of filters, dealing in particular with non-monotonic collections. However, their general framework concerning the class of estimators requires in return additional regularity assumptions which we intend to relax in this paper. We focus on a specific class of projection estimators that extends the spectral cut-off to nonmonotonic collections of filters. Precisely, we consider the collection of unrestricted binary filters λ i ∈ {0, 1}. The computation of the estimator relies on the choice of a proper set of coefficients m ⊂ {1, ..., n}, which considerably increases the number of possibilities compared to the spectral cut-off procedure. We show this method to satisfy a non-asymptotic exact oracle inequality, when the oracle is computed in the class of binary filters. Moreover, we show our estimator to nearly achieve the rate of convergence of the best linear estimator in the maximal class of filters, i.e. when no restriction is made on λ i .

It many actual situations, the operator A n is not known precisely and only an approximation of it is available. Regularization of inverse problems with approximate operator is studied in [START_REF] Cavalier | Adaptive estimation for inverse problems with noisy operators[END_REF], [START_REF] Efromovich | On inverse problems with unknown operators[END_REF] and [START_REF] Hoffmann | Nonlinear estimation for linear inverse problems with error in the operator[END_REF]. In this paper, we tackle the problem of estimating x 0 in the situation where we observe independently a noisy version bi of each eigenvalue b i . We consider a new framework where the observations bi are made once and for all, and are seen as non-random. We provide a bound on the conditional risk of the estimator, given the values of bi , in the form of a conditional oracle inequality.

The paper is organized as follows. We introduce the problem in Section 2. We define our estimator in Section 3, and provide two types of oracle inequalities. Section 4 is devoted to an application of the method to inverse problems with noisy operators. The proofs of the results are postponed to the Appendix.

Problem setting

Let (X , . ) be a Hilbert space and A n : X → R n (n > 2) a linear operator. We want to recover an unknown signal x 0 ∈ X based on the indirect observations

y = A n x 0 + ε, (1)
where ε is a random noise vector. We assume that ε is centered with covariance matrix σ 2 I, where I denotes the identity matrix. We endow R n with the scalar product u, v n = n -1 n i=1 u i v i and the associated norm . n and we note A * n : R n → X the adjoint of A n . Let K n be the kernel of A n and K ⊥ n its orthogonal in X which we assume to be of dimension n. The surjectivity of A n ensures that the observation y provides information in all directions. If this condition is not met, one may simply reduce the dimension of the image in order to make A n surjective.

The efficiency of the estimator relies first of all on the accuracy of the discrete operator A n and how "close" it is to the true value A. The convergence of the estimator towards x 0 is subject to the condition that the distance of x 0 to the set K ⊥ n tends to 0, which is reflected by a proper asymptotic behavior of the design t 1 , ..., t n . This aspect is not discussed here, we consider a framework where we have no control over the design t 1 , ..., t n and we focus on the convergence of the estimator towards the projection x † . Let {b i ; φ i , ψ i } i=1,...,n be a singular system for the linear operator A n , that is,

A n φ i = b i ψ i and A * n ψ i = b i φ i and b 2 1 ≥ ... ≥ b 2
n > 0 are the ordered non-zero eigenvalues of the self-adjoint operator A * n A n . The φ i 's (resp. ψ i 's) form an orthonormal system of K ⊥ n (resp.

R n ).
In this framework, the available information on x 0 consists in a noisy version of A n x 0 . As a result, estimating the part of x 0 lying in K n is impossible, based only on the observations. The best approximation of x 0 one can get without prior information is the orthogonal projection of x 0 onto K ⊥ n . This projection, noted x † , is called best approximate solution and is obtained as the image of A n x 0 through the generalized Moore-Penrose inverse operator

A † n = (A * n A n ) † A * n , where (A * n A n ) † denotes the inverse of A * n A n , restricted to K ⊥ n .
By construction, the generalized Moore-Penrose inverse A † n can also be defined as the operator for which {b -1 i ; ψ i , φ i } i=1,...,n is a singular system. We refer to [START_REF] Engl | of Mathematics and its Applications[END_REF] for further details.

Searching for a solution in the subspace K ⊥ n allows to reduce the number of regressors to n. Then, estimating x † can be made using a classical linear regression framework where the number of regressors is equal to the dimension of the observation. Decomposing the observation in the singular basis {ψ i } i=1,...,n leads to the following model

y i = b i x i + ε i , i = 1, ..., n,
where we set y i = y, ψ i n , ε i = ε, ψ i n and x i = x 0 , φ i . It now suffices to divide each term by the known singular value b i to observe the coefficient x i , up to a noise term η i := b -1 i ε i . Equivalently, this is obtained by applying the Moore-Penrose inverse A † n in the model [START_REF] Besov | skiȋ. Integral representations of functions and imbedding theorems[END_REF]. We thus consider the function y † = A † n y ∈ K ⊥ n , defined as the inverse image of y through A n with minimal norm. Identifying y † with the vector of its coefficients y † i = b -1 i y i in the basis {φ i } i=1,...,n , we obtain

y † i = x i + η i , i = 1, ..., n. (2)
The covariance matrix of the noise η = (η 1 , ..., η n ) ′ is diagonal in this model, as we have

E(η i η j ) = n -1 b -1 i b -1 j σ 2 ψ i , ψ j n which is null for all i = j and equal to σ 2 i := σ 2 b -2 i /n if i = j.
Thus, the model can be interpreted as a linear regression model with heteroscedastic noises, the variances σ 2 i being inversely proportional to the eigenvalues b 2 i . In the case where ε in the original model ( 1) is Gaussian with distribution N (0, σ 2 I), the noises η i remain Gaussian in (2). This representation points out the effect of the decay of the singular values b i on the noise level, making the problem ill-posed. To control the noise with a too large variance σ 2 i , a solution is to consider weighted versions of y † . For some filter λ = (λ 1 , ..., λ n ) ′ , note x(λ) ∈ K ⊥ n the function defined by x(λ), φ i = λ i y † i for i = 1, ..., n. Filter-based methods aim to cancel out the high frequency noises by allocating low weights to the components y † i corresponding to small singular values. A widely used example is the Tikhonov regularization, with weights of the form λ i = (1 + τ σ 2 i ) -1 for some τ > 0. The Tikhonov solution can be expressed as the minimizer of the functional

y -A n x 2 + τ x 2 , x ∈ X ,
which makes the method particularly convenient in cases where the SVD of A * n A n or the coefficients y † i are not easily computable. We refer to [START_REF] Cavalier | Nonparametric statistical inverse problems[END_REF] and [START_REF] Tikhonov | Solutions of ill-posed problems[END_REF] for further details. Another common filter-based method is the truncated singular value decomposition or spectral cut-off studied in [START_REF] Bissantz | Convergence rates of general regularization methods for statistical inverse problems and applications[END_REF], [START_REF] Engl | of Mathematics and its Applications[END_REF] and [START_REF] Per | The truncated SVD as a method for regularization[END_REF]. An estimator of x 0 is obtained as a truncated version of y † , where all coefficient y † i corresponding to arbitrarily small singular values are replaced by 0. This approach can be viewed as a principal component analysis, where only the highly explanatory directions are selected. The spectral cut-off estimator is associated to filter factors of the form λ i = 1{i ≤ k}, where 1{.} denotes the indicator function and k is a bandwidth to be determined. Data-driven methods for selecting suitable values of k are discussed in [START_REF] Cavalier | Nonparametric statistical inverse problems[END_REF], [START_REF] Cavalier | Risk hull method and regularization by projections of ill-posed inverse problems[END_REF], [START_REF] Per | The truncated SVD as a method for regularization[END_REF], [START_REF] Varah | On the numerical solution of ill-conditioned linear systems with applications to ill-posed problems[END_REF] and [START_REF] Varah | A practical examination of some numerical methods for linear discrete ill-posed problems[END_REF].

A natural way to generalize the spectral cut-off procedure is to enlarge the class of estimators by considering non-ordered truncated versions of y † , as made in [START_REF] Loubes | l 1 penalty for ill-posed inverse problems[END_REF], [START_REF] Loubes | Adaptive complexity regularization for linear inverse problems[END_REF] or [START_REF] Loubes | Penalized estimators for non linear inverse problems[END_REF] (see also Examples 1 and 2 in [START_REF] Cavalier | Oracle inequalities for inverse problems[END_REF]). This approach reduces to a model selection issue where each model is identified with a set of indices m ⊂ {1, ..., n}. Precisely, for m a given model, define xm ∈ K ⊥ n as the orthogonal projection of

y † onto X m := span{φ i , i ∈ m}, that is, xm satisfies xm , φ i = y † i if i ∈ m, 0 otherwise.
The objective is to find a model m that makes the expected risk E xmx 0 2 small. The computation of the estimator no longer relies on the choice of one parameter k ∈ {1, ..., n} as for spectral cut-off, but on the choice of a set of indices m ⊂ {1, ..., n}, which increases the number of possibilities. In particular, this approach allows non-monotonic collections of filters that may perform better than decreasing sequences obtained by spectral cut-off. To see this, write the bias-variance decomposition of the estimator xm for a deterministic model m:

E xm -x 0 2 = E x 0 -x † 2 + i / ∈m x 2 i + i∈m σ 2 i .
In these settings, it appears that in order to minimize the risk, best is to select indices i for which the component x 2 i is larger than the noise level σ 2 i . A proper choice of filter should depend on both the variance σ 2 i and the coefficient x 2 i . Consequently, the resulting sequence {λ i } i=1,...,n has no reason of being a decreasing function of σ 2 i if some coefficients x 2 i are large enough to compensate for a large variance.

3 Non-ordered variable selection

Threshold regularization

The construction of the estimator by non-ordered variable selection reduces to finding a proper set m. Following the discrepancy principle, an optimal value for m (minimizing the risk) is obtained by keeping small simultaneously the bias term i / ∈m x 2 i and the variance term i∈m σ 2 i in the expression of the risk E xmx 0 2 . Following the previous argument, a minimizer of the risk E xmx 0 2 is obtained by selecting only the indices i for which the coefficient x 2 i is larger than the noise level σ 2

i . An optimal model is thus given by m * := {i :

x 2 i ≥ σ 2 i }.
The coefficients x i being unknown to the practitioner, the optimal set m * can not be computed in practical cases. For this reason it will be referred to as an oracle.

We shall now provide a model m constructed from the available information, that mimics the oracle m * . Fixing a threshold on the coefficients x i being impossible, we propose to use a threshold on the coefficients y † i . Precisely, consider the set

m = i : y †2 i ≥ 4σ 2 i µ i ,
for {µ i } i=1,...,n a sequence of positive parameters to be chosen. Obviously, the behavior of the resulting estimator x m relies on the choice of the sequence {µ i } i=1,...,n : the larger the µ i 's, the more sparse is x m . It must be chosen so that the resulting set m contains only the indices i for which the noise level is small compared to the actual value of x i . Although, the only knowledge of the observations y † i and the variances σ 2 i makes it a difficult task.

There exist general filter-based methods that can be applied to arbitrary classes of filter estimators. One example is the unbiased risk estimation discussed in [START_REF] Cavalier | Oracle inequalities for inverse problems[END_REF], which defines an estimator of x 0 via the minimization of an unbiased estimation of the risk, over an arbitrary set of filters. When restricted to the class of binary filters λ i ∈ {0, 1}, unbiased risk estimation reduces to minimizing over M the criterion

m → y † -xm 2 + 2 i∈m σ 2 i .
The minimum can be shown to be reached for the set m = {i : y †2 i ≥ 2σ 2 i }, which corresponds to taking µ i = 1/2 in our method. This choice is shown to be asymptotically efficient in Proposition 2 in [START_REF] Cavalier | Oracle inequalities for inverse problems[END_REF], although additional restrictions are made on the λ i 's which we intend to relax here in an asymptotic framework. If these conditions are not met, the accuracy of the choice µ i = 1/2 is not clear. We investigate in the next section a different choice for µ i which turns out to be nearly optimal in a general framework.

In a general point of view, the estimator x m can be obtained via a minimization procedure, using a BIC-type criterion for heteroscedastic models,

x m = argmin x∈X y † -x 2 + 4 n i=1 σ 2 i µ i 1{ x, φ i = 0} .
In a certain way, this can be seen as a hard-thresholding version of the estimator considered in [START_REF] Loubes | Penalized estimators for non linear inverse problems[END_REF], obtained with a ℓ 1 penalty. However, expressing the estimator as the solution to a minimization equation does not ease the computation. The method requires in any case calculation of the SVD of A * n A n and the coefficients y † i , which may be computationally expensive. On the other hand, the computation of the estimator is simple once the decomposition of y † in the SVD of A * n A n is known, as it suffices to compare each coefficient y † i to the threshold 4σ 2 i µ i .

Oracle inequalities

In the definition of m, the choice of the parameters µ i is crucial. Too large values of µ i will result in an under-adjustment, keeping too few relevant components y † i to estimate x 0 . On the contrary, a small value of µ i increases the probability of selecting a component y † i that is highly affected with noise. Thus, it is essential to find a good balance between these two types of errors. In the next theorem, we provide a nearly optimal choice for the parameters µ i , under the condition that ε has finite exponential moments.

For i = 1, ..., n, note γ i := η 2 i /σ 2 i = nε 2 i /σ 2 . We make the following assumption.

A1. There exist K, β > 0 such that ∀t > 0, ∀i = 1, ..., n, P(γ i > t) ≤ Ke -t/β .

In a Gaussian model, the γ i 's have χ 2 distribution with one degree of freedom. The condition A1 holds for any β > 2, taking K = 1 -2/β.

Theorem 3.1 Assume that the condition A1 holds. Set

µ i = max{β log(n 2 σ 2 i ), 0}, the estimator x m satisfies E x m -x † 2 ≤ E xm * -x † 2 + (K 1 log n + K 2 ) i∈m * σ 2 i + K 3 n , with K 1 = 12β, K 2 = 2 + β log x † 2 and K 3 = 2Kβ.
Remark 1. This theorem establishes a non-asymptotic oracle inequality with exact constant. The residual term is similar to that in Corollary 1 in [START_REF] Loubes | l 1 penalty for ill-posed inverse problems[END_REF]. The fact that the term x † depends on n is not problematic here as it can in any case be bounded by the norm of x 0 . Remark 2. The method requires knowledge of the operator A n , the variance σ 2 and the constant β in the condition A1. Note however that knowing the constant K is not necessary to build the estimator. Remark 3. The set m contains all indices i for which σ 2 i ≤ 1/n 2 , as we have in this case µ i = 0. This suggests that the error caused by wrongfully selecting indices i for which the variance is smaller than 1/n 2 is negligible, regardless of the value of y † i .

Remark 4. In an asymptotic concern, the accuracy of the result stated in Theorem 3.1 relies on the convergence rate of the residual term to zero, compared to the risk of the oracle. The residual term i∈m * σ 2 i is actually the variance term in the bias-variance decomposition of xm * , and therefore, it is bounded by the risk of the oracle. As a result, the estimator x m is shown to reach at least the same rate of convergence as the oracle, up to a logarithmic term, which warrants good adaptivity properties. The logarithmic term vanishes in the convergence rate if the bias term i / ∈m * x 2 i dominates in the risk of the xm * . Precisely, the oracle inequality is asymptotically exact as soon as the residual term log n i∈m * σ 2 i is negligible compared to the bias term i / ∈m * x 2 i . In this case, it follows from Theorem 3.1 that

E x m -x † 2 = (1 + o(1)) E xm * -x † 2 .
Of course, this condition is hard to verify in practice and assuming it is true reduces to make strong regularity assumptions on the asymptotic behavior of x 0 and A n . In a non-asymptotic framework, the theorem warrants that the estimator x m is close to the oracle as soon as the variance term i∈m * σ 2 i is small compared to the bias term i / ∈m * x 2 i in the bias-variance decomposition of the oracle.

The estimator x m being built using binary filters λ i ∈ {0, 1}, it is natural to measure its efficiency by comparing its risk to that of the best linear estimator in this class. Nevertheless, we see in the next corollary that a similar oracle inequality holds if we consider the oracle in the maximal class of filters, that is, allowing the λ i 's to take any real value. Corollary 3.2 Assume that the condition A1 holds, the estimator x m of Theorem 3.1 satisfies

E x m -x † 2 ≤ K 4 log n inf λ∈R n E x(λ) -x † 2 + K 5 n ,
for some constants K 4 , K 5 independent of n.

This result is a straightforward consequence of Lemma 5.2 in the Appendix, where it is shown that the oracle in the class of binary filters λ i ∈ {0, 1} achieves the same rate of convergence up to a factor 2, as the best filter estimator obtained with non-random values of λ. The class of unrestricted binary filters leads to a simple solution while it induces a slight loss of efficiency compared to the maximal class.

Interest of oracles lies in the fact that the best estimator in a given class will often reach the optimal rate of convergence. In many situations, comparing the risk of the estimator to that of an oracle might be sufficient to deduce optimality results, as well as adaptivity properties, as discussed in [START_REF] Cavalier | Nonparametric statistical inverse problems[END_REF]. In the literature of inverse problems, rates of convergence of oracles are obtained under regularity conditions on the map x 0 and the spectrum of A n . These conditions can be gathered into a single assumption, generally referred to as source condition, relating the behavior of x 0 to the regularity of the operator A n (see for instance [START_REF] Bissantz | Convergence rates of general regularization methods for statistical inverse problems and applications[END_REF], [START_REF] Engl | of Mathematics and its Applications[END_REF] or [START_REF] Fermín | A statistical view of iterative methods for linear inverse problems[END_REF]). Another point of view widely adopted in the literature is the minimax approach (see [START_REF] Cavalier | Nonparametric statistical inverse problems[END_REF]), aiming to determine the behavior the worst possible value of x 0 in a given class of functions. Typically, the condition can be a polynomial decay of the coefficients x i , which reduces to assuming that x 0 lies in the unit ball in a proper Besov space. For rates of convergence with a minimax approach, we refer to [START_REF] Besov | skiȋ. Integral representations of functions and imbedding theorems[END_REF], [START_REF] Donoho | Neo-classical minimax problems, thresholding and adaptive function estimation[END_REF] and [START_REF] Loubes | Adaptive estimation with soft thresholding penalties[END_REF]. In our framework, rates of convergence for x m can be deduced from Theorem 2 in [START_REF] Loubes | l 1 penalty for ill-posed inverse problems[END_REF], under a polynomial decay of the coefficients x i and the eigenvalues b i .

Regularization with unknown operator

In many actual situations, the operator A n is not precisely known. In this section, we consider the framework where the operator A n is observed independently from y. This situation is treated in [START_REF] Cavalier | Adaptive estimation for inverse problems with noisy operators[END_REF], [START_REF] Efromovich | On inverse problems with unknown operators[END_REF] or [START_REF] Hoffmann | Nonlinear estimation for linear inverse problems with error in the operator[END_REF]. The method discussed in the previous section does not apply for such problems since it requires complete knowledge of the operator A n . As in [START_REF] Cavalier | Adaptive estimation for inverse problems with noisy operators[END_REF], we assume that the eigenvectors φ i and ψ i are known. This seemingly strong assumption is actually met in many situations, for instance if the problem involves convolution or differential operators which can be decomposed in Fourier basis (see also the examples in [START_REF] Cavalier | Nonparametric statistical inverse problems[END_REF]). Thus, only the eigenvalues b i are unknown and we assume they are observed independently of y, with a centered noise ξ i with known variance

s 2 > 0: bi = b i + ξ i , i = 1, ..., n.
The method discussed in this paper is different according to whether the eigenvalues are known exactly or observed with a noise. Thus, we need to assume here that s is positive and the known operator framework can not be seen as a particular case. Moreover, we assume the ξ i 's are independent and satisfy the two following conditions. A2. There exist K ′ , β ′ > 0 such that ∀t > 0, ∀i = 1, ..., n, P(ξ 2 i /s 2 > t) ≤ K ′ e -t/β ′ .

A3. There exist C, α > 0 such that ∀i = 1, ..., n, min{P(ξ i < -αs), P(ξ i > αs)} ≥ C.

As discussed previously, the condition A2 means that that the ξ i 's have finite exponential moments. The condition A3 is hardly restrictive, and is fulfilled for instance as soon as the ξ i 's are identically distributed. As we shall see in the sequel, the method requires knowledge of the constant α (or at least an upper bound for it), but no information on the constants β ′ , K ′ or C is needed to build the estimator.

Knowing the eigenvectors of A * n A n allows us to write the model in the form

y i = b i x i + ε i , i = 1, ..., n.
In our framework where the actual eigenvalues b i are unknown, a natural estimator of each component x i is obtained by ỹi = b-1 i y i , provided that bi = 0. However, it is clear that this estimate is not satisfactory if bi is far from the true value (consider for instance the extreme case where bi = 0 or if bi and b i are of opposite signs). Actually, the naive estimator b-1 i can not be used efficiently to estimate b -1 i because it may have an infinite variance. In [START_REF] Cavalier | Adaptive estimation for inverse problems with noisy operators[END_REF], the authors fix a threshold w the estimate can not exceed and consider an estimator of b -1

i equal to b-1 i if | bi | > 1/
w and null otherwise. As we will see below, we use the same idea here, although the threshold fixed on the bi 's is implicitly part of the variable selection process.

We can reasonably assume that null values of bi do not provide any relevant information and can not be used to estimate x 0 . Thus, to avoid considering trivial situations, we assume that all bi are non-zero. In all generality, the ỹi 's can be viewed as noisy observations of x i by writing ỹi = x i + ηi , i = 1, ..., n,

with ỹi = b-1 i y, ψ i n and ηi = b-1 i (ε i -ξ i x i ),
where we recall ε i = ε, ψ i n . As in the previous section, we propose a threshold procedure to filter out the observations ỹi that are potentially highly contaminated with noise. Here, the noise ηi is more difficult to deal with because it depends on the unknown coefficient x i .

Our objective is to find an optimal variable selection criterion conditionally to the bi 's. In order to do so, we consider a framework where the bi 's are observed once and for all, and are treated as non-random. Thus, we define as an oracle, a model m * ξ minimizing the conditional risk E ξ xmx † 2 , where E ξ (.) denotes the expectation knowing ξ = (ξ 1 , ..., ξ n ) ′ . Following a similar argument as in the previous section, a model minimizing the conditional risk contains only the indices i for which the coefficient x 2 i is larger than the noise level. Hence, we may define m * ξ = {i :

x 2 i > E ξ (η 2 i )}.
A notable difference here is that the noise ηi actually depends on the value

x i . Let σ2 i = n -1 b-2 i σ 2 , we can calculate the conditional expectation of η2 i , given by E ξ (η 2 i ) = σ2 i + b-2 i ξ 2 i x 2 i .
After simplifications, it appears that the optimal model conditionally to the ξ i 's can be expressed in the two following explicit forms

m * ξ = i : 2| bi | > σ 2 n|b i |x 2 i + |b i | = i : x 2 i > σ 2 n( b2 i -ξ 2 i ) , | bi | > |b i | 2 .
In the first expression, we see that the oracle selects indices i for which the observation bi exceeds a certain value depending on both x i and b i . Interestingly, components ỹi corresponding to observations bi smaller than half the true eigenvalue b i are not selected in the oracle, regardless of the coefficient x i . Here again, the optimal model m * ξ can not be used in practical cases since it involves the unknown values x i and ξ i . We can only try to mimic the optimal threshold, based on the observations ỹi and bi . Consider the set

m ξ = i : ỹ2 i > 8σ 2 i ν i , | bi | > αs ,
where {ν i } i=1,...,n are parameters to be chosen and α is the constant defined in A3. With this definition, only the indices for which the observation bi is larger than a certain value, namely αs, are selected. This conveys the idea discussed in [START_REF] Cavalier | Adaptive estimation for inverse problems with noisy operators[END_REF], that when b i is small compared to the noise level, the observation bi is potentially mainly noise. Remark however that in [START_REF] Cavalier | Adaptive estimation for inverse problems with noisy operators[END_REF], the lower limit for the observed eigenvalues is s log 2 (1/s), while in our method, it is chosen of the same order as the standard deviation s.

Define the set M = {i : |b i | < 2αs}.

Theorem 4.1 Assume that the condition A1 holds. The threshold estimator obtained with ν i = max{β log(n 2 σ2 i ), 0} satisfies,

E ξ x m ξ -x † 2 ≤ (K ′ 1 log n + K ′ 2 ) E ξ xm * ξ -x † 2 + i∈M x 2 i + κ(ξ), with K ′ 1 = max{18β, 4α -2 β ′ }, K ′ 2 = max{9(β log x † 2 + 1), 1}, and
κ(ξ) = 4Kβ n + 4 i / ∈m * ξ ξ 2 i x 2 i α 2 s 2 1{ξ 2 i > s 2 β ′ log n}. Moreover, if A2 holds, E(κ(ξ)) = O(n -1 log n).
The main interest of this result lies in the fact that it provides an oracle inequality, conditionally to the bi 's. In particular, the conditional oracle xm * ξ is more efficient than the estimator obtained by minimizing the expected risk m → E xmx † 2 , since the optimal set m * ξ is allowed to depend on the ξ i 's. We see that the estimator x m ξ performs almost as well as the conditional oracle. Indeed, the residual term κ(ξ) is independent from ξ with high probability, and its expectation is negligible under A2 as pointed out in the theorem. The non-random term i∈M x 2 i is small if the eigenvalues b i are observed with a good precision, i.e. if the variance s 2 is small. Moreover, this term can be shown to be of the same order as the risk under the condition A3.

E x m ξ -x † 2 ≤ K ′ 4 log n E xm * ξ -x † 2 + K ′ 5 log n n ,
for some constants K ′ 4 and K ′ 5 independent from n and s 2 .

With a noisy operator, we manage to provide an estimator that achieves the rate of convergence of the conditional oracle, regardless of the precision of the approximation of the spectrum of A n . Indeed, the constants K ′ 4 and K ′ 5 in Corollary 4.2 do not involve the variance s 2 of ξ. Actually, the variance only plays a role in the accuracy of the oracle. The result is non-asymptotic and requires no assumption on s 2 .

Appendix

Technical lemmas

Lemma 5.1 Assume the condition A1 holds. We have

• E ((η 2 i -x 2 i )1{i ∈ m}) ≤ 2Kβσ 2 i e -µ i /β . • E ((x 2 i -η 2 i )1{i / ∈ m}) ≤ σ 2 i (6µ i + 2).
Proof. Using the inequality (a + b) 2 ≤ 2a 2 + 2b 2 , we find that η

2 i -x 2 i ≤ 2η 2 i -y †2 i /2. By definition of m, we get (η 2 i -x 2 i )1{i ∈ m} ≤ 2σ 2 i (γ i -µ i )1{i ∈ m} ≤ 2σ 2 i (γ i -µ i )1{γ i ≥ µ i },
where we used that X ≤ X1{X ≥ 0}. We finally obtain for all i / ∈ m * ,

E (η 2 i -x 2 i )1{i ∈ m} ≤ 2σ 2 i ∞ 0 P(γ i ≥ t + µ i ) dt ≤ 2Kβσ 2 i e -µ i /β ,
as a consequence of A1. For the second part of the lemma, write

x 2 i -η 2 i = y †2 i -2η i y † i which is bounded by 3y †2 i /2 + 2η 2 i , using the inequality 2ab ≤ 2a 2 + b 2 /2. This leads to E (x 2 i -η 2 i )1{i / ∈ m} ≤ σ 2 i (6µ i + 2). Lemma 5.2 inf m∈M E xm -x † 2 ≤ 2 inf λ∈R n E x(λ) -x † 2 .
Proof. The minimal values of the expected risks can be calculated explicitly in the two classes considered here. Minimizing over R n the function λ → E x(λ)x † 2 , we find that the optimal value of λ i is reached for λ

* i = x 2 i /(x 2 i + σ 2 i ).
On the other hand, we know that m → E xmx † 2 reaches its minimum at m * = {i :

x 2 i ≥ σ 2 i }, yielding inf λ∈R n E x(λ) -x † 2 = n i=1 x 2 i σ 2 i x 2 i + σ 2 i and inf m∈M E xm -x † 2 = i∈m * σ 2 i + i∈m * x 2 i . By definition, if i ∈ m * , 2x 2 i /(x 2 i + σ 2 i ) ≥ 1.
In the same way, 2σ 2 i /(x 2 i + σ 2 i ) ≥ 1, for all i / ∈ m * . We conclude by summing all the terms. Lemma 5.3 Assume the condition A1 holds. We have, for all i = 1, ..., n,

• E ξ ((η 2 i -x 2 i )1{i ∈ m ξ }) ≤ 4Kβ σ2 i e -ν i /β + 4ξ 2 i x 2 i α 2 s 2 . • E ξ ((x 2 i -η2 i )1{i / ∈ m ξ }) ≤ 9σ 2 i ν i + 8E ξ (η 2 i ) + x 2 i 1{| bi | ≤ αs}. Proof. Remark that η2 i = b-2 i (ε i -ξ i x i ) 2 ≤ 2 b-2 i ε 2 i + 2 b-2 i ξ 2 i x 2 i . Using that x 2 i ≥ ỹ2 i /2 -η2 i , we deduce η2 i -x 2 i ≤ 4 b-2 i ε 2 i + 4 b-2 i ξ 2 i x 2 i - ỹ2 i 2 . Writing m ξ = {ỹ 2 i > 8σ 2 i ν i } ∩ {| bi | > αs}, we find (η 2 i -x 2 i )1{i ∈ m ξ } ≤ 4σ 2 i (γ i -ν i )1{γ i ≥ ν i } + 4 b-2 i ξ 2 i x 2 i 1{| bi | > αs},
where we recall that

γ i = nε 2 i /σ 2 . Clearly, b-2 i 1{| bi | > αs} < α -2 s -2
and the result follows using the condition A1. For the second part of the lemma, remark that the complement of

m ξ is {ỹ 2 i ≤ 8σ 2 i ν i , | bi | > αs}∪{| bi | ≤ αs}. Using the inequality x 2 i -η2 i ≤ (1+θ -1 )ỹ 2 i +θη 2 i for θ = 8, we get (x 2 i -η2 i )1{i / ∈ m ξ } ≤ 9σ 2 i ν i + 8η 2 i + x 2 i 1{| bi | ≤ αs}.
Lemma 5.4 If A2 holds, we have

ξ 2 i ≤ s 2 β ′ log n + ξ 2 i 1{ξ 2 i > s 2 β ′ log n}, with E (ξ 2 i 1{ξ 2 i > s 2 β ′ log n}) = O(n -1 log n). Proof. Write ξ 2 i ≤ s 2 β ′ log n 1{ξ 2 i ≤ s 2 β ′ log n} + ξ 2 i 1{ξ 2 i > s 2 β ′ log n}.
To bound the first term, we use the crude inequality 1{ξ 2 i ≤ s 2 β ′ log n} ≤ 1. For the second term, we have as a consequence of A2,

E ξ 2 i 1{ξ 2 i > s 2 β ′ log n} = ∞ 0 P ξ 2 i 1{ξ 2 i /s 2 > β ′ log n} > t dt = s 2 β ′ log n P(ξ 2 i /s 2 > β ′ log n) + s 2 ∞ β ′ log n P(ξ 2 i /s 2 > t) dt ≤ K ′ β ′ s 2 (1 + log n) n .

Proofs

Proof of Theorem 3.1. Write

x m -x 0 2 = xm * -x 0 2 + i / ∈m * (η 2 i -x 2 i )1{i ∈ m} + i∈m * (x 2 i -η 2 i )1{i / ∈ m}.
The objective is to bound the terms E((η 2 ix 2 i )1{i ∈ m}) and E((x 2 iη 2 i )1{i / ∈ m}). First, assume that σ 2 i > 1/n 2 , i.e. µ i = β log (n 2 σ 2 i ). By Lemma 5.1, we know that

E (η 2 i -x 2 i )1{i ∈ m} ≤ 2Kβσ 2 i e -µ i /β ≤ 2Kβ n 2 .
The same bound holds if σ 2 i ≤ 1/n 2 with µ i = 0, as a straight-forward consequence of Lemma 5.1. On the other hand, note that if i / ∈ m, then µ i = β log(n 2 σ 2 i ). Lemma 5.1 warrants E (x 2 iη 2 i )1{i / ∈ m} ≤ σ 2 i 6β log(n 2 σ 2 i ) + 2 . Since i ∈ m * , log(n 2 σ 2 i ) ≤ 2 log n + log x † 2 . We conclude by summing all the terms.

Proof of Theorem 4.1. The proof starts as in Theorem 3.1. We have

x m ξ -x † 2 = xm * ξ -x † 2 + i / ∈m * ξ (η 2 i -x 2 i )1{i ∈ m ξ } + i∈m * ξ (x 2 i -η2 i )1{i / ∈ m ξ },
and the objective is to bound the conditional expectation of each term separately. Using successively Lemma 5.3 and Lemma 5.4, we get The result follows by summing all the term, using that the risk of the oracle x m ξ is

E ξ (η 2 i -x 2 i )1{i ∈ m ξ } ≤ 4Kβ n 2 + 4α
E ξ xm * ξ -x † 2 = i / ∈m * ξ x 2 i + i∈m * ξ E ξ (η 2 i ).
Proof of Corollary 4.2. It suffices to show that the term i∈M x 2 i is of the same order as the risk of the oracle. Write For all i ∈ M, the probability P(| bi | ≤ |b i |/2) is greater than C as a consequence of A3. We deduce i∈M x 2 i ≤ C -1 E xm * ξx † 2 .

Corollary 4 . 2

 42 If the conditions A1, A2 and A3 hold, the threshold estimator defined in Theorem 4.1 satisfies

  | bi | ≤ |b i |/2).

  -2 s -2 ξ 2 By Lemma 5.4, we know that κ(ξ) = i / On the other hand, Lemma 5.3 gives, for θ = 8,E ξ (x 2 i -η2 i )1{i / ∈ m ξ } ≤ 9σ 2 i ν i + 8E ξ (η 2 i ) + x 2 i 1{| bi | ≤ αs}. For all i ∈ m * ξ , we know that | bi | ≥ |b i |/2. Thus, if i ∈ m * ξ , 1{| bi | ≤ αs} ≤ 1{i ∈ M}, where we recall M = {i : |b i | < 2αs}. We know also that, if i ∈ m * ξ , then σ2 i ≤ x 2 i . Thus, ν i = β log(n 2 σ2 i ) ≤ 2β log n + β log x † 2 . Noticing that σ2 i ≤ E ξ (η 2 i ), we find E ξ (x 2 i -η2 i )1{i / ∈ m ξ } ≤ (18β log n + 9β log x † 2 + 8)E ξ (η 2 i ) + x 2 i 1{i ∈ M}.

					∈m * ξ	κ i (ξ) is such that
		E(κ(ξ)) ≤	4(Kβ + 2α -2 K ′ β ′ x † 2 log n) n	= O	log n n	.
						i x 2 i ≤	4β ′ log n α 2	x 2 i + κ i (ξ),
	with	κ i (ξ) =	4Kβ n 2 +	4ξ 2 i x 2

i α 2 s 2 1{ξ 2 i > s 2 β ′ log n}.


est destinée au dépôt

Modelling the transport of particulate suspensions and formation damage during the deep injection of carbon dioxide

Mohammed Adil Sbai, Mohamed Azaroual

To cite this version: ! "

Mohammed
# $ % & ! " ! # $ ! ! ! % & '! ' & ! " & ! ! ( ) ( ) ! ( ) ! ( ) " % ! ! * % ( )! ! ! + , ! ( )! - ! '&( ) * + , *

  Adil Sbai, Mohamed Azaroual. Modelling the transport of particulate suspensions and formation damage during the deep injection of carbon dioxide. 1st French-German Symposium on Geological Storage of CO2, Jun 2007, Potsdam, Germany. hal-00589938


Joseph Onderi
email: joseph.orero@lip6.fr

Orero �a

Florent Levillain

Marc Damez-Fontaine

Maria Rifqi

Bernadette Bouchon-Meunier

ASSESSING GAMEPLAY EMOTIONS FROM PHYSIOLOGICAL SIGNALS A FUZZY DECISION TREES BASED MODEL

Keywords: Emotion Recognition, Video Games, Physiological Signals, Fuzzy Sets

As video games become a widespread form of entertainment, there is need to develop new evaluative methodologies for acknowledging the various aspects of the player's subjective experience, and especially the emotional aspect. Video game developers could benefit from being aware of how the player reacts emotionally to specific game parameters. In this study, we addressed the possibility to record physiological measures on players involved in an action game, with the main objective of developing adequate models to describe emotional states. Our goal was to estimate the emotional state of the player from physiological signals so as to relate these variations of the autonomic nervous system to the specific game narratives. To achieve this, we developed a fuzzy set theory based model to recognize various episodes of the game from the user's physiological signals. We used fuzzy decision trees to generate the rules that map these signals to game episodes characterized by a variation of challenge at stake. A specific advantage to our approach is that we automatically recognize game episodes from physiological signals with explicitly defined rules relating the signals to episodes in a continuous scale. We compare our results with the actual game statistics information associated with the game episodes.

INTRODUCTION

In the recent years video games have enjoyed a dramatic increase in popularity, the growing market being echoed by a genuine interest in the academic field, with an attempt to justify the existence of video game theory in academia [START_REF] Wolf | Introduction[END_REF]. In this context of flourishing technological and theoretical efforts, it is regrettable that effective methods of evaluation of the gaming technology still lag behind. Until now, the evaluation of video games has mainly relied on traditional methods through subjective self-reports such as questionnaires [START_REF] Nel | Automated video games evaluation based on the template formalism[END_REF], interviews and focus groups [START_REF] Grassioulet | A cognitive ergonomics approach to the process of game design and development[END_REF]. These subjective reports have limitations, one of them being that they only generate data at the moment the question is asked, and not through a continuous process. Also, subjects are known to be poor in self reporting their behaviors during game situations [START_REF] Pagulayan | User-centered design in games[END_REF]. Since the main goal of video games is to entertain through a continuous renewal of the user's interest, controlling the emergence of certain affective states is a key determinant in achieving a truly immersive experience. Therefore, in evaluating video games, it would be more appropriate to assess dynamically the emotional responses in relation to variations in challenge difficulty and to specific tasks at hand in the game.

Secondly, it has a particular interest regarding the field of developing affective human-computer interaction applications. The video games industry appears as a laboratory where we can develop innovative devices of interaction that could be used as a standard for designing affective-oriented interfaces keeping track of users's emotional states.

It is thus important to develop algorithms to automatically compute the affective state during the play process. In this study, we addressed the possibility of recording physiological measures on players involved in an action game, with the main objective of developing adequate models to describe emotional states. Specifically, our main goal was to estimate the emotional state of the player from physiological signals so as to: (i) relate these variations of the autonomic nervous system to the specific game narratives and (ii) differentiate these game narratives according to their level of challenge from physiological signature independent of the player.

The rest of the paper is organized as follows: In Section 2 we give an overview of affective computing. In Section 3, we outline the methods of classification and physiological measures used. We then describe the experimental setup in Section 4 and give our results in Section 5. Finally we give conclusions and future work in Section 6.

AFFECTIVE COMPUTING

Ever since Picard's highlight [START_REF] Picard | Affective Computing[END_REF] on the subject, affective computing has received considerable interest from many researchers in Human Computer Interaction (HCI). In particular the recognitions of emotions from physiological signals has become a major research focus in the recent past [START_REF] Fairclough | Fundamentals of physiological computing[END_REF]. Among a vast range of possible ways to access a user's emotional responses such as facial gesture or voice recognition through video and audio recording, subjective self report measures, physiological measures stand out. They grant an access to non conscious and non reportable processes [START_REF] Cacioppo | Psychophysiological science: Interdisciplinary approaches to classic questions about the mind[END_REF] and may to a certain extent be unobtrusively monitored [START_REF] Scheirer | Frustrating the user on purpose: A step toward building an affective computer[END_REF]. In addition, physiological measures may provide a fine-grained account of the rise, surge or extinction of an emotional response.

Considerable research progress has been made in the use of physiological measures to reliably predict specific emotions [START_REF] Picard | Toward machine emotional intelligence: Analysis of affective physiological state[END_REF][START_REF] Wagner | From physiological signals to emotions: Implementing and comparing selected methods for feature extraction and classification[END_REF][START_REF] Kim | Emotion recognition based on physiological changes in music listening[END_REF]. However, recognition of emotions through physiological signals is not yet a clear systematic task. There is lack of standardized methodology in the analysis of physiological data for affective system in terms of models, measures, features and algorithms to use [START_REF] Arroyo-Palacios | Towards a standardization in the use of physiological signals for affective recognition systems[END_REF]. The difficulty to uniquely map physiological patterns onto specific emotion types is due to the nature of these bio data. Physiological data tends to vary considerably from one person to another and may even display considerable differences within individuals on different occasions [START_REF] Picard | Toward machine emotional intelligence: Analysis of affective physiological state[END_REF][START_REF] Wagner | From physiological signals to emotions: Implementing and comparing selected methods for feature extraction and classification[END_REF][START_REF] Kim | Emotion recognition based on physiological changes in music listening[END_REF]. The aim of this work was to address some of these difficulties.

First, one of the main difficulty in emotion recognition is to ensure that people are really in the expected emotions. Researchers have employed various methods of eliciting emotions such as guided imagery technique, movie clips, math problems and music songs [START_REF] Picard | Toward machine emotional intelligence: Analysis of affective physiological state[END_REF][START_REF] Nasoz | Emotion recognition from physiological signals for presence technologies[END_REF][START_REF] Kim | Emotion recognition based on physiological changes in music listening[END_REF]. However, use of video games to elicit emotions is more advantageous in many ways. Particulary, due to their interactive nature, video games tend to promote a natural sense of immersion, the player being concerned by the consequences of its own actions, especially when mastering a skill is at stake. In this work, we took advantage of this propensity to immerse the user with the hope that we could drive specific physiological reactions corresponding to variations in the challenge proposed by the game at different points.

Secondly, emotions are better represented in fuzzy terms rather than discrete categories. As Goulev [START_REF] Goulev | An Investigation into the Use of AffectiveWare in Interactive Computer Applications[END_REF] already pointed out, even human ability to identify emotions from physical appearance is rarely hundred percent and as such, it is necessary to express in fuzzy terms mapping emotions from physiological data. Moreover, the physiological data from sensors is itself imperfect such that it is difficult to express the results in crisp terms [START_REF] Bouchon-Meunier | Aggregation and Fusion of Imperfect Information[END_REF]. Also, one could express more than one emotional state category at a particular time. A fuzzy set theory based approach could better represent both the multi-state complexity and uncertainties as well as imperfections in emotion recognition.

CLASSIFICATION AND PHYSIOLOGICAL MEASURES

Classification Methods

In emotion recognition from physiological signals, various methods of classification have been used in the past such as linear discriminant analysis, k-nearest-neighbor (KNN), multilayer perceptron network and decision trees (DT) [START_REF] Picard | Toward machine emotional intelligence: Analysis of affective physiological state[END_REF][START_REF] Nasoz | Emotion recognition from physiological signals for presence technologies[END_REF][START_REF] Wagner | From physiological signals to emotions: Implementing and comparing selected methods for feature extraction and classification[END_REF][START_REF] Rani | Anxiety-based affective communication for implicit human machine interaction[END_REF][START_REF] Kim | Emotion recognition based on physiological changes in music listening[END_REF]. A comparison of various methods using optimal features [START_REF] Wagner | From physiological signals to emotions: Implementing and comparing selected methods for feature extraction and classification[END_REF] seems to suggest that results depend on the method chosen and the nature of the experiment. In particular, decision trees have been found to perform comparably well [START_REF] Rani | Anxiety-based affective communication for implicit human machine interaction[END_REF]. A major superiority characteristic of decision trees is that they give explicitly defined rules used in classification. This is important in our case to explicitly know the relationship between the physiological signals attributes and the emotional states.

Secondly, as we have already noted, it is preferable to use a fuzzy sets theory based approach in classification. In fuzzy set theory [START_REF] Zadeh | Fuzzy sets[END_REF], a fuzzy set is represented by a membership function, µ � : A → [�� 1], indicating the degree to which an element belongs to a given set A. This is a contrast to ��� 1} in a crisp set, in which an element can only belong to a given set (membership value of 1) or not (membership value of 0). It is interesting for this kind of recognition to express the emotional states in a continuous scale of values [�� 1]. Therefore, we choose to use fuzzy decision trees (FDT) in this work. Fuzzy approach has been used in the past to model emotions from physiological signals [START_REF] Goulev | An Investigation into the Use of AffectiveWare in Interactive Computer Applications[END_REF][START_REF] Mandryk | A fuzzy physiological approach for continuously modeling emotion during interaction with play technologies[END_REF]. Particularly, Mandryk and Atkins [START_REF] Mandryk | A fuzzy physiological approach for continuously modeling emotion during interaction with play technologies[END_REF], developed a fuzzy logic model of rules that were used to transform arousal and valence based on psychophysiology literature. The aim of our experiment was to automatically generated these rules that define the psychophysiology relations from the fuzzy decision trees based on the players' physiological data. We used Marsala's Salammbô Fuzzy Decision Tree [START_REF] Marsala | Fuzzy decision trees to help flexible querying[END_REF] and compared our results with Quilan's C4.5 decision trees [START_REF] Quinlan | C4.5: Programs for Machine Learning[END_REF] and KNN [START_REF] Cover | Nearest neighbor pattern classficiation[END_REF] .

Physiological Signals and Features

Based on previous literature, we chose to collect galvanic skin response (GSR) which is a measure of the conductivity of the skin. GSR is considered as an effective correlate to arousal [START_REF] Lang | The emotion probe studies of motivation and attention[END_REF][START_REF] Bradley | Affective picture processing[END_REF][START_REF] Detenber | Roll' em� the effects of picture motion on emotional responses[END_REF] and has been extensively used in the domain of affective computing [START_REF] Scheirer | Frustrating the user on purpose: A step toward building an affective computer[END_REF][START_REF] Ward | Physiological responses to different web page designs[END_REF][START_REF] Mandryk | A fuzzy physiological approach for continuously modeling emotion during interaction with play technologies[END_REF]. We also collected heart rate (HR) through a measure of cardiovascular activity. HR may differentiate between positive and negative emotions [START_REF] Winton | Facial and autonomic manifestations of the dimensional structure of emotion[END_REF][START_REF] Papillo | Principles of psychophysiology: Physical, social, and inferential elements[END_REF][START_REF] Bradley | Affective picture processing[END_REF] although the possible correlation between heart rate and valence remains debated [START_REF] Neumann | Similar patterns of cardiovascular response during emotional activation as a function of affective valence and arousal and gender[END_REF][START_REF] Ritz | Airway response of healthy individuals to affective picture series[END_REF].

After collecting the signals, the raw signals were expressed in values [�� 1] for each subject's data. Then as in [START_REF] Picard | Toward machine emotional intelligence: Analysis of affective physiological state[END_REF], we calculated from each of these signals, six statistical features (the mean, the standard deviation, and the mean of the first and the second absolute difference of the raw signal and the normalized signal) for each game episode segment. However, in our case, for each feature, we subtracted the subject's baseline value. The baseline value was calculated by taking the value of each player's signal from the two minutes period preceding the beginning of the game. This was to account for the variations from subject to subject in signal values and also to inform us whether the attribute value increases or decreases on playing the particular game episode.

Feature Selection

A key research question is to determine which set of measures and features of these measures provide the optimal combination for discriminating the different emotional categories [START_REF] Picard | Affective Computing[END_REF]. There is no agreement on a particular feature selection method [START_REF] Arroyo-Palacios | Towards a standardization in the use of physiological signals for affective recognition systems[END_REF], it highly depends on the nature of the data and method used in classification. In this work we tested Sequential Backward Search (SBS), Sequential Forward Search (SFS) and Sequential Floating Forward Search (SFFS) which have been found to perform well in this kind of data [START_REF] Picard | Toward machine emotional intelligence: Analysis of affective physiological state[END_REF][START_REF] Wagner | From physiological signals to emotions: Implementing and comparing selected methods for feature extraction and classification[END_REF][START_REF] Yang | Emotion feature selection from physiological signal based on bpso[END_REF][START_REF] Kim | Emotion recognition based on physiological changes in music listening[END_REF]. However, there was no difference in results between them, the optimal feature set was the same for all the three. This could be expected since, unlike other classifiers, decision trees also perform feature selection by selecting the best attribute to discriminate the classes in each node of the tree. Nevertheless, to be able to rank features in order of their relevance in classification, we utilized SFS.

THE EXPERIMENTAL SETUP

Game Episodes

In this study we tested a game belonging to a popular genre in the game industry, Halo3, which is a First-Person Shooter (FPS). This game appears as one of the most immersive genre in that it propels the player at the heart of action through a first-person perspective. To segment the game session in different episodes we used specific events scripted inside the game. Scripted events are triggered by certain actions from the player, they are mandatory and are typically used by game designers either to mark transitions between different periods of the game, to inform the player about a certain state of the world or to draw her attention toward key information. They manifest through dialogs, visual and/or sound effects (e.g. the rumor of a battleship afar indicating the proximity of the enemy). By taking advantage of these scripted events we could delimitate three different combat episodes. Each episode is characterized by distinct geographical cues and corresponds approximatively to a game level in the traditional acceptation where boundaries of a level are defined by the completion of a goal (in our case the complete "cleaning" of a zone).

Before playing the game, participants were told to rest for two minutes during which a movie was presented, depicting abstract forms evolving smoothly and slowly so as to help the participant to relax. This period was considered as a baseline for the subsequent physiological recording. The first combat sequence (Combat1) is characterized by a succession of short skirmishes. The overall difficulty in this episode is low as it is almost possible to progress without shooting while fellow soldiers take care of the enemies. This is followed by the second combat session (Combat2) which is globally a bit more difficult than Combat1. The last combat (Combat3) begins when the player arrives at an enemy site. This is the most intense episode with enemies outnumbering and taking cover. In order to contrast combat episodes with more relaxing periods, we distinguished two additional sessions consisting in transitory periods (REST), triggered by the death of the last enemy in a level, wherein the player is able to rest, collect ammunitions and explore the environment.

Participants and Setting

The experiment was conducted at the LUTIN (Laboratoire des Usages en Technologies d' Information Numérique) in the Cité des Sciences et de l'Industrie, Paris, France. Twelve male participants aged between 18 to 40 were recruited from visitors of the Cité des Sciences et de l'Industrie to participate in the experiment. No specific expertise in the field of video games was required, although we selected participants able to manipulate a gamepad and to orient themselves in a virtual environment. Halo3 was played on a Microsoft Xbox 360 on a 32-inch LCD television. A camera captured the TV screen. Participants were seated at approximatively one meter from the screen, they were explicitly told not to move and keep the game pad onto their laps in order to avoid any muscular artefact in the physiological recording. All participants played to the game Halo3, at the easiest to difficult level. Participants were allowed to play at their own pace, until the experimenter told them to quit the game. The game ended when participants reached a specific point in the game (end of Combat3).

To collect the physiological measures we used the Biopac MP35 acquisition unit and the software BSLPro to visualize the data. To record the HR, we measured the electrocardiography (ECG) through a Einthoven derivation II placing pre-gelled surface electrodes on the ankles and on the wrist. We recorded GSR using surface electrodes attached with Velcro Straps that were placed on two fingers of the left hand. The fingers wearing the electrodes remained wedged under the gamepad. In order to synchronize the video with the physiological data, we used sound markers emitted at the beginning and at the end of the baseline video which were sent to the acquisition unit. ECG and GSR data were collected at 200 Hz. As noisy ECG data may produce failures in computing the HR, we inspected the HR data and corrected manually every erroneous samples. As we have already noted, our main goal in this experiment was to automatically discriminate episodes of the game characterized by high challenge from those characterized by low challenge based on physiological signals. The three different combat episodes could be characterized in term of their level of difficulty, based on certain information provided by the player activity. Specifically, we took advantage of a warning signal triggered each time the player is in danger, i.e. when his energy bar is fully depleted. This kind of information is a fair indicator of the level of intensity of an episode, as it reflects the amount of stress placed upon the player in a circumstance of likely death. As can be seen in Figure 1, there is a global and significant increase of the number of warning signals in Combat3 when compared to Combat1. However, there is only a slight increase of warning signals in Combat2 when compared to Combat1. As a consequence, we chose to contrast between Combat1 and Combat3.

RESULTS

�

� � � �� � � � � � � � � � �� �� �� ������������ �������������� ������� ������� �������
As earlier explained, we calculated 6 features from each physiological signals and utilized these attribute values as input to the classifiers. Although the output from the fuzzy decision tree was in continuous values [�� 1], we deffuzzified the values into ��� 1} for quantitative analysis and to compare the results with other classifiers. In the first attempt, a classifier was trained by leaving out one subject's data (test sample) against the rest of the subjects' data (training sample) using the aggregate of the entire episode (3 to 5 minutes). This was to give us a global recognition model of various game episodes. But in a related way, we considered a second alternative, in which we subdivided each game episode into 20 segments (approximately 10 seconds length). This was to validate if the model could apply on a continuous scale at various times during the game episodes. Half (six) of the subjects' samples was used for training and the other half for testing. In the analysis, we give the results of both alternatives.

First, we performed machine learning to differentiate between game combats of high challenge from those of low challenge with reasonable success as shown in Table 1 andTable 2. In discriminating these episodes, the GSR (mean amplitude and the mean of the second difference of the normalized) signal was the most relevant when compared across the classifiers using both the aggregate of the entire episode and the 10 seconds segments samples. Secondly, we tried to discriminate between episodes wherein the player is involved in combat and transitory periods (REST). As already noted, the transitory periods (REST) came after the death of the last enemy in the combat episode and thus index the moments the players feels relieved to have accomplished a goal. The player is not involved in any shooting while moving forward to start the next combat episode. During our machine learning process, we managed to differentiate these episodes (combat and REST) with 68.75% and 70.4% success as shown in Table 3 and Table 4 respectively. It turned out that, in addition to GSR (standard deviation) signal, the HR (the mean of the first and second difference) signal was very relevant in discriminating the REST/combat categories across the classifiers using both the aggregate of the entire episode and the 10 seconds segments samples. into discrete categories. However, fuzzy decision trees will be more applicable as they can be used on a continuous scale.

CONCLUSIONS

In this study we tested a model to automatically recognize specific episodes in an action game. Our aim was to relate a set of physiological signals features to specific game episodes. By taking advantage of the game narratives and the corresponding variations in the challenge proposed by the game through its different levels, we trained fuzzy decision trees to predict the episodes and generate rules that map the features to these periods of the game.

Although we couldn't differentiate clearly between combat episodes and transitory episodes where the player is not involved in a fight, we managed to identify with considerable success different combat episodes characterized by a variation in the challenge at stake. Our results point to galvanic skin response as a very relevant measure of the level of arousal of a player since it best discriminated between very intense episodes and not so intense episodes. Similarly, we found that, in addition to galvanic skin response, heart rate was also very relevant measure in discriminating between shooting and the subsequent rest episodes. We anticipate that variations of these signals could index the moments the players feels relieved when a goal has been accomplished. However, several issues inhibited us from obtaining optimal results in this experiment. One issue regards to the players level of expertise. As we couldn't control properly the level of expertise of our participants, we may have failed to obtain homogeneous samples. Confronted to the same level of difficulty, players with sensibly different playing background may react very differently. For instance a well trained player would react positively to an especially intense episode, as it may represent a fair challenge to him. Inversely, a beginner would consider that the same episode is outrageously punishing and thus, frustrating. Therefore, we would hope to obtain an even better level of recognition by testing people with the same amount of experience, and possibly the same kind of reaction confronted to an increase in intensity. Also, the indexation method we used to differentiate the game episodes may have failed to be the most accurate in differentiating clearly between combat episodes and transitory episodes where the player is not involved in a fight. It is possible that psychological boundaries of an episode (i.e. when a period of the game is subjectively felt as starting and ending) only partially coincide with objective markers such as scripted events. In consequence, it would be useful to consider the possibility to index events according to psychological markers.

Much is still to be done before getting access to the structure of the player's emotional pro-cesses. In particular, we need to improve our methods in order to define episodes with respect to psychological states and their variations. As a next step toward recognizing levels of affective involvement in a game, we would first need to consider more physiological measures. The present study did not exhaust all the possible signals available and improvements can be expected from the consideration of a large amount of physiological features. Additionally, we would need to differentiate physiological signatures not only from a level of intensity we assume they elicit, but from the players' direct assessment. We could therefore evaluate episodes based on various subjective scales (e.g. how frustrated I am, how concentrated I am, how fun is my experience, etc.) so as to define the episodes on several dimensions. Altogether, the road map for the forthcoming investigation of affective states in video game will get through a clear definition of the most relevant dimensions to account for the emotional response we target, as well as a thorough examination of the best physiological features. This way, we hope for a truly systematic affective recognition procedure to be incorporated to the games evaluation routines.

Figure 1 :

 1 Figure 1: Numbers of Red Warnings during Combat1, Combat2 and Combat3 for each participant

Table 1 :

 1 Low/High Challenge Episodes Pre-labeled Vs Predicted Results (Aggregate of Entire Episode)

	-	Low High Salammbô FDT Quilan C4.5 DT KNN�k=1)
	Low	8	4	66.67%	50.00%	66.67%
	High	1	11	91.67%	83.33%	58.33%
	Average			79.17%	66.67%	62.50%

Table 2 :

 2 Low/High Challenge Episodes Pre-labeled Vs Predicted Results (10 Seconds Segments)

	-	Low High Salammbô FDT Quilan C4.5 DT KNN�k=10)
	Low	202	38	84.17%	85.00 %	80.00%
	High	49	191	79.58%	77.08%	86.67%
	Average			81.88%	81.04%	83.33%

Table 3 :

 3 Combat/REST Episodes Pre-labeled Vs Predicted Results (Aggregate of Entire Episode) Although using the aggregate of the entire episode gives us a global classification model, it drastically reduces the sample size for good machine learning. The performance of classifiers is also comparable especially when episodes are segmented, if judged by classification

	-	Combat REST Salammbô FDT Quilan C4.5 DT KNN�k=1)
	Combat 19	5	79.17%	41.70%	70.08 %
	REST	10	14	58.33%	54.20%	70.08 %
	Average			68.75%	47.95%	70.08 %
	Generally, segmentation of episodes seems to have increased the classification rate across all
	the classifiers.

Table 4 :

 4 Combat/REST Episodes Pre-labeled Vs Predicted Results (10 Seconds Segments)

	-	Combat REST Salammbô FDT Quilan C4.5 DT KNN�k=10)
	Combat 377	103	78.54%	76.0%	77.3 %
	REST	181	299	62.29%	66.5%	69.8 %
	Average			70.42%	71.25%	73.55 %