Valdério Anselmo Reisen
email: valderioanselmoreisen@gmail.com

Céline Lévy-Leduc

Higor Cotta

Pascal Bondon

Márton Ispány

Higor Henrique Aranda Cotta

Marton Ispany

Paulo Roberto Prezotti Filho

An Overview of Robust Spectral Estimators

Keywords: Time series, M -estimation, Q N -estimation, long-memory, periodic processes, robustness

come   L'archive ouverte pluridisciplinaire

Introduction

It is well known that outlying observations may completely destroy most of the standard estimators and several authors developed robust approaches in order to mitigate the impact of additive outliers, specially in time series models which is the process considered in this paper. However, most of the work is devoted to the robust estimation of the location, scale and other statistical tools. In this direction, the classical periodogram is the natural estimator of the spectral density of a time series and recent studies indicate that the periodogram is highly sensitive to the presence of outliers, and, thus, it becomes useless in any sub-sequential analysis. As a viable approach to attenuate this issue, the M -regression method applied to build alternative spectral estimators given in [START_REF] Fajardo | M-periodogram for the analysis of long-range-dependent time series[END_REF] and [START_REF] Reisen | An M-estimator for the long-memory parameter[END_REF]) and the Q N -periodogram introduced in [START_REF] Fajardo | Robust estimation in long-memory processes under additive outliers[END_REF] are some methodologies proposed recently in the literature of time series to handle additive outliers.

The M -periodogram is discussed in [START_REF] Fajardo | M-periodogram for the analysis of long-range-dependent time series[END_REF] and [START_REF] Reisen | An M-estimator for the long-memory parameter[END_REF] for the long-memory time series. The short-range process was still an open problem and is one main contribution of this paper. The asymptotic property of the M -periodogram is derived for the process which is identified to have short-memory property such as an ARMA model (Theorem 1). As a second contribution of this paper, the recent results given [START_REF] Fajardo | M-periodogram for the analysis of long-range-dependent time series[END_REF] and [START_REF] Reisen | An M-estimator for the long-memory parameter[END_REF], for long-memory model, are summarized and these methods are compared empirically with Q N -periodogram and the classical periodogram which is widely used in modelling time series data. Here, these methods are empirically studied and compared in time series with and without additive outliers with the aim to verify their finite sample size robustness properties, that is, to verify their capacity to accommodate the additive outlier's effect.

The use of M -and Q N -periodograms in periodic ARMA (PARMA) models is also discussed here in the context of handling atypical or aberrant observations (additive outliers). This becomes the third contribution of this paper.

This paper is organized as follows: Section 2 discusses robust periodograms based on M -regression method and Q N function for short and long-memory time series. Section 3 presents some simulation results for the methods discussed in Section 2. Section 4 gives some applications of the alternative periodograms in short and long-memory and periodic processes.

Robust periodograms

Let {Y t } t∈Z be a second order stationary process. Since this paper deals with short and longmemory processes, additional assumptions on the process {Y t } t∈Z will be given in the sequel of the paper. For a sample {Y 1 , Y 2 , ..., Y N }, the classical periodogram function, at the Fourier frequency λ j = 2πj/N, j = 1, . . . , [N/2], is defined as

I N (λ j ) = 1 2πN N k=1 Y k exp(ikλ j ) 2 .
(1)

Next subsections deal with alternative periodogram functions which present similar performance (from theoretical and empirical meaning) to I N (λ), λ ∈ (-π, π), but with robustness property against additive outliers and asymmetric and heavy-tail distributions.

M -periodogram

One alternative way to derive the periodogram function I N (λ j ) is based on the Least Square (LS) estimates of a two-dimensional vector β = (β (1) , β (2) ) in the linear regression model

Y i = c N i β + ε i = β (1) cos(iλ j ) + β (2) sin(iλ j ) + ε i , 1 ≤ i ≤ N, β ∈ R 2 , (2)
where ε i denotes the deviation of Y i from c N i β and E(ε i ) = 0 and E(ε 2 i ) < ∞. In the sequel (ε i ) is assumed to be a function of a stationary Gaussian process, see (10) for a precise definition. Then,

βLS N (λ j ) = Arg min β∈R 2 N i=1 (Y i -c N i (λ j )β) 2 , (3)
where

c N i (λ j ) = (cos(iλ j ) sin(iλ j )) . (4)
The solution of (3) is βLS

N (λ j ) = (C C) -1 C Y , (5)
where Y = (Y 1 , . . . , Y N ) , C and C C are defined by

C =      cos(λ j ) sin(λ j ) cos(2λ j ) sin(2λ j ) . . . . . . cos(N λ j ) sin(N λ j )      (6)
and

C C = N k=1 cos(kλ j ) 2 N k=1 cos(kλ j ) sin(kλ j ) N k=1 cos(kλ j ) sin(kλ j ) N k=1 sin(kλ j ) 2 = N 2 Id 2 (7)
where Id 2 is the identity matrix 2 by 2. Hence,

βLS N (λ j ) = 2 N C Y = 2 N N k=1 Y k cos(kλ j ) N k=1 Y k sin(kλ j ) = ( βLS,(1) N (λ j ), βLS,(2) N (λ j )) . (8)
In view of (1),

I N (λ j ) = N 8π βLS N (λ j ) 2 = N 8π ( βLS,(1) N (λ j )) 2 + ( βLS,(2) N (λ j )) 2 =: I LS N (λ j ) , (9)
where • denotes the classical Euclidean norm and βLS N (λ j ) = ( βLS,( 1)

N (λ j ), βLS, (2)
N (λ j )) is the least square estimates of β = (β (1) , β (2) ) see, for example, [START_REF] Fajardo | M-periodogram for the analysis of long-range-dependent time series[END_REF] and [START_REF] Reisen | An M-estimator for the long-memory parameter[END_REF] and references therein. Note that I N (λ j ) (9) can be derived for different choices of ε i , i = 1, . . . , N .

It is supposed here that

ε i = G(η i ) . (10)
In ( 10), G is a non null real-valued and skew symmetric measurable function (i.e. G(-x) = -G(x), for all x) and (η i ) i≥1 is a stationary Gaussian process with zero mean and unit variance. Additional assumptions of (η i ) i≥1 will be given in the sequel of the paper. Let ψ(.) be a function satisfying the following assumptions.

(A1) 0 < E[ψ 2 (ε 1 )] < ∞ .
(A2) The function ψ is absolutely continuous with its almost everywhere derivative ψ satisfying

E[|ψ (ε 1 )|] < ∞ and such that the function z → E[|ψ (ε 1 -z) -ψ (ε 1 )|] is continuous at zero. (A3) ψ is nondecreasing, E[ψ (ε 1 )] > 0 and E[ψ (ε 1 ) 2 ] < ∞.
(A4) ψ is skew symmetric, i.e. ψ(-x) = -ψ(x), for all x.

It is now introduced the M -periodogram based on the M -estimator βM N of the parameter β defined in Equation ( 2). The M -estimator βM

N = ( β(1) N , β(2) N ) is defined as the solution (t 1 , t 2 ) of N i=1 cos(iλ j ) ψ(Y i -cos(iλ j )t 1 ) = 0 and N i=1 sin(iλ j ) ψ(Y i -sin(iλ j )t 2 ) = 0. (11)
β(1) N and β(2) N can be also seen as the minimizers with respect to t 1 and t 2 , respectively, of

N i=1 cos(iλ j ) ψ(Y i -cos(iλ j )t 1 ) and N i=1 sin(iλ j ) ψ(Y i -sin(iλ j )t 2 ) , (12)
where ψ satisfies the same assumptions as in [START_REF] Koul | Second order behavior of M-estimators in linear regression with long-memory errors[END_REF]. By analogy to (9), the robust periodogram I M N (λ j ) at λ j = 2πj/N, j = 1, . . . , [N/2], is defined by

I M N (λ j ) = N 8π βM N (λ j ) 2 = N 8π ( β(1) N (λ j )) 2 + ( β(2) N (λ j )) 2 . ( 13
)
2.1.1 M -Periodogram in short-memory processes In this subsection the asymptotic properties of βM N are established in the short-range dependence framework. For this, the following assumptions are introduced. This result helps to establish the theoretical properties of the robust periodogram I M N given in Corollary 1. (A5) Let η t , t ∈ Z, be i.i.d. standard Gaussian random variables and let a j be real numbers such that j≥0 |a j |< ∞ and a 0 = 1. Then,

ε i = j≥0 a j η i-j .
(A6) ψ is the Huber function that is ψ(x) = max[min(x, c), -c], for all x in R, where c is a positive constant.

Theorem 1. Assume that (A5) and (A6) hold and that β = 0 in (2) so that Y i = ε i . Then, for any fixed j, βM N defined by (12) satisfies

N 2 (F (c) -F (-c)) βM N (λ j ) d -→ N 0, ∆ (j) , N → ∞ ,
where F is the c.d.f. of ε 1 and

∆ (j) = k∈Z E{ψ(ε 0 )ψ(ε k )} cos(kλ j ) sin(kλ j )
-sin(kλ j ) cos(kλ j ) .

Theorem 1 is proved in Section 5.

Corollary 1. Under the assumptions of Theorem 1, I M N (λ j ) defined in (13) satisfies for any fixed j,

I M N (λ j ) d -→ X 2 + Y 2 4π(F (c) -F (-c)) 2 , as N → ∞ , where X ∼ N 0, k∈Z E{ψ(ε 0 )ψ(ε k )} cos(kλ j ) , Y ∼ N 0, k∈Z E{ψ(ε 0 )ψ(ε k )} cos(kλ j ) and Cov(X, Y ) = k∈Z E{ψ(ε 0 )ψ(ε k )} sin(kλ j ).
The proof of Corollary 1 is a straightforward consequence of Theorem 1 and (13).

2.1.2 M -periodogram for long-memory processes Now, consider the following assumption for (η i ) i≥1 in the case of long-memory process. The results in this subsection are well detailed in [START_REF] Fajardo | M-periodogram for the analysis of long-range-dependent time series[END_REF].

(A7) (η i ) i≥1 is a stationary zero-mean Gaussian process with covariances ρ(k) = E(η 1 η k+1 ) satisfying:

ρ(0) = 1 and ρ(k) = k -D L(k), 0 < D < 1 ,
where the function L is slowly varying at infinity and is positive for large k. Recall that a slowly varying function L(x), x > 0 is such that L(xt)/L(x) → 1, as x → ∞ for any t > 0. Constants and logarithms are example of slowly varying functions. Moreover, the spectral density f of (η i ) i≥1 can be expressed as:

f (λ) = |1 -exp(-iλ)| -2d f * (λ) , (14)
where d ∈ (0, 1/2) and f * is an even, positive, continuous function on (-π, π], bounded above and bounded away from zero.

Note that

D = 1 -2d , ( 15
)
where D is defined in Assumption (A7) and d is the standard long-memory parameter notation given in the literature of long-memory models. The fact that (η i ) i≥1 is required to satisfy (A7) essentially means that both L(x), x ≥ 1 and f * (λ), λ in (-π, π] satisfy some smoothness properties.

Theorem 2. Assume that (A7), (A1), (A2), ( A3) and (A4) hold and that β = 0 in (2) so that Y i = ε i . Then, for any fixed j, βM N (λ j ) defined by ( 12) satisfies

N 2 βM N (λ j ) = J 1 E[ψ (ε 1 )] 2 N N i=1 cos(iλ j ) sin(iλ j ) η i + o p (N (1-D)/2 ) , as N → ∞ , (16)
where

J 1 = E[ψ(G(η)
)η] = 0, η being a standard Gaussian random variable and D = 1 -2d. Moreover,

N D/2 βM N (λ j ) d -→ N 0, J 2 1 (E[ψ (ε 1 )]) 2 Γ , N → ∞ , (17)
where

Γ = lim N →∞ 4 N 2-D 1≤k, ≤N c N k (λ j )c T N (λ j )ρ(k -) (18) = 8π × (2πj) -2d f * (0) L 1 0 0 L 2 . ( 19
)
In Relation (18), the vector c N k (λ j ) is defined in (4),

L 1 = 1 π R sin 2 (λ/2) (2πj -λ) 2 λ 2πj -2d dλ - 1 π R sin 2 (λ/2) (2πj -λ)(2πj + λ) λ 2πj -2d dλ , (20)
and

L 2 = 1 π R sin 2 (λ/2) (2πj -λ) 2 λ 2πj -2d dλ + 1 π R sin 2 (λ/2) (2πj -λ)(2πj + λ) λ 2πj -2d dλ . (21)
Corollary 2. Under the assumptions of Theorem 2, the periodogram I M N defined in (13) satisfies

N D-1 I M N (λ j ) d -→ (Z 2 1 + Z 2 2 ) , as N → ∞ , (22)
where (Z 1 , Z 2 ) is a zero-mean uncorrelated Gaussian vector with covariance matrix equal to

J 2 1 8π(E[ψ (ε 1 )]) 2 Γ , (23)
with Γ defined in (18).

Theorem 2 and Corollary 2 are proved in [START_REF] Fajardo | M-periodogram for the analysis of long-range-dependent time series[END_REF].

Q N -periodogram

Another possible approach to obtain the classical periodogram ( 1) is to write it in terms of the sample autocovariance function

I N (λ j ) = 1 2π N -1 h=-(N -1) γ(h) cos(hλ j ), (24)
where λ j = 2πj/N, j = 1, . . . , [N/2] and γ(h) is the classical sample autocovariance function for a sample {Y 1 , ..., Y N } .

A straightforward approach to robustify 24 is to plug in a robust autocovariance function replacing the classical one. This methodology is now addressed.

For a sample x 1 , ..., x N [START_REF] Rousseeuw | Alternatives to the median absolute deviation[END_REF] proposed a robust scale estimator function Q N (•) which is based on the τ th order statistic of N 2 distances {|x j -x k |, j < k}, and can be written as

Q N (x) = κ × {|x j -x k |; j < k} (τ ) , ( 25
)
where κ is a constant used to guarantee consistency (κ = 2.2191 for the Gaussian distribution) and τ = ( N 2 + 2)/4 + 1. The above function can be evaluated using the algorithm proposed by [START_REF] Croux | Time-efficient algorithms for two highly robust estimators of scale[END_REF], which is computationally efficient.

Based on Q N (•), [START_REF] Ma | Highly robust estimation of the autocovariance function[END_REF] proposed a highly robust estimator for the autocovariance function:

γ Q N (h) = 1 4 Q 2 N -h (u + v) -Q 2 N -h (u -v) , (26)
where u and v are vectors containing the initial N -h and the final N -h observations of x 1 , ..., x N , respectively. The robust estimator for the autocorrelation function is

ρ Q N (h) = Q 2 N -h (u + v) -Q 2 N -h (u -v) Q 2 N -h (u + v) + Q 2 N -h (u -v) . ( 27
)
It can be shown that | ρ Q N (h)|≤ 1 for all h. Now, returning to (24), the robust Q N -periodogram for a sample {Y 1 , ..., Y N } is defined by

I Q N N (λ j ) = 1 2π N -1 h=-(N -1) γ Q N (h) cos(hλ j ), ( 28
)
where

λ j = 2πj/N, j = 1, . . . , [N/2].
The theoretical properties of I Q N N are still under study. Therefore, in the sequel, the asymptotic properties of γQ N are summarized for short and long memory processes. These are well detailed in [START_REF] Lévy-Leduc | Robust estimation of the scale and the autocovariance function of Gaussian short and long-range dependent processes[END_REF]. Theorem 3. Assume that (A8) holds and let h be a non negative integer. Then, the autocovariance estimator γ Q N (h) satisfies the following Central Limit Theorem:

Main asymptotic results for short memory process

√ N ( γ Q N (h) -γ(h)) d -→ N (0, σ2 h ), N → ∞, where σ2 (h) = E[ζ 2 (Y 1 , Y 1+h )] + 2 k≥1 E[ζ(Y 1 , X 1+h )ζ(Y k+1 , Y k+1+h )] , (29)
and the function ζ is defined by

ζ : (x, y) → (γ(0) + γ(h)) IF x + y 2(γ(0) + γ(h)) , Q, Φ -(γ(0) -γ(h)) IF x -y 2(γ(0) -γ(h)) , Q, Φ . ( 30
)
where IF is defined by

IF(x, Q, Φ) = κ 1/4 -Φ(x + 1/κ) + Φ(x -1/κ) R φ(y)φ(y + 1/κ)dy , ( 31
)
where Φ and φ denote the c.d.f. and p.d.f. of a standard Gaussian random variable, respectively with κ defined in (25).

Theorem 3 is proved in Lévy-Leduc et al. ( 2011).

Main asymptotic results for long-memory process

The following results concern the robust autocovariance function for long-memory process see, also, Lévy-Leduc et al. ( 2011).

(A9) (Y i ) i≥1 is a stationary zero-mean Gaussian process with autocovariance γ(h

) = E(Y 1 Y h+1 ) satisfying: γ(h) = h -D L(h), 0 < D < 1 ,
where L is slowly varying at infinity and is positive for large h. Note that, as previously stated, D = 1 -2d.

Theorem 4. Assume that (A9) holds and that L has three continuous derivatives. Assume also that 1), for some in (0, D), as x tends to infinity, for all i = 0, 1, 2, 3, where L (i) denotes the ith derivative of L. Let h be a non negative integer. Then, γ Q N (h) satisfies the following limit theorems as N tends to infinity.

L i (x) = x i L (i) (x) satisfy: L i (x)/x = O(
(i) If D > 1/2, √ N ( γ Q N (h) -γ(h)) d -→ N (0, σ2 (h)) , where σ2 (h) = E[ζ 2 (Y 1 , Y 1+h )] + 2 k≥1 E[ζ(Y 1 , Y 1+h )ζ(Y k+1 , Y k+1+h )] , ζ being defined in (30). (ii) If D < 1/2, β(D) N D L(N ) ( γ Q N (h) -γ(h)) d -→ γ(0) + γ(h) 2 (Z 2,D (1) -Z 1,D (1) 2 )
where β(D) = B((1 -D)/2, D), B denotes the Beta function, the processes Z 1,D (•) and Z 2,D (•) are defined as follows:

Z 1,D (t) = R t 0 (u -x) -(D+1)/2 + du dB(x), 0 < D < 1 , (32)
Z 2,D (t) = R 2 t 0 (u -x) -(D+1)/2 + (u -y) -(D+1)/2 + du dB(x)dB(y), 0 < D < 1/2 , ( 33
)
and

L(N ) = 2L(N ) + L(N + h)(1 + h/N ) -D + L(N -h)(1 -h/N ) -D , ( 34
)
where B is the standard Brownian motion. The symbol means that the domain of integration excludes the diagonal.

Theorem 4 is proved in Lévy- [START_REF] Lévy-Leduc | Robust estimation of the scale and the autocovariance function of Gaussian short and long-range dependent processes[END_REF].

Monte Carlo simulation

In this section, small sample size experiments are conducted with the aim to clarify the empirical performance of the spectral estimates discussed previously in a different context such as time series with additive outliers. Based on this, some standard questions, such as ( 1) what is the best method to be used in a real application? (2) which method ( if any) should be considered when dealing with outliers? (3) Does the large observation ( if any) make similar outlier´s effect on the statistical time series modelling functions, that is, on the ACF and periodogram functions? among others, are expected to be answered or, at least, clarified.

Let {X t } t=1,...,N be a sample from a Gaussian second order stationary process and let {Y t } t=1,...,N be a sample of the process defined by

Y t = X t + ωW t ( 35
)
where the parameter ω represents the magnitude of the outlier, and W t is a random variable with probability distribution

P (W t = -1) = P (W t = 1) = δ/2 and P (W t = 0) = 1 -δ ,
where 35) is based on the parametric models proposed by [START_REF] Fox | Outliers in time series[END_REF]. W t is the product of Bernoulli(δ) and Rademacher random variables; the latter equals 1 or -1, both with probability 1/2. X t and W t are independent random variables. Note that, if ω = 0.0 {Y t } is an outlier free time series.

E[W t ] = 0 and E[W 2 t ] = Var(W t ) = δ. Note that (
In order to compare the performance of M -and Q N -periodogram, a Monte Carlo investigation was carried out under different contamination scenarios. For the simulations, the number of replications was 5000, the samples {X t } of size N = 500 were generated according to a model autocorrelation structure, which is given in what follows, and the contaminated data Y t were generated from ( 35) with δ = 0.01 for magnitudes ω = 0 (no outliers) and 10.

The comparison between the methods is performed by estimating α in the linear regression log(I(λ j )) const + α log(λ j ) + E j , j = 1, . . . , N 0.7 , where I(.) is either I N (.), I M N (.) or I Q N N (.). The data were generated based on

X t = (1 -B) -d Z t = j≥0 Γ(j + d) Γ(j + 1)Γ(d) t-j , (36)
where t is an AR(1) model, that is, t = φ t-1 + η t , where η t , t = 1, ..., N , are i.i.d. standard Gaussian random variables.

In the finite sample size investigation, the model correlation structures are divided in two cases:

1. An AR(1) model with φ = 0.6 and d = 0.

2. An ARFIMA(0, d, 0) model with d = 0.3.

Figure 1 displays the plots of the empirical densities of αI N , αI M N and αI Q N N for the case of AR(1) models without contamination (ω = 0). Although, αI M N has a slight better performance than αI Q N N , that is, the first method and the classical periodogram presented very close densities, all the methods provided similar results showing that, even for small sample sizes, the empirical density is very close which corroborate the theoretical results discussed previously. Based on the asymptotic theory and the empirical results all three methods can be used to estimate the spectral density of a time series when there is no contamination of additive outliers. This opens an important contribution in the context that alternative spectral estimators such as I M N and I Q N N can be used instead of the classical periodogram I N in the step procedure for modelling time series data. For example, these estimators can be an alternative tools to be used in the Whittle function to obtain the parameter estimates. This will be also discussed in what follows. Note that, the disadvantage of

I Q N N over I M
N and I N is that the ACF using Q N (.) does not have the positive definite property. When the data is contaminated with additive outliers the scenario changes significantly. As well known, the periodogram, which depends on the classical autocovariance, is corrupted by the outliers. Therefore, the alternative methods are almost unaffected. This is displayed in Figure 2 in which ω = 10 and δ = 0.01. The empirical density of αI N is shifted to the right side which is an expected result since the variance increases with outliers. The empirical densities of αI M N and αI Q N N remain almost unchangeable. In the case of long-memory process, the empirical density plots are given in Figures 3 and4 are drawn. That is, in the uncontaminated scenarios, all three methods displayed similar densities although the method M and the classical one (periodogram) are very close. In the contaminated case, the classical one is totally affected by the additive outliers. Reinforcing that the ACF using Q N does not have the positiveness property.

- 4 Applications of M and Q N -periodograms

Robust estimation of the fractional parameter

Based on the theoretical results discussed previously, this section introduces some applications related to the use of M -regression and Q N estimation functions. The application is divided in two cases: (a) Estimation of the fractional parameter d in long-memory processes; (b) Estimation in periodic AR (PAR) processes. Some finite sample size investigation is also addressed in the context of time series with and without outliers.

(a) Estimation of the fractional parameter in long-memory process

The estimation methods of the fractional parameter d discussed here are derived from the well-known semi-parametric regression method (GPH) originally proposed by [START_REF] Geweke | The estimation and application of long memory time series model[END_REF]. The regression estimation methods based on I M N and I Q N N were previously introduced in [START_REF] Reisen | An M-estimator for the long-memory parameter[END_REF] and [START_REF] Fajardo | Robust estimation in long-memory processes under additive outliers[END_REF], respectively, papers where the reader will find more details related to theoretical and empirical results of these estimation methodologies. (A10) (ε i ) i≥1 is a stationary mean-zero Gaussian process with spectral density given in Assumption (A7).

For estimating the fractional parameter d of long-memory processes having their spectral density satisfying ( 14), it is usual to use the standard GPH [START_REF] Geweke | The estimation and application of long memory time series model[END_REF]) estimator defined in the following. This estimator is motivated heuristically by starting from 37) where X j = log|2 sin(λ j /2)| and and, by (37),

log(f (λ j )) = -2d log(|2 sin(λ j /2)|) + log(f * (λ j )) = -2dX j + log(f * (λ j )) = log(f 0 ) -2dX j + log(f j /f 0 ), (
f j = f (λ j ). If ε R j = log I N (λ j ) f (λ j ) , ( 38
) then log(I N (λ j )) = ε R j + log(f (λ j )),
log(I N (λ j )) = log(f 0 ) -2dX j + log(f j /f 0 ) + ε R j . ( 39
)
The GPH estimator is given by dGPH = -0.5

m N j=1 (X j -X) log(I LS N (λ j )) m N k=1 (X k -X) 2 , ( 40
)
where

X j = log|2 sin(λ j /2)|, X = m N j=1 X j /m N , I LS N (λ j ) is defined in (9) and m N is a function of N .
Based on the above discussion, one way to define a M -regression estimator of d consists in replacing I LS N in ( 40) by I M N defined in ( 13):

dM = -0.5 m N j=1 (X j -X) log(I M N (λ j )) m N k=1 (X k -X) 2 , ( 41
)
where

X j = log|2 sin(λ j /2)|, X = m N j=1 X j /m N and m N is a function of N which is specified in Theorem 5.
The theoretical properties of dM are established under the following assumptions. The random process (ε j ) is obtained through a moving average process:

ε j = k≤j a j-k ζ k , a j = L(j)j -(1+D)/2 , j ≥ 1 , ( 42
)
for some D in (0, 1), where L(•) is a positive slowly varying function at infinity and where the random variables ζ k are i.i.d. with zero mean and variance 1. It is assumed that the distribution of

ζ 0 satisfies E(e iuζ0 ) ≤ C(1 + |u|) -δ , u ∈ R . ( 43
)
where C < ∞ and δ > 0 are constants. Note that, Conditions ( 42) and ( 43) imply that the cumulative distribution function F ε0 of ε 0 is infinitely boundedly differentiable, see [START_REF] Koul | Second order behavior of M-estimators in linear regression with long-memory errors[END_REF].

Theorem 5. Let Y i = ε i , for all i in {1, . . . , N }, where ε i satisfy ( 42) and (A10). Assume that 1/D is not an integer and that β = 0 in (2). Assume moreover that 42) and satisfies (43), ν 1 = 0, ν 2 = 0 and ν 3 = 0, where the ν k are defined by

E(ζ 4∨2k 0 ) < ∞, where k = [1/D], ζ 0 is defined in (
ν k = ∞ 0 ψ(y) 1 -(-1) k f (k) (y)dy, for all integer k ≥ 0 , ( 44
)
where ψ is the Huber function. Then, if

1/3 < D < 1, √ m N ( dM -d) d -→ N (0, π 2 /24), as N → ∞, ( 45
)
where dM is defined in (41) and

m N = N β with 0 < β < (1 -D)/3.
This result is proved in [START_REF] Reisen | An M-estimator for the long-memory parameter[END_REF].

Another way of defining a robust estimator of d is to consider:

dQ N = -0.5 m N j=1 (X j -X) log(I Q N N (λ j )) m N k=1 (X k -X) 2 , ( 46
)
where

X j = log|2 sin(λ j /2)|, X = m N j=1 X j /m N , I Q N N (λ j
) is defined in (28) and m N is a function of N . For further information, see [START_REF] Fajardo | Robust estimation in long-memory processes under additive outliers[END_REF]. The asymptotic property of dQ N is still an open problem, however, the empirical results given in [START_REF] Fajardo | Robust estimation in long-memory processes under additive outliers[END_REF] support the use of this method under time series with and without outliers. The performance of fractional estimators dGP H , dM and dQ N is the motivation of the next subsection for long-memory time series with and without additive outliers.

4.1.1-Finite sample size investigation

In this subsection, the numerical experiments were carried out in accordance with the model of Section 3. For the simulations, N = 500, ω = 10 and δ = 0.01 for 5000 replications. The results are displayed in Figures 5,6 and Table 1. Since there is not short-memory component in the model m N was fixed at N 0.7 for all tree methods.

Figure 5 presents the boxplots with the results of dGP H , dM and dQ N estimators for the uncontaminated scenario. dM and dQ N seem to present positive bias and, surprisingly, dQ N displays smaller deviation. However, in general, all methods perform similarly, i.e., all estimation methods leaded to comparable estimates close to the real values of d.

Figure 6 displays the boxplots of dGP H , dM and dQ N when the series has outliers. As can be perceived from the boxplots, the GPH estimator is clearly affected by additive outliers while the robust ones keep almost the same picture as the one of the non-contaminated scenario, except that the bias of dQ N becomes negative, that is, this estimator tends to overestimate the true parameter.

The empirical mean, bias and mean square root are displayed in Table 1. This numerically corroborates the results discussed based on Figures 5,6, that is, the estimators have similar performance in the absence of outliers in the data. While the performance of dGP H changes dramatically in the presence of outliers, the estimates from dQ N and dM keep almost unchangeable. As a general conclusion, the empirical result suggests that all the methods can be used to estimate the parameter d when there is not a suspicion of additive or abrupt observation. However, in the existence of a single atypical observation, the methods dQ N and dM should be preferred. Similar conclusions are given in [START_REF] Fajardo | Robust estimation in long-memory processes under additive outliers[END_REF] and [START_REF] Reisen | An M-estimator for the long-memory parameter[END_REF] for dQ N and dM , respectively. where for each season ν ( 1 ≤ ν ≤ S) where S is the period, p ν and q ν are the AR and MA orders, respectively, φ ν,1 , . . . , φ ν,pν and θ ν,1 , . . . , θ ν,qν are the AR and MA coefficients, respectively, and φ ν,0 = θ ν,0 = 1. The sequence {ε t } t∈Z is zero-mean and uncorrelated, and has periodic variances with period S, i.e. E(ε 2 rS+ν ) = σ 2 ν for ν = 1, . . . , S. In the following, p = max ν p ν , q = max ν q ν , φ ν,j = 0 for j > p ν , θ ν,k = 0 for k > q ν , and ( 47) is referred as the PARMA(p, q) S model (see, for example, [START_REF] Basawa | Large sample properties of parameter estimates for periodic ARMA models[END_REF] and [START_REF] Sarnaglia | Periodic ARMA models: Application to particulate matter concentrations[END_REF]).

To deal with outliers effect in the estimation of PAR model, [START_REF] Sarnaglia | Robust estimation of periodic autoregressive processes in the presence of additive outliers[END_REF] proposed the use of the Q N (.) function in this model. Following the same lines of the linear time series model described previously, the Q N (.) function is used to compute an estimator of the periodic autocovariance function γ (ν) (h) at lag h and this sample ACF based on Q N (.) estimator, denoted here as γ (ν) Q (h), replaces the classical periodic ACF γ (ν) (h) in the Yule-Walker periodic equations (see, for example, [START_REF] Mcleod | Diagnostic checking periodic autoregression models with application[END_REF] and [START_REF] Sarnaglia | Robust estimation of periodic autoregressive processes in the presence of additive outliers[END_REF]) to derive an alternative parameter estimator method for a periodic AR model. The authors derived some asymptotic and empirical properties of the proposed estimator. They showed that the method well accommodate the effect of additive outliers, that is, it presented robustness against these type of observations in the finite sample size series as well as in a real data set.

Let now Z 1 , ...Z N , where N = nS, be a sample from PAR process which is a particular case of the model definition in (47) with q ν = 0 and let now Q N (.) for PAR process be defined as

Q (ν) N (Z) = Q N ({Z rS+ν } 0≤r≤N ) . ( 48
) Based on Q (ν)
N (Z), the authors derived the sample ACF for periodic stationary processes γ(ν) Q (h). Under some model assumptions, they proved the following main results.

1. For a fixed lag h, γ(ν) Q (h) satisfies the following central limit theorem:

As N -→ ∞, √ N γ(ν) Q (h) -γ (ν) (h) D -→ N (0, σ2 h ) ,
where γ (ν) (h) is the periodic ACF function and σ2 h is the variance, more details are given in [START_REF] Sarnaglia | Robust estimation of periodic autoregressive processes in the presence of additive outliers[END_REF]. 2. The Q (ν) N Yule-Walker estimators ( φν,i ) 1≤i≤pν ,ν=1,...,S satisfy φν,i -φ ν,i = O P (N -1/2 ) for all i = 1, . . . , p ν and ν in {1, . . . , S}.

Recently, [START_REF] Solci | Empirical study of robust estimation methods for PAR models with application to PM 10 data[END_REF] compared the Yule-Walker estimator (YWE), the robust least squares estimator [START_REF] Shao | Robust estimation for periodic autoregressive time series[END_REF]) and the ACF Q n estimator (γ (ν) Q (h), denoted here RYWE, in the context of estimating the parameters in PAR models with and without outliers. Their main conclusion is similar to the cases discussed previously, that is, for the case of ARFIMA model γ(ν) Q (h) displayed good performance in estimating the parameters in PAR models, periodic samples with and without outliers. As expected, the YWE estimator performed very poorly with the presence of outliers in the data. One of their simulation results is reproduced in the table below (Table 2) in which n = 100, 400 (cycles), S = 4, t is a Gaussian white noise process and δ = 0.01 (outlier's probability) and magnitude ω = 10. The results correspond to the mean of 5000 replications. 0 N (0, 1) 100 0.9 -0.007 0.077 -0.003 0.103 0.8 -0.002 0.065 0.004 0.084 0.7 0.000 0.063 -0.001 0.083 0.6 -0.005 0.066 -0.003 0.083 400 0.9 -0.001 0.037 -0.001 0.047 0.8 -0.001 0.031 0.000 0.038 0.7 -0.001 0.032 0.001 0.038 0.6 0.000 0.032 0.000 0.039 7 N (0, 1) 100 0.9 -0.181 0.247 0.014 0.120 0.8 -0.118 0.176 0.012 0.096 0.7 -0.105 0.157 0.015 0.091 0.6 -0.097 0.151 0.012 0.091 400 0.9 -0.183 0.203 0.017 0.055 0.8 -0.129 0.144 0.012 0.046 0.7 -0.108 0.124 0.013 0.044 0.6 -0.103 0.119 0.014 0.043

As an alternative estimator of φν,i , [START_REF] Sarnaglia | A robust estimation approach for fitting a PARMA model to real data[END_REF] proposed the use of M -periodogram function to obtain estimates of the parameters in PARMA models. The estimator is based on the approximated Whittle function suggested in [START_REF] Sarnaglia | Periodic ARMA models: Application to particulate matter concentrations[END_REF]. Basically, the Whittle M -estimator of PARMA parameters is derived by the ordinary Fourier transform with the nonlinear M -regression estimator for periodic processes in the harmonic regression equation that leads to the classical periodogram. The empirical simulation investigation in [START_REF] Sarnaglia | A robust estimation approach for fitting a PARMA model to real data[END_REF] considered the scenarios of periodic time series with presence and absence of additive outliers. Their small sample size investigation leaded to a very promising estimation method under the context of modelling periodic time series with additive outliers and heavy-tailed distributions. The theoretical justification of the proposed estimator is still an open problem and it is now a current research theme of the authors.

Table 3 displays results of a simple simulation example to show the empirical performance of the Whittle M -estimator with the Huber function ψ(x) [START_REF] Huber | Robust estimation of a location parameter[END_REF]) compared to the maximum Gaussian and Whittle likelihood estimators to estimate a PAR(2) model with parameters φ 1,1 = -0.2, φ 2,1 = -0.5, σ 2 1,1 = 1.0 and σ 2 2,1 = 1.0. The sample sizes are N = nS = 300, 800 (n = 150, 400, respectively) and the Huber function was used with constant equal to 1.345, which ensure that the M -estimator is 95% as efficient as the least squares estimator for univariate multiple linear models with independent and identically distributed Gaussian white noise. The sample root mean square error (RMSE) was computed over 5000 replications. The PAR(2) model with additive outliers was generated with outlier's probability δ = 0.01 and magnitude ω = 10. The values with " * " refer to the RMSE for the contaminated series. In the absence of outliers, in general, all estimators present similar behaviour. Relating to the estimation of the variance of the innovations, the MLE and WLE seem to be more precise which is an expected result since the data is Gaussian with zero-mean and these two methods are asymptotically equivalents. The RMSE of the estimators decreases as the sample size increases. When the simulated data has outliers, as an expected result the MLE and WLE estimates are totally corrupted by the atypical observations while the RWLE estimator presents generally accurate estimates. This simple example of simulation leads to the same conclusions of the models discussed previously in which M -regression method was also considered.

The methods discussed above give strong motivation to use the methodology in practical situations in which periodically correlated time series contain additive outliers. For example, [START_REF] Sarnaglia | Robust estimation of periodic autoregressive processes in the presence of additive outliers[END_REF] applied the robust ACF estimator γ(ν) Q (h) to fit a model for the quarterly Fraser River data. [START_REF] Sarnaglia | A robust estimation approach for fitting a PARMA model to real data[END_REF] and [START_REF] Solci | Empirical study of robust estimation methods for PAR models with application to PM 10 data[END_REF] analysed air pollution variables using the robust methodologies discussed in these papers. In the first paper, the authors considered the daily average SO 2 concentrations and, in the second one, it was analysed the daily average PM 10 concentrations. Both data set were collected at Automatic Air Quality Monitoring Network (RAMQAr) in the Great Vitória Region GVR-ES, Brazil, which is composed by nine monitoring stations placed in strategic locations and accounts for the measuring of several atmospheric pollutants and meteorological variables in the area. In general, the models well fitted the series and all these applied examples revealed outliers effects on the estimates.

  In the short-memory scenario, the process under study (Y i ) i≥1 satisfies the following assumption (see, also, Lévy-Leduc et al. (2011)): (A8) (Y i ) i≥1 is a stationary zero-mean Gaussian process with autocovariance sequence γ(h) = E(Y 1 Y h+1 ) satisfying: h≥1 |γ(h)|< ∞ .

Fig. 2 .

 2 Fig. 1. Densities of αI N , αI M N and αI Q N N

Fig. 3 .Fig. 4 .

 34 Fig. 3. Densities of αI N , αI M N and αI Q N N

Fig. 5 .Fig. 6 .

 56 Fig. 5. Boxplots of dGP H , dM and dQ N when δ = 0.

4. 2

 2 Q n and M -estimators in PARMA models One of the most popular periodic causal process is the PARMA model which generalizes the ARMA model. {Z t } t∈Z is said to be a PARMA model if it satisfies the difference equationpν j=0 φ ν,j Z rS+ν-j = qν k=0 θ ν,k ε rS+ν-k , r ∈ Z(47)

Table 1 .

 1 Empirical Mean, Bias and RMSE of dGP H , dM and dQ N when ω = 10 and δ = 0, 0.01, 0.05.

	d δ	MEAN dGP H dM	dQ N dGP H	BIAS dM	RMSE dQ N dGP H dM	dQ N
		0.0 0.3029 0.2950 0.2933 0.0029 -0.0049 -0.0066 0.0601 0.0596 0.0558
	0.3	0.01 0.2226 0.2899 0.3052 -0.0773 -0.0101 0.0052 0.0972 0.0581 0.0584
		0.05 0.1225 0.2681 0.3236 -0.1775 -0.0318 0.0236 0.1873 0.0689 0.0682

Table 2 .

 2 Bias and RMSE for Model 1 and outliers with probability δ = 0.01.

				YWE	RYWE
	ω	t	n φν,1	Bias RMSE Bias RMSE

Table 3 .

 3 Empirical RMSE results for estimating an PAR(2) model. MLE 800 0.048; 0.101 * 0.079; 1.122 * 0.046; 0.239 * 0.074; 1.253 * 300 0.068; 0.121 * 0.117; 1.368 * 0.079; 0.252 * 0.111; 1.364 * WLE 800 0.048; 0.101 * 0.079; 1.122 * 0.046; 0.239 * 0.074; 1.253 * 300 0.067; 0.067 * 0.147; 0.179 * 0.083; 0.089 * 0.147; 0.189 * RWLE 800 0.051; 0.054 * 0.118; 0.149 * 0.051; 0.058 * 0.108; 0.152 *

	Method N 300 0.067; 0.121 φ1,1	σ 2 1,1	φ2,1	σ 2 2,1

* 0.117; 1.366 * 0.079; 0.252 * 0.111; 1.363 *

Acknowledgements

V. A. Reisen gratefully acknowledges partial financial support from FAPES/ES, CAPES/Brazil and CNPq/Brazil and CentraleSupélec. Màrton Ispàny was supported by the EFOP-3.6.1-16-2016-00022 project. The project is cofinanced by the European Union and the European Social Fund. Paulo Roberto Prezotti Filho and Higor Cotta are Ph.D students under supervision of V. A. Reisen and P. Bondon. The authors would like to thank the referee for the valuable suggestions.

Proof of Theorem 1

By Propositions 1 and 4 and Example 1 of [START_REF] Wu | M-estimation of linear models with dependent errors[END_REF] the assumptions of Theorem 1 of [START_REF] Wu | M-estimation of linear models with dependent errors[END_REF] hold. Thus,

where

Observe that

.

By observing that

tends to zero as N tends to infinity and that the same holds for N -1 N -|k| =1 sin(2 + k), this concludes the proof.