Rosa M Crujeiras
email: rcrujeir@usc.es.

Rubén Fernández-Casal

Wenceslao González-Manteiga

An L 2 -test for comparing spatial spectral densities

Keywords: hypothesis testing, spectral density, spatial process

Consider Z 1 and Z 2 stationary spatial processes. We study the asymptotic behaviour of a Cramer-von-Mises type test statistic for testing the hypothesis H 0 : f 1 = f 2 , where each f l denotes the spectral density of each observed process, for l = 1, 2. Asymptotic distribution theory under H 0 and under local alternatives is given. A simulation study showing the performance of the test is also provided.

Introduction

An important problem in spatial statistic is modeling the dependence structure of spatial data, both from parametric and nonparametric approaches. This problem can be focused from the spatial domain, where the variogram or the covariogram are the functions that describe the dependence. From this point of view, on goodness-of-fit testing for dependence structures, [START_REF] Diblasi | On the use of the variogram in checking for independence in spatial data[END_REF] propose a test for independence and [START_REF] Maglione | Exploring a valid model for the variogram of an isotropic spatial process[END_REF] extend the former technique for choosing a valid model for a variogram, based on smoothed versions of the observed variables.

Spectral techniques constitute an alternative way for studying dependent data and this methodology has been broadly used in time series analysis (e.g. [START_REF] Priestley | Spectral analysis and time series[END_REF]). Despite the extension of these techniques to multidimensional settings is not straightforward, this approach is gaining acceptance in spatial data analysis. From this point of view, the target function is no longer the variogram or the covariogram, but the spatial spectral density. The classical nonparametric estimator for the spectral density is the periodogram. Its extension to the spatial context has been studied by [START_REF] Fuentes | Spectral methods for nonstationary spatial processes[END_REF].

Under this spectral scheme, [START_REF] Crujeiras | Goodness-of-fit tests for the spatial spectral density[END_REF] provide two testing techniques for goodness-of-fit testing, using distances on the spectral and on the log-spectral domain. These test statistics take advantage of the representation of the spatial periodogram as the response variable in a multiplicative regression model. By a logarithmic transform, the spatial log-periodogram can be written as the exogenous variable in a regression model, where the regression function is the log-spectral density.

In the goodness-of-fit for regression models literature, [START_REF] King | Testing the equality of two regression curves using linear smoothers[END_REF] study the problem of comparing two regression curves, using linear smoothers, under independendent and Gaussian errors. The general case of comparing L ≥ 2 regression curves is studied in [START_REF] Dette | Nonparametric analysis of covariance[END_REF], under heterocesdastic errors. In [START_REF] Vilar-Fernández | Nonparametric comparison of curves with dependent errors[END_REF], the authors provide a procedure for testing the equality of regression curves, under fixed design and dependent errors. Based on the ideas in [START_REF] Vilar-Fernández | Nonparametric comparison of curves with dependent errors[END_REF], the goal of our work is to provide a test for testing the hypothesis that the spectral densities of L spatial random processes are equal, without specifying a parametric model. For that purpose, we consider a Cramer-von-Mises type functional, as in [START_REF] Dette | Nonparametric analysis of covariance[END_REF] and [START_REF] Vilar-Fernández | Nonparametric comparison of curves with dependent errors[END_REF]. This paper is organized as follows. In Section 2 we provide some background on spatial spectral methods and nonparametric regression. In Section 3, we study the asymptotic distribution of the test under the null hypothesis and under local alternatives. A brief simulation study is given in Section 4. Section 5 is devoted to technical details.

2 Some background on spectral techniques and nonparametric regression.

Let Z l be a zero mean second-order stationary spatial process, observed on a regular grid D l , for l = 1, . . . , L. That is, {Z l (s), s ∈ D l = a l + D}, with D = {1, . . . , d 1 } × {1, . . . , d 2 }. The case a 1 = . . . = a L implies that the processes are observed on the same grid of locations. Denote by N d = d 1 d 2 the number of points in any of the grids D l , with l = 1, . . . , L. The covariance function of the processes are defined by:

C l (u) = E(Z l (s), Z l (s + u)), s, u ∈ Z 2 . (1)
Assuming that u |C l (u)| < ∞, by Khinchin's theorem (e.g. [START_REF] Yaglom | Correlation theory of stationary and related random functions[END_REF]), the covariance function of a stationary random process can be written, for l = 1, . . . , L as:

C l (u) = Π 2 e -iu T λ f l (λ)dλ, Π 2 = [-π, π] × [-π, π] (2)
where f l , the spectral density, is bounded and continuous for all l and T denotes the transpose operator.

The classical nonparametric estimator of the spectral density is the periodogram, which is given by:

I l (λ k ) = 1 (2π) 2 N d s∈D l Z l (s)e -is T λ k 2 , (3)
where s T λ k denotes the scalar product in R 2 . The periodogram is usually computed at the set of bidimensional Fourier frequencies, λ T k = (λ k1 , λ k2 ):

λ kj = 2πk j d j , k j = 0, ±1, . . . , ±n j = ⌊ d j -1 2 ⌋, j = 1, 2 (4)
and denote by N = (2n 1 + 1)(2n 2 + 1) the number of Fourier frequencies. The periodogram (3) can be also written in terms of the sample covariances as:

I l (λ k ) = 1 (2π) 2 u∈U Ĉl (u)e -iu T λ k , l = 1, . . . , L (5)
where U = {u = (u 1 , u 2 ); u j = 1 -d j , . . . , d j -1, j = 1, 2} and the sample covariances, for Z l with l = 1, . . . , L, are given by:

Ĉl (v) = 1 N d s∈D l (v) Z l (s)Z l (s + v), D l (v) = {s ∈ D l ; s + v ∈ D l }. (6)
Assume that the spatial processes Z l , l = 1, . . . , L can be represented as:

Z l (s) = ∞ j=-∞ ∞ k=-∞ ψ l jk ε l (s 1 -j, s 2 -k), (7)
where the error variables ε l are independent (and independent among themselves) and identically distributed (i.i.d.) as N (0, σ 2 ε l ), for l = 1, . . . , L. Note that any Gaussian stationary process can be represented as in (7). Then, the corresponding spectral density f l can be written as:

f l (λ) = |A l (λ)| 2 f ε l (λ), λ ∈ Π 2 (8)
where

f ε l (λ) = σ 2 ε l (2π) 2 and A l (λ) = ∞ j=-∞ ∞ l=-∞ ψ l jk e -i(j,k)λ , (j, k)λ = jλ 1 + kλ 2 .
In this case, the periodogram for each process Z l , with l = 1, . . . , L, admits the following representation:

I l (λ k ) = f l (λ k )V l k + R l N (λ k ), (9)
where the variables V l k are i.i.d. standard exponential distributed, and V l k and V l ′ k , with l = l ′ are also independent. The residual term R l N (λ k ) is uniformly bounded (see [START_REF] Crujeiras | Goodness-of-fit tests for the spatial spectral density[END_REF]). Applying logarithms in (9) we have:

Y l k = m l (λ k ) + z l k + r l k , l = 1, . . . , L (10)
where m l = log f l is the log-spectral density, the variables z l k = log V l k are i.i.d. with density function h(x) = e x-e x , and the residual term r l k is given by:

r l k = log 1 + R l N (λ k ) f l (λ k )V l k .
Several nonparametric estimators of the spatial log-spectral density could be obtained considering a smoothed combination of log-periodogram values, that is:

ml (λ k ) = i W l i (λ k )Y l i . (11)
The weights W l i can be defined as Gasser-Muller weights, for instance:

W l i (λ) = (2π) 2 |H| 1/2 B i K(H -1/2 (λ -ν))dν, (12)
where K is a bidimensional kernel function, H is a bidimensional bandwidth matrix and the integration region is given by:

B i = [a i1-1 , a i1 ] × [a i2-1 , a i2 ], λ i ∈ B i , ∪ i B i = B, B i ∩ B j = ⊘, for i = j.
The sets B i in the partition of B must be Jordan measurable and max i µ(B i ) = O(N -1 ) (see [START_REF] Müller | Nonparametric regression analysis of longitudinal data[END_REF]). Other options are Priestley-Chao weights:

W l i (λ) = (2π) 2 N K H (λ -λ i ) = (2π) 2 N |H| 1/2 K(H -1/2 (λ -λ i )), (13)
or Nadaraya-Watson weights. Another alternative consists of considering a local-linear estimator for the spatial log-spectral density. All these weights, under a fixed design setting, are asymptotically equivalent.

Asymptotic Analysis

Consider {Z l (s), s ∈ D l }, with l = 1, . . . , L, L realizations of a spatial stochastic process (for instance, realizations taken on L time moments) or L realizations of different spatial processes. Our main purpose is to test whether the dependence structure of {Z l , l = 1, . . . , L} is the same. In terms of the log-spectral densities m l , the testing problem can be written as H 0 : m 1 = . . . = m L vs. H a : m l = m j , for some l = j. In this context, the comparison can be made by considering nonparametric estimators of the spatial log-spectral densities.

Consider the following test statistic, based on a L 2 -distance:

Q = L l=2   l-1 j=1 Π 2 ( m l (λ) -m j (λ)) 2 ω(λ)dλ   , ( 14
)
where ω is a positive, bounded weight function with support Π 2 . This weight function is usually chosen to avoid edge-effects. In the spectral context, this function can be chosen in order to filter frequencies where the periodogram presents higher variability, as the origin or those frequencies with π-valued components.

A1 The spatial processes Z l can be represented as in ( 7) and

i j |i| 1/2 |j||ψ l ij | ≤ ∞ and i j |i||j| 1/2 |ψ l ij | ≤ ∞.
A2 The spectral densities are non-vanishing:

inf λ∈Π 2 f l (λ) > 0, for l = 1, . . . , L.
A3 We consider Gasser-Muller type weights, given by ( 12), or Priestley-Chao weights, given by ( 13). Besides,

W 1 i = . . . = W L i . A4 The bidimensional kernel function K is continuously differentiable, with compact sup- port and K 2 (u)du < ∞.
A5 The bidimensional bandwidth matrix,

H satisfies N |H| 1/2 → ∞, as N → ∞, with n 1 , n 2 → ∞ and n 1 /n 2 → c, for some constant c.
Consider the hypothesis testing problem:

H 0 : m 1 = m 2 , H a : m 1 = m 2
and assume that both Z 1 and Z 2 have been observed on grids with the same design. This implies that the corresponding Fourier frequencies are the same in both cases. By Riemann's approximation, Q can be approximated by:

Q = (2π) 2 N k ( m 1 (λ k ) -m 2 (λ k )) 2 ω(λ k ). ( 15
)
Theorem 1. Assume conditions (A1)-(A5) hold. Then, under the null hypothesis that H 0 : m 1 = m 2 , we have that, as N → ∞:

N 2 |H| 1/2 Q - (2π) 4 12N |H| 1/2 C K I ω → N (0, σ 2 Q ), ( 16
)
in distribution, with

C K = K 2 (u)du, I ω = ω(v)
dv and the asymptotic variance is

σ 2 Q = (2π) 8 72 (K * K) 2 (u)du ω 2 (v)dv,
where * denotes the convolution operator.

Also in this context of comparing two dependence structures, consider that the null hypothesis is false and assume that:

m 1 (λ) -m 2 (λ) = C N p(λ), (17)
where p(λ) is a non-zero function. We will see that the test statistic Q allows for detecting local alternatives at a distance of order N -1/2 |H| -1/8 .

Theorem 2. Assume conditions (A1)-(A5) hold. Then, if (17) holds and C 2 N = (N 2 |H| 1/2 ) -1/2 , we have that, as N → ∞:

N 2 |H| 1/2 Q - (2π) 4 12N |H| 1/2 C K I ω → N p 2 (v)ω(v)dv, σ 2 Q ,
in distribution, with C K , I ω and σ 2 Q as in Theorem 1.

Theorems 1 and 2 can be generalized for stationary random fields on R d , under a similar asymptotic framework, assuming that the sampling grid increases at the same rate in all directions. A d-variate kernel function K satisfying condition A4 and a d-dimensional bandwidth matrix H, satisfying condition A5 must be considered. The corresponding asymptotic mean and variance in ( 16) are given by:

4 3 π 2d 1 N |H| 1/2 C K I ω , and σ 2 Q,d = 2 2d+1 9 π 4+2d (K * K) 2 (u)du ω 2 (v)dv,
where the weight function ω is now defined on Π d = [-π, π] d . Thus, in the particular case of d = 1, we provide a testing technique for comparing spectral densities in time series context. In this case, we have a scalar bandwidth parameter h, which plays the role of H 1/2 in the general dimension setting.

As a particular case, this technique can be adapted to detect changes on the dependence structure of a process observed at different time moments. If this is the case, denote Z l (s) = Z(s, t l ) the observation of the process Z at location s ∈ D (taking a l = 0, for l = 1, . . . , L) and time moment t l . This capacity makes the technique relevant when studying spatial-temporal processes. Invariance of the spatial dependence along time makes feasible the use of stationary spatial-temporal dependence models (e.g. [START_REF] Fernández-Casal | Flexible spatio-temporal stationary variogram models[END_REF]).

If the spatial process Z l are observed on regular grids with different sizes, then the corresponding frequency spectrum is not the same. The asymptotic behaviour of Q could be determined following similar arguments to those in (Vilar-Fernández and González-Manteiga (2004), Theorem 3), under some conditions on the asymptotic rates of the samples.

In the regression context, [START_REF] Zhang | A power comparison between nonparametric regression tests[END_REF] consider three major types of nonparametric regression test, among them, the Cramer-von-Mises test we consider in this work. The L 2 -test result the most powerful, from a local asymptotic point of view.

We provide in this section some simulation results. One thousand simulations of Z 1 and Z 2 have been drawn in a 51 × 51 and a 101 × 101 regular grids. These processes have been considered as bidimensional autoregressive, introduced by [START_REF] Martin | A subclass of lattice processes applied to a problem in planar sampling[END_REF], with autoregression parameters (β 1 , β 2 ) (from now on BAR(β 1 , β 2 )). For the simulation study, we consider that both processes present the same model BAR(0.5, 0.5), under the null hypothesis; then, we set the following alternatives: H a : Z 1 = BAR(0.0, 0.0) and Z 2 = BAR(0.05, 0.05) and also H a : Z 1 = BAR(0.0, 0.0) and Z 2 = BAR(0.1, 0.1).

If we denote by (i, j) the points on an observation grid, a BAR(β 1 , β 2 ) satisfies the following expression:

Z(i, j) = β 1 Z(i -1, j) + β 2 Z(i, j -1) -β 1 β 2 Z(i -1, j -1) + ε(i, j),
where ε(i, j) is a white noise process with variance σ 2 (which has been set to one in the simulations). The corresponding spectral density is given by: In order to reproduce as exactly as possible the behaviour of the test statistic, we have considered Priestley-Chao weights (see Appendix), although a local linear estimator could be used in practice, in order to avoid edge-effects. Weight function ω in ( 14) has been chosen such that frequency (0, 0) and frequencies with components next to ±π have been filtered (two of them in the smaller grid and four of them for the largest one). We have considered diagonal bandwidth matrix, with elements proportional to the spacing between frequencies:

f (λ) = σ 2 (2π) 2 (1 + β 1 -2β 1 cos(λ 1 )) -1 (1 + β 2 -2β 2 cos(λ 2 )) -1 .
H = rh • diag 2π d 1 , 2π d 2 ,
and a range of rh has been explored. In Tables 1 and2, we show the percentage of rejections in the simulation schemes described above. In general, the performance of the test agrees with the theoretical results, although with a bit conservative behaviour, showing the influence of the smoothness parameter, denoted by rH. In addition, a bias correction on the test statistic would reduce the influence of this parameter. As it is expected, the results are better as the size of the grid increases, both in terms of size and power.

Appendix

Let's introduce the following notation. Consider the following regression model:

Y l * k = m l (λ k ) + z l * k , l = 1, 2. ( 18
)
where

Y l * k = Y l k -C 0 -r l k , and z l * k = z l k -C 0 , where C 0 = E(z l k ) is the Euler constant. Denote by m l
* the nonparametric estimator of m l as in ( 11) for ( 18) equation and denote by

B l k = j W j (λ k )r l j , l = 1, 2.
Lemma 1. The test statistic Q can be decomposed in three addends:

Q = Q1 + Q2 + Q3 , where Q1 = (2π) 2 N k ( m 1 * (λ k ) -m 2 * (λ k )) 2 ω(λ k ), Q2 = (2π) 2 N k (B 1 k -B 2 k ) 2 ω(λ k ), Q3 = 2 (2π) 2 N k ( m 1 * (λ k ) -m 2 * (λ k ))(B 1 k -B 2 k )ω(λ k ).
Proof. It is straightforward from the definitions of the non parametric estimator in regression model ( 18) and B l k , for l = 1, 2.

Lemma 2. Under conditions (A1)-(A5) and under H 0 , we have that:

Q2 = O P log 2 N N 2 |H| 1/2 and Q3 = O P log 2 N N .
Proof. Q2 can be decomposed as:

Q2 = (2π) 2 N k j W 2 j (λ k )(r 1 j -r 2 j ) 2 ω(λ k ) + (2π) 2 N k i j =i W j (λ k )W i (λ k )(r 1 j -r 2 j )(r 1 i -r 2 i )ω(λ k ) = Q2,1 + Q2,2 .
By a Taylor expansion on the residual part r l j around 0, for l = 1, 2:

r l j = - R l N (λ j ) f l (λ j )V j - 1 2(1 + x j ) 2 R l N (λ j ) f l (λ j )V j 2 , where x j ∈ 0, R l N (λ j ) f l (λ j )V j . Since, max j |R l N (λ j )| = O P (N -1/2 log N ),
for l = 1, 2, just following [START_REF] Kooperberg | Rate of convergence for logspline spectral density estimation[END_REF], the Lagrange remainder in the Taylor expansion, denoted by LR l j can be uniformly bounded by O P log 2 N N

. Then, under the null hypothesis, Q2,1 can be decomposed as Q2,

1 = Q1 2,1 + Q2 2,1 + Q3 2,1 , where Q1 2,1 = (2π) 2 N k j W 2 j (λ k ) R 2 N (λ j ) f (λ j )V - R 1 N (λ j ) f (λ j )V 2 ω(λ k ), (19)
Q2 2,1 = (2π) 2 N k j W 2 j (λ k )(LR 2 j -LR 1 j ) 2 ω(λ k ), ( 20
) and Q3 2,1 = 2 (2π) 2 N k j W 2 j (λ k ) R 2 N (λ j ) f (λ j )V - R 1 N (λ j ) f (λ j )V (LR 2 j -LR 1 j )ω(λ k ). (21)
Let's find a bound for Q2 2,1 in (20), the addend involving the Lagrange remainders.

Q2 2,1 = (2π) 2 N k j W 2 j (λ k )(LR 1 j -LR 2 j ) 2 ω(λ k ) = O P log 4 N N 2 (2π) 2 N k j 1 N 2 |H| K 2 (H -1/2 (λ k -λ j ))ω(λ k ) = O P log 4 N N 2 (2π) 2 N   k j =k (2π) 4 N 2 |H| K 2 (H -1/2 (λ k -λ j ))ω(λ k ) + k 1 N |H| K 2 (H -1/2 0)ω(λ k )   = O P log 4 N N 4 |H| 1/2 + O P log 4 N N 2 (2π) 6 N 3 |H| k j =k K 2 (H -1/2 (λ k -λ j ))ω(λ k ) ≤ O P log 4 N N 4 |H| 1/2 + O P log 4 N N 2 (2π) 6 N 2 |H| max k j =k K 2 (H -1/2 (λ k -λ j )) = O P log 4 N N 3 |H| 1/2 ,
where the last inequality follows from max k ω(λ k ) ≤ c, for some constant c. Following similar arguments as above, we will found bounds for Q2,1 and Q2,3 . For Q1 2,1 in (19):

Q1 2,1 = (2π) 2 N k j W 2 j (λ k ) R 2 N (λ j ) f (λ j )V - R 1 N (λ j ) f (λ j )V 2 ω(λ k ) = Q1,1 2,1 + Q1,2 2,1 + Q1,3 2,1
where

Q1,1 2,1 = (2π) 2 N k j W 2 j (λ k ) R 1 N (λ j ) f (λ j )V 2 ω(λ k ), Q1,3 2,1 = 2 (2π) 2 N k j W 2 j (λ k ) R 1 N (λ j ) f (λ j )V R 2 N (λ j ) f (λ j )V ω(λ k )
and Q1,2 2,1 is the similar to Q1,1 2,1 , but replacing each R 1 N (λ j ) for R 2 N (λ j ). We will find a bound for Q1,1 2,1 .

Q1,1 2,1 ≤ (2π) 2 N N max j R 1 N (λ j ) f (λ j )V 2 j k W 2 j (λ k )ω(λ k ) = O P log 2 N N j k 1 N 2 |H| K 2 (H -1/2 (λ k -λ j ))ω(λ k ) = O P log 2 N N 2 |H| 1/2 .
Similar computations lead to the same bound for the other addends. Let's find a bound for the third addend, Q3

2,1 in (21):

Q3 2,1 = 2 (2π) 2 N k j W 2 j (λ k ) R 1 N (λ j ) f (λ j )V - R 2 N (λ j ) f (λ j )V (LR 2 j -LR 1 j )ω(λ k ) = O P log N N 1/2 O P log 2 N N (2π) 6 N j k 1 N 2 |H| K 2 (H -1/2 (λ k -λ j ))ω(λ k ) = O P log 3 N N 5/2 |H| 1/2 . Therefore, Q2 ≤ O P log 2 N N 2 |H| 1/2 .
Finally, Q3 can be written as:

Q3 = 2 (2π) 2 N k ( m 1 (λ k ) -m 2 (λ k )) (B 1 k -B 2 k )ω(λ k ) ≤ 2 max k |B 1 k -B 2 k | (2π) 2 N k j W j (λ k )(Y 1 j -Y 2 j )ω(λ k ) ≤ O P log N N 1/2 max j |Y 1 j -Y 2 j | (2π) 4 N k j 1 N |H| 1/2 K(H -1/2 (λ k -λ j ))ω(λ k ) = O P log 2 N N since H 0 : m 1 = m 2 implies that Y 1 j -Y 2 j = r 1 j -r 2 j
, for every Fourier frequency. Lemma 3. Under conditions (A1)-(A5) and under H 0 , we have that

N 2 |H| 1/2 Q1 - (2π) 4 12N |H| 1/2 C K I ω → N (0, σ 2 Q ),
in distribution, where C K , I ω and σ 2 Q as in Theorem 1.

Proof. Define the following random variables Λ

k = z 1 * k -z 2 * k , with E(Λ k ) = 0, E(Λ 2 k ) = π 2 3
and Cov(Λ k , Λ j ) = 0 for j = k. The statistic Q1 can be decomposed in two addends as Q1 = Q1,1 + Q1,2 , where :

Q1,1 = (2π) 2 N k j W 2 j (λ k )Λ 2 j ω(λ k ), Q1,2 = (2π) 2 N k j i =j W j (λ k )W i (λ k )Λ j Λ i ω(λ k ). Define b i,j = (2π) 2 N k W i (λ k )W j (λ k )ω(λ k ). Then: Q1,1 = j b j,j Λ 2 j , and Q1,2 = i =j b i,j Λ i Λ j .
First, we will study the behaviour of Q1,1 . For simplicity, consider Priestley-Chao weights. Taking expectations and using Riemann's approximation, it is easy to see that:

E( Q1,1 ) ≈ (2π) 4 12 1 N |H| 1/2 K 2 (u)du ω(v)dv. (22)
Let's check the order of the variance of Q1,1 . Denote by c 2 = V ar(Λ 2 j ) This variance can be computed taking into account that:

V ar( Q1,1 ) = (2π) 12 N 6 |H| 2 c 2 j α 2 j , α j = k K 2 (H -1/2 (λ k -λ j ))ω(λ k ).
Then,

V ar( Q1,1 ) = c 2 (2π) 12 N 6 |H| 2 k ω(λ k ) k ′ ω(λ k ′ ) j K 2 (H -1/2 (λ k -λ j ))K 2 (H -1/2 (λ k ′ -λ j ))
square matrix with size N . The new coefficients, with one dimensional indexes, are given by the following matrix: A = (a ij ) , A ∈ M N ×N , and each entry of this matrix is defined by a ij = π √ 3 b ij and a ii = 0, where the bidimensional indexes are given by: i = (i 1 , i 2 ) = (k, k 0 ), if (k -1)κ 2 ≤ i ≤ kκ 2 and i = (k -1)κ 2 + k 0 , (24) j = (j 1 , j 2 ) = (l, l 0 ), if (l -1)κ 2 ≤ j ≤ lκ 2 and j = (l -1)κ 2 + l 0 .

(25)

Now, define the variables X i = √ 3 π Λ i , where i and i satisfy: i = (i 1 -1)κ 2 + i 2 . Q1,2 can be written as a quadratic form with one-dimensional indexes:

Q1,2 = i,j a i,j X i X j , a ij = π √ 3 (2π) 2 N k W i (λ k )W j (λ k )ω(λ k ),
where i and j are determined by ( 24) and ( 25), respectively. Asymptotic normality is proved if the following conditions are satisfied:

(i) There exists a sequence of real numbers k(N ) such that k(N ) 4 σ 2 Q max i j a 2 ij → 0.

(ii) The random variables X i satisfy max i E X 2 i 1 {|Xi|>k(N )} → 0, N → ∞.

(iii) The eigenvalues of the matrix A = (a ij ) are negligible:

σ 2 Q max i µ 2 i → 0, N → ∞.
In order to check (i) -(iii) note that:

a 2 ij = π 2 (2π) 4 3N 2 k W i (λ k )W j (λ k )ω(λ k ) 2 = π 2 (2π) 4 3N 2 k W 2 i (λ k )W 2 j (λ k )ω 2 (λ k ) + π 2 (2π) 4 3N 2 k W i (λ k )W j (λ k )ω(λ k ) k =k ′ W i (λ k ′ )W j (λ k ′ )ω(λ k ′ ) = a 2A ij + a 2B ij .
Besides, max i j

a 2A ij ≈ C K π 2 (2π) 8 3N 4 |H| max i (ω 2 * K 2 )(λ i ) = O N -4 |H| . ( 26
)
Following similar arguments, the same rate is obtained for max i j a 1B ij . Therefore,

k 4 (N )σ 2 Q max i j a ij = O k 4 (N ) N 6 |H| , (27)
which tends to zero if the sequence k(N ) → ∞ satisfies that k 4 (N ) N 6 |H| → 0.

Condition (ii) holds since the variables X i are i.i.d. with second order moment E(X 2 i ) = 1. It remains to show that condition (iii) also holds. Since max i j |a ij | = O(N -1 ) and

  51 × 51 H0 : (0.5, 0.5) Ha : (0.0, 0.0), (0.05, 0.05) Ha : (0.0, 0.0), (0.1, 0.1) rh α = 0.01 α = 0.05 α = 0.10 α = 0.01 α = 0.05 α = 0.10 α = 0.01 α = 0.05 α = 0.10 5 0

  01 α = 0.05 α = 0.10 α = 0.01 α = 0.05 α = 0.10 α = 0.01 α = 0.05 α = 0.10

Table 1 :

 1 Percentage of rejections for 51 × 51 regular grid. α: significance level.

		.042	0.094	0.139	0.177	0.265	0.340	0.648	0.775	0.837
	10	0.039	0.088	0.123	0.229	0.323	0.397	0.792	0.880	0.907
	15	0.045	0.079	0.105	0.280	0.388	0.456	0.703	0.820	0.857
	20	0.051	0.078	0.099	0.149	0.215	0.249	0.484	0.606	0.682
	25	0.053	0.079	0.101	0.112	0.169	0.197	0.333	0.419	0.471
	30	0.051	0.082	0.102	0.089	0.136	0.175	0.215	0.296	0.345
	35	0.043	0.070	0.098	0.066	0.106	0.136	0.154	0.154	0.251
	40	0.035	0.059	0.087	0.047	0.080	0.108	0.090	0.150	0.190
	45	0.029	0.046	0.066	0.025	0.061	0.086	0.057	0.091	0.132
	50	0.019	0.039	0.053	0.019	0.032	0.063	0.020	0.064	0.082

Table 2 :

 2 Percentage of rejections for 101 × 101 regular grid. α: significance level.

	5.00	0.031	0.074	0.109	0.550	0.711	0.785	0.998	1.000	1.000
	12.78	0.045	0.089	0.127	0.978	0.992	0.996	1.000	1.000	1.000
	20.55	0.052	0.089	0.114	0.983	0.990	0.994	1.000	1.000	1.000
	28.33	0.057	0.096	0.118	0.976	0.983	0.988	1.000	1.000	1.000
	36.10	0.055	0.093	0.111	0.958	0.976	0.983	1.000	1.000	1.000
	43.88	0.055	0.092	0.115	0.940	0.961	0.971	1.000	1.000	1.000
	51.65	0.054	0.088	0.118	0.890	0.941	0.959	1.000	1.000	1.000
	59.43	0.055	0.082	0.116	0.787	0.884	0.916	1.000	1.000	1.000
	67.20	0.050	0.076	0.108	0.647	0.755	0.817	1.000	1.000	1.000
	74.98	0.041	0.073	0.094	0.525	0.618	0.685	0.979	0.994	0.999

Acknowledgements

The authors acknowledge the financial support of the Spanish Ministry of Education and Science and FEDER, project MTM2005-00820. Research of Rosa Crujeiras was partially supported by grant BES2003-0581 and Xunta de Galicia project PGIDIT06PXIB207009PR. In addition, research of Rubén Fernández-Casal was partially supported by Xunta de Galicia project PGIDIT05TIC00701CT.

A C C E P T E D M A N U S C R I P T ACCEPTED MANUSCRIPT

which can be approximated, using a changes of variable and Riemann's sums, by:

Therefore, applying Markov's inequality, it follows that:

Besides, E( Q1,2 ) = 0 since Λ i and Λ j are uncorrelated, for i = j and its variance is given by:

V ar( Q1,2 ) =

The same applies for u = j. For E(Λ i Λ j Λ u Λ v ) to be different from zero, one of these two conditions must hold: i = u and j = v or i = v and j = u. Then:

Consider the following approximation for the product of two b ij coefficients:

Then,

Therefore, we obtain the following expression for σ 2 Q :

In order to prove the asymptotic normal distribution of Q1,2 , we will apply Theorem 5.2 in de Jong (1987). For that purpose, we must write Q1,2 as a quadratic form, namely Q1,2 = i,j a i,j X i X j , where i and j are one-dimensional indexes and X i are i.i.d. random variables with zero mean and unit variance.

First, define a new subindex for the Fourier frequencies λ k , with k = (k 1 , k 2 ) and k l = 0, ±1, . . . , ±n l , for l = 1, 2. Consider

ACCEPTED MANUSCRIPT

taking into account that the spectral ratio of a matrix (max i |ν i |) is bounded by its supremun norm,

Then, the asymptotic convergence to a normal distribution is proved.

Proof of Theorem 1. The theorem is proved combining the results from Lemmas 1-3.

Proof of Theorem 2. Consider the decomposition of the test statistic given in Lemma 1: Q = Q1a + Q2a + Q3a . The sketch of the proof is as follows: find bounds for Q2a and Q3a and Q1a is decomposed in three addends Q1a1 , Q1a2 and Q1 . The asymptotic normality of Q1 is proved in Theorem 1. Besides, Q1a2 can be also bounded and

where the approximation holds fof C 2 N = N -1 |H| -1/4 .


Uwe K Üchler

Stefan Tappe

ON THE SHAPES OF BILATERAL GAMMA DENSITIES

Keywords: bilateral Gamma distributions, selfdecomposability, unimodality, asymptotic behaviour, density shapes

We investigate the four parameter family of bilateral Gamma distributions. The goal of this paper is to provide a thorough treatment of the shapes of their densities, which is of importance for assessing their fitting properties to sets of real data. This includes appropriate representations of the densities, analyzing their smoothness, unimodality and asymptotic behaviour.

Introduction

In many fields of applications it is important to find appropriate classes of distributions for fitting observed data. For this issue, normal distributions often provide only a poor fit. Specific examples are given by the logarithmic returns of stock prices, because their empirical densities typically possess heavier tails and much higher located modes than normal distributions.

Thus, several authors have looked for other appropriate classes of distributions. We mention the generalized hyperbolic distributions [START_REF] Barndorff-Nielsen | Exponentially decreasing distributions for the logarithm of particle size[END_REF] and their subclasses, which have been applied to finance in [START_REF] Eberlein | Hyperbolic distributions in finance[END_REF], the Variance Gamma distributions [START_REF] Madan | Purely discontinuous asset pricing processes[END_REF] and CGMYdistributions [START_REF] Carr | The fine structure of asset returns: an empirical investigation[END_REF].

Recently, another family of distributions was proposed in [START_REF] Küchler | Bilateral Gamma distributions and processes in financial mathematics[END_REF]: Bilateral Gamma distributions. In the mentioned article, bilateral Gamma distributions are fitted to observed stock prices and compared to other classes of distributions considered in the literature.

In order to provide a general overview about their fitting properties -also in view of other applications than finance -we present a thorough treatment of the shapes of their densities. After recalling the basic properties of bilateral Gamma distributions in Section 2, we provide suitable representations of the densities in Section 3, which we can use in order to obtain density plots with a computer program. Afterwards, the investigation of the shapes of bilateral Gamma distributions starts: Section 4 concerns the smoothness of the densities, Section 5 the unimodality and Section 6 is devoted to the asymptotic behaviour of the densities near zero and for x → ±∞. In Section 7 we characterize typical shapes of the densities and draw implications concerning the fitting properties of bilateral Gamma distributions.

Date: February 25, 2008.

Bilateral Gamma distributions

In this section, we define bilateral Gamma distributions and review some of their properties. For details and more informations, we refer to [START_REF] Küchler | Bilateral Gamma distributions and processes in financial mathematics[END_REF].

A bilateral Gamma distribution with parameters α + , λ + , α -, λ -> 0 is defined as the distribution of X -Y , where X and Y are independent, X ∼ Γ(α + , λ + ) and Y ∼ Γ(α -, λ -).

The characteristic function of a bilateral Gamma distribution is

ϕ(z) = λ + λ + -iz α + λ - λ -+ iz α - , z ∈ R (2.1)
where the powers stem from the main branch of the complex logarithm.

If X is bilateral Gamma distributed with parameters (α + , λ + ; α -, λ -), then for any c > 0 the random variable cX has, by (2.1), again a bilateral Gamma distribution with parameters (α + , λ + c ; α -, λ - c ). Note that, also by (2.1), the sum of two independent bilateral Gamma random variables with parameters (α + 1 , λ + ; α - 1 , λ -) and (α + 2 , λ + ; α - 2 , λ -) has again a bilateral Gamma distribution with parameters (α

+ 1 + α + 2 , λ + ; α - 1 + α - 2 , λ -).
In particular, bilateral Gamma distributions are stable under convolution, and they are infinitely divisible. It follows from [START_REF] Sato | Lévy processes and infinitely divisible distributions[END_REF]Ex. 8.10] that both, the drift and the Gaussian part in the Lévy-Khintchine formula (with truncation function h = 0), are equal to zero, and that the Lévy measure is given by

F (dx) = α + x e -λ + x 1 (0,∞) (x) + α - |x| e -λ -|x| 1 (-∞,0) (x) dx.
(2.2) Thus, we can also express the characteristic function ϕ as

ϕ(z) = exp R e izx -1 k(x) x dx , z ∈ R (2.3)
where k : R → R is the function

k(x) = α + e -λ + x 1 (0,∞) (x) -α -e -λ -|x| 1 (-∞,0) (x), x ∈ R (2.4)
which is decreasing on each of (-∞, 0) and (0, ∞). It is an immediate consequence of [START_REF] Sato | Lévy processes and infinitely divisible distributions[END_REF]Cor. 15.11] that bilateral Gamma distributions are selfdecomposable, and hence of class L in the sense of [START_REF] Sato | On distribution functions of class L. Zeit[END_REF] and [START_REF] Sato | On higher derivatives of distribution functions of class L[END_REF]. This is a key property for analyzing their densities, which is exploited in Sections 4, 5 and 6.

Using the characteristic function (2.1), we can specify the following quantities.

Mean:

α + λ + -α - λ -, Variance: α + (λ + ) 2 + α - (λ -) 2 , Skewness: 2 α + (λ + ) 3 -α - (λ -) 3 α + (λ + ) 2 + α - (λ -) 2 3/2 , Kurtosis: 3 + 6 α + (λ + ) 4 + α - (λ -) 4 α + (λ + ) 2 + α - (λ -) 2 2 .

Representations of the densities

Bilateral Gamma distributions are absolutely continuous with respect to the Lebesgue measure, because they are the convolution of two Gamma distributions. Since the densities satisfy the symmetry relation

f (x; α + , λ + , α -, λ -) = f (-x; α -, λ -, α + , λ + ), x ∈ R \ {0} (3.1)
it is sufficient to analyze the density functions on the positive real line. As the convolution of two Gamma densities, they are for x ∈ (0, ∞) given by

f (x) = (λ + ) α + (λ -) α - (λ + + λ -) α -Γ(α + )Γ(α -) e -λ + x ∞ 0 v α --1 x + v λ + + λ - α + -1 e -v dv. (3.2)
We can express the density f by means of the Whittaker function W λ,µ (z) [5, p. 1014]. According to [5, p. 1015], the Whittaker function has the representation

W λ,µ (z) = z λ e -z 2 Γ(µ -λ + 1 2 ) ∞ 0 t µ-λ-1 2 e -t 1 + t z µ+λ-1 2 dt for µ -λ > - 1 2 . (3.3) From (3.2) and (3.3) we obtain for x > 0 f (x) = (λ + ) α + (λ -) α - (λ + + λ -) 1 2 (α + +α -) Γ(α + ) x 1 2 (α + +α -)-1 e -x 2 (λ + -λ -) (3.4)
× W 1 2 (α + -α -), 1 2 (α + +α --1) (x(λ + + λ -))
. By [5, p. 1014], we can express the Whittaker function W λ,µ (z) by the Whittaker functions M λ,µ (z), namely it holds

W λ,µ (z) = Γ(-2µ) Γ( 1 2 -µ -λ) M λ,µ (z) + Γ(2µ) Γ( 1 2 + µ -λ) M λ,-µ (z).
For the Whittaker function M λ,µ (z) the identity [5, p. 1014]

M λ,µ (z) = z µ+ 1 2 e -z 2 Φ(µ -λ + 1 2 , 2µ + 1; z) is valid, with Φ(α, γ; z) denoting the confluent hypergeometric function [5, p. 1013] Φ(α, γ; z) = 1 + α γ z 1! + α(α + 1) γ(γ + 1)
z 2 2! + α(α + 1)(α + 2) γ(γ + 1)(γ + 2)
z 3 3! + . . . (3.5)
Because of the series representation (3.5) of Φ(α, γ; z), we can use (3.4) in order to obtain density plots with a computer program, which is done in Section 7.

If one of α + , α -is an integer, the representation becomes more convenient at one half of the real axis.

3.1. Proposition. Assume α + ∈ N = {1, 2, . . .}. Then it holds for each x ∈ (0, ∞)

f (x) = (λ + ) α + (λ -) α - (λ + + λ -) α -(α + -1)! α + -1 k=0 a k x k e -λ + x ,
where the coefficients (a k ) k=0,...,α + -1 are given by

a k = α + -1 k 1 (λ + + λ -) α + -1-k α + -2-k l=0 (α -+ l).
Proof. Since α + is an integer, we can compute the integral appearing in (3.2) by using the binomial expansion formula. The calculations are obvious.

The symmetry relation (3.1) and the identity [5, p. 1017]

W 0,µ (z) = z π K µ z 2 ,
where K µ (z) denotes the Bessel function of the third kind, imply that in the case α + = α -=: α the density (3.4) is of the form

f (x) = 1 Γ(α) λ + λ - λ + + λ - α |x| α-1 e -x 2 (λ + +λ -) |x|(λ + + λ -) π K α-1 2 |x| 2 (λ + + λ -) (3.6) for x ∈ R \ {0}. The density of a Variance Gamma distribution V G(µ, σ 2 , ν) is, according to [7, Sec. 6.1.5], given by h(x) = 2 exp( µx σ 2 ) ν 1/ν √ 2πσΓ( 1 ν ) x 2 2σ 2 ν + µ 2 1 2ν -1 4 K 1 ν -1 2 1 σ 2 x 2 2σ 2 ν + µ 2 . (3.7) Inserting the parametrization (µ, σ 2 , ν) := α λ + - α λ -, 2α λ + λ -, 1 α (3.8)
into (3.7), we obtain the density (3.6), showing that bilateral Gamma distributions with α + = α -=: α are Variance Gamma with parameters given by (3.8).

Conversely, for a bilateral Gamma distribution which is Variance Gamma it necessarily holds α + = α -, see [START_REF] Küchler | Bilateral Gamma distributions and processes in financial mathematics[END_REF]Thm. 3.3].

Smoothness

As we have pointed out in Section 2, bilateral Gamma distributions are selfdecomposable. Therefore, we may use the results of [START_REF] Sato | On distribution functions of class L. Zeit[END_REF] and [START_REF] Sato | On higher derivatives of distribution functions of class L[END_REF] in the sequel.

The smoothness of the density depends on the parameters α + and α -. Let N ∈ N 0 = {0, 1, 2, . . .} be the unique nonnegative integer satisfying

N < α + + α -≤ N + 1. 4.1. Theorem. It holds f ∈ C N (R \ {0}) and f ∈ C N -1 (R) \ C N (R).
Proof. This is an immediate consequence of [9, Thm. 1.2].

Thus, the N -th order derivative of the density f is not continuous. The only point of discontinuity is zero. In Section 6, we will explore the behaviour of f (N ) near zero.

The densities of bilateral Gamma distributions satisfy the following integro-differential equation.

Proposition. f satisfies for

x ∈ R \ {0} xf (x) = (α + + α --1)f (x) -α + λ + ∞ 0 f (x -u)e -λ + u du -α -λ - ∞ 0 f (x + u)e -λ -u du.
Proof. The assertion follows from Cor. 2.1 of [START_REF] Sato | On distribution functions of class L. Zeit[END_REF].

Unimodality

Bilateral Gamma distributions are strictly unimodal, which is the content of the next result.

5.1.

Theorem. There exists a point x 0 ∈ R such that f is strictly increasing on (-∞, x 0 ) and strictly decreasing on (x 0 , ∞).

Proof. The existence of the mode x 0 is a direct consequence of [9, Thm. 1.4], because neither the distribution function of a bilateral Gamma distribution nor its reflection is of type I 4 in the sense of [START_REF] Sato | On distribution functions of class L. Zeit[END_REF].

We emphasize that the mode x 0 from Theorem 5.1 can, in general, not be determined explicitly. However, we get the following result, which narrows the location of the mode. 5.2. Proposition. If α + , α -≤ 1, then x 0 = 0. Presumed α + > 1 and α -≤ 1, it holds x 0 ∈ (0, α + -1 λ + ). In the case α + , α -> 1 we have x 0 ∈ (-α --1 λ -, α + -1 λ + ), and it holds

x 0 = 0 if and only if λ -α + -λ + α -= λ --λ + , x 0 > 0 if and only if λ -α + -λ + α -> λ --λ + , x 0 < 0 if and only if λ -α + -λ + α -< λ --λ + .
Proof. The first statement is a consequence of parts (viii) and (ix) of [START_REF] Sato | On distribution functions of class L. Zeit[END_REF]Thm. 1.3].

Since the mode of a Γ(α, λ)-distribution with α > 1 is given by α-1 λ , parts (ii) and (iii) of [START_REF] Sato | On distribution functions of class L. Zeit[END_REF]Thm. 4.1] yield the second assertion.

In the case α + , α -> 1, part (iv) of [START_REF] Sato | On distribution functions of class L. Zeit[END_REF]Thm. 4.1] shows that x 0 ∈ (-α --1 λ -, α + -1 λ + ). According to Theorem 4.1, the density f is continuously differentiable. Using the representation (3.2) and Lebesgue's dominated convergence theorem, we obtain the first derivative for x ∈ (0, ∞)

f (x) = (λ + ) α + (λ -) α - (λ + + λ -) α -Γ(α + )Γ(α -) -λ + e -λ + x ∞ 0 v α --1 x + v λ + + λ - α + -1 e -v dv + (α + -1)e -λ + x ∞ 0 v α --1 x + v λ + + λ - α + -2 e -v dv .
Applying Lebesgue's dominated convergence theorem again, by the continuity of f and the fact Γ(x + 1) = xΓ(x), x > 0 we get

f (0) = (λ + ) α + (λ -) α - (λ + + λ -) α + +α --2 Γ(α + + α --2) Γ(α + -1)Γ(α -) 1 - λ + λ + + λ -• α + + α --2 α + -1 ,
which yields the remaining statement of the proposition.

A particular consequence of Proposition 5.2 is that for λ + α + and λ - α -the mode x 0 is necessarily close to zero.

Asymptotic behaviour

We have seen in Section 4 that for N := α + + α --1 the N -th order derivative of the density f is not continuous. The only point of discontinuity is zero. We will now explore the behaviour of f (N ) near zero. For the proof of the upcoming result, Theorem 6.1, we need the following properties of the Exponential Integral [1, Chap. 5]

E 1 (x) := ∞ 1 e -xt t dt, x > 0.
The Exponential Integral has the series expansion

E 1 (x) = -γ -ln x - ∞ n=1 (-1) n n • n! x n , (6.1)
where γ denotes Euler's constant

γ = lim n→∞ n k=1 1 k -ln(n) .
The derivative of the Exponential Integral is given by

∂ ∂x E 1 (x) = - e -x x . (6.2)
Due to symmetry relation (3.1) it is, concerning the behaviour of f (N ) near zero, sufficient to treat the case x ↓ 0.

Theorem. Let

N := α + + α --1. (1) lim x↓0 f (N ) (x) is finite if and only if α + ∈ N = {1, 2, . . .}. (2) If α + / ∈ N and α + + α -/ ∈ N, then f (N ) (x) ∼ C 1 x α as x ↓ 0 for constants C 1 = 0, α ∈ (0, 1). (3) Let α + / ∈ N be such that α + + α -∈ N. Then f (N ) (x) ∼ M (x) as x → 0, where M is a slowly varying function as x → 0 satisfying lim x→0 M (x) = ∞. Moreover, it holds lim x↓0 (f (N ) (x) -f (N ) (-x)) = C 2 ∈ R.
The constants in Theorem 6.1 are given by

α = N + 1 -α + -α -, C 1 = (λ + ) α + (λ -) α -sin(α + π) Γ(α + + α --N ) sin((α + + α -)π) , C 2 = (λ + ) α + (λ -) α - 2 (-1) N +1 cos(α + π) + cos(α -π) .
Proof. For α + ∈ N we conclude the finiteness of the limit in the first statement from [10, Thm. 3], since for each β ∈ (0, 1) (recall that the function k was defined in (2.4))

lim u↓0 u β-1 (α + -k(u)) = lim u↓0 u β-1 α + (1 -e -λu ) = 0.
In order to prove the rest of the theorem, we evaluate expressions (1.8)-(1.10) in [START_REF] Sato | On distribution functions of class L. Zeit[END_REF], and then we apply [START_REF] Sato | On distribution functions of class L. Zeit[END_REF]Thm. 1.7]. The constant c in [9, eqn. (1.8)] is in the present situation

c = exp (α + + α -) 1 0 e -u -1 u du + (α + + α -) ∞ 1 e -u u du - ∞ 1 α + e -λ + u + α -e -λ -u u du . (6.3)
The first integral appearing in (6.3) is by (6.2) and the series expansion (6.1)

1 0 e -u -1 u du = lim x↓0 -E 1 (u) -ln u u=1 u=x = -E 1 (1) -γ.
Using (6.2) and the fact lim

x→∞ E 1 (x) = 0, see [1, Chap. 5], for each constant λ > 0 the identity ∞ 1 e -λu u du = lim x→∞ -E 1 (λu) u=x u=1 = E 1 (λ)
is valid. Thus, we obtain c = e -(α + +α -)γ-α + E 1 (λ + )-α -E 1 (λ -) . (6.4)

The function K(x) in [9, eqn. (1.9)] is in the present situation

K(x) = exp 1 |x| α + + α --α + e -λ + u -α -e -λ -u u du .
Since by (6.2)

1 |x| e -λu u du = E 1 (λ|x|) -E 1 (λ) for λ > 0,
we obtain

K(x) = e α + E 1 (λ + )+α -E 1 (λ -) |x| -(α + +α -) e -α + E 1 (λ + |x|)-α -E 1 (λ -|x|) . (6.5)
Using the series expansion (6.1), we get (6.6) showing that for the slowly varying function

lim x→0 K(x) = (λ + ) α + (λ -) α -e (α + +α -)γ+α + E 1 (λ + )+α -E 1 (λ -) ,
L(x) =
W λ,µ (z) ∼ e -z 2 z λ H(z)
with H denoting the function

H(z) = 1 + ∞ k=1 µ 2 -(λ -1 2 ) 2 µ 2 -(λ -3 2 ) 2 • • • µ 2 -(λ -k + 1 2 ) 2 k!z k .
Obviously, it holds H(z) ∼ 1 for z → ∞. Taking into account (3.1) and (3.4), for x → ±∞ the density has the asymptotic behaviour

f (x) ∼ C 3 x α + -1 e -λ + x as x → ∞, f (x) ∼ C 4 |x| α --1 e -λ -|x| as x → -∞,
where the constants C 3 , C 4 > 0 are given by

C 3 = (λ + ) α + (λ -) α - (λ + + λ -) α -Γ(α + ) and C 4 = (λ + ) α + (λ -) α - (λ + + λ -) α + Γ(α -) .
As a consequence, we obtain for the logarithmic density function ln f

lim x→∞ ln f (x) x = -λ + and lim x→-∞ ln f (x) x = λ -.
In particular, the density of a bilateral Gamma distribution is semiheavy tailed.

Shapes of the densities

The shapes of bilateral Gamma distributions can have remarkable differences. Using the results of the previous sections, we characterize typical shapes of their densities. We infer that the density has a pole at the mode x 0 = 0. Notice that in the special case α + + α -= 1 the difference f (x) -f (-x) tends to a finite value as x ↓ 0 by the third statement of Theorem 6.1. We observe that densities with (7.1) are appropriate for fitting data sets with many observations accumulating closely around zero.

If 1 < α + + α -≤ 2, then, by Theorem 4.1, f is continuous on R, but its derivative is not continuous at zero. Let us have a closer look at the behaviour of f near zero in this case.

• If α + , α -∈ (0, 1) and α + + α -∈ (1, 2), it holds, according to Theorem 6.1,

lim x↑0 f (x) = ∞ and lim x↓0 f (x) = -∞.
Hence, we have a steep mode of the density at zero with exploding slope from the left and from the right. This shape is also suitable for sets of data with many observations being close to zero. • If α -< 1 < α + , applying Theorem 6.1 yields

lim x↑0 f (x) = ∞ and lim x↓0 f (x) = ∞.
Hence, the mode x 0 is located at the positive half axis and f has infinite slope at zero. We remark that in the special case α + + α -= 2 the difference f (x) -f (-x) tends to a finite value as x ↓ 0 by the third statement of Theorem 6.1.

• If α + , α -= 1, we have a two-sided exponential distribution, which is in particular Variance Gamma, as we have shown at the end of Section 3. We obtain lim with finite constants C -, C + ∈ (0, ∞). Consequently, we have a peak mode of the density at zero with finite slope from both sides. If α + + α -> 2, then the density is smooth, that is at least of type C 1 (R) by Theorem 4.1. Choosing λ + α + and λ - α -, the mode x 0 is necessarily close to zero by Proposition 5.2. Such shapes are in particular applicable for observations of financial data. We refer to [START_REF] Küchler | Bilateral Gamma distributions and processes in financial mathematics[END_REF]Sec. 9], where for a specific data set of stock returns the maximum likelihood estimation α + = 1.55, λ + = 133.96, α -= 0.94, λ -= 88.92 provided a good fit.

Summarizing the preceding results, Figure 1 provides an overview about typical shapes of bilateral Gamma densities.

Figure 1 .

 1 Figure 1. The shapes of f for λ + = λ -. Different choices of λ + and λ -may shift the mode and change the skewness.

  If α + + α -≤ 1, then f is not continuous at zero by Theorem 4.1. According to Theorem 6.1, it holds lim x↑0 f (x) = ∞ and lim x↓0 f (x) = ∞. (7.1)

x↑0f

  (x) = C -and lim x↓0 f (x) = -C +


Wojciech Sarnowski
email: wojciech.sarnowski@pwr.wroc.pl

Krzysztof Szajowski
email: krzysztof.szajowski@pwr.wroc.pl

On-line detection of a part of a sequence with unspecified distribution

Keywords: Disorder problem, sequential detection, optimal stopping, Markov process, change point. 1991 MSC: Primar 60G40, 60K99; Secondary 90D60

We register a Markov process. At two random moments θ 1 , θ 2 , where θ 1 < θ 2 , the distribution of observed sequence changes. It is known before θ 1 and after θ 2 . Between these instants is unknown, chosen randomly from a set of distributions. The optimal stopping rule which stops observation of the sequence between disorders θ 1 and θ 2 is identified.

Introduction

In this paper we present a generalization of the double disorder problem considered by Yoshida [START_REF] Yoshida | Probability maximizing approach for a quickest detection problem with complicated Markov chain[END_REF] and Szajowski [START_REF] Szajowski | Optimal on-line detection of outside observation[END_REF]. The considerations are inspired by the problem regarding how can we protect ourselves against a second fault in a technological system after the occurrence of an initial fault. The proposed procedure assumes that the distribution of observations between disorders is unknown in advance and it is chosen randomly by "nature" (see Bojdecki & Hosza [START_REF] Bojdecki | On a generalized disorder problem[END_REF]). In order to find the optimal stopping rule we use the approach of maximizing the probability of stopping between these two faults. The formulation of the problem can be found in section 2. The main result is presented in section 3.

Formulation of the problem

Let (X n , n ∈ N) be an observable sequence of random variables defined on the space (Ω, F, P ) with values in (E, B), where E is a subset of R. On (E, B) there are σ-additive measures {µ x } x∈E . Space (Ω, F, P ) supports variables θ 1 , θ 2 . They are F-measurable variables with values in N. We assume the following distributions:

P (θ 1 = j) = p j-1 1 q 1 , P (θ 2 = k | θ 1 = j) = p k-j-1 2 q 2 ; k > j, j = 1, 2, ... (1)
Moreover we introduce ε, which is F-measurable, independent of θ 1 , θ 2 . ε has distribution:

P (ε = i) = e i ; i = 1, 2, ..., d; d < ∞; d i=1 e i = 1. Additionally we consider d + 2 Markov processes on (Ω, F, P ): (X 1 n , G 1 n , µ 1 x ), (X 2,i n , G 2,i n , µ 2,i x ), i = 1, .., d, (X 3 n , G 3 n , µ 3
x ) defined in the following way:

X n = X 1 n • I {θ 1 >n} + X 2,i n • I {θ 1 ≤n<θ 2 , ε=i} + X 3 n • I {θ 2 ≤n} . (2)
σ-fields G 1 n , G 2,i n , G 2
n are the smallest σ-fields for which (X1 ), (X2,i ), (X 3 ) are adapted (respectively).

We make inference only based on observable sequence (X n , n ∈ N) 1 . Thus for further consideration we define filtration {F n } n∈N where F n = σ(X 0 , X 1 , ..., X n ). Variables θ 1 , θ 2 and ε are not stopping time with respect to F n and σ-fields G • n . Moreover we assume that θ 1 , θ 2 are independent of (X 1 n , n ∈ N), (X 2,i n , n ∈ N), i = 1, ..., d and (X 3 n , n ∈ N). Measures µ satisfy the relations:

µ 1 x (dy) = 1 • µ x (dy), µ 2,i
x (dy) = f i x (y)µ x (dy), i = 1, ..., d, µ 3 x (dy) = g x (y)µ x (dy). Where f i x (.) = 1, g x (.) = 1. Furthermore the functions f 1 x (.), ..., f d x (.),g x (.) are different and g x (y)/f i x (y) < ∞ for all i = 1, ..., d and all x, y ∈ E. We assume that the measures µ 1

x , µ 2,1 x , ..., µ 2,d x , µ 3 x , x ∈ E are known in advance and we have that

P (X 1 1 ∈ A | X 0 = x) = µ 1 x (A) for every A ∈ B.
The model presented has the following heuristic justification: two disorders take place in the observed sequence (X n ). They affect distributions by changing their parameters. Disorders occur at two random moments of time θ 1 and θ 2 . They split the sequence of observations into three segments. The first segment is described by (X 1 n ), the second one -for θ 1 ≤ n < θ 2 -by (X 2 n ). The third is given by (X 3 n ) and is observed when n ≥ θ 2 . The second epoch of observation is the important one. We do not know the distribution of observations in advance. After the first disorder takes the place there is a "switch" from the initial distribution to one of d possible distributions (given by densities f i x ). The choice of this distribution is random. We model this using an unobservable variable ε that indicates which of them is drawn (if ε = j, then the µ j

x distribution is chosen with density f j x with respect to measure µ x ). Next, at the random time θ 2 the distribution of observations becomes µ 3

x . We assume that the variables ε, θ 1 , θ 2 are unobservable.

Our aim is to stop the observed sequence between the two disorders.This can be interpreted as a strategy for protecting against a second failure when the first has already happened. In the mathematical approach -putting T as the set of all stopping times with respect to the filtration {F n } n∈N -we want to control the probability V (τ ) = P (τ < ∞, θ 1 ≤ τ < θ 2 ) by choosing the stopping time τ * ∈ T for which

V (τ * ) = sup τ ∈T V (τ ) (3)
3 Solution of the problem For X 0 = x let us define:

Z n = P x (θ 1 ≤ n < θ 2 | F n ) for n = 0, 1, 2, . . ., Y n = esssup {τ ∈T , τ ≥n} P x (θ 1 ≤ τ < θ 2 | F n ) for n = 0, 1, 2, . . . and τ 0 = inf{n : Z n = Y n } (4) Notice that, if Z ∞ = 0, then Z τ = P x (θ 1 ≤ τ < θ 2 | F τ ) for τ ∈ T . Since F n ⊆ F τ (when n ≤ τ ) we have Y n = ess sup τ ≥n P x (θ 1 ≤ τ < θ 2 | F n ) = ess sup τ ≥n E x (I {θ 1 ≤τ <θ 2 } | F n ) = ess sup τ ≥n E x (Z τ | F n ) Lemma 1
The stopping time τ 0 defined by formula ( 4) is the solution of problem (3).

Proof. From the theorems presented in [START_REF] Bojdecki | Probability maximizing approach to optimal stopping and its application to a disorder problem[END_REF] it is enough to show that lim n→∞ Z n = 0. For all natural numbers n, k, where n ≥ k for each x ∈ E we have:

Z n = E x (I {θ 1 ≤n<θ 2 } | F n ) ≤ E x (sup j≥k I {θ 1 ≤j<θ 2 } | F n ) From Levy's theorem lim sup n→∞ Z n ≤ E x (sup j≥k I {θ 1 ≤j<θ 2 } | F ∞ ) where F ∞ = σ ( ∞ n=1 F n ).
I {θ 1 ≤j<θ 2 } | F ∞ ) = 0 a.s.
what ends the proof of the Lemma.

For further considerations let us introduce the following notation: let u , v = d i=1 u i v i . Thus . , . stands for inner product. Let us now define some necessary processes and functions: Π

1 n = P x (θ 1 ≤ n | F n ), Π 2 n = P x (θ 2 ≤ n | F n ), for n = 0, 1, 2, . . ., H(y, z, α, β) = (1 -α)p 1 + [(1 -α)q 1 + (α -β)p 2 ] e , f y (z) + [(α -β)q 2 + β] g y (z) (5) Π 1 (y, z, α, β) = 1 -(1 -α)p 1 H -1 (y, z, α, β) (6) Π 2 (y, z, α, β) = [(α -β)q 2 + β]g y (z)H -1 (y, z, α, β) (7)
Using these definitions we formulate and prove a technical Lemma which will be useful when proving the main theorem.

Lemma 2 For each x ∈ E and each Borel function u : R -→ R the equations given below are true

Π 1 n+1 = Π 1 (X n , X n+1 , Π 1 n , Π 2 n ) (8) Π 2 n+1 = Π 2 (X n , X n+1 , Π 1 n , Π 2 n ), (9)
with boundary condition:

Π 1 0 = Π 2 0 = 0 E x u(X n+1 )(1 -Π 1 n+1 ) | F n = (1 -Π 1 n )p 1 E u(y)µ Xn (dy), (10)
E x u(X n+1 )(Π 1 n+1 -Π 2 n+1 ) | F n = (1 -Π 1 n )q 1 + (Π 1 n -Π 2 n )p 2 E u(y) e , f Xn (y) µ Xn (dy), (11)
E x u(X n+1 )Π 2 n+1 ) | F n = q 2 (Π 1 n -Π 2 n ) + Π 2 n E u(y)g Xn (y)µ Xn (dy) (12)
Proof. It is obvious that Π 1 0 = Π 2 0 = 0. Notice that since θ 1 < θ 2 we obtain the inclusion {ω : θ 1 ≤ n} ⊇ {ω : θ 2 ≤ n} and:

Z n = P x (θ 1 ≤ n < θ 2 | F n ) = Π 1 n -Π 2 n ( 13
)
Next suppose that A i ∈ F i , i ≤ n + 1 and let us consider the probability

P x (θ 1 > n + 1 | X i ∈ A i , i ≤ n + 1) = P x (θ 1 > n + 1, X n+1 ∈ A n+1 | X i ∈ A i , i ≤ n) P x (X n+1 ∈ A n+1 | X i ∈ A i , i ≤ n) E x (u(X n+1 )Π 2 n+1 | F n ) = E x (E x (u(X n+1 )1 {θ 2 ≤n+1} | F n ) | F n ) = E x E u(y)P x (dy | F n , θ 2 ≤ n + 1)P x (θ 2 ≤ n + 1 | F n ) | F n = E u(y)g Xn (y)µ Xn (dy)P x (θ 2 ≤ n + 1 | F n ) = E u(y)g Xn (y)µ Xn (dy)(P x (θ 2 = n + 1, θ 1 ≤ n < θ 2 | F n ) + P x (θ 2 ≤ n | F n )) = (Π 1 n -Π 2 n )q 2 + Π 2 n E u(y)g Xn (y)µ Xn (dy)
We used the properties of conditional expectation here. Similar transformations give us equations ( 10) and ( 11). This proves Lemma 2.

Notice that Lemma 2 enables us to calculate E x (u(X n+1 ) | F n ) in easy way. Combining ( 10)-( 12) we have

E x (u(X n+1 ) | F n ) = E u(y)H(X n , y, Π 1 n , Π 2 n )µ Xn (dy) (14) Lemma 3 System X x = {X x n }, where X x n = (X n-1 , X n , Π 1 n , Π 2
n ) forms a family of random Markov functions.

Proof. a function:

ϕ(x 1 , x 2 , α 1 , α 2 ; z) = (x 2 , z, Π 1 (x 2 , z, α 1 , α 2 ), Π 1 (x 2 , z, α 1 , α 2 )) (15)
Observe that

X x n = ϕ(X n-2 , X n-1 , Π 1 n-1 , Π 2 n-1 ; X n ) = ϕ(X x n-1 ; X n )
Hence X x n can be interpreted as function of previous state X x n-1 and random variable X n . Moreover, applying (14), we get that conditional distribution of X n given σ-field F n-1 depends only on X x n-1 . According to [START_REF] Shiryaev | Optimal Stopping Rules[END_REF] (pp. 102-103) system X x is a family of random Markov functions. This fact implies that we can reduce initial problem (3) to the problem of optimal stopping fourdimensional process (

X n-1 , X n , Π 1 n , Π 2 n ) with reward h(x 1 , x 2 , α 1 , α 2 ) = α 1 -α 2 (16)
The reward function results from equation ( 13). Thanks to Lemma 3 we construct the solution using standard tools of optimal stopping theory (c.f [START_REF] Shiryaev | Optimal Stopping Rules[END_REF] ).

For any Borel function v : E 2 × [0, 1] d+1 -→ [0, 1] let us define two operators:

T x v(y, z, α, β) = E x (v(X n , X n+1 , Π 1 n+1 , Π 2 n+1 ) | X n-1 = y, X n = z, Π 1 n = α, Π 2 n = β) Q x v(y, z, α, β) = max{v(y, z, α, β), T x v(y, z, α, β)}
From well known theorems of optimal stopping theory ( [START_REF] Shiryaev | Optimal Stopping Rules[END_REF]), we infer that the solution of the problem posed is the Markov time τ 0 :

τ 0 = inf{h(X n , X n+1 , Π 1 n+1 , Π 2 n+1 ) ≥ h * (X n , X n+1 , Π 1 n+1 Π 2 n+1 )}
where:

h * (y, z, α, β) = lim k→∞ Q k x h(y, z, α, β) Of course Q k x v(y, z, α, β) = max{Q k-1 x v, T x Q k-1 x v} = max{v, T x Q k-1 x v}
To obtain a clearer formula for 0 , we formulate the basic result of this paper:

Theorem 1 (a) The solution of problem ( 3) is given by:

τ * = inf{n : (X n , X n+1 , Π 1 n+1 , Π 2 n+1 ) ∈ B * } (17)
Set B * is of the form:

B * = {(y, z, α, β) : (α -β) ≥ (1 -α) × p 1 E R * (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β))µ y (du) + q 1 E S * (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β)) e , f y (u) µ y (du)
+ (α -β)p 2 E S * (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β)) e , f y (u) µ y (du)
Where:

R * (y, z, α, β) = lim k→∞ R k (y, z, α, β) , S * (y, z, α, β) = lim k→∞ S k (y, z, α, β)
Functions R k and S k are defined recursively:

R 1 (y, z, α, β) = 0 , S 1 (y, z, α, β) = 1 R k+1 (y, z, α, β) = (1 -I R k (y, z, α, β)) (18)
× p 1 E R k (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β))µ y (du) +q 1 E S k (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β)) e , f y (u) µ y (du) , S k+1 (y, z, α, β) = I R k (y, z, α, β) + (1 -I R k (y, z, α, β)) (19) ×p 2 E S k (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β)) e , f y (u) µ y (du)
Where the set R k is:

A C C E P T E D M A N U S C R I P T ACCEPTED MANUSCRIPT R k = (y, z, α, β) : h(y, z, α, β) ≥ T x Q k-1 x h(y, z, α, β) = {(y, z, α, β) : (α -β) ≥ (1 -α) × p 1 E R k (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β))µ y (du) + q 1 E S k (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β)) e f y (u) µ y (du) + (α -β)p 2 E S k (y, u, Π 1 (y, u, α, β), Π 2 (y, u, α, β)) e , f y (u) µ y (du) (b)
The value problem. The optimal value for (3) is given by the formula

V (τ * ) = p 1 E R * (x, u, ϕ x (u), 0)µ x (du) +q 1 E S * (x, u, ϕ x (u), 0) e , f x (u) µ x (du)
where: ϕ x (u) = 1 -

p 1 p 1 +q 1 e ,f x (u)
Proof. Part (a). Part (a) results from Lemma 2 -the problem reduces to the problem of optimal stopping of the Markov process (X n-1 , X n , Π 1 n , Π 2 n ) with payoff function h(y, z, α, β) = α -β. Given (11) with the function u equal to unity we get:

T x h(y, z, α, β) = E Π 1 n+1 -Π 2 n+1 | F n | X n-1 =y,Xn=z,Π 1 n =α,Π 2 n =β = (Π 1 n -Π 2 n )p 2 E e , f Xn (u) µ Xn (du) + (1 -Π 1 n )q 1 | X n-1 =y,Xn=z,Π 1 n =α,Π 2 n =β = (1 -α)q 1 + (α -β)p 2
From the definition of R 1 and S 1 it is clear that h(y, z, α, β) = α -β = (1 -α)R 1 (y, z, α, β) + (α -β)S 1 (y, z, α, β)

Also R 1 = {(y, z, α, β) : h(y, z, α, β) ≥ T x h(y, z, α, β)}. From the definition of Q x and the facts above we obtain:

Q x h(y, z, α, β) = (1 -α)R 2 (y, z, α, β) + (α -β)S 2 (y, z, α, β)
where: R 2 (y, z, α, β) = q 1 (1 -I R 1 (y, z, α, β)) and S 2 (y, z, α, β) = p 2 + (1p 2 )I R 1 (y, z, α, β)). Suppose the following induction hypothesis holds:

Q k-1 x h(y, z, α, β) = (1 -α)R k (y, z, α, β) + (α -β)S k (y, z, α, β)
where R k and S k are given by equations ( 18), (19) respectively. We will show that:

Q k x h(y, z, α, β) = (1 -α)R k+1 (y, z, α, β) + (α -β)S k+1 (y, z, α, β)
From the induction assumption and equations ( 10), (11) we obtain:

A C C E P T E D M A N U S C R I P T ACCEPTED MANUSCRIPT T x Q k-1 x h(y, z, α, β) = T x (1 -α)R k (y, z, α, β) (20) +T x (α -β)S k (y, z, α, β) = (1 -α)p 1 E R k (y, u, α, β)µ y (du) + [(1 -α)q 1 + (α -β)p 2 ] E S k (y, u, α, β) e , f y (u) µ y (du) = (1 -α) p 1 E R k (y, α, β)µ y (du) + q 1 E S k (y, u, α, β) × e , f y (u) µ y (du) + (α -β)p 2 E S k (y, u, α, β) e , f y (u) µ y (du) Notice that (1 -α)R k+1 (y, z, α, β) + (α -β)S k+1 (y, z, α, β) is equal α -β = h(y, z, α, β) = Q k x h(y, z, α, β) for (y, z, α, β) ∈ R k and, taking into account (20), it is equal T x Q k-1 x h(y, z, α, β) = Q k x h(y, z, α, β) for (y, z, α, β) / ∈ R k ,
where R k is given by (20). Finally we get:

Q k x h(y, z, α, β) = (1 -α)R k+1 (y, z, α, β) + (α -β)S k+1 (y, z, α, β)
This proves ( 18) and ( 19). Using the monotone convergence theorem and theorems of optimal stopping theory ( [START_REF] Shiryaev | Optimal Stopping Rules[END_REF]) we conclude that the optimal stopping time τ * is given by (17).

Proof. Part (b). First, notice that Π 1 1 = ϕ x (u) = 1 -p 1 p 1 +q 1 e ,f x (u) . Under the assumption τ * < ∞ a.s. we get:

P x (τ * < ∞, θ 1 ≤ τ * < θ 2 ) = sup τ EZ τ = E max{h(x, X 1 , Π 1 1 , Π 2 1 ), T x h * (x, X 1 , Π 1 1 , Π 2 1 )} = E lim k→∞ Q k x h(x, X 1 , Π 1 1 , Π 2 1 ) = E (1 -Π 1 1 )R * (x, X 1 , Π 1 1 , Π 2 1 ) + (Π 1 1 -Π 2 1 )S * (x, X 1 , Π 1 1 , Π 2 1 )
= p 1 E R * (x, u, ϕ x (u), 0)µ x (du) + q 1 E S * (x, u, ϕ x (u), 0) e , f x (u) µ x (du)

We used Lemma 2 here and simple calculations for Π 1 1 , Π 2 1 . This ends the proof.

Remarks

It is notable that the solution of formulated problem depends only on twodimensional vector of posterior processes. None additional process containing knowledge about ε is necessary. This is a consequence of considered payoff function [START_REF] Brodsky | Nonparametric Methods in Change-Point Problems[END_REF]. The presence of variable ε and random choice of the distribution on the second epoch are reflected in weighted average e , f y (u) . What implies that the problem can be reduced to the case when we just fix the distribution on the middle segment of observations.

  It is true that: lim sup j≥k, k→∞ I {θ 1 ≤j<θ 2 } = 0 a.s. and by the dominated convergence theorem we get lim

k→∞ E x (sup j≥k

It should be emphasized that the sequence (X n , n ∈ N) is not markovian one under admitted assumption as it has been written in[START_REF] Szajowski | Optimal on-line detection of outside observation[END_REF],[START_REF] Yakir | Optimal detection of a change in distribution when the observations form a Markov chain with a finite state space[END_REF] and[START_REF] Dube | A framework for quickest detection of traffic anomalies in networks[END_REF]. However the sequence has Markov property given θ 1 and θ

(see[START_REF] Szajowski | A two-disorder detection problem[END_REF] and[START_REF] Moustakides | Quickest detection of abrupt changes for a class of random processes[END_REF]).

The formulas obtained are very general and for this reason -quite complicated. We simplify the model by assuming that P (θ 1 > 0) = 1 and P (θ 2 > θ 1 ) = 1. However, it seems that some further simplifications can be made in special cases. Further research should be carried out in this direction. From a practical point of view, computer are necessary to construct B * -the set in which we stop our observable sequence. In the paper an extension of the sequential approach proposed by [6] is applied. There variety of models related to disorder in random sequences are presented in the monography by [3].

ACCEPTED MANUSCRIPT

This follows from Bayes' formula. Let us transform the probability appearing in the numerator:

Now we split the probability in the denominator into three parts

The first part is:

The second one:

And the third:

Thus, taking into account ( 5) we obtain for Π 6), it can be seen that ( 8) is satisfied. Applying similar reasoning and transformations to the process Π 2 n+1 we get:

which leads to formula (9). In the case of the conditional expectations we show that equation (12) holds. To do this we need to first define σ-field F n = σ(ε, θ 1 , θ 2 , X 0 , ..., X n ). Notice that F n ⊂ F n . We have: