Elena Di Bernardino
email: elena.dibernardino@cnam.fr

Anne Estrade
email: anne.estrade@parisdescartes.fr

José R León
email: jose.leon@ciens.ucv.ve

A test of Gaussianity based on the Euler characteristic of excursion sets

Keywords: MSC 2010 subject classifications: Primary 62G10 ; secondary 60G10, 60G15, 60G60 Test of Gaussianity, Gaussian fields, Excursion sets, Level sets, Euler characteristic, Crossings

In the present paper, we deal with a stationary isotropic random field X : R d → R and we assume it is partially observed through some level functionals. We aim at providing a methodology for a test of Gaussianity based on this information. More precisely, the level functionals are given by the Euler characteristic of the excursion sets above a finite number of levels. On the one hand, we study the properties of these level functionals under the hypothesis that the random field X is Gaussian. In particular, we focus on the mapping that associates to any level u the expected Euler characteristic of the excursion set above level u. On the other hand, we study the same level functionals under alternative distributions of X, such as chi-square, harmonic oscillator and shot noise. In order to validate our methodology, a part of the work consists in numerical experimentations. We generate Monte-Carlo samples of Gaussian and non-Gaussian random fields and compare, from a statistical point of view, their level functionals. Goodnessof-fit p-values are displayed for both cases. Simulations are performed in one dimensional case (d = 1) and in two dimensional case (d = 2), using R.

Introduction

The question of the Gaussianity of a phenomenon is a historical and fundamental problem in statistical literature. This type of information can be crucial in many application problems: oceanography and waves behavior, hydrology and climatology, agronomy, neurology and spike behavior, insurance and finance, astrophysics. For instance, in this latter application, during the last decade a large number of researchers joined efforts to decide whether the Cosmic Microwave Background (CMB) temperature is Gaussian or not. The problem of determining whether an i.i.d random sample comes from a Gaussian distribution has been studied extensively. In the case where the mean and the variance of the random variable are known, one can use a classical goodness-of-fit test. However, if these parameters need to be estimated the test of Lillifors and a variant of the Cramer-Von Mises test, with estimated parameters, are well adapted. These tests are no more distribution-free and depend on the true value of the parameters. The p-values must be obtained by simulations. The situation becomes more complicated when the sample comes from a stationary stochastic process satisfying some mixing conditions. For this type of problem, some tests have been designed to determine whether the one dimensional marginals are Gaussian. We can cite, by way of illustration: the Eps test [START_REF] Epps | Testing that a stationary time series is Gaussian[END_REF] based on the empirical characteristic function; the test built by Lobato and Velasco [START_REF] Lobato | A simple test of normality for time series[END_REF] that uses a test of symmetry and kurtosis; the test of Subba and Gabr [START_REF] Subba Rao | A test for linearity of stationary time series[END_REF] where the bi-spectral density is the basis for the test. A remarkable exception is the test built in Cuesta-Albertos et al. [START_REF] Cuesta-Albertos | A random-projection based procedure to test if a stationary process is Gaussian[END_REF]. There, the method for constructing the test uses a one dimensional random projection, after which the projected sample is subjected to a test that infers whether the one dimensional marginal is Gaussian. Actually, the random projection procedure allows one to test whether the whole distribution of the process is Gaussian or not, and hence is not limited to marginal distributions.

In the present article, we deal with a real valued stationary isotropic random field and we use the information given by level functionals of a single realisation of the process to build a test of Gaussianity. The level functional is the Euler characteristic (EC) of the excursion sets above some levels. Our first motivation comes from the article [START_REF] Adler | Some New Random Field Tools for Spatial Analysis[END_REF] Section 7 (Model identification). In the aforementioned paper, Adler suggests to use the EC of the excursion sets as a way to determine what is the actual distribution of the observed process. His words can describe better than ours the main goal "Suppose that we are given a real valued random field f on R 3 , which we observe on the unit cube [0, 1] 3 . Our problem is that we are not certain what the distribution of the field is. For simplicity, we shall assume that the choices are Gaussian and χ 2 , where both are assumed to be stationary and isotropic. Then one way to choose between the two is to calculate, from the data, an empirical Euler characteristic curve, ....." c.f. [START_REF] Adler | Some New Random Field Tools for Spatial Analysis[END_REF]. If the data is Gaussian, then the theoretical curve is a very precise one (see Figures 1 and9, left panels), depending on the second spectral moment of the process and some other invariant quantities. Otherwise, if it is a χ 2 -process, Kramer oscillator process or a shot noise process, a completely different curve appears (see Figures 6,[START_REF] Cabaña | Affine processes: a test of isotropy based on level sets[END_REF] left panel, 10 left panel, and 11). In what follows, we offer a methodology that tries to implement these ideas.

The idea of observing level functionals of a random field in order to infer some information on the distribution of the field is not new. In [START_REF] Lindgren | Spectral Moment Estimation by Means of Level Crossings[END_REF], Lindgren provides estimators for the second spectral moment of an univariate Gaussian process that are based on the number of upcrossings at various levels. In [START_REF] Cabaña | Affine processes: a test of isotropy based on level sets[END_REF], Cabaña builds a test of isotropy, for a two dimensional random field, that is based on the perimeter and the area of excursion sets. In [START_REF] Taheriyoun | Testing the covariance function of stationary Gaussian random fields[END_REF], the covariance function of a bivariate Gaussian field is inferred from the excursion sets Euler characteristic. The same idea has inspired many precursors working in materials science, see for instance [START_REF] Roberts | Chord-distribution functions of three-dimensional random media: approximate first-passage times of Gaussian processes[END_REF], [START_REF] Nott | Multi-phase image modelling with excursion sets[END_REF]. In those papers, the modelling of images or slices of a two-phases materials (even more complicated) is achieved by using a two dimensional stationary Gaussian field that has been thresholded at a certain level. The observed data are the lengths of the respective two phases along any line extracted from the image. The aim of these studies is to identify the Gaussian covariance function. Let us also mention [START_REF] Adler | Estimating thresholding levels for random fields via Euler characteristics[END_REF] where the authors start from the observation of a neurological space-time signal at some moderate levels and deduce some parameters that help in estimating the probability of exceeding very high values. Not far from this thema, one can find the question of exceedances or the study of extreme values, when considering high levels (see, e.g., the seminal work of Rice [START_REF] Rice | Mathematical Analysis of Random Noise[END_REF] and [START_REF] Worsley | Estimating the number of peaks in a random field using the Hadwiger characteristic of excursion sets, with applications to medical images[END_REF]). We will not go further in that direction and, at the opposite, stay with the observation of moderate levels. In all the mentioned papers, the field that is under study is assumed to be Gaussian. On the contrary, in the present paper, Gaussianity is not assumed but has to be inferred without knowing the spatial correlation.

To be more precise, we aim at proving that the function that associates u to the mean excursions EC at level u provides a kind of signature of the distribution of the random field under study. Although a so complex information as the knowledge of the whole distribution cannot be summarize in a single function, our guess is that the shape of its graph could be enough to discriminate between Gaussianity and non Gaussianity. Our main tool will be a Central Limit Theorem for the EC of an excursion set of a stationary isotropic Gaussian random field. This asymptotic normality takes place when the domain grows to R d . The result is proved in [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF] with the help of a Breuer-Major Theorem [START_REF] Nourdin | Quantitative Breuer-Major theorems[END_REF]. We will also need some generalisations, extensions and explicit computations of this result. In particular, in the present work, we extend the results in [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF] by showing that the random variables

Z i k = |T i | -1/2 (φ(X, T i , u k ) -E[φ(X, T i , u k )]
) for i = 1, . . . , m; k = 1, . . . , p;

with φ(X, T i , u k ) the EC of a standard Gaussian field X associated to the level u k in the domain T i and E[φ(X, T i , u k )] the theoretical Gaussian mean of excursions EC at level u k , are asymptotically jointly distributed as Gaussian when the disjoint domains T i grow up and satisfy asymptotic independence (see Proposition 5 for further details). Furthermore, we provide a tractable expression of the associated asymptotic variances (see Propositions 2 and 8).

The aforementioned results are used to build a Gaussianity test for a standard random field X. Indeed, I. if the H0 hypothesis: "the random field X is Gaussian" holds true, then the test statistic based on the sum of the scaled Z i k 's follows an appropriate chi-squared distribution. This case is for instance, illustrated by the pretty QQ-plots alignments and high p-values of goodness-of-fit tests in Sections 4.1.3 and 5.1. Conversely, II. if the underlying standard random field X is not Gaussian with a given second spectral moment (in particular if X is a χ 2 or a Kramer oscillator process), we deliberately center again, as in the I. case, the Z i k 's variables by using the (wrong) theoretical Gaussian mean of excursions EC at level u k (see Equation ( 6)) with the same second spectral moment. In this case, we obtain a very small goodness-of-fit p-values for the chi-squared distribution associated to the considered test statistics based on Z i k 's (see Sections 4.2.1 and 4.2.2). Then we are able to reject the H0 hypothesis under considered alternatives. The crucial point is that the theoretical mean of excursions EC has a very different shape in the Gaussian case with respect to the considered alternatives: χ 2 -process or a Kramer oscillator process. This evident difference allows us to easily discriminate between Gaussian and non Gaussian setting.

As mentioned above, we consider two alternative hypothesis versus the Gaussianity: χ 2 -process and Kramer oscillator process. These processes are chosen for two types of reasons. First, the mean excursions EC curve for the considered processes can be analytically known. Secondly, they are typical examples of models with specific properties: positive asymmetry or oscillation. This variety of behaviors can cover real-life situations in the applications. For the same types of reasons, we also consider a discontinuous shot noise process as an alternative. However, due to the non-smoothness property of this process, we separate its study in a specific section (see Section 6). A further discussion in this sense is postponed to the conclusion section.

Main contributions of the paper. Under the assumption of stationarity and isotropy for the Gaussian field X, we give a new explicit formula for the second moment of the excursions EC (see Proposition 1 and Proposition 2). We extend the results in [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF] by showing a Central Limit Theorem for the joint excursions EC concerned with different levels and disjoint domains (see Proposition 5). We also show finiteness of the asymptotic variance of the excursions EC of a chi-square field (see Proposition 7). A numerically tractable formula for the asymptotic variance in the univariate Gaussian case is given in Proposition 8. We propose a statistical methodology to test the Gaussianity distribution and we implement that on simulated data-sets to evaluate the finite-sample performance of our strategy.

Outline of the paper. Section 1 contains the general setting and the definition of the observation tool, namely the Euler characteristic of excursion sets. In Section 2, we focus on the Gaussian hypothesis. Assuming that the field under study X is a stationary isotropic Gaussian field, we give explicit formulas for the first two moments of the excursions EC. We also recall the Central Limit Theorem satisfied by the excursions EC when the domain tends to R d (see Theorem 4) and give an extension for the joint excursions EC concerned with different levels and disjoint domains. It allows us to describe the situation in term of a statistical model with observations whose distribution is known under the null hypothesis H0. Section 3 is concerned with the study of two alternative distributions of the observed random field: χ 2 and Kramer oscillator. In both cases, we give an explicit formula for the mean EC of excursion sets. Section 4 is devoted to numerical illustrations for univariate processes. We generate trajectory samples of stationary processes, Gaussian and non Gaussian, and compare the theoretical mean function of the excursions EC and the empirical one. We also build some chi-square statistics and associated goodness-of-fit tests in order to quantify the deviation between Gaussian and non Gaussian case. At last, in Section 5, we go further in the numerical study by considering two dimensional random fields. We generate Gaussian and χ 2 samples and compare with the theoretical situation. All the random generations and numerical computations are performed with R. In the last Section 6, we consider a shot noise process as a potential alternative, giving a preliminar study that could lead to a more general test of Gaussianity for non smooth processes. In order to keep the methodological spirit of the present paper, we have reported the technical proofs to the appendix section.

Setting

Let us start by introducing some definitions and terminology. Let X be a real valued random field defined on R d . We will sometimes call it a multivariate process. We say that it is stationary if its distribution is invariant under translations in R d and isotropic if its distribution is invariant under rotations in R d . Hence, the distribution of a stationary isotropic random field is invariant under the action of the isometry group of R d . Assuming X is smooth enough, we write X ′ for its derivative, which is a random vector field from R d to R d . We also write X ′′ for its second derivative. For any t ∈ R d , we denote by ∇ 2 X(t) the 1 2 d(d + 1) random vector that contains the upper coefficients of the symmetric Hessian matrix X ′′ (t) and by X(t) the d + 1 2 d(d + 1) + 1 random vector (X ′ (t), ∇ 2 X(t), X(t)). Let us assume now that X is stationary, isotropic and centered. If almost every realisation of X is of class C 1 , then the covariance matrix of X ′ (0) has the following form

E(X ′ (0)X ′ (0) t ) = λ I d , (1)
for some λ > 0 usually named as second spectral moment.

All over the paper, we consider a real valued random field X defined on R d that satisfies the following assumption.

Assumption (A):

The random field X is stationary, isotropic, E(X(0)) = 0, Var(X(0)) = 1 and almost all realisations belong to C 3 (R d ). For any fixed t in R d , the covariance matrix of the random vector X(t) has full rank.

At last, the covariance function r is such that,

ψ(t) → 0 when ∥t∥ → +∞ and ψ ∈ L 1 (R d ) ,
where ψ(t) = max

( ∂ k r ∂t k (t) ; k = (k 1 , . . . , k d ) ∈ N d , k 1 + . . . + k d ≤ 4
) .

Remark: Assumption (A) could appear too strong for the task described in this paper. Actually, we are not interesting in working under the weakest assumption, our aim is rather to provide a methodology that applies for any dimension d. Since the theoretical tools that we use are mainly derived from [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF], we chose to work with the same assumption setting as in this paper.

Notations.

• for any u ∈ R and any compact T ⊂ R d , we call "excursion set of X above the level u within the domain T " the following set {t ∈ T : X(t) ≥ u},

• p Z (.) denotes the probability density function of any random vector Z (assuming it exists),

• | • | denotes without any ambiguity, either the absolute value, or the d-dimensional Lebesgue measure.

Euler characteristic.

The Euler characteristic of a compact domain K in R d can be heuristically defined in the case d = 1 as the number of connected components of K, or in the case d = 2 as the number of connected components minus the number of holes in K. In the case where K is an excursion set {t ∈ T : X(t) ≥ u}, with T a rectangle in R d and u a real number, there exists a rather tractable formula that uses the theory of Morse functions (see [START_REF] Adler | Random Fields and Geometry[END_REF] Chapter 9, for instance). It states that the Euler characteristic of {t ∈ T : X(t) ≥ u} is equal to a sum of two terms. The first one only depends on the restriction of X to the interior of T , it is given by the quantity φ(X, T, u) defined in Equation (2) below. The second one exclusively depends on the behaviour of X on the l-dimensional faces of T , with 0 ≤ l < d. From now on, we focus on the term φ(X, T, u), named as "modified Euler characteristic" in [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF], and we still call it Euler characteristic (EC). It is defined by the following,

φ(X, T, u) = d ∑ k=0 (-1) k µ k (T, u), where (2)
µ k (T, u) = #{t ∈ T : X(t) ≥ u, X ′ (t) = 0, index(X ′′ (t)) = d -k},
with T the interior of T and the "index" stands for the number of negative eigenvalues.

Special case 1: dimension one. When d = 1, writing [0, T ] instead of T for a while, Equation (2) becomes

φ(X, [0, T ], u) = #{local maxima of X above u in (0, T )} (3)
-#{local minima of X above u in (0, T )}.

Morse's theorem says that this quantity is linked with the number of up-crossings,

U (X, [0, T ], u) = #{t ∈ [0, T ] : X(t) = u, X ′ (t) ≥ 0} , by the relation φ(X, [0, T ], u) + 1 {X(0)>u,X ′ (0)<0} + 1 {X(T )>u,X ′ (T )>0} = U (X, [0, T ], u) + 1 {X(0)>u} .
Taking expectation in both expressions and using stationarity yield the next formula that we will use in the forthcoming sections, namely Sections 3 and 4,

E[φ(X, [0, T ], u)] = E[U (X, [0, T ], u)]. ( 4
)
Special case 2: dimension two. When d = 2, Equation ( 2) can be rewritten in the following way.

With the notations introduced within this equation, µ 0 (T, u) denotes the number of local maxima above u, µ 2 (T, u) denotes the number of local minima above u and µ 1 (T, u) the number of local saddle points above u. Hence,

φ(X, T, u) = #{local extrema of X above u in T } ( 5
)
-#{local saddle points of X above u in T }.

Under Gaussian hypothesis

In this section, we assume that X is Gaussian and satisfies Assumption (A) described in Section 1.

First two moments of the Euler characteristic of an excursion set

Let T be a cube in R d . This section is devoted to explicit formulas for the first two moments of φ(X, T, u). They are based on the decomposition in (2) and on Rice formulas for the factorial moments of µ k (T, u) (see for instance [START_REF] Adler | Random Fields and Geometry[END_REF] Chapter 11 or [START_REF] Azaïs | Level Sets and Extrema of Random Processes and Fields[END_REF] Chapter 6).

In particular, using the stationarity of X, the expectation can be computed as follows

E[φ(X, T, u)] = (-1) d ∫ T E[1 [u,∞) (X(t)) det(X ′′ (t)) | X ′ (t) = 0] p X ′ (t) (0) dt = (-1) d |T | (2πλ) -d/2 E[1 [u,∞) (X(0)) det(X ′′ (0))] ,
where we recall that λ is the second spectral moment of X, see [START_REF] Adler | Some New Random Field Tools for Spatial Analysis[END_REF]. Moreover, it is proved in [START_REF] Adler | Random Fields and Geometry[END_REF] Lemma 11.7.1, through a regression and due to Wick's formula, that

E[1 [u,∞) (X(0)) det(X ′′ (0))] = (-1) d (2π) -1/2 λ d H d-1 (u) e -u 2 /2 ,
where H k stands for the Hermite polynomial of order k. Hence, the next formula holds

E[φ(X, T, u)] = |T | (2π) -(d+1)/2 λ d/2 H d-1 (u) e -u 2 /2 . ( 6
)
In what follows, we will be particularly interested in the next function

C(u) = (2π) -(d+1)/2 λ d/2 H d-1 (u) e -u 2 /2 , ( 7
)
that yields E[φ(X, T, u)] = |T | C(u). Equation [START_REF] Biermé | Level total curvature integral: Euler characteristic and 2D random fields[END_REF] shows that C(u) implicitly depends on X through its dimension parameter d and its second spectral moment λ. Whenever necessary in the next sections, we will emphasize this dependence by writing C(u) = C(u, λ).

For the second moment, a so nice formula as [START_REF] Biermé | A Fourier approach for the level crossings of Shot Noise processes with jumps[END_REF] seems to be out of reach. Nevertheless, in the next proposition, we express the second moment as an integral that can be numerically evaluated. Let us mention the paper [START_REF] Taheriyoun | A note on the higher moments of the Euler characteristic of the excursion sets of random fields[END_REF] where a similar formula is established but concerning another modified Euler characteristic, namely the differential topology (DT) characteristic. In this reference, the authors also give the second moment of the excursions DT characteristic as an integral of conditional distribution functionals and they propose a Monte-Carlo simulation method to evaluate these functions. We will use the following functions, defined for u ∈ R and t ∈ R d ,

g(u) = E[ 1 [u,∞) (X(0)) | det(X ′′ (0))| ] D(t) = (2π) 2d det(λ 2 I d -r ′′ (t) 2 ) G(u, t) = E[1 [u,∞) (X(0)) 1 [u,∞) (X(t)) det(X ′′ (0)) det(X ′′ (t)) | X ′ (0) = X ′ (t) = 0]      (8)
In the one dimensional case, explicit formulas for the expectation functions g(u) and G(u, t) are given in the numerical Section 4 (see Proposition 8).

Proposition 1. Assume that X is Gaussian and satisfies Assumption (A).

Then, for any u ∈ R, the map t → G(u, t) D(t) -1/2 is integrable on any compact set in R d and

E[φ(X, T, u) 2 ] = ∫ R d |T ∩ (T -t)| G(u, t) D(t) -1/2 dt + |T | (2πλ) -d/2 g(u) .
Proof. Integrability comes from [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF] Proposition 1.1 since X ∈ C 3 . In order to compute the expectation of φ(X, T, u) 2 , let us start with (2). It yields

φ(X, T, u) 2 = ∑ 0≤k≤d µ k (T, u) + ∑ 0≤k≤d µ k (T, u)(µ k (T, u) -1) + ∑ 0≤k̸ =l≤d (-1) k+l µ k (T, u)µ l (T, u).
The expectation of the first term is equal to

E[ ∑ 0≤k≤d µ k (T, u)] = E[ #{t ∈ T : X(t) ≥ u, X ′ (t) = 0} ] = ∫ T E[ 1 [u,∞) (X(t)) | det X ′′ (t)| | X ′ (t) = 0] p X ′ (t) (0) dt = |T | g(u) p X ′ (0) (0) ,
where we have used Rice Formula to get the second line and stationarity as well as the independence between X ′ (0) and (X(0), X ′′ (0)) to get the third one.

For the second and third terms, we introduce, for k, l = 1, . . . , d and s, t ∈ R d ,

F k,l u (s, t) = E[1 [u,∞) (X(s)) | det X ′′ (s)| 1 D k (X ′′ (s)) ×1 [u,∞) (X(t)) | det X ′′ (t)| 1 D l (X ′′ (t)) | X ′ (s) = X ′ (t) = 0],
where D k denotes the set of symmetric matrices having index equal to d -k. Hence, Rice Formula for the second factorial moment allows us to obtain

E[ ∑ 0≤k≤d µ k (T, u)(µ k (T, u) -1)] = ∫ T ×T ( ∑ 0≤k≤d F k,k u (s, t)) p X ′ (s),X ′ (t) (0, 0) dsdt.
Given that D k ∩ D l = ∅ for k ̸ = l, one can adapt the proof of the second moment Rice formula to get

E[ ∑ 0≤k̸ =l≤d (-1) k+l µ k (T, u)µ l (T, u)] = ∫ T ×T ( ∑ 0≤k̸ =l≤d (-1) k+l F k,l u (s, t)) p X ′ (s),X ′ (t) (0, 0) dsdt.
Let us remark that

∑ 0≤k≤d F k,k u (s, t) + ∑ 0≤k̸ =l≤d (-1) k+l F k,l u (s, t) = G(u, t -s) .
It yields

E[φ(X, T, u) 2 ] = ∫ R d |T ∩ (T -t)| G(u, t) p X ′ (0),X ′ (t) (0, 0) dt + |T | g(u) p X ′ (0) (0).
It remains to compute the probability density function of (X ′ (0), X ′ (t)). The covariance matrix of this vector is equal to

( λI d -r ′′ (t) -r ′′ (t) λI d )
and so p X ′ (0),X ′ (t) (0, 0) = D(t) -1/2 . Hence the result.

Asymptotic variance

In the next proposition, we let the cube T grow to R d and we give a formula for the asymptotic variance of φ(X, T, u). Actually, we consider

T (N ) = {N t : t ∈ T }
the image of a fixed cube T by the dilation t → N t and we let N grow to +∞.

Proposition 2. Assume that X is Gaussian and satisfies Assumption (A) and let T be a cube in R d . Then for any u in R,

lim N →+∞ Var[|T (N ) | -1/2 φ(X, T (N ) , u)] = V (u) < +∞ with V (u) = ∫ R d (G(u, t) D(t) -1/2 -C(u) 2 ) dt + (2πλ) -d/2 g(u) , ( 9
)
where C(u), g(u), G(u, t), D(t) have been defined in [START_REF] Biermé | Level total curvature integral: Euler characteristic and 2D random fields[END_REF] and [START_REF] Cabaña | Affine processes: a test of isotropy based on level sets[END_REF].

Proof. From Proposition 1, for a fixed cube T , we have

Var[|T | -1/2 φ(X, T, u)] = ∫ R d |T ∩ (T -t)| |T | G(u, t) D(t) -1/2 dt + (2πλ) -d/2 g(u) -|T | -1 (Eφ(X, T, u)) 2 = ∫ R d |T ∩ (T -t)| |T | (G(u, t) D(t) -1/2 -C(u) 2 ) dt + (2πλ) -d/2 g(u) ,
where we have used the relation

|T | 2 = ∫ R d |T ∩ (T -t)
| dt to get the last line. Hence, the asymptotic formula can easily be derived using Lebesgue dominated convergence theorem conditionally to the fact that t → G(u, t) D(t) -1/2 -C(u) 2 belongs to L 1 (R d ). This point is the matter of the next lemma.

Lemma 3. The map t → G(u, t) D(t)

-1/2 -C(u) 2 belongs to L 1 (R d ).
The proof of Lemma 3 is postponed in the appendix section.

Beyond the existence of a finite asymptotic variance as stated in the previous proposition, φ(X, T (N ) , u) satisfies a Central Limit Theorem.

Theorem 4 (Theorem 2.6 in [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF]). Assume that X is Gaussian and satisfies Assumption (A) and let T be a cube in R d . Then for any u in R, the next convergence holds in distribution

|T (N ) | -1/2 (φ(X, T (N ) , u) -E[φ(X, T (N ) , u)]) L -→ N →∞ N (0, V (u)) , ( 10
)
where N (0, V (u)) stands for the centered Gaussian distribution with variance V (u) (see Equation ( 9)).

Disjoint domains and various levels

We now consider two domains T 1 and T 2 that are disjoint and two levels u 1 and u 2 that can be equal or not.

Proposition 5. Assume that X is Gaussian and satisfies Assumption (A). (a)

Let T 1 and

T 2 be two cubes in R d such that |T 1 | = |T 2 | and dist(T 1 , T 2 ) > 0 and let u 1 and u 2 belong to R. For any integer N > 0, we introduce Z (N ) i = |T (N ) i | -1/2 (φ(X, T (N ) i , u i ) -E[φ(X, T (N ) i , u i )]) for i = 1, 2.
As N → +∞,

( Z (N ) 1 , Z (N ) 2
) converges in distribution to a centered Gaussian vector with diagonal covariance matrix

( V (u 1 ) 0 0 V (u 2
)

)
where V (u i ) is prescribed by [START_REF] Cammarota | A Quantitative Central Limit Theorem for the Euler-Poincare Characteristic of Random Spherical Eigenfunctions[END_REF].

(b)

Let T be a cube in R d and let u 1 and u 2 belong to R. For any integer N > 0, we introduce

ζ (N ) i = |T (N ) | -1/2 (φ(X, T (N ) , u i ) -E[φ(X, T (N ) , u i )]) for i = 1, 2.
As N → +∞,

( ζ (N ) 1 , ζ (N ) 2
) converges in distribution to a centered Gaussian vector with covariance matrix

( V (u 1 ) V (u 1 , u 2 ) V (u 1 , u 2 ) V (u 2 )
) , where V (u i ) is prescribed by [START_REF] Cammarota | A Quantitative Central Limit Theorem for the Euler-Poincare Characteristic of Random Spherical Eigenfunctions[END_REF] and V (u 1 , u 2 ) is a finite real number.

Comment. In [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF] Theorem 2.5, the covariance V (u 1 , u 2 ) is prescribed by a series expansion that, although convergent, is so awkward that it cannot be used to evaluate it in practise.

Proof of Proposition 5.

Concerning item (a), we first have to establish that the fields Z (N ) 1

and

Z (N ) 2
are asymptotically independent. Intuitively, it comes from the fact that the distance between the cubes T goes to infinity and the covariance function of X has a sufficient rate of decay at infinity due to Assumption (A). Precisely, the asymptotic decorrelation is given in the next lemma.

Lemma 6. Under assumptions of Proposition 5 and using the same notation, it holds that

Cov ( Z (N ) 1 , Z (N ) 2 ) → N →+∞ 0.
Lemma 6 is proved in the appendix section. By using this result and Theorem 4, we know that the covariance matrix of the random vector (Z

(N ) 1 , Z (N ) 2
) tends to

( V (u 1 ) 0 0 V (u 2 )
)

. One can use the same arguments as those in [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF] (which were inspired by the Breuer-Major Theorem of [START_REF] Nourdin | Quantitative Breuer-Major theorems[END_REF]) to establish that any linear combination xZ

(N ) 1 + yZ (N ) 2
has a Gaussian limit in distribution. Item (b) is proved in [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF], Theorem 2.5.

Statistical model

We are now able to build a test for the H0 hypothesis: "the random field X is Gaussian". Actually, we assume that X is observed through the family (Y i k ) 1≤i≤m;1≤k≤p where

Y i k = φ(X, T i , u k ) |T i | , i = 1, . . . , m ; k = 1, . . . , p ,
for disjoint domains T 1 , . . . , T m that have the same large volume and are at large distance one from each other, and various levels u 1 ≤ . . . ≤ u p . Under these conditions, by using Proposition 5 and Equation ( 6), one can write the following approximation valid under H0,

Y i k ≈ C(u k , λ) + ϵ i k , i = 1, . . . , m ; k = 1, . . . , p , (11)
where

(ϵ i k ) i,k is a mp-dimensional centered Gaussian vector with covariance Cov(ϵ i k , ϵ j l ) = Cov(Y i k , Y j l ) = δ ij V (u k , u l )/|T i | . ( 12
)
The deterministic value C(u k , λ) is given by Equation [START_REF] Biermé | Level total curvature integral: Euler characteristic and 2D random fields[END_REF]. It depends on the level u k and on the standard field X only through its second spectral moment λ.

For k = l, V (u k , u l ) = V (u k
) is given by ( 9), and for k ̸ = l, since no explicit formula is available for V (u k , u l ) in practise, it has to be estimated. Using the statistical model in ( 11)-( 12), we will focus on two particular cases.

(a) Diagonal case. We take m disjoint domains T 1 , . . . , T m and m levels u 1 , . . . , u m , such that each level u k is associated to a single domain T k . In this setting we have

m observations Y 1 1 , . . . , Y m m and a diagonal covariance matrix equal to diag(V (u 1 )/|T 1 |, . . . , V (u m )/|T m |). (b) Crossed case.
We take m disjoint domains T 1 , . . . , T m and p levels u 1 , . . . , u p , such that different levels u k are associated to the same domain T i . In this setting we have mp observations (Y i k ) i,k and their covariance matrix is given by (12).

Under alternative hypothesis

χ 2 hypothesis

In this section, we deal with χ 2 distributions instead of Gaussian ones. Let us fix s, a non negative integer, as the degrees of freedom. We start with {X i (.)} s i=1 , an independent sample of standard stationary Gaussian fields on R d that satisfy Assumption (A) of Section 1. We denote by r X their covariance function and recall that r X (0) = 1. Consider now the following stationary fields

χ 2 s (•) = s ∑ i=1 (X i (•)) 2 and Z (s) (•) = 1 √ 2s ( s ∑ i=1 (X i (•)) 2 -s) .
Note that for any t ∈ R d , χ 2 s (t) is a chi-square random variable with s degrees of freedom. One get readily that

E[Z (s) (0)] = 0, Var[Z (s) (0)] = 1, E[Z (s) (t)Z (s) (0)] = r X (t) 2 .
Therefore, Z (s) also satisfies Assumption (A). Moreover, its second spectral moment is equal to

λ with λ = -((r X ) 2 ) ′′ (0) = -2r ′′ X (0).
We are now interested in the expectation of the Euler characteristic of

{t ∈ T : Z (s) (t) ≥ u} = {t ∈ T : χ 2 s (t) ≥ s + u √ 2s} ,
for a fixed cube T ⊂ R d and a fixed level u in R. A formula for the mean EC of excursion sets of χ 2 fields is given in [START_REF] Worsley | Local maxima and the expected Euler characteristic of excursion sets of χ 2 , F and t fields[END_REF], Theorem 3.5 (see also Theorem 15.10.1 in [START_REF] Adler | Random Fields and Geometry[END_REF]). It applies in our context, although we have to handle carefully with the second spectral moment of the X i 's, which is equal to λ/2. It yields

E[φ(Z (s) , T, u)] (13) = |T | ( λ π ) d/2 e -(s+u √ 2s)/2 (s + u √ 2s) (s-d)/2 2 (s-2+2d)/2 Γ(s/2) P d,s (s + u √ 2s) 1 [0,∞) (s + u √ 2s) ,
where P d,s (•) is a polynomial of degree d -1 with integer coefficients (see Section 3.3 in [START_REF] Worsley | Local maxima and the expected Euler characteristic of excursion sets of χ 2 , F and t fields[END_REF]). In particular, we stress that P 1,s (u) = 1 and

P 2,s (u) = u -s + 1.
Let us recall that, in dimension one, the mean Euler characteristic of the excursion above the level u is equal to the mean number of upcrossings at level u, see Equation ( 4). With this point of view, Formula (13) can also be found in [START_REF] Sharpe | Some properties of the crossings process generated by a stationary χ 2 process[END_REF] for instance. Next proposition is concerned with the second moment in dimension one, in the same spirit as Proposition 2.

Proposition 7. Let d = 1. Let us assume that the X ′ i s are one dimensional i.i.d. Gaussian processes that satisfy Assumption (A) for i = 1, . . . , s. Then, for any u in R and T in [0, +∞), φ(Z (s) , [0, T ], u) admits a finite second moment. Morever, there exists v s (u) ∈ [0, +∞) such that lim

T →+∞ 1 T Var[φ(Z (s) , [0, T ], u)] = v s (u) .
The proof of Proposition 7 can be found in the appendix section.

Kramer oscillator hypothesis

We work in dimension d = 1. Let us consider the following system of stochastic differential equations, well known as Kramer oscillator system,

dQ(t) = P (t)dt, dP (t) = σdW (t) -(cP (t) + V ′ (Q(t))dt, } (14)
where V (q) = a 0 q 4 -a 1 q 2 for positive constants a 0 and a 1 , σ and c are also positive constants, and W is a Brownian motion. The asymptotic properties of such a system have been studied for instance in [START_REF] Wu | Large and moderate deviations and exponential convergence for stochastic damping Hamiltonian systems[END_REF]. Besides, it is well known that the Markov process (Q, P ) has an invariant measure µ that can be written (up to a numerical constant factor) as

dµ(p, q) = exp ( - 2c σ 2 ( p 2 2 + V (q))
) dp dq.

From now on, we assume that (Q(0), P (0))'s distribution is proportional to µ, so that (Q, P ) is stationary. Here, we are interested in the process Q. It is stationary, centered, but certainly not Gaussian since its distribution is proportional to exp

( -2c σ 2 V (q)
) dq. Nevertheless, the distribution of its derivative process, Q ′ = P , is actually Gaussian with zero mean and variance equal to σ 2

2c . An application of Rice formula gives the mean number of upcrossings of Q in [0, T ] at any level u and hence, using (4), we get

E[φ(Q, [0, T ], u)] = E[U (Q, [0, T ], u)] = T ( 1 √ 2π(σ 2 /2c) ∫ ∞ 0 p exp -( c σ 2 p 2 )dp ) p Q (u) = T σ 2 √ πc p Q (u) with p Q (u) = exp(-2c σ 2 V (u)) ∫ R exp(-2c σ 2 V (q))dq . ( 15
)
Moreover, let us remark that a suitable choice of the parameters σ, c, a 0 , a 1 allows us to get Var(Q(0)) = 1 and Var(Q ′ (0)) = λ, so that Q satisfies the same moments constraints as the generic process X in Section 1. Actually, it is sufficient to prescribe σ 2 2c = λ and to solve the following non linear equation

∫ R exp ( - 2c σ 2 (a 0 q 4 -a 1 q 2 ) ) dq = ∫ R q 2 exp ( -( 2c σ 2 (a 0 q 4 -a 1 q 2 ) ) dq. ( 16
)

Univariate numerical illustrations

In this section, we focus on the one dimensional case and hence only deal with univariate processes. Both for Gaussian and alternative distributions, we first compare the theoretical formulas for the moments of excursions EC with the empirical ones, obtained by Monte Carlo simulations. Secondly, we perform the test statistics and calculate associated goodness-of-fit p-values. A graphical illustration of our test is also provided by drawing QQ-plots.

Gaussian process

First and second moments: E[φ(X, T, u)] and V (u)

We rewrite Equations ( 6) and ( 9) of Section 2 in the case d = 1 under H0 hypothesis, i.e. when X is a stationary isotropic standard Gaussian process with covariance function r and second spectral moment λ = -r ′′ (0):

E[φ(X, T, u)] = |T | C(u, λ) with C(u, λ) = (2π) -1 λ 1/2 e -u 2 /2 , ( 17
)
and Var[φ(X,

T, u)] ∼ |T | V (u) as |T | → +∞, with V (u) = ∫ R (G(u, t) D(t) -1/2 -C(u, λ) 2 ) dt + (2πλ) -1/2 g(u) , ( 18
)
where D(t) = (2π) 2 (λ 2 -r ′′ (t) 2 ) and g(u), G(u, t) are given by ( 8). In the one dimensional case, it yields

g(u) = E[ 1 [u,∞) (X(0)) |X ′′ (0)| ] G(u, t) = E[1 [u,∞) (X(0)) 1 [u,∞) (X(t)) X ′′ (0) X ′′ (t) | X ′ (0) = X ′ (t) = 0].
Our first task is now to provide formulas for g(u) and G(u, t) that can be numerically evaluated. They will obviously include Gaussian integrals and, therefore, we need to describe the Gaussian distributions that are involved.

On the one hand, introducing γ = r (4) (0), the covariance matrix of (X(0), X ′′ (0)) is

( 1 -λ -λ γ
) .

On the other hand, the conditional distribution

L(X(0), X(t), X ′′ (0), X ′′ (t) | X ′ (0) = X ′ (t) = 0)
is a 4-dimensional centered Gaussian distribution with covariance matrix given by

Γ(t) = C 11 (t) -C 12 (t)C 22 (t) -1 C 12 (t) T ,
where C 11 (t) is the covariance matrix of the vector (X(0), X(t), X ′′ (0), X ′′ (t)), C 22 (t) is the covariance matrix of the vector (X ′ (0), X ′ (t)) and C 12 (t) is the matrix of the covariances between those two vectors.

We still need some extra notations:

• ϕ denotes the standard Gaussian density, Φ the standard Gaussian distribution, and

Φ(•) = 1 -Φ(•); • for any (u, α, m) ∈ R 3 , we introduce        Q (0) (u, α, m) = ∫ ∞ u Φ(αx + m)ϕ(x)dx ; Q (1) (u, α, m) = - ∫ ∞ u ϕ(αx + m)ϕ ′ (x)dx Q (2) (u, α, m) = - ∫ ∞ u ϕ ′ (αx + m)ϕ(x)dx ; Q (3) (u, α, m) = ∫ ∞ u Φ(αx + m)ϕ ′′ (x)dx.
In Proposition 8 below, we finally give explicit formulas for the functions g and G. They will be useful to give a numerical evaluation of V (u) in [START_REF] Lachièze-Rey | Two-dimensional Kac-Rice formula. Application to shot noise processes excursions[END_REF] for various values of u.

Proposition 8. Let X be a Gaussian process that satisfies Assumption (A).

1. For any u in R,

g(u) = λϕ(u) ( 2Φ( λu (γ -λ 2 ) 1/2 ) -1 ) + ( 2γ π
) 1/2 Φ ( ( γ (γ -λ 2 ) ) 1/2 u ) .
2. Let u ∈ R and t ∈ R be fixed. We denote by (l ij ) 1≤i≥j≤4 the coefficients of any lower triangular matrix

L(t) such that L(t)L(t) T = Γ(t). Then G(u, t) = (l 41 l 31 + l 42 l 32 + l 43 l 33 ) Q (0) ( u l 11 , - l 21 l 22 , u l 22 ) +(l 42 l 31 + l 41 l 32 ) Q (1) ( u l 11 , - l 21 l 22 , u l 22
)

+l 42 l 32 Q (2) ( u l 11 , - l 21 l 22 , u l 22 ) + l 41 l 31 Q (3) ( u l 11 , - l 21 l 22 , u l 22
) .

The proof of this proposition is postponed in the appendix section.

In order to illustrate ( 17) and ( 18), we generate a 300 Monte-Carlo sample of a stationary standard Gaussian process with covariance function r(t) = e -t 2 . Note that it implies the second spectral moment λ = 2. By using the Cholesky function chol in R programm we are able to numerically evaluate L(t) matrix, for all t. These quantities are necessary to build g, G and numerically approximate the integral that gives the asymptotic variance V (u) (see Equation ( 18)). In order to evaluate φ(X, T, u) on each realisation of X, we use Equation ( 3) for various values of u. Comparison between theoretical formulas and empirical counterparts are shown in Figure 1. ) with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 . In this case λ = 2.

Estimation of the second spectral moment

It has been already quoted that the second spectral moment of X, denoted by λ, plays an important role in Equation [START_REF] Kauderer | Modes in n-dimensional first-order systems[END_REF]. Actually, all the influence of X in E[φ(X, T, u)] is summarised through this parameter. From our statistical point of view, since we aim at inferring the distribution of X from observations of the excursions of X, parameter λ is unknown a priori. It has to be estimated from the observed level functionals X.

We estimate γ = λ 1/2 by using an unbiased estimator introduced by Lingren [START_REF] Lindgren | Spectral Moment Estimation by Means of Level Crossings[END_REF] for stationary zeromean Gaussian processes. It is based on p different levels u 1 < u 2 < . . . < u p by the following prescription,

γ = p ∑ k=1 c k γ u k , with γ u k = 2 π T -1 e u 2 k /2 φ(X, [0, T ], u k ) . ( 19
)
Actually, in the paper of Lingren, the number of up-crossings of level u k in the interval [0, T ], namely U (X, [0, T ], u k ), is used instead of φ(X, [0, T ], u k ) (recall the analogy described in Section 1). As a general rule, considering

p = 3 , (u 1 , u 2 , u 3 ) = (-u, 0, u), with u = 2 3 √ V ar(X(t)) , c 1 = c 2 = c 3 = 1 3
seems to be an acceptable choice (see discussion in [START_REF] Lindgren | Spectral Moment Estimation by Means of Level Crossings[END_REF]). In this case

λ 1/2 = γ = 1 3 ( γ -u + γ 0 + γ -u
). An illustration of this estimation procedure is given in Figure 2. This estimation of λ 1/2 is also used in Figure 1 (left), dashed line. Boxplot, based on 300 Monte Carlo simulations, of the ratio between theoretical value λ 1/2 and estimated value γ given by [START_REF] Lelievre | Free Energy Computations[END_REF]. X(•) is a Gaussian univariate process (d = 1) with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 . In this case λ = 2. We choose p = 3, c 1 = c 2 = c 3 = 1 3 and u = (-2/3, 0, 2/3).

The same estimation strategy could be developed in the alternative hypothesis framework as soon as we have at our disposal a closed formula linking the mean EC with the second spectral moment (see for instance Equation [START_REF] Lindner | Analytical approach to the stochastic FitzHugh-Nagumo system and coherence resonance[END_REF] in the univariate chi-square case).

Chi-square statistics

In the following we will consider particular sub-cases of the two cases presented in Section 2.4, diagonal case (a) and crossed case (b). In particular, we focus on the four next models. We will now associate chi-square statistics with each of the four aforementioned models. For the diagonal case (a.1), we consider

F a1 := m ∑ i=1 ( φ(X, T i , u) -E[φ(X, T i , u)] √ |T i | V (u) ) 2 F a1 := m ∑ i=1   φ(X, T i , u) -E[φ(X, T i , u)] √ V ar(φ(X, T i , u))   2 where E[φ(X, T i , u)] = |T i | C(u, λ
) is given by ( 17) and V (u) is given by [START_REF] Lachièze-Rey | Two-dimensional Kac-Rice formula. Application to shot noise processes excursions[END_REF]. Furthermore, E and V ar in F a1 are respectively the empirical mean and the empirical variance on considered Monte-Carlo sample generations. Under H0 hypothesis, a consequence of ( 11)-( 12) is that both random variables F a1 and F a1 are approximately χ 2 m distributed, i.e. central chi-square with m degrees of freedom. We evaluate F a1 and F a1 on 300 Monte Carlo simulations. We choose m = 3 and u = 1.2. The QQ-plot comparison between the obtained empirical quantiles with the theoretical quantiles of a centered χ 2 m distribution is gathered in Figure 3 (first and second panels).

1 = -1.2, u 2 = 0, u 3 = 1.2. X(•) is a Gaussian univariate process (d = 1)
with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 with λ = 2. See associated goodness-of-fit test p-values in Table 1.

Consider now the diagonal case (a.2). Let

F a2 := m ∑ i=1 ( φ(X, T i , u i ) -E[φ(X, T i , u i )] √ |T i |V (u i ) ) 2 F a2 := m ∑ i=1   φ(X, T i , u i ) -E[φ(X, T i , u i )] √ V ar(φ(X, T i , u i ))   2
where E and V ar are respectively the empirical mean and the empirical variance on considered Monte Carlo sample generations. By using again the statistical model in ( 11)-( 12), both F a2 and F a2 are approximately χ 2 m distributed. An illustration is presented in Figure 3 (third and fourth panels), by choosing 300 Monte Carlo simulations, m = 3, |T i | = 200 and different levels

u 1 = -1.2, u 2 = 0, u 3 = GOF test p-
F b1 := ||Λ -1 (Z -E[Z])|| 2 , F b1 := || Λ -1 (Z -E[Z])|| 2 ,
where Z is the p-dimensional Gaussian vector given by

Z = ( φ(X, T, u k ) |T | 1/2 ) 1≤k≤p , E is

Consider now the crossed case (b.2).

In this setting we have m disjoint domains T 1 , . . . , T m and p different levels u 1 < . . . < u p . Let

F b2 := m ∑ i=1 ||Λ -1 (Z i -E[Z i ])|| 2 , F b2 := m ∑ i=1 || Λ (i) -1 (Z i -E[Z i ])|| 2 ,
where for any i ∈ {1, . . . , m}, Z i is the p-dimensional Gaussian vector given by

Z i = ( φ(X, T i , u k ) |T i | 1/2 ) 1≤k≤p
and E is the empirical mean on the Monte-Carlo simulations. Let Γ = (V (u k , u l )) 1≤kl,≤p be the (theoretical) covariance matrix of Z i . Let Λ stand for any square root of Γ. Similarly, let Γ (i) be the empirical covariance matrix of Z i evaluated on the Monte-Carlo sample generations, and let Λ (i) be any of its square root matrix. Hence,

||Λ -1 (Z i -E[Z i ])|| 2 and || Λ (i) -1 (Z i -E[Z i ])|| 2 are both approximately χ 2 p distributed.
Moreover, since for 1 ≤ i ̸ = j ≤ m, the Gaussian vectors Z i and Z j are independent, F b2 and F b2 are still centered χ 2 distributed with now mp degrees of freedom.

An illustration of the behaviour of F b2 is presented in Figure 4 (second panel), by choosing 300 Monte Carlo simulations, m = 3 disjoint domains T 1 , T 2 , T 3 with |T i | = 200 and p = 3 different levels u 1 = -1.5, u 2 = 0, u 3 = 1.5. Furthermore, in Table 2 we display the goodness-of-fit test p-values for F b1 and F b2 in Figure 4. Obtained p-values allow us to statistically accept the chi-squared distribution for these test statistics. 2.

GOF test p-values for χ 2 m p

Cramér-von Mises test

Anderson-Darling test

F b1 (m = 1, p = 3) 0.7725 0.8548 F b2 (m = 3, p = 3) 0.9745 0.9886

Table 2 Goodness-of-fit test p-values associated to test statistics F b1 and F b2 in Figure 4. The considered process X(•) is a

Gaussian univariate process (d = 1) with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 .

Alternative processes

In this section, we simulate trajectories of an alternative standard process Z. We compare it with a standard Gaussian process X having the same second spectral moment.

Similarly to Section 4.1.3, we now consider a test statistic. For sake of brevity, we only focus on the diagonal case (a.1) with a unique level u and m disjoint domains T 1 , . . . , T m , and we introduce

D := m ∑ i=1   φ(Z, T i , u) -E[φ(X, T i , u)] √ V ar(φ(Z, T i , u i ))   2 , ( 20
)
where E[φ(X, T i , u)] is the mean Euler characteristic of the excursions of X given by ( 6) with d = 1 and V ar(φ(Z, T i , u i )) is the empirical variance on considered Monte-Carlo sample generations. In order to use this empirical variance, one has to be sure that the variance of the Euler characteristic of the excursions of Z is finite. We recall that this is the case when Z is chi-square distributed (see Proposition 7). Under H0 hypothesis, the D statistic is asymptotically distributed as a chi-square with m degrees of freedom.

First alternative: χ 2 process

In conformity with Section 3.1, we consider

Z (s) (•) = 1 √ 2s (χ 2 s (•)) -s) with χ 2 s (•) = ∑ 1≤i≤s X i (•) 2 ,
where the X i 's are independent copies of a stationary standard Gaussian process with covariance function r(t) = e -t 2 /2 . Hence, we get λ = -2r ′′ (0) = 2 and then the obtained process Z (s) has the same variance and the same second spectral moment λ as the previous Gaussian one presented in Section 4.1.

In the first and third panels of Figure 5, we display the boxplot for the ratio between the empirical mean of 300 Monte Carlo values of φ(Z (s) , T, u) and the theoretical mean given by (13) in the case d = 1, i.e.,

E[φ(Z (s) , T, u)] = |T | ( λ π ) 1/2 (s + u √ 2s) (s-1)/2 2 s/2 Γ(s/2) e -(s+u √ 2s)/2 1 [0,∞) (s + u √ 2s). (21)
This expectation has to be compared with the expectation given by ( 17) when X is a stationary standard Gaussian process with second spectral moment λ = 2. Therefore, in the second and fourth panels of Figure 5, we display the boxplot for the ratio between the empirical 300 Monte Carlo mean value of φ(Z (s) , T, u) and the Gaussian theoretical expectation given in [START_REF] Kauderer | Modes in n-dimensional first-order systems[END_REF].

Remark. When the degrees of freedom s tend to infinity, the CLT implies that Z (s) tends in distribution to a stationary centered Gaussian process with covariance function equal to t → r(t) 2 , which implies a variance equal to 1 and a second spectral moment equal to 2 (exactly as the Gaussian process X considered in Section 4.1). On the other hand, using Stirling approximation for the Γ function, one can prove that the right-hand side of (21) tends to the right-hand side of (17), i.e.,

E[φ(Z (s) , T, u)] → |T | λ 1/2 (2π) -1 e -u 2 /2 , when s → ∞.
See also Formula (3.4) in [START_REF] Sharpe | Some properties of the crossings process generated by a stationary χ 2 process[END_REF] for the same remark. The comparison between the second and fourth panels in Figure 5, as well as Figure 6, illustrate this convergence.

In Figure 7 below, we show the QQ-plot of the test statistic D (see [START_REF] Lindgren | Spectral Moment Estimation by Means of Level Crossings[END_REF]) versus the quantiles of χ 2 m . Since the considered chi-square process Z (s) is not Gaussian, a deviation can be observed. In order and the theoretical mean given by Equation [START_REF] Lindner | Analytical approach to the stochastic FitzHugh-Nagumo system and coherence resonance[END_REF]. Second and fourth panels: Boxplot for the ratio between the same empirical values and the theoretical mean in the Gaussian case given by Equation [START_REF] Kauderer | Modes in n-dimensional first-order systems[END_REF]. In both cases λ = 2. The degrees of freedom s is chosen equal to 2 and 10 respectively.

Boxplot using Chi-square and Gaussian Expectations, s=2

Levels u -0.9 -0.75 -0.5 0 0.5 0.75 0.9 1.5 to quantify how we are really able to reject the null hypothesis H0 under the alternative, in Table 3 we display the generated p-values for the goodness-of-fit test of D. One can easily see that this test allows us to reject the fact that the considered underlying process Z (s) is Gaussian. This is particularly evident if one compares these p-values with those in Tables 1 and2.

Figure 7 shows a huge deviation from the bisector of the first orthant in the case s = 2 whereas this deviation is less significant in the case s = 10. It emphasizes the convergence in distribution of Z (s) towards a Gaussian as s tends to infinity. 7 versus a χ 2 m distribution with m = 6. The considered process Z (s) (•) is a chi-square univariate process (d = 1) with E(Z (s) (0)) = 0, Var(Z (s) (0)) = 1 and λ = 2.

Second alternative: Kramer oscillator process

In this section we generate a 300 Monte-Carlo sample of a Kramer oscillator process as defined in Section 3.2. In order to obtain a process Q with zero mean, unit variance and second spectral moment equal to 2, we solve Equation ( 16) and choose

σ = 2 , c = 1 , a 0 = 1 , a 1 = 2.3373 , C = 4.886.
The generation procedure is the following. We define a discretized schema to simulate the solution (Q(t), P (t)) of the stochastic differential system in [START_REF] Gabriel | Estimating and testing zones of abrupt change for spatial data[END_REF]. Actually, we use the Metropolization of the Euler-Verlet schema with a sufficient small discretization step. The interested reader is referred to Algorithm 2.11 (Generalized Hybrid Monte-Carlo) in [START_REF] Lelievre | Free Energy Computations[END_REF].

The comparison between the (theoretical) expectation given by ( 15) and the empirical one is shown in left panel of Figure 8 below. Furthermore, as we did in previous Section 4.2.1, we also consider the D test statistic. The QQ-plot of D versus χ 2 m distribution is shown in right panel of Figure 8. Since the considered process Q is not Gaussian, a strong deviation is observed. It allows us to reject the H0 hypothesis. To quantify this deviation we present the generated p-values for the goodness-of-fit test associated to D in Table 4. 8 (right) versus a χ 2 m distribution with m = 6. The considered process Q(•) is a Kramer oscillator process as defined in Section 3.2 with E(Q(0)) = 0, Var(Q(0)) = 1 and second spectral moment λ = 2.

Bivariate numerical illustration

In this section, we focus on two dimensional random fields. For both Gaussian and chi-square distributions, we apply the same methodology as in the previous section. In particular, we use the test statistic D introduced in (20), its definition being not specific to dimension one.

Under H0 hypothesis

We consider a stationary centered Gaussian random field X = {X(t) : t ∈ R 2 }. Its restriction to a finite regular grid included can be seen as a model for a grey level image. The modified Euler characteristic of an excursion set in the rectangle domain T above level u is given by Equation ( 5). On the other hand, Equation ( 6) in dimension d = 2 gives its expectation,

E[φ(X, T, u)] = |T | C 2 (u, λ) with C 2 (u, λ) = (2π) -3/2 λ u e -u 2 /2 . ( 22
)
In what follows, we generate a 300 Monte-Carlo sample of a bivariate stationary centered Gaussian random field X with covariance function r(t) = e -||t|| 2 , t ∈ R 2 . In that case, E(X(0)) = 0, Var(X(0)) = 1 and the second spectral moment λ is equal to 2. We use [START_REF] Azaïs | Level Sets and Extrema of Random Processes and Fields[END_REF] in order to compute φ(X, T, u) for a fixed cube T and various values of u. The local extremum points of X are given by the R function extrema2dC in the EMD package. To identify the saddle points, we find all the stationary points of X in the considered domain T (i.e., all the points with an associated null gradient function) and we exclude the points previously identified as local extremum points.

In the left panel of Figure 9, the comparison between the theoretical mean given by ( 22) and the empirical mean based on the simulations is illustrated.

In the right panel, we consider the chi-square statistic D introduced in ( 20) for a unique level u and m disjoint domains T 1 , . . . , T m . The considered process is now the Gaussian bivariate process X. Furthermore, in In this case λ = 2.

Alternative: χ 2 bivariate process

As an alternative to hypothesis H0, let us consider

Z (s) (•) = 1 √ 2s (χ 2 s (•)) -s) with χ 2 s (•) = ∑ 1≤i≤s X i (•) 2 ,
where the X i 's are independent copies of a centered stationary Gaussian two dimensional field with covariance function r(t) = e -||t|| 2 /2 . As described in Section 3.1, for any fixed t ∈ R 2 , χ 2 s (t) has a central χ 2 distribution with s degrees of freedom. Furthermore, the field Z (s) is centered, stationary and its covariance function is given by r(t) = e -||t|| 2 . Hence its variance is equal to 1 and its second spectral moment is equal to λ = -2r ′′ (0) = 2. The expectation of φ(Z (s) , T, u) is given by ( 13) with

d = 2.
We generate a sample of 300 realisations of the bivariate process Z (s) for s = 2 on a fixed cube T . We use Equation ( 5) to empirically compute φ(Z (2) , T, u) for various values of u. On the other hand, Equation ( 13) yields the following in the case

d = s = 2, E[φ(Z (2) , T, u)] = |T | ( λ 4π ) (2u + 1) e -(u+1) 1 [0,∞) (u + 1). (23)
In the left panel of Figure 10, we compare Equation ( 23) and its empirical counterpart based on the Monte Carlo simulations of Z (2) . In the right panel, we show the QQ-plot of the test statistic D versus the χ 2 m distribution, where D is computed with the realisations of Z (2) . Since the considered chi-square process Z (2) is not Gaussian, a large deviation can be observed. Associated generated p-value for the goodness-of-fit test of D is given in Table 6. This quantifies how we are really able to reject the null hypothesis H0 under the alternative. 10 (right) versus a χ 2 m distribution with m = 3. The considered processZ (2) (•) is a chi-square bivariate process (d = 2) with covariance function r(t) = e -||t|| 2 . In this case E(Z (2) )(0) = 1, Var(Z (2) (0)) = 1 and λ = 2.

Non smooth processes: a case study

In this section, we sketch a preliminary study concerning non smooth processes. For this purpose, the shot noise process is a convenient toy model. We will only consider the one dimensional case pointing out two recent references that deal with dimension two, [START_REF] Biermé | Level total curvature integral: Euler characteristic and 2D random fields[END_REF] and [START_REF] Lachièze-Rey | Two-dimensional Kac-Rice formula. Application to shot noise processes excursions[END_REF]. Let d = 1 and let us introduce the following shot noise process S,

S(t) =   ∑ ξ∈Ψ 1 [0,a] (t -ξ)   -νa , t ∈ R , ( 24
)
where a > 0 and Ψ is a homogeneous Poisson point process on R with intensity ν > 0. The process S is clearly stationary with zero mean and variance equal to νa. Moreover, its values almost surely belong to the discrete set {k -νa, k ∈ N}.

Z (2) (•) is a chi-square bivariate process (d = 2) with covariance function r(t) = e -||t|| 2 . In this case E(Z (2) )(0) = 1, Var(Z (2) (0)) = 1 and λ = 2.
See associated goodness-of-fit test p-value in Table 6.

In what follows, we still aim at testing the Gaussianity hypothesis. If a whole trajectory of S could be observed, its discrete shape would clearly indicate the right choice between Gaussian and shot noise. We now assume that only excursion sets of S above a few levels are observed and we consider the number of "upcrossings". Note that the notion of upcrossings has to be properly defined in this context, since S is not continuous. Actually, it is proved in [START_REF] Biermé | A Fourier approach for the level crossings of Shot Noise processes with jumps[END_REF] (see also [START_REF] Hsing | On the intensity of crossings by a shot noise process[END_REF] for a different approach) that for any level u ∈ R \ {k -νa, k ∈ N} and any interval T ,

E[U (S, T, u)] = |T | ν e -νa ∑ k≥0 (νa) k k! 1 {k<u+νa<k+1} . ( 25
)
This formula, which obviously differs from [START_REF] Biermé | A Fourier approach for the level crossings of Shot Noise processes with jumps[END_REF], is an indicator that could enable to discriminate between a Gaussian or a shot noise model.

In the following, we generate 300 trajectories of such a process on a fixed interval T with a = ν = 1. Since aν = 1, then Var(S(0)) = 1. A comparison between the theoretical expectation and the Monte-Carlo empirical mean is presented in Figure 11. One can observe the evident different shape of the mean EC curve compared to the Gaussian case (see Figure 1, left).

Conclusion

In the present paper, we have proposed a methodology to test whether a random field defined on R d is Gaussian. The test statistics are computed from the observation of a single realisation, precisely the Euler characteristic of excursion sets at moderate levels is concerned. Our approach requires two theoretical ingredients that have to be valid under the Gaussian hypothesis: a Central Limit Theorem satisfied by the Euler characteristic of excursion sets when the domain tends to R d , and a closed formula that gives the mean Euler characteristic of the excursion sets. We have established both of them. We consider this methodological work as a first step, there are many open questions and possible future researches extending and going further this work.

There are two major domains where the Gaussianity hypothesis is really relevant and where our results could certainly be applied: neurology and cosmology. The neurological signals that are collected when studying brain activity are commonly observed through their level sets (see [START_REF] Worsley | Estimating the number of peaks in a random field using the Hadwiger characteristic of excursion sets, with applications to medical images[END_REF] and [START_REF] Lindner | Analytical approach to the stochastic FitzHugh-Nagumo system and coherence resonance[END_REF]). A natural question concerns the type of the signal distribution: Gaussian or Poisson? Gaussian or Oscillator?

The test that we have built, together with the alternatives that we have considered, could be an appropriate tool for studying this problem. On the other hand, in the analysis of the Cosmic Microwave Background (CMB) radiation, the question of Gaussianity has been tackled in a huge amount of publications (see [START_REF]Isotropy and statistics of the CMB[END_REF] for a recent overview). In that context, the random field under study is defined on the two dimensional celestial sphere. Hence, in order to be applied, our methodology should be extended beyond the Euclidean case. Let us mention that recent studies, like [START_REF] Cheng | Excursion probability of Gaussian random fields on sphere[END_REF], [START_REF] Marinucci | High-Frequency Asymptotics for Lipschitz-Killing Curvatures of Excursion Sets on the Sphere[END_REF], [START_REF] Cammarota | A Quantitative Central Limit Theorem for the Euler-Poincare Characteristic of Random Spherical Eigenfunctions[END_REF] for instance, contain theoretical results on high energy behavior of spherical random fields that could provide the required background for a test of Gaussianity in the same spirit as the one presented here.

Furthermore, it is not difficult to imagine that the same question of Gaussianity is actually asked in other real-life situations, and that the potential alternatives take various shapes, like a chi-square when the process under study is always positive, or an oscillator when the process under study presents an almost periodic structure (sea waves) or a shot noise process when sudden peaks are observed. Our methodology could also serve as a goodness-of-fit test in the opposite way. For instance, in the area of geostatistics, the authors of [START_REF] Gabriel | Estimating and testing zones of abrupt change for spatial data[END_REF] observe the Euler characteristic of excursion sets to detect abrupt changes on soil data sets. They use as an a priori model for the data a chi-square random field. A natural extension to this study could be a test of chi-square distribution versus not chi-square.

In order to go further in this direction, it would be necessary to establish Central Limit Theorems for the Euler characteristic of excursion sets for new classes of processes or random fields. Beside the present work, we implemented simulations that lead us to postulate that such a Central Limit Theorem could exist, in particular for the two alternative processes that we have considered.

At last, another appealing question concerns the hypothesis that the process under study is expected to satisfy, such as smoothness of the realisations and fast decay of the covariance function. As already mentioned at the beginning of the paper, Assumption (A) seems to strong at least in the one dimensional case. Indeed, let us remark that in that case, two important results of Section 2 remain valid under weaker assumptions. Actually, on the one hand, the Rice formula gives the mean number of upcrossings under the only assumption of a finite second spectral moment. On the other hand, it is well known that a sufficient condition for the number of upcrossings to have a finite second moment is the so-called Geman condition, which only requires that r ′′ , the second derivative of the covariance function of X, exists and satisfies some integrability condition near 0. Such considerations allows us to believe that the theoretical results that we have presented deserve to be studied under weaker assumptions, even in dimension greater than one. Another extension could concern non continuous random fields, like the shot noise processes.

Appendix: Proofs

A.1. Proof of Lemma 3.

The level u is fixed, so we do not refer to it in the following lines. For t ∈ R d \ {0}, we define h(t) = G(u, t) D(t) -1/2 -C(u) 2 . We have to establish that h is integrable on R d . We already know that h is integrable on any compact set of R d , so it remains to study the behavior of h(t) for large ∥t∥.

We write K = 1 2 d(d + 1) + 1 and we introduce the function

(z, x) ∈ R K = R 1 2 d(d+1) × R → f (z, x) = det(z) 1 [u,∞) (x) ,
where det(z) stands for the determinant of the d × d symmetric matrix whose upper coefficient are given by the 1 2 d(d + 1) dimensional vector z. Recall that we denote by X(.) the stationary and isotropic random vector field (X ′ (t), ∇ 2 X(t), X(t)), which has dimension

D = d + 1 2 d(d + 1) + 1 = d + K. We have h(t) = ∫ R 2K f (ẍ 1 , x 1 )f (ẍ 2 , x 2 ) ( p X(0),X(t) (0, ẍ1 , x 1 , 0, ẍ2 , x 2 ) -p X(0) (0, ẍ1 , x 1 )p X(t) (0, ẍ2 , x 2 ) ) dx 1 dẍ 1 dx 2 dẍ 2 .
Let Σ be the covariance matrix of X(0) and, for any t ∈ R d , let Γ X (t) be the covariance matrix of (X(0), X(t)). The latter can be written by blocks, each of one having dimension D × D, as follows

Γ X (t) = ( Σ Γ 12 (t) Γ 12 (t) T Σ
) .

Moreover we know that all the terms of matrix Γ X (t) are uniformly bounded in absolute value by ψ(t). Using Theorem 3.1 of [START_REF] Withers | Expansions for the multivariate normal[END_REF] we get, for any

z 1 , z 2 in R D , p X(0),X(t) (z 1 , z 2 ) -p X(0) (z 1 )p X(t) (z 2 ) = p X(0) (z 1 )p X(0) (z 2 )   ∑ J∈N D×D ,|J|>0 H A J (z 1 , Σ)H B J (z 2 , Σ) Γ 12 (t) J J!   ,
where

• the sum runs over all D × D matrices J = (J jk ) 1≤j,k≤D with non negative integer entries and

|J| = ∑ 1≤j,k≤D J jk > 0, • A J and B J are D-dimensional indices defined by (A J ) j = ∑ 1≤k≤D J jk and (B J ) k = ∑ 1≤j≤D J jk , • H k (z, Σ) = ϕ Σ (z) -1 (-∂ ∂z ) k ϕ Σ (z)
for any D-dimensional index k and any z ∈ R D , with ϕ Σ the probability density function of any D-dimensional Gaussian vector N (0, Σ),

• M J = ∏ 1≤j,k≤D (M jk ) J jk for any matrix M = (M jk ) 1≤j,k≤D , • and J! = ∏ 1≤j,k≤D (J jk )!. Thus h(t) = ∫ R 2K f (x 1 ) f (x 2 ) p X(0) (0, x 1 ) p X(0) (0, x 2 )   ∑ J,|J|>0 H A J ((0, x 1 ), Σ) H B J ((0, x 2 ), Σ) Γ 12 (t) J J!   dx 1 dx 2 .
We can bound |h(t)| in the following form

|h(t)| ≤ ∑ J,|J|>0 ψ(t) |J| J! D(A J ) D(B J ), (26)
where

D(A J ) = ∫ R K |f (x 1 )H A J ((0, x 1 ), Σ)| p X(0) (0, x 1 ) dx 1 ≤ E ( | det(∇ 2 X(0)) H A J ((0, ∇ 2 X(0), X(0)), Σ)| ) ≤ E( det(∇ 2 X(0) 2 ) 1/2 E ( H A J ((0, ∇ 2 X(0), X(0)), Σ) 2 ) 1/2 ≤ ψ(0) d E ( H A J ((0, ∇ 2 X(0), X(0)), Σ) 2 ) 1/2 , ( 27
)
where the last inequality is a consequence of expanding det(∇ 2 X(0)) as a multivariate polynomial function of degree d evaluated at the coordinates of the matrix X ′′ (0). We now concentrate on H A J ((0, ∇ 2 X(0), X(0)), Σ). Following [START_REF] Withers | Expansions for the multivariate normal[END_REF], we have

H A J ((0, x), Σ) = E [ ( Σ -1 ((0, x) + iZ) ) A J ] ,
where Z is any D-dimensional centered Gaussian random vector with covariance matrix Σ. Remember that Σ is the covariance matrix of X(0) = (X ′ (0), ∇ 2 X(0), X(0)), so it can be factorized as Σ = ( λI d 0 0 Σ 1

)
where Σ 1 has size

K × K. Hence Σ -1 = ( λ -1 I d 0 0 Σ -1 1 )
and Z can be expanded as Z = (Z 0 , Z 1 ) where Z 0 is a d-dimensional centered Gaussian vector with covariance matrix λI d and Z 1 is a K-dimensional centered Gaussian vector with covariance matrix Σ 1 , Z 0 and Z 1 being independent. Then Σ -1 ((0, x)

+ iZ) = ( iλ -1 Z 0 , Σ -1 1 (x + iZ 1 )
) ,

where the two blocks of coordinates are independent. So

H A J ((0, x), Σ) = H A (0) J (0, λI d ) H A (1) J (x, Σ 1 ),
where the D-dimensional vector A J is equal to (A (0)

J , A (1)
J ) with blocks of respective size d and K. For the first term, we note that H A (0)

J (0, λI d ) = λ -1 2 |A (0) J | H A (0) J (0)
where H k (.) is the usual multidimensional Hermite polynomial of multi-order k. Concerning the second term, we compute E

( H A (1) J ((∇ 2 X(0), X(0)), Σ 1 ) 2 )
thanks to a formula for multivariate Gaussian integrals that can be found in [START_REF] Kauderer | Modes in n-dimensional first-order systems[END_REF] for instance (see equation [START_REF] Worsley | Estimating the number of peaks in a random field using the Hadwiger characteristic of excursion sets, with applications to medical images[END_REF] in this reference),

E ( H A (1) J ((∇ 2 X(0), X(0)), Σ 1 ) 2 ) = ∫ R K H A (1) J (x, Σ 1 ) 2 p ∇ 2 X(0),X(0) (x) dx = 2 |A (1) J |-D (A (1)
J )! diag(Σ -1 1 ) A (1)

J

where diag(B) denotes the vector containing the diagonal coefficients of any square matrix B.

Finally, since H A (0)

J (0) 2 ≤ (A (0) J )!, E ( H A J ((0, ∇ 2 X(0), X(0)), Σ) 2 ) ≤ (2∥Σ -1 ∥) |J| (A J )!
where we have used that |A (k)! is a multinomial coefficient). Hence (A J ) j ! ≤ D ∑ 1≤k≤D J jk ∏ 1≤k≤D J jk ! and so

(A J )! = ∏ 1≤j≤D (A J ) j ! ≤ D |J| J!
Coming back to inequality [START_REF] Rice | Mathematical Analysis of Random Noise[END_REF] yields

D(A J ) ≤ ψ(0) d (2D∥Σ -1 ∥) |J|/2 √ J! .
The same inequality can be established for D(B J ). Hence, from (26) we get

|h(t)| ≤ ψ(0) 2d ψ(t) ∑ J,|J|>0 (2D∥Σ -1 ∥ ψ(t)) |J|-1 .
The series in the right hand is convergent for ∥t∥ large enough. Indeed we can choose A such that for ∥t∥ > A, we have 2D∥Σ -1 ∥ ψ(t) ≤ η < 1 and so the series is bounded by (1 -η) -1 . The result follows since on the other hand, by assumption (A), ∫ ∥t∥>A ψ(t)dt < +∞.

A.2. Proof of Lemma 6.

Without loss of generality, we can assume that |T . Hence, by stationarity

|Cov(Z (N ) 1 , Z (N ) 2 )| = N -d | ∑ s,s ′ ∈IN E ( φ(X, [0, 1) d , u 1 ) φ(θ N (b-a)+s ′ -s X, [0, 1) d , u 2 ) ) | ≤ ∑ s∈J N |E ( φ(X, [0, 1) d , u 1 ) φ(θ N (b-a)+s X, [0, 1) d , u 2 ) ) | , with J N = [[-N + 1, N -1]] d .
It is proved in [START_REF] Estrade | A Central Limit Theorem for the Euler characteristic of a Gaussian excursion set[END_REF] Proposition 1.3 that, for a cube T in R d and a level u in R, the following expansion holds both a.s. and in L 2 (Ω).

φ(X, T, u) = (-1) d ∞ ∑ q=1 ∑ n∈N D ;|n|=q a(n, u) ∫ T H n (Y (t)) dt (28)
where D = d + 1 2 d(d + 1) + 1, H n denotes the D dimensional Hermite coefficient of multi-order n ∈ N D and Y (t) = Λ -1 X(t) with Λ any square root of the covariance matrix of X(0). Moreover the Hermite coefficients a(n, u) are such that for any q ≥ 1,

∑ n∈N D ;|n|=q a(n, u) 2 n! ≤ C q d ( 29
)
where C is some positive constant that only depends on d, Λ and u. Hence, thanks to [START_REF] Roberts | Chord-distribution functions of three-dimensional random media: approximate first-passage times of Gaussian processes[END_REF] and to the orthogonality of Hermite polynomials, we have for any

τ ∈ R d , E ( φ(X, [0, 1) d , u 1 ) φ(θ τ X, [0, 1) d , u 2 ) ) = ∞ ∑ q=1 ∫ [0,1] d ∫ [0,1] d E (F q (Y (t 1 ), u 1 )F q (Y (τ + t 2 ), u 2 )) dt 1 dt 2 ,
where

F q (y, u) = ∑ n∈N D ;|n|=q a(n, u)H n (y) , y ∈ R D .
One can adapt Arcones inequality, for instance following the proof of [START_REF] Arcones | Limit theorems for nonlinear functionals of a stationary Gaussian sequence of vectors[END_REF] Lemma 1 step by step, to establish that there exists a constant κ > 0 depending on the covariance function of X such that

|E (F q (Y (t 1 ), u 1 )F q (Y (τ + t 2 ), u 2 ))| ≤ κ q ψ(τ + t 2 -t 1 ) q ∑ n∈N D ; |n|=q 1 2 (a(n, u 1 ) 2 + a(n, u 2 ) 2 )n!.
Let now η > 0 be fixed. Thanks to Assumption (A), we can choose N large enough such that κψ(τ ) ≤ 1/2 as soon as ∥τ ∥ ≥ N δ, ∫ ∥t∥≥N ψ(t) dt < η and N δ > 2. Hence, using (29), we deduce that for ∥τ ∥ ≥ N δ,

E ( φ(X, [0, 1) d , u 1 ) φ(θ τ X, [0, 1) d , u 2 ) ) ≤ C ∞ ∑ q=1 q d κ q ∫ [0,1] d ∫ [0,1] d ψ(τ + t 2 -t 1 ) q dt 1 dt 2 ≤ C ( ∞ ∑ q=1 q d (1/2) q-1 ) ∫ [-1,1] d ψ(τ + t) dt .
Note that C stands here for a positive constant that may change from a line to another. Coming back to the covariance of Z (N ) 1

and

Z (N ) 2 yields |Cov(Z (N ) 1 , Z (N ) 2 )| ≤ C ∑ s∈JN ∫ [-1,1] d ψ(N (b -a) + s + t) dt ≤ 2C ∫ ∥t∥≥N ψ(t) dt ≤ 2C η.
Thus, Lemma 6 is proved.

A.3. Proof of Proposition 7.

Firstly, let us remark that, since Z (s) (.) = 1 √ 2s (χ 2 s (.) -s), it is sufficient to prove the statement of Proposition 7 for the process χ 2 s instead of Z (s) . Secondly, let us recall the relation ( 4), so that we can replace φ(Z (s) , [0, T ], u) by U (χ 2 s , [0, T ], u) in the proposition.

The finiteness of V s (T, u)

:= E[U (χ 2 s , [0, T ], u)(U (χ 2 s , [0, T ], u) -1)
] is proved in [START_REF] Sharpe | Some properties of the crossings process generated by a stationary χ 2 process[END_REF] under the following conditions,

r X (t) = 1 -λ 2 t 2 2 + λ 4 t 4 4 + o(t 4 ) as t → 0 and r X (t) = o(t -α ) as t → ∞ (30)
for some positive finite constants λ 2 , λ 4 , α (see (4.1) and (4.2) in the cited paper). They are clearly satisfied in our case since Assumption (A) is in force. This establishes that U (χ 2 s , [0, T ], u) admits a finite variance.

We now focus on the asymptotic variance as T goes to infinity. We will actually prove that lim

T →+∞ 1 T ( V s (T, u) -(E[U (χ 2 s , [0, T ], u)]) 2 ) < + ∞.
First, in a similar way as (13), let us remark that the expectation of

U (χ 2 s , [0, T ], u) is proportional to T , K s (u) := 1 T E[U (χ 2 s , [0, T ], u)] = ( λ π ) 1/2 e -u/2 u (s-1)/2 2 s/2 Γ(s/2) . ( 31
)
On the other hand, introducing the function G(x) = x1 [0,∞) (x), Rice formula for the above factorial moment gives

V s (T, u) = 2 ∫ T 0 (T -t) I(t) dt,
where

I(t) = E[G((χ 2 s ) ′ (0)) G((χ 2 s ) ′ (t)) | χ 2 s (0) = χ 2 s (t) = u] p χ 2 s (0),χ 2
s (t) (u, u). For computing this conditional expectation, we begin with the following one that is easier to handle with. Let us define ⃗ X(t) = (X 1 (t), . . . , X s (t)) and, for z

= (z 1 , • • • , z s ) and w = (w 1 , • • • , w s ) in R s , J(t, z, w) := E[G((χ 2 s ) ′ (0)) G((χ 2 s ) ′ (t)) | ⃗ X(0) = z, ⃗ X(t) = w, χ 2 s (0) = χ 2 s (t) = u].
Let us fix z and w such that ||z|| 2 = ||w|| 2 = u. In that case, we have

J(t, z, w) = 4 E[G( s ∑ j=1 X j (0)X ′ j (0)) G( s ∑ j=1 X j (t)X ′ j (t)) | ⃗ X(0) = z, ⃗ X(t) = w].
We will use, as is classical, the following regression model. For each j ∈ {1, . . . , s}, we write

X ′ j (0) = α 1 (t)X j (0) + β 1 (t)X j (t) + ε 1j (t) X ′ j (t) = α 2 (t)X j (0) + β 2 (t)X j (t) + ε 2j (t).
It is important to point out that, because of the independence, the coefficients above do not depend on j. We readily obtain that

α 1 (t) = r ′ X (t)r X (t) 1 -r 2 X (t) , β 1 (t) = -r ′ X (t) 1 -r 2 X (t) , β 2 (t) = -α 1 (t), α 2 (t) = -β 1 (t). Thus J(t, z, w) = 4E[G( s ∑ j=1 z j (α 1 (t)z j + β 1 (t)w j + ε 1j (t))) G( s ∑ j=1 w j (β 1 (t)z j + α 1 (t)w j + ε 2j (t)))] = 4E[G(α 1 (t)u + β 1 (t)⟨w, z⟩ + √ uσ ε W (t)) G(α 1 (t)u + β 1 (t)⟨w, z⟩ + √ uσ ε Z(t))] = 4E[G(d u (t) + √ uσ ε (t)W (t)) G(d u (t) + √ uσ ε (t)Z(t))],
where

σ 2 ε (t) = -r ′′ X (0) -(α 2 1 (t) + β 2 1 (t) + 2α 1 (t)β 1 (t)r X (t)) is the variance of each ε 1,j (t) or ε 2,j (t), W (t) = 1 √ uσ ε (t) s ∑ i=1 w j ε 1,j (t) and Z(t) = 1 √ uσ ε (t) s ∑ i=1 z j ε 2,j (t)
are independent N (0, 1) random variables, and d u (t) = α 1 (t)u + β 1 (t)⟨w, z⟩.

In particular, the latter shows that J(t, z, w) only depends on w and z through their scalar product. So, from now on, for ease of notation, we will write J(t, ⟨z, w⟩).

We have

β 1 (t) → ∞ when t → 0 and d u (t) = (α 1 (t) + β 1 (t))u -β 1 (t)u(1 -1 u ⟨w, z⟩) → t→0 -∞, for all w ̸ = z. The latter convergence is equivalent to -1 t . Moreover, d u (t) → 0 whenever t → ∞,

because of Assumption (A).

Then we can write

J(t, ⟨z, w⟩) = 4 u σ 2 ε (t) E[G( d u (t) √ uσ ε (t) + W (t)) G( d u (t) √ uσ ε (t) + Z(t))] = 4 u σ 2 ε (t) ∞ ∑ k=0 d 2 uk (t)k! ( Corr(W (t), Z(t))) k , ( 32
)
where the last equality is a consequence of Mehler's formula.

) 33
Thus, for 0 < t < δ, we have J(t, ⟨z, w⟩) ≤ C u σ 2 ε (t). Moreover, as this bound only depends on u, it yields immediately that, I(t) ≤ C u σ 2 ε (t) , for 0 < t < δ.

We now focus on the behaviour of I(t) for t close to +∞. We have ) sin s-2 ψ dψ.

I(t) = ∫ R s ×R s J(t,
Following [START_REF] Taheriyoun | A note on the higher moments of the Euler characteristic of the excursion sets of random fields[END_REF], we now expand J(t, x) as J(t, x) := 4 u σ 2 ε (t)

∑ ∞ k=0 Q k (t, x), so that 1 T V s (u, T ) = 2 ∫ T 0 (1 - t T ) I(t) dt = ∫ T 0 (1 - t T ) 4u s-1 σ 2 ε (t) (2π) s (1 -r 2 X (t)) s/2 σ s-1 (S s-1 ) 2 W s-2 × ∫ π 0 ∞ ∑ k=0 Q k (t, u cos ψ 1 )e -u( 1-r X (t) cos ψ 1-r 2 X (t)
)

sin s-2 ψ dψ dt.

In first place, let us study the limit when t → ∞ of the (k = 0)-term in the integrand. We introduce the following notation,

I 0 (t) = 4u s-1 σ 2 ε (t) (2π) s (1 -r 2 X (t)) s/2 σ s-1 (S s-1 ) 2 W s-2 ∫ π 0 Q 0 (t, u cos ψ 1 )e -u( 1-r X (t) cos ψ 1-r 2 X (t)
) sin s-2 ψ dψ.

We have

I 0 (t) → t→∞ 4u s-1 (-r ′′ X (0)) (2π) s+1
e -u σ s-1 (S s-1 ) 2 = K s (u) and using [START_REF] Taheriyoun | Testing the covariance function of stationary Gaussian random fields[END_REF]. Therefore, recalling the identity

) 34
This last equality holds true if we show that the difference appearing into the integrand is bounded by an L 1 (R) function, outside of a compact interval. Actually, under (A), it is easy to prove that for t large enough, |I 0 (t) -K s (u) 2 | ≤ C| du(t) √ uσε(t) | ≤ Cψ(t). Hence, the limit (34) is established.

In the sequel we are going to study the asymptotic behaviour of the remaining terms. Let us introduce

R(t) := 4u s-1 σ 2 ε (t) (2π) s (1 -r 2 (t)) s/2 σ s-1 (S s-1 ) 2 2W s-2 × ∫ π 0 ∞ ∑ k=1 Q k (t, u cos ψ 1 )e -u(

Computation of g(u)

. By Equation ( 8), we have We make the following change of variable, ẍ+λx (γ-λ 2 ) 1/2 = y, to get

g(u) = E[ 1 [u,∞) (X(0)) |X ′′ (0)| ] = 1 2π(γ -λ 2 ) 1/2 ∫ ∞ u ∫ R |ẍ|e - 1 2(γ-λ 2 ) (γx 2 +2λxẍ+ẍ 2 ) dẍdx = 1 2π(γ -λ 2 ) 1/2 ∫ ∞ u e -1
I 1 (x, λ, γ) = 1 √ 2π ∫ ∞ λx (γ-λ 2 ) 1/2 ((γ -λ 2 ) 1/2 y -λx)e -1 2 y 2 dy = (γ -λ 2 ) 1/2 ϕ( λx (γ -λ 2 ) 1/2 ) -λx Φ( λx (γ -λ 2 ) 1/2 ) I 2 (x, λ, γ) = (γ -λ 2 ) 1/2 ϕ( λx (γ -λ 2 ) 1/2 ) + λx Φ(- λx (γ -λ 2 ) 1/2 ).
But it holds that Φ(x) = Φ(-x). Thus I 1 (x, λ, γ) + I 2 (x, λ, γ) = 2(γ -λ 2 ) 1/2 ϕ( λx (γ -λ 2 ) 1/2 ) + λx (2Φ( λx (γ -λ 2 ) 1/2 ) -1).

Summing up, we have

g(u) = 1 √ 2π ∫ ∞ u e -1 2 x 2 [ 2(γ -λ 2 ) 1/2 ϕ( λx (γ -λ 2 ) 1/2 ) + λx (2Φ( λx (γ -λ 2 ) 1/2 ) -1) ] dx.
Integrating by parts the second integral, we get that it is equal to

λϕ(u)(2Φ( λu (γ -λ 2 ) 1/2 ) -1) + 2λ 2 √ 2π(γ -λ 2 ) 1/2 ∫ ∞ u e -1 2 x 2 ϕ( λx (γ -λ 2 ) 1/2 )dx, yielding g(u) = λϕ(u)(2Φ( λu (γ -λ 2 ) 1/2 ) -1) + 2γ √ 2π(γ -λ 2 ) 1/2 ∫ ∞ u e -1 2 x 2 ϕ( λx (γ -λ 2 ) 1/2 )dx.
The second term is equal to

2γ 2π(γ -λ 2 ) 1/2 ∫ ∞ u e -1 2 x 2 e -λ 2 x 2 2(γ-λ 2 ) dx = γ π(γ -λ 2 ) 1/2 ∫ ∞ u e -γx 2 2(γ-λ 2 ) dx = ( 2γ π ) 1/2 Φ(( γ (γ -λ 2 )
) 1/2 u).

Finally,

g(u) = λϕ(u)(2Φ( λu (γ -λ 2 ) 1/2 ) -1) + ( 2γ π ) 1/2 Φ(( γ (γ -λ 2 )
) 1/2 u).

Computation of G(u, t).

Since Γ(t) is the covariance matrix of (X(0), X(t), X"(0), X"(t)) conditioned to (X ′ (0) = X ′ (t) = 0), one can write

G(u, t) = E(1 Z1>u 1 Z2>u Z 3 Z 4 ) where Z = (Z 1 , Z 2 , Z 3 , Z 4 ) L = N (0, Γ(t)). (36)
From now on, we remove the dependence on t. We begin by writing Z = L Y, where L is a lower triangular matrix such that LL T = Γ(t). Denoting by l ij for i ≥ j the elements of L, we have

Z 1 = l 11 Y 1 Z 2 = l 21 Y 1 + l 22 Y 2 Z 3 = 3 ∑ i=1 l 3i Y i Z 4 = 4 ∑ i=1 l 4i Y i ,
and the expectation (36) can be written as

G(u, t) = ∫ R 4 1 {l11y1>u} 1 {l21y1+l22y2>u} ( 3 ∑ i=1 l 3i y i )( 4 ∑ i=1 l 4i y i )ϕ(y)dy.
By expanding the second sum, the integral can be written as the sum of the following terms (starting with index i = 4 term)

I 4 = l 44 ∫ R 4

Fig 1 .

 1 Fig 1. Left: Theoretical u → E[φ(X, T, u)] from Equation (17) for |T | = 200 (full line). We also display, for different levels u, the empirical counterpart E[φ(X, T, u)] (red dots) based on 300 Monte Carlo simulations. We also display the associated empirical intervals. The map u → |T | C(u, λ) from Equation (17), by using the estimated spectral moment λ as in Section 4.1.2, is represented in dashed line. Right: Theoretical u → V (u) from Equation (18) for various values of u (black triangles) and empirical variance of |T | -1/2 φ(X, T, u) (red dots) based on 300 Monte Carlo simulations with |T | = 200 and the same values of u. X(•) is a Gaussian univariate process (d = 1) with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 . In this case λ = 2.

  Fig 2.Boxplot, based on 300 Monte Carlo simulations, of the ratio between theoretical value λ 1/2 and estimated value γ given by[START_REF] Lelievre | Free Energy Computations[END_REF]. X(•) is a Gaussian univariate process (d = 1) with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 . In this case λ = 2. We choose p = 3, c 1 = c 2 = c 3 = 1 3 and u = (-2/3, 0, 2/3).

( a ). 1 )(a. 2 ). 1 )(b. 2 )

 a1212 Diagonal case. Each level u k is associated to a single domain T k . We have m observations Y 1 1 , . . . , Y m m in the statistical model (11). (aCase with same level u, i.e. u 1 = . . . = u m = u, and m disjoint domains T 1 , . . . , T m . Case with different levels u 1 < . . . < u m and m disjoint domains T 1 , . . . , T m . (b) Crossed case. Different levels u k can be associated to the same domain T i . We have mp observations (Y i k ) 1≤i≤m,1≤k≤p in (11). (bCase with one single domain T (m = 1) and p different levels u 1 < . . . < u p . Case with m disjoint domains T 1 , . . . , T m and p different levels u 1 < . . . < u p .

Fig 3 .

 3 Fig 3. QQ-plot based on 300 Monte Carlo simulations. First and second panels: quantiles of F a1 , F a1 versus quantiles of the χ 2 m distribution for m = 3 disjoint domains T 1 , T 2 , T 3 with |T i | = 200 and a single level u = 1.2. Third and fourth panels: quantiles of F a2 , F a2 versus quantiles of the χ 2 m distribution for m = 3 disjoint domains T 1 , T 2 , T 3 with |T i | = 200 and different levels u1 = -1.2, u 2 = 0, u 3 = 1.2. X(•) is a Gaussian univariate process (d = 1)with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 with λ = 2. See associated goodness-of-fit test p-values in Table1.

  For the crossed case (b.1), let us now consider different levels u 1 < . . . < u p and one single domain T , i.e. m = 1. Let us define

  the empirical mean on the Monte-Carlo simulations and the matrices Λ and Λ are defined below. Let Γ = (V (u k , u l )) 1≤kl,≤p be the (theoretical) covariance matrix of Z and let Λ stand for any square root of Γ. Similarly, let Γ be the empirical covariance matrix of Z evaluated on the Monte-Carlo sample generations, and let Λ be any of its square root matrix. Hence, F b1 and F b1 are both approximately χ 2 m p distributed, with m = 1. An illustration of the behaviour of F b1 is presented in Figure4(first panel), by choosing 300 Monte-Carlo simulations, m = 1, |T | = 200 and p = 3 different levels u 1 = -1.5, u 2 = 0, u 3 = 1.5.

5 )Fig 4 .

 54 Fig 4. QQ-plot based on 300 Monte Carlo simulations. First panel: quantiles of F b1 versus quantiles of the χ 2 p distribution for m = 1, with |T | = 200 and p = 3 different levels u 1 = -1.5, u 2 = 0, u 3 = 1.5. Second panel: quantiles of F b2 versus quantiles of the χ 2 mp distribution for m = 3 disjoint domains T 1 , T 2 , T 3 with |T i | = 200 and p = 3 different levels u 1 = -1.5, u 2 = 0, u 3 = 1.5. X(•) is a Gaussian univariate process (d = 1) with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 with λ = 2. See associated goodness-of-fit test p-values in Table2.

2 Fig 5 .

 25 Fig 5. First and third panels: Boxplot for the ratio between the empirical 300 Monte Carlo values of φ(Z (s) , T, u) and the theoretical mean given by Equation[START_REF] Lindner | Analytical approach to the stochastic FitzHugh-Nagumo system and coherence resonance[END_REF]. Second and fourth panels: Boxplot for the ratio between the same empirical values and the theoretical mean in the Gaussian case given by Equation[START_REF] Kauderer | Modes in n-dimensional first-order systems[END_REF]. In both cases λ = 2. The degrees of freedom s is chosen equal to 2 and 10 respectively.

Fig 6 .

 6 Fig 6. Boxplots of the empirical 300 Monte Carlo values of φ(Z (s) , T, u) for the considered chi-square univariate process Z (s) and different values of u. Red points represent the theoretical means given by Equation (21) for the same values of u; blue ones are the Gaussian means given by Equation (17). In this case λ = 2, |T | = 200, s = 2 (left panel) and s = 10 (right panel).

Fig 7 . 2 m

 72 Fig 7. QQ-plot: quantiles of D versus quantiles of the χ 2 m distribution for fixed level u = 1.2 and m = 6 disjoint domains. Here |T i | = 200 for i ∈ {1, . . . , m}. The dashed red line represents the bisector of the first orthant. We take 300 Monte Carlo simulations. The considered process Z (s) (•) is a chi-square univariate process (d = 1) with s = 2 (first panel), s = 10 (second panel), E(Z (s) (0)) = 0, Var(Z (s) (0)) = 1 and λ = 2. See associated goodness-of-fit test p-values in Table3.

2 Fig 8 . 2 m

 282 Fig 8. Left: Theoretical u → E[φ(Q, T, u)] from Equation (15) for |T | = 200 (full line). We also display, for different levels of u, the empirical counterpart E[φ(Q, T, u)] (red dots), with associated empirical intervals, based on 300 Monte Carlo simulations. Right: QQ-plot: quantiles of D versus quantiles of the χ 2 m distribution for fixed level u = 1.2 and m = 6 disjoint domains. The dashed red line represents the bisector of the first orthant. Here |T i | = 200 for i ∈ {1, . . . , m}. We take 300 Monte Carlo simulations. Q(•) is a Kramer oscillator process as defined in Section 3.2 with E(Q(0)) = 0, Var(Q(0)) = 1 and second spectral moment λ = 2. See associated goodness-of-fit test p-value in Table4.

2 Fig 9 . 2 m 5

 2925 Fig 9. Left: Theoretical u → E[φ(X, T, u)] from Equation[START_REF] Lobato | A simple test of normality for time series[END_REF] for |T | = 196 (full line). We also display, for different levels u, the empirical counterpart E[φ(X, T, u)] (red dots), with associated empirical intervals, based on 300 Monte Carlo simulations. Right: QQ-plot: quantiles of D versus quantiles of the χ 2 m distribution for fixed level u = 1.2 and m = 3 disjoint domains. Here |T i | = 196 for i ∈ {1, . . . , m}. We take 300 Monte Carlo simulations. X(•) is a Gaussian bivariate process (d = 2) with covariance function r(t) = e -||t|| 2 , t ∈ R 2 . In this case E(X(0)) = 1, Var(X(0)) = 1 and λ = 2. See associated goodness-of-fit test p-value in Table5.

GOF test p-values for χ 2 m 9

 29 Cramér-von Mises test Anderson-Darling testD

Fig 10 .

 10 Fig 10. Left: Theoretical u → E[φ(Z(2) , T, u)] from Equation[START_REF] Marinucci | High-Frequency Asymptotics for Lipschitz-Killing Curvatures of Excursion Sets on the Sphere[END_REF] for |T | = 196 (full line). We also display, for different levels u, the empirical counterpart E[φ(Z(2) , T, u)] (red dots), with associated empirical intervals, based on 300 Monte Carlo simulations. Right: QQ-plot: quantiles of D versus quantiles of the χ 2 m distribution for fixed level u = 3 and m = 3 disjoint domains. Here |T i | = 196 for i ∈ {1, . . . , m}. We take 300 Monte Carlo simulations. Z(2) (•) is a chi-square bivariate process (d = 2) with covariance function r(t) = e -||t|| 2 . In this case E(Z(2) )(0) = 1, Var(Z(2) (0)) = 1 and λ = 2. See associated goodness-of-fit test p-value in Table6.

Fig 11 .

 11 Fig 11. Theoretical u → E[U (S, T, u)] from Equation (25) for |T | = 200 (full line). We also display, for different levels u, the empirical counterpart E[U (S, T, u)] (red dots) with associated empirical intervals based on 300 Monte Carlo simulations. S(•) is a shot noise univariate process (d = 1) defined by Equation (24) with a = ν = 1.

J

  )! = (A J )!. Note that for any D-dimensional index k = (k 1 , . . . , k D ), we have |k|! ≤ D |k| (k)! (actually, |k|!

1 = N -d/ 2 ∑

 12 1 | = |T 2 | = 1 and that there exists a and b in R d such that T 1 = a + [0, 1) d and T 2 = b + [0, 1) d . Moreover, since the distance between T 1 and T 2 is non negative, we have ∥b -a∥ ∞ > 1. In the following lines, we denote δ = ∥b -a∥ ∞ -1. We use the additivity of the Euler characteristic to expand Z s∈IN φ(θ N a+s X, [0, 1) d , u 1 ) , where φ(X, T, u) = φ(X, T, u) -E[φ(X, T, u)], θ s stands for the shift operator and I N for the set of integer valued indices [[0, N -1]] d . The same holds for Z (N ) 2

2 = 2 ∫ ∞ 0 (

 220 (t)dt -T K s (u) I 0 (t) -K s (u) 2 ) dt. (

)

  sin s-2 ψ dψ. So, by using[START_REF] Withers | Expansions for the multivariate normal[END_REF], we obtain that for t large enough,| ∞ ∑ k=1 Q k (t, u cos ψ 1 )| ≤ | ( Corr(W (t), Z(t)) ) | ∞ ∑ k=1 d 2 uk (t)k! ≤ Cψ(t).

( 2 ∫ ∞ 0 ( 2 ∫ ∞ 0 RA. 4 .

 20204 V s (T, u) -(E[U (χ 2 s , [0, T ], u)]) 2 ) = I 0 (t) -K s (u) 2 )dt + (t)dt < ∞. Proof of Proposition 8.

2 dẍ 2 )

 22 Let us split the inner integral into ẍ > 0 and ẍ < 0, i.e.1 √ 2π(γ -λ 2 ) (ẍ+λx) dẍ := I 1 (x, λ, γ) + I 2 (x, λ, γ),whereI 1 (x, λ, γ) = 1 √ 2π(γ -λ 2 ) -λ 2 ) (ẍ-λx) 2 dẍ = I 1 (x, -λ, γ).

1 {l11y1>u} 1 {l21y1+l22y2>u} ( 3 ∑l 4 1 {l11y1>u} 1 {l21y1+l22y2>u} ( 3 ∑l= l 43 l 33 ∫ R 2 1 {l11y1>u} 1 3 )

 34333213 i=1 3i y i )y 4 ϕ(y)dy = 0,I 3 = l 43 ∫ R i=1 3i y i )y 3 ϕ(y)dy {l21y1+l22y2>u} ϕ(y 1 )ϕ(y 2 )dy 1 dy 2 = l 43 l 33 Q (0) (

Table 1

 1 values for χ 2 Goodness-of-fit test p-values associated to test statistics F a1 , F a1 , F a2 and F a2 in Figure3. The considered processX(•) is a Gaussian univariate process (d = 1) with E(X(0)) = 0, Var(X(0)) = 1 and covariance function r(t) = e -t 2 .1.2. Furthermore, in Table1we display the goodness-of-fit test p-values for F a1 , F a1 , F a2 and F a2 in Figure3. Obtained p-values allow us to statistically accept the χ 2 m distribution for these test statistics. Moreover, we highlight how p-values increase when one evaluates the EC curve for several levels u i (see F a2 and F a2 cases).

	m	Cramér-von Mises test Anderson-Darling test
	F a1	0.7162	0.7519
	F a1	0.6379	0.7585
	F a2	0.8933	0.9576
	F a2	0.7518	0.7731

Table 3

 3 Goodness-of-fit test p-values associated to test statistics D in Figure

Table 4

 4 Goodness-of-fit test p-value associated to test statistic D in Figure

Table 5 ,

 5 we display the goodness-of-fit test p-values for D. The obtained high p-values allow us to statistically accept the χ 2 m distribution for the test statistic D.

Table 6

 6 Goodness-of-fit test p-value associated to test statistic D in Figure

  The coefficients d uk (t) are given by ′′X (t) -2α 1 (t)β 1 (t) -(α 2 1 + β 2 1 (t))r X (t)).Let us continue our study with the behaviour of J(t, ⟨z, w⟩) for t close to 0. Let us start with the series. It holds

		d uk (t) = =	1 k! 1 k!	∫ ∫ ∞ R G( -du (t) d u (t) √ uσ ε (t) √ uσε (t) ( d u (t) + y)H k (y)ϕ(y)dy √ uσ ε (t) + y)H k (y)ϕ(y)dy,
	with ϕ the standard Gaussian density, and
	Corr(W (t), Z(t)) = (-r ∞ ⟨z, w⟩ uσ 2 ε (t) ∑ k=0 d 2 uk (t)k! = ∫ ∞ -du(t) √ ( d u (t) √ uσ ε (t) uσε (t) ≤ 2(( ∫ ∞ d u (t) √ ) 2 uσ ε (t) -du (t) + y) 2 ϕ(y)dy, ϕ(y)dy + √ uσε (t)	∫ ∞ -du(t) √ uσε (t)
	k=0	d 2 uk (t)k! ≤ 2(2(-	d u (t) √ uσ ε (t)	)ϕ(	d u (t) √ uσ ε (t)	) + (-	√ d u (t) uσ ε (t)	)ϕ(	d u (t) √ uσ ε (t)	)) = O(1).	(

y 2 ϕ(y)dy).

For the last integral, we use the upper bound

∫ ∞ a y 2 ϕ(y)dy ≤ aϕ(a) + 1 a ϕ(a),

that yields, as t → 0, ∞ ∑

  (τ 1 )dσ s-1 (τ 2 ), where σ s-1 stands for the surface Lebesgue measure on S s-1 . A consequence of Funk-Hecke formula says that, for any reasonable function F , ∫S s-1 ×S s-1 F (⟨τ 1 , τ 2 ⟩)dσ s-1 (τ 1 )dσ s-1 (τ 2 ) = σ s-1 (S s-1 ) 2 2W s-2

			⟨z, w⟩) p ⃗ X(0), ⃗ X(t),χ 2 s (0),χ 2 s (t) (z, w, u, u) dz dw
	=	u s-2 (2π) s (1 -r 2 X (t)) s/2	S s-1 ×S s-1 ∫	J(t, u⟨τ 1 , τ 2 ⟩) e	-u(	1-r X (t)⟨τ 1 ,τ 2 ⟩ 1-r 2 X (t)
							∫ π
							F (cos ψ) sin s-2 ψ dψ,
							0
	with W s the Wallis integral, i.e. W s =	∫ π/2 0	sin s (ψ)dψ. It implies
	I(t) =	(2π)	4 u s-2

)

dσ s-1 s (1 -r 2 X (t)) s/2 σ s-1 (S s-1 ) 2 2W s-2 ∫ π 0 J(t, u cos ψ) e -u( 1-r X (t) cos ψ 1-r 2 X (t)

Acknowledgment. The authors sincerely express their gratitude to the referees for their valuable comments that help to motivate the study and make the paper clearer.

* Partially supported by the INRIA International Chairs program

1 {l11y1>u} 1 {l21y1+l22y2>u} (l 31 y 1 y 2 + l 32 [(y 2 2 -1) + 1])ϕ(y 1 )ϕ(y 2 )dy 1 dy 2 ).

In this manner, we get the result.