Randal Douc
email: randal.douc@telecom-sudparis.eu

François Roueff
email: roueff@telecom-paristech.fr

Tepmony Sim
email: sim@telecom-paristech.fr

Keywords: 62F12; Secondary: 60J05 consistency, ergodicity, maximum likelihood, observation-driven models, time series of counts

Handy sufficient conditions for the convergence of the maximum likelihood estimator in observation-driven models

Observation-driven time series models have been widely used in various disciplines such as in economics, finance, epidemiology, population dynamics, etc. These models have been introduced by [START_REF] Cox | Statistical analysis of time-series: some recent developments[END_REF] and later considered by [START_REF] Streett | Some observation driven models for time series of counts[END_REF], [START_REF] Davis | Observation-driven models for Poisson counts[END_REF], [START_REF] Fokianos | Poisson autoregression[END_REF], [START_REF] Michael | Absolute regularity and ergodicity of Poisson count processes[END_REF], [START_REF] Doukhan | On weak dependence conditions for Poisson autoregressions[END_REF], [START_REF] Davis | Theory and inference for a class of observationdriven models with application to time series of counts[END_REF] and [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF]. The celebrated GARCH(1, 1) model, see [START_REF] Bollerslev | Generalized autoregressive conditional heteroskedasticity[END_REF], as well as most of the models derived from this one, see [START_REF] Bollerslev | Glossary to arch (garch)[END_REF] for a list of some of them, are typical examples of observation-driven models. Observation-driven models have the nice feature that the associated (conditional) likelihood and its derivatives are easy to compute and the prediction is straightforward. The consistency of the maximum likelihood estimator (in short, MLE) for the class of these models can be cumbersome, except when it can be derived using computations specific to the studied model (the GARCH(1,1) case being one of the most celebrated example). When the observed variable is discrete, general consistency results have been obtained only recently in [START_REF] Davis | Theory and inference for a class of observationdriven models with application to time series of counts[END_REF] or [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF] (see also in [START_REF] Henderson | Stationarity of generalized autoregressive moving average models[END_REF] for the existence of stationary and ergodic solutions to some observation-driven time series models). However, the consistency result of [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF] applies to some restricted class of models and does not cover the case where the distribution of the observations given the hidden variable also depends on an unknown parameter. We now introduce three simple examples, to which the results of [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF] can not be directly applied. The first one is the negative binomial integer-valued GARCH (NBIN-GARCH) model, which was first introduced by [START_REF] Zhu | A negative binomial integer-valued GARCH model[END_REF] as a generalization of the Poisson IN-GARCH model.

The NBIN-GARCH model belongs to the class of integer-valued GARCH models that account for overdispersion (i.e., variability is larger than mean) and potential heavy tails in the high values. In [START_REF] Zhu | A negative binomial integer-valued GARCH model[END_REF], the author applied this model to treat the data of counts of poliomyelitis cases in the USA from 1970 to 1983 reported by the Centres for Disease Control, where data overdispersion was detected. The estimation result showed that NBIN-GARCH(1, 1) outperformed among some commonly used models such as Poisson and Double Poisson models. The NBIN-GARCH(1, 1) model is formally defined as follows.

Example 1 (NBIN-GARCH(1, 1) model). Consider the following recursion.

X k+1 = ω + aX k + bY k , Y k+1 |X 0:k+1 , Y 0:k ∼ N B r, X k+1 1 + X k+1 , (1)
where X k takes values in X = R + , Y k takes values in Z + and θ = (ω, a, b, r) ∈ (0, ∞) 4 is an unknown parameter. In [START_REF] Alexander | Normal mixture garch (1, 1): Applications to exchange rate modelling[END_REF], N B(r, p) denotes the negative binomial distribution with parameters r > 0 and p ∈ (0, 1), that is: if Y ∼ N B(r, p), then P(Y = k) = Γ(k+r) k!Γ(r) (1p) r p k for all k ≥ 0, where Γ stands for the Gamma function. Though substantial analysis on this model has been carried out in the literature, to the best of our knowledge, the consistency of the MLE has not been treated, see the end of the discussions of Section 6 in [START_REF] Zhu | A negative binomial integer-valued GARCH model[END_REF].

The second example is the univariate normal mixture GARCH model (NM-GARCH) proposed by [START_REF] Haas | Mixed normal conditional heteroskedasticity[END_REF] and later considered by [START_REF] Alexander | Normal mixture garch (1, 1): Applications to exchange rate modelling[END_REF]. The NM-GARCH model is another natural extension of GARCH processes, where the usual Gaussian conditional distribution of the observations given the hidden volatility variable is replaced by a mixture of Gaussian distributions given a hidden vector volatility variable. The NM-GARCH model has the ability of capturing time variation in both conditional skewness and kurtosis, while the classical GARCH cannot. In [START_REF] Alexander | Normal mixture garch (1, 1): Applications to exchange rate modelling[END_REF], the NM-GARCH(1, 1) model was applied to examine the data of exchange rates consisting of daily prices in US dollars of three different currencies (British pound, euro and Japanese yen) from 2 January 1989 to 31 December 2002. The empirical evidence suggested the best performance of NM(2)-GARCH(1, 1) when compared to the classical GARCH [START_REF] Alexander | Normal mixture garch (1, 1): Applications to exchange rate modelling[END_REF][START_REF] Alexander | Normal mixture garch (1, 1): Applications to exchange rate modelling[END_REF], standardized symmetric and skewed t-GARCH(1, 1) models applied to this same data. The definition of this model is formally stated as follows.

Example 2 (NM(d)-GARCH(1, 1) model). Let d ∈ N \ {0} and consider the following recursion.

X k+1 = ω + AX k + Y 2 k b , Y k+1 |X 0:k+1 , Y 0:k ∼ G θ (X k+1 ; •) , (2)
dG θ (x; •) dν (y) = d ℓ=1 γ ℓ e -y 2 /2x ℓ (2πx ℓ ) 1/2 , x ∈ (0, ∞) d , y ∈ R ,
where ν is the Lebesgue measure on R,

X k = [X 1,k . . . X d,k ] T takes values in X = R d + ; γ = [γ 1 . . . γ d ]
T a d-dimensional vector of mixture coefficients belonging to the d-dimensional simplex

P d = γ ∈ R d + : d ℓ=1 γ ℓ = 1 , (3)
ω, b are d-dimensional vector parameters with positive and non-negative entries, respectively and A is a d × d matrix parameter with non-negative entries. Here we have θ = (γ, ω, A, b). Note that G θ depends on θ only through the mixture coefficients γ 1 , . . . , γ d . If d = 1, we obtain the usual conditionally Gaussian GARCH(1,1) process. In such a case, since γ = γ 1 = 1, G θ no longer depends on θ. Up to our knowledge, the usual consistency proof of the MLE for the GARCH cannot be directly adapted to this model.

Finally, we consider the following new example, where a threshold is added to the usual INGARCH model in the conditional distribution.

Example 3 (Threshold INGARCH model). Consider the following recursion.

X k+1 = ω + aX k + bY k , Y k+1 |X 0:k+1 , Y 0:k ∼ P (X k+1 ∧ τ ) , (4)
where X k takes values in X = (0, ∞), Y k takes values in Z + and θ = (ω, a, b, τ ) ∈ (0, ∞) 4 is an unknown parameter. Comparing with the usual INGARCH model, a threshold τ has been added in the conditional observation distribution. This corresponds to the practical case where the hidden variable has an influence on the observation up to this threshold.

For a well-specified model, a classical approach to establish the consistency of the MLE generally involves two main steps: first the maximum likelihood estimator (MLE) converges to the maximizing set Θ ⋆ of a limit criterion, and second the maximizing set indeed reduces to the true parameter θ ⋆ , which is usually referred to as solving the identifiability problem. In this paper, we are interested in solving the problem involved in the first step, that is, the convergence of MLE. We extend the convergence result of MLE obtained in [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF], which is valid for a restricted class of models, to a larger class of models in which the three examples introduced above are embedded. More precisely, we show the convergence of MLE in observation-driven models where the probability distributions of observations explicitly depend on the unknown parameters. Moreover, we provide very simple conditions that are easy to check, as shown by the three illustrating examples.

The paper is organized as follows. Specific definitions and notation are introduced in Section 2. Then, Section 3 contains the main contribution of the paper, that is, sufficient conditions for the existence of ergodic solutions and for the consistency of the MLE. These results are then applied in Section 4 to the three examples introduced above. Numerical experiments for the NBIN-GARCH(1, 1) model are given in Section 5. Finally, Section 6 provides the proofs of the main results, mainly inspired from [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF].

Definitions and notation

Consider a bivariate stochastic process {(X k , Y k ) : k ∈ Z + } on X × Y, where (X, d) is a complete and separable metric space endowed with the associated Borel σ-field X and (Y, Y) is a Borel space. Let (Θ, ∆), the set of parameters, be a compact metric space, {G θ : θ ∈ Θ} be a family of probability kernels on X × Y and {(x, y) → ψ θ y (x) : θ ∈ Θ} be a family of measurable functions from (X × Y, X ⊗ Y) to (X, X ). The observation-driven time series model can be formally defined as follows.

Definition 1. A time series {Y k : k ∈ Z + } valued in Y is said to be distributed according to an observation-driven model with parameter θ ∈ Θ if there is a bivariate Markov chain {(X k , Y k ) : k ∈ Z + } on X × Y whose transition kernel K θ satisfies K θ ((x, y); dx ′ dy ′ ) = δ ψ θ y (x) (dx ′ ) G θ (x ′ ; dy ′ ) , (5)
where δ a denotes the Dirac mass at point a. Moreover, we will say that the observation-driven time series model is dominated by some σ-finite measure ν on (Y, Y) if for all x ∈ X, the probability kernel G θ (x; •) is dominated by ν.

In this case we denote by g θ (x; •) its Radon-Nikodym derivative, g θ (x; y) = dG θ (x;•) dν (y), and we always assume that for all (x, y) ∈ X × Y and for all θ ∈ Θ, g θ (x; y) > 0 .

A dominated parametric observation-driven model is thus characterized by the collection {(g θ , ψ θ ) : θ ∈ Θ}. The class of observation-driven time series models is a particular case of partially-observed Markov chains since only Y k 's are observed, whereas X k 's are hidden variables. Note that our notation for observation-driven models is slightly different from that of [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF] where their sequence {Y k } corresponds to our sequence {Y k-1 }. Note also that the process {X k : k ≥ 1} by itself is a Markov chain with transition kernel defined by

R θ (x; A) = ½ A (ψ θ y (x)) G θ (x; dy), x ∈ X, A ∈ X . (6)
However, observation-driven time series models do not belong to the class of hidden Markov models. This can be seen in the following recursive relation, which holds for all k ≥ 0,

X k+1 = ψ θ Y k (X k ) , Y k+1 | F k ∼ G θ (X k+1 ; •) , (7)
where F k = σ (X ℓ , X ℓ+1 , Y ℓ : ℓ ≤ k, ℓ ∈ Z + ) and which can be represented graphically as below. The most popular example is the GARCH(1,1) process, where G θ (x; •) is a centered (say Gaussian) distribution with variance x and ψ θ y (x) is an affine function of x and y 2 . One can readily check that Examples 1 and 2 are other instances of dominated observation-driven models.

X k X k+1 Y k+1 Y k X k+2 Y k+2 ψ θ G θ
The inference about model parameter is carried out by relying on the conditional likelihood of the observations (Y 1 , . . . , Y n ) given X 1 = x for an arbitrary x ∈ X. The corresponding conditional density function with respect to ν ⊗n is, under parameter θ, for all x ∈ X,

y 1:n → n k=1 g θ ψ θ y 1:k-1 (x); y k , (8)
where, for any vector y 1:p = (y 1 , . . . , y p ) ∈ Y p , ψ θ y 1:p is the X → X function obtained as the successive composition of ψ θ y 1 , ψ θ y 2 , ..., and ψ θ yp ,

ψ θ y 1:p = ψ θ yp • ψ θ y p-1 • • • • • ψ θ y 1 , (9)
with the convention ψ θ y s:t (x) = x for s > t. Then, the corresponding (conditional) Maximum Likelihood Estimator (MLE) θx,n of the parameter θ, is defined by θx,n ∈ argmax

θ∈Θ L θ x,n Y 1:n , (10)
where

L θ x,n y 1:n := n -1 n k=1 ln g θ ψ θ y 1:k-1 (x); y k . (11)
In this contribution, we study the convergence of θx,n as n → ∞ for some well-chosen value of x under the assumption that the model is well specified and the observations are in a steady state. This means that we assume that the observations {Y k : k ∈ Z + } are distributed according to Pθ⋆ with θ ⋆ ∈ Θ, where, for all θ ∈ Θ, Pθ denotes the stationary distribution of the observation-driven time series corresponding to the parameter θ. However whether such a distribution is well defined is not always obvious. We will use the following ergodicity assumption.

(A-1) For all θ ∈ Θ, the transition kernel K θ of the complete chain admits a unique stationary distribution π θ on X × Y.

With this assumption, we can now define Pθ . The following notation and definitions will be used throughout the paper.

Definition 2. For any probability distribution µ on X × Y, we denote by P θ µ the distribution of the Markov chain {(X k , Y k ), k ≥ 0} with kernel K θ and initial probability mesure µ. Under Assumption (A-1), we denote by π θ 1 and π θ 2 the marginal distributions of π θ on X and Y, respectively and by P θ and Pθ the probability distributions defined respectively as follows. a) P θ denotes the extension of P θ π θ on the whole line (X × Y) Z . b) Pθ is the corresponding projection on the component Y Z .

The probability distributions P θ and Pθ are more formally defined by setting, for all m ∈ Z and B ∈ Y

⊗(m+Z * + ) , Pθ Y m+Z -× B = P θ X Z × Y m+Z -× B = P θ π θ X m+Z * + × B , (12)
or equivalently, using the canonical functions

Y k , k ∈ Z, Pθ (Y m+1:∞ ∈ B) = P θ (Y m+1:∞ ∈ B) = P θ π θ (Y m+1:∞ ∈ B) . (13)
Here and in what follows, we abusively use the same notation Y k both for the canonical projection defined on Y Z and for the one defined on (X × Y) Z + . We also use the symbols E θ and Ẽθ to denote the expectations corresponding to P θ and Pθ , respectively.

3 Main results

Preliminaries

In this section, we follow the same lines as in [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF] to derive the convergence of the MLE θx,n for a general class of observation-driven models. The approach is to establish that, as the number of observations n → ∞, there exists a [START_REF] Fokianos | Poisson autoregression[END_REF], for some appropriate value of x, can be approximated by

(Y Z , Y ⊗Z ) → (R, B(R)) measurable function p θ (•|•) such that the normalized log-likelihood L θ x,n Y 1:n defined in
n -1 n k=1 ln p θ (Y k |Y -∞:k-1 ) .
To define p θ (•|•), we set, for all y -∞:1 ∈ Y Z -, whenever the following limit is well defined,

p θ (y 1 | y -∞:0 ) =    lim m→∞ g θ ψ θ y -m:0 (x); y 1 if the limit exists, ∞ otherwise. (14)
By (A-1), the process Y is ergodic under Pθ⋆ and provided that Ẽθ⋆ ln

+ p θ (Y 1 |Y -∞:0 ) < ∞ , it follows that lim n→∞ L θ x,n Y 1:n = Ẽθ⋆ ln p θ (Y 1 |Y -∞:0 ) , Pθ⋆ -a.s.
In this paper we show that with probability tending to one, the MLE θx,n eventually lies in a neighborhood of the set

Θ ⋆ = argmax θ∈Θ Ẽθ⋆ ln p θ (Y 1 |Y -∞:0 ) , (15)
which only depends on θ ⋆ . In this contribution, we provide easy-to-check sufficient conditions implying

lim n→∞ ∆( θx,n , Θ ⋆ ) = 0, Pθ⋆ -a.s., (16)
but, for the sake of brevity, we do not precisely determine the set Θ ⋆ . Many approaches have been proposed to investigate this problem, which is often referred to as the identifiability problem. In particular cases, one can prove that Θ ⋆ = {θ ⋆ }, in which case the strong consistency of the MLE follows from ( 16). We will mention a general result which precises how the set Θ ⋆ is related to the true parameter θ ⋆ in Remark 3. For the moment, let us mention that we have

θ ⋆ ∈ Θ ⋆ , (17)
provided that the following assumption holds:

(B-1) For all θ, θ ⋆ ∈ Θ, we have (i) If θ = θ ⋆ , y → p θ (y|Y -∞:0 ) is a density function Pθ⋆ -a.s. (ii) Under Pθ⋆ , the function y → p θ⋆ (y|Y -∞:0 ) is the conditional density function of Y 1 given Y -∞:0 .
Indeed, [START_REF] Michael | Absolute regularity and ergodicity of Poisson count processes[END_REF] follows by writing for all θ ∈ Θ,

Ẽθ⋆ ln p θ⋆ (Y 1 |Y -∞:0 ) -ln p θ (Y 1 |Y -∞:0 ) = Ẽθ⋆ ln p θ⋆ (Y 1 |Y -∞:0 ) p θ (Y 1 |Y -∞:0 ) = Ẽθ⋆ Ẽθ⋆ ln p θ⋆ (Y 1 |Y -∞:0 ) p θ (Y 1 |Y -∞:0 ) Y -∞:0 ,
which is nonnegative under (B-1) since it is the expectation of a conditional Kullback-Leibler divergence.

Convergence of the MLE

In this part, we always assume that (A-1) holds. The following is a list of additional assumptions on which our convergence result relies.

(A-2) There exists a function V :

X → R + such that, for all θ ∈ Θ, π θ 1 ( V ) < ∞.
Remark 1. Assumption (A-2) is usually obtained as a byproduct of the proof of Assumption (A-1), see Section 3.3. It is here stated as an assumption for convenience.

The following set of conditions can readily be checked on g θ and ψ θ .

(B-2) For all y ∈ Y, the function (θ, x) → g θ (x; y) is continuous on Θ × X. (B-3) For all y ∈ Y, the function (θ, x) → ψ θ y (x) is continuous on Θ × X.
The function V appearing in (B-4)(viii) below is the same one as in Assumption (A-2). Moreover, in this condition and throughout the paper we write f V for a real-valued function f and a nonnegative function V defined on the same space X, whenever there exists a positive constant c such that |f (x)| ≤ cV (x) for all x ∈ X.

(B-4) There exist x 1 ∈ X, a closed set X 1 ⊆ X, ̺ ∈ (0, 1), C ≥ 0 and measurable functions ψ :

X 1 → R + , H : R + → R + and φ : Y → R + such that the following assertions hold. (i) For all θ ∈ Θ and (x, y) ∈ X × Y, ψ θ y (x) ∈ X 1 . (ii) sup (θ,x,y)∈Θ×X 1 ×Y g θ (x; y) < ∞.
(iii) For all θ ∈ Θ, n ∈ Z + , x ∈ X, and

y 1:n ∈ Y n , d ψ θ y 1:n (x 1 ), ψ θ y 1:n (x) ≤ ̺ n ψ(x) , (18)
(iv) ψ is locally bounded.

(v) For all θ ∈ Θ and y ∈ Y, ψ(ψ θ y (x 1 )) ≤ φ(y). (vi) For all θ ∈ Θ and (x, x ′ , y) ∈ X 1 × X 1 × Y, ln g θ (x; y) g θ (x ′ ; y) ≤ H(d(x, x ′ )) e C (d(x 1 ,x)∨d(x 1 ,x ′ )) φ(y) , (19)
(vii) H(u) = O(u) as u → 0.

(viii) If C = 0, then, for all θ ∈ Θ,

G θ ln + φ V , (20)
otherwise, for all θ ∈ Θ,

G θ φ V . (21)
Let us now state our main result as follows.

Theorem 3. Assume that (A-1), (A-2), (B-2), (B-3) and (B-4) hold. Then, letting x 1 ∈ X as in (B-4), the function p θ (•|•) defined by ( 14) with x = x 1 satisfies (B-1) and the convergence ( 16) of the MLE holds with the set Θ ⋆ defined by [START_REF] Meyn | Markov Chains and Stochastic Stability[END_REF].

For convenience, the proof is postponed to Section 6.1.

Remark 2. As noticed in [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF], the techniques used to prove Theorem 3 also apply in the misspecified case, where Y is not distributed according to Pθ⋆ . We do not pursue in this direction in this contribution.

The consistency of the MLE then follows from Theorem 3 by the following remark.

Remark 3. In many specific cases, one can show that Θ ⋆ defined by [START_REF] Meyn | Markov Chains and Stochastic Stability[END_REF] is the singleton {θ ⋆ }. However this task appears to be quite difficult in some cases such as Example 3. Instead one can use [START_REF] Douc | The maximizing set of the asymptotic normalized log-likelihood for partially observed markov chains[END_REF]Section 4.2], where it is shown that the assumptions of Theorem 3 imply that Θ ⋆ is exactly the set of parameters θ such that Pθ = Pθ⋆ . Thus we can conclude that the MLE converges to the equivalence class of the true parameter. This type of consistency has been introduced by [START_REF] Leroux | Maximum-likelihood estimation for hidden Markov models[END_REF] in the context of hidden Markov models in order to disentangle the proof of the consistency from the problem of identifiability. Recall that the model is identifiable if and only if the equivalent classes {θ : Pθ = Pθ⋆ } reduce to singletons {θ ⋆ } for all θ ⋆ ∈ Θ.

Ergodicity

In this section, the observation-driven model is studied to prove the condition (A-1). Since this is a "for all θ (...)" condition, to save space and alleviate the notational burden, we will drop the superscript θ from, for example, G θ , R θ and ψ θ and respectively write G, R and ψ, instead.

Ergodicity of Markov chains are usually studied using ψ-irreducibility. This approach is well known to be quite efficient when dealing with fully dominated models, see [START_REF] Meyn | Markov Chains and Stochastic Stability[END_REF]. It is not at all the same picture for observationdriven models, where other tools need to be invoked, see [START_REF] Fokianos | Log-linear poisson autoregression[END_REF][START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF]. Since the ergodicity is studied for a given parameter θ, the ergodicity results of [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF] directly apply, even though observation-driven models are restricted to the case where g does not depend on the unknown parameter θ in this reference. Our main contribution here is to focus on an easy-to-check list of assumptions yielding the ergodicity conditions (A-1) and (A-2). We also provide a lemma (Lemma 5) which gives the construction of the instrumental functions α and φ used in the list of assumptions.

(A-3) The measurable space (X, d) is a locally compact, complete and separable metric space and its associated σ-field X is the Borel σ-field.

(A-4) There exist (λ, β) ∈ (0, 1) × R + and a measurable function V : X → R + such that RV ≤ λV + β and {V ≤ M } is compact for any M > 0.

(A-5) The Markov kernel R is weak Feller, that is, for any continuous and bounded function f defined on X, Rf is continuous and bounded on X.

(A-6) The Markov kernel R has a reachable point, that is, there exists x 0 ∈ X such that, for any x ∈ X and any neighborhood N of x 0 , R m (x; N ) > 0 for at least one positive integer m.

(A-7) We have sup

(x,x ′ ,y)∈X 2 ×Y x =x ′ d(ψ y (x), ψ y (x ′ )) d(x, x ′ ) < 1.
(A-8) There exist a measurable function α from X 2 to [0, 1], a measurable function φ : X 2 → X and a measurable function W : X 2 → [1, ∞) such that the following assertions hold.

(i) For all (x, x ′ ) ∈ X 2 and y ∈ Y,

min g(x; y), g(x ′ ; y) ≥ α(x, x ′ )g φ(x, x ′ ); y . ( 22
)
(ii) For all x ∈ X, W (x, •) is finitely bounded in a neighborhood of x, that is, there exists γ x > 0 such that sup

x ′ ∈B(x,γx) W (x, x ′ ) < ∞. (iii) For all (x, x ′ ) ∈ X 2 , 1 -α(x, x ′ ) ≤ d(x, x ′ )W (x, x ′ ). (iv) sup Y W (ψ y (x), ψ y (x ′ )) G(φ(x, x ′ ); dy) -W (x, x ′ ) < ∞,
where the sup is taken over all (x, x ′ ) ∈ X 2 .

We can now state the main ergodicity result.

Theorem 4. Conditions (A-3), (A-4), (A-5), (A-6), (A-7) and (A-8) imply that K admits a unique stationary distribution

π on X × Y. Moreover π 1 V < ∞ for every V : X → R + such that V V .
The proof of Theorem 4 is postponed to Section 6.2 for convenience. The first conclusion of Theorem 4 can directly be applied for all θ ∈ Θ to check (A-1). The second conclusion can be used to check (A-2). In doing so, one must take care of the fact that although V may depend on θ, V does not.

Assumptions (A-4), (A-5) and (A-6) have to be checked directly on the Markov kernel R defined by [START_REF] Davis | Theory and inference for a class of observationdriven models with application to time series of counts[END_REF]. To this end it can be useful to define, for any given x ∈ X, the distribution

Px := P δx⊗G(x;•) (23)
on (X × Y) Z + , where P µ is defined for any distribution µ on X × Y as in Definition 2. Then the first component process {X k , k ∈ Z + } associated to Px is a Markov chain with Markov kernel R and initial distribution δ x .

We now provide a general framework for constructing α and φ that appear in (A-8).

Lemma 5. Suppose that X = C S for some measurable space (S, S) and C ⊆ R. Thus for all x ∈ X, we write x = (x s ) s∈S , where x s ∈ C for all s ∈ S. Suppose moreover that for all x = (x s ) s∈S ∈ X, we can express the conditional density g(x; •) as a mixture of densities of the form j(x s )h(x s ; •) over s ∈ S. This means that for all t ∈ C, y → j(t)h(t; y) is a density with respect to ν and there exists a probability measure µ on (S, S) such that g(x; y) = S j(x s )h(x s ; y)µ(ds) , y ∈ Y .

(

) 24
We moreover assume that h takes non-negative values and that one of the two following assumptions holds.

(F-1) For all y ∈ Y, the function h(•; y) : t → h(t; y) is non-decreasing.

(F-2) For all y ∈ Y, the function h(•; y) : t → h(t; y) is non-increasing.

For all (x, x ′ ) ∈ X 2 , denoting x ∧ x ′ := (min{x s , x ′ s }) s∈S and x ∨ x ′ := (max{x s , x ′ s }) s∈S , we define α(x, x ′ ) and φ(x, x ′ ) as

       α(x, x ′ ) = inf s∈S j(x s ∨ x ′ s ) j(x s ∧ x ′ s ) and φ(x, x ′ ) = x ∧ x ′ under (F-1) ; α(x, x ′ ) = inf s∈S j(x s ∧ x ′ s ) j(x s ∨ x ′ s ) and φ(x, x ′ ) = x ∨ x ′ under (F-2) .
Then α and φ defined above satisfy (A-8)(i).

Proof. We only prove this result under Condition (F-1). The proof is similar under (F-2). Since for all t ∈ C, y → j(t)h(t; y) is a density with respect to ν, we have

j(t) = h(t; y)ν(dy) -1 > 0 .
Thus j is non-increasing on C. Clearly, the defined α takes values on [0, 1] and φ defines a function from X 2 to X. For all (x, x ′ ) ∈ X 2 and y ∈ Y, we have

g(x; y) = S j(x s )h(x s ; y)µ(ds) ≥ S j(x s ∨ x ′ s )h(x s ∧ x ′ s ; y)µ(ds) ≥ S j(x s ∨ x ′ s ) j(x s ∧ x ′ s ) j(x s ∧ x ′ s )h(x s ∧ x ′ s ; y)µ(ds) ≥ S inf s∈S j(x s ∨ x ′ s ) j(x s ∧ x ′ s ) j(x s ∧ x ′ s )h(x s ∧ x ′ s ; y)µ(ds) = α(x, x ′ )g(φ(x, x ′ ); y) .
By symmetry of α and φ, we get (22) and thus (A-8)(i) holds.

Examples

Let us now apply these results to prove the convergence of MLE of Examples 1, 2 and 3.

NBIN-GARCH model

Example 1 is a specific case of Definition 1 where ν is the counting measure on Y = N,

ψ θ y (x) = ω + ax + by , (25)
g θ (x; y) = Γ(y + r) y!Γ(r) 1 1 + x r x 1 + x y , (26)
with θ = (ω, a, b, r) in a compact subset Θ of (0, ∞) 4 and X = (0, ∞).

In [20, Theorem 1], the equation satisfied by the mean of the observations

µ k = E[Y k ] is derived and is shown to admit a constant solution if and only if rb + a < 1 . (27)
This clearly implies that this condition is necessary to have a stationary solution {Y k } with finite mean. However it does not imply the existence of such a solution. In fact, the following result shows that ( 27) is indeed a necessary and sufficient condition to have a stationary solution {Y k } with finite mean. It also shows that all the assumptions of Theorem 3 hold, which, with Remark 3, provides the consistency of the MLE θx 1 ,n for any x 1 ∈ X.

Theorem 6. Suppose that all θ = (ω, a, b, r) in Θ satisfy Condition (27).

Then Assumptions (A-1), (A-2), (B-2), (B-3) and (B-4) hold with V being defined as the identity function on X and with any x 1 ∈ X.

Proof. For convenience, we divide the proof into two steps.

Step 1. We first prove Assumptions (A-1) and (A-2) by applying Theorem 4. We set V (x) = V (x) = x and thus we only need to check (A-3), (A-4), (A-5), (A-6), (A-7) and (A-8). Condition (A-3) holds. We have for all θ ∈ Θ,

RV (x) = ω + (a + br)x = (a + br)V (x) + ω,
which yields (A-4). The fact that the kernel R is weak Feller easily follows by observing that, as p → p ′ , N B(r, p) converges weakly to N B(r, p ′ ), so (A-5) holds. We now prove (A-6). Let x ∞ = ω/(1-a). Let x ∈ R and define recursively the sequence x 0 = x, x k = ω + ax k-1 for all positive integers k. Since 0 < a < 1, this sequence converges to the fixed point x ∞ . Therefore, defining P x as in (23), for any neighborhood N of x ∞ , there exists some n such that x n ∈ N and we have

R n (x; N ) = Px (X n ∈ N ) ≥ Px (X k = x k for all k = 1, . . . , n) = Px (Y 0 = . . . = Y n-1 = 0) > 0.
So (A-6) holds. Assumption (A-7) holds since we have for all (x, x ′ , y) ∈

X 2 × Y with x = x ′ , |ψ y (x) -ψ y (x ′ )| |x -x ′ | = a < 1 .
To prove (A-8), we apply Lemma 5 with C = X, S = {1} (so µ boils down to the Dirac measure on {1}). For all (x, y) ∈ X × Y, let j(x) = . Indeed, h satisfies (F-1). Thus by Lemma 5, for all (x, x ′ ) ∈ X 2 and y ∈ Y, we get that

α(x, x ′ ) = 1 + x ∧ x ′ 1 + x ∨ x ′ r ∈ (0, 1] and φ(x, x ′ ) = x ∧ x ′ satisfy (A-8)(i).
For any given r > 0, let a function W :

X 2 → [1, ∞) be defined by, for all (x, x ′ ) ∈ X 2 , W (x, x ′ ) = 1 ∨ r.
By definition of W , as a constant function, (A-8)(ii) and (A-8)(iv) clearly hold. Moreover, (A-8)(iii) holds since for all (x, x ′ ) ∈ X 2 , we have that

1 -α(x, x ′ ) ≤ (1 ∨ r)|x -x ′ | = W (x, x ′ )|x -x ′ | .
Therefore, (A-8) holds, which completes Step 1.

Step 2. We now prove (B-2), (B-3) and (B-4). By assumption on Θ, then there exists (ω, ω, b, b, r, r, α, ᾱ) ∈ (0, ∞) 6 × (0, 1) 2 such that

ω ≤ ω ≤ ω, ≤ b ≤ b, r ≤ r ≤ r, α ≤ a + br ≤ ᾱ .
Clearly, (B-2) and (B-3) hold by definitions of ψ θ y (x) and g θ (x; y). It remains to check (B-4) for a well-chosen closed subset X 1 and any x 1 ∈ X. Let X 1 = [ω, ∞) ⊂ X so that (B-4)(i) holds. By noting that for all (θ, x, y) ∈ Θ×X×Y, g θ (x; y) ≤ 1, we have (B-4)(ii). From ( 9) and (25), we have for all s ≤ t, y s:t ∈ Y t-s+1 , x ∈ X and θ ∈ Θ,

ψ θ y s:t (x) = ω 1 -a t-s+1 1 -a + a t-s+1 x + b t-s j=0 a j y t-j . (28)
Using (28), we have, for all θ ∈ Θ, x ∈ X and y 1:n ∈ Y n ,

ψ θ y 1:n (x 1 ) -ψ θ y 1:n (x) = a n |x 1 -x| ≤ ᾱn |x 1 -x| .
This gives (B-4)(iii) and (B-4)(iv) by setting ̺ = ᾱ < 1 and ψ(x) = |x 1 -x|. Next we set φ, H and C to meet Conditions (B-4)(v) and (B-4)(vi) and (B-4)(vii). Let us write, for all θ ∈ Θ and y ∈ Y,

x 1 -ψ θ y (x 1 ) ≤ ω + (1 + a)x 1 + by ≤ ω + (1 + ᾱ)x 1 + by and, for all (x, x ′ ) ∈ X 2 1 = [ω, ∞) 2 , ln g θ (x; y) g θ (x ′ ; y) = (r + y) ln(1 + x ′ ) -ln(1 + x) + y ln x -ln x ′ ≤ (r + y)(1 + ω) -1 + y ω -1 |x -x ′ | ≤ r + y (1 + ω -1 ) |x -x ′ | . Setting φ(y) = ω∨r+(1+ ᾱ)x 1 + b ∨ (1 + ω -1
) y, H(x) = x and C = 0 then yield Conditions (B-4)(v), (B-4)(vi) and (B-4)(vii). Now (B-4)(viii) follows from ln

+ y G θ (x, dy) ≤ y G θ (x, dy) = rx ≤ r V (x) .
This concludes the proof.

NM-GARCH model

The NM(d)-GARCH(1, 1) of Example 2 is a specific case of Definition 1 where X = R d + and ν is the Lebesgue measure on Y = R,

ψ θ y (x) = ω + Ax + y 2 b , (29)
g θ (x; y) = d ℓ=1 γ ℓ e -y 2 /2x ℓ (2πx ℓ ) 1/2 , (x, y) ∈ X × Y , (30)
and θ = (γ, ω, A, b) ∈ Θ, a compact subset of P d × (0, ∞) d × R d×d + × R d + , with P d defined by [START_REF] Bollerslev | Glossary to arch (garch)[END_REF].

In [START_REF] Haas | Mixed normal conditional heteroskedasticity[END_REF], it is shown that the equation satisfied by the variance of a univariate NM(d)-GARCH(1, 1) process admits a constant solution if and only if

|λ| max (A + bγ T ) < 1 , (31)
where, for any square matrix M, |λ| max (M) denotes the spectral radius of M. It follows that the existence of a weakly stationary solution implies (31) but it does not say anything about the existence of stationary or weakly stationary solution. The result below shows that (31) is indeed a sufficient condition for the existence of a stationary solution with finite variance. It moreover provides with Theorem 3 and Remark 3 the consistency of the MLE θx 1 ,n for any x 1 ∈ X.

Theorem 7. Suppose that all θ = (γ, ω, A, b) in Θ satisfy Condition (31).

Then Assumptions (A-1), (A-2), (B-2), (B-3) and (B-4) hold with V being defined as any norm on X.

Proof. In this proof section, we set

V (x) = |x| = d ℓ=1 |x ℓ | , (32)
for all x = (x ℓ ) ∈ X. As in Theorem 6, we divide the proof into two steps.

Step 1. We first show that Assumptions (A-1) and (A-2) hold with the above V by applying Theorem 4. Define V on X by setting

V (x) = (1 + x 0 ) T x ,
where 1 is the vector of X with all entries equal to 1 and x 0 is defined by

1 + x 0 = (I -A + bγ T T ) -1 1 .
We indeed note that by Condition (31) the above inversion is well defined and moreover

(I -(A + bγ T ) T ) -1 = I + k≥1 A T + γb T k ,
and, since A, b, γ all have non-negative entries, it follows that x 0 has nonnegative entries. Thus, for all x = (x ℓ ) ∈ X,

V (x) = 1 T x ≤ V (x) ,
so that V V . Hence by Theorem 4, we thus only need to check (A-3), (A-4), (A-5), (A-6), (A-7) and (A-8) with V defined as above for a given θ = (γ, ω, A, b) ∈ Θ (so we drop θ in the notation in the remaining of Step 1). Condition (A-3) holds for any metric d associated to a norm on the finite dimensional space X. (The precise choice of d is postponed to the verification of (A-7).) We have

RV (x) = V (ω + Ax + y 2 b) G(x, dy) = (1 + x 0 ) T ω + (1 + x 0 ) T A + bγ T x = V (ω) + 1 T (I -A + bγ T ) -1 A + bγ T -I + I x = V (ω) + x T 0 x ≤ V (ω) + λV (x) ,
where we set λ = max ℓ {x 0,ℓ /(1 + x 0,ℓ )} < 1. Hence (A-4) holds. Condition (A-5) easily follows from the continuity of the Gaussian distribution with respect to its variance parameter. We now prove (A-6). From ( 9) and (29), we have for all n ≥ 1, y 0:n-1 ∈ Y n and x ∈ X,

ψ θ y 0:n-1 (x) = A n x + n-1 j=0 A j (ω + y 2 n-1-j b) . ( 33
)
Let us use the norm

M = max j i |M i,j | = sup |x|≤1 |Mx| on d × d matrices.
Note that by (31), there exists δ ∈ (0, 1) and c > 0 such that, for any k ≥ 1,

A + bγ T k ≤ c δ k . ( 34
)
Using that A, b, γ all have nonnegative entries, we have

A k ≤ A + bγ T k . ( 35
)
Hence (I-A) -1 = I+ k≥1 A k is well defined and we set x ∞ = (I -A) -1 ω so that, with (29), we have

ψ θ y 0:n-1 (x) -x ∞ = A n x + j≥n A j ω + n-1 j=0 y 2 n-1-j A j b .
Then, using definition (23), we get that, Px -a.s., for all n ≥ 1,

|X n -x ∞ | = |ψ Y 0:n-1 (x) -x ∞ | ≤ |A n (x -x ∞ )| + j≥n A j ω + max 0≤j≤n-1 Y 2 j n-1 j=0 A j b .
With (34) and ( 35), this implies

Px |X n -x ∞ | ≤ c δ n |x -x ∞ | + |ω| 1 -δ + |b| 1 -δ max 0≤j≤n-1 Y 2 j = 1 .
To obtain (A-6), it is sufficient to observe that, since g takes positive values in (30), for any positive ǫ, x ∈ X and any n ≥ 1,

Px max

0≤j≤n-1

Y 2 j < ǫ > 0 .
Next we prove (A-7). We have

ψ y (x) -ψ y (x ′ ) = A(x -x ′ ) .
Since (34) and ( 35) imply that |λ| max (A) < 1, there exists a vector norm which makes A strictly contracting. Choosing the metric d on X as the one derived from this norm, we get (A-7). To show (A-8), we again rely on Lemma 5. Let us set C = (0, ∞) and S = {1, . . . , d} and define the probability measure µ on S by µ({s}) = γ s , for all s ∈ S. For all (t, y) ∈ C×Y, let j(t) = 1 (2πt) 1/2 and h(t; y) = exp -y 2 /2t . Obviously, Relation (24) holds and h satisfies (F-1). Hence, Lemma 5 implies that α and φ defined respectively for all x = (x 1 , . . . , x d ),

x ′ = (x ′ 1 , . . . , x ′ d ) ∈ X by α(x, x ′ ) = min 1≤ℓ≤d x ℓ ∧ x ′ ℓ x ℓ ∨ x ′ ℓ 1 2 ∈ (0, 1] and φ(x, x ′ ) = (x 1 ∧x ′ 1 , . . . , x d ∧x ′ d ), satisfy (A-8)(i). For x = (x 1 , . . . , x d ), x ′ = (x ′ 1 , . . . , x ′ d ) ∈ X, we have 1 -α(x, x ′ ) = 1 -min 1≤ℓ≤d 1 - |x ℓ -x ′ ℓ | x ℓ ∨ x ′ ℓ 1 2 ≤ max 1≤ℓ≤d |x ℓ -x ′ ℓ | x ℓ ∨ x ′ ℓ ≤ min 1≤ℓ≤d (x -1 ℓ ∧ x ′-1 ℓ ) |x -x ′ | ≤ W (x, x ′ ) d(x, x ′ ) ,
where d is the metric previously defined and W is defined by

W (x, x ′ ) = 1 ∨ c d min 1≤ℓ≤d (x -1 ℓ ∧ x ′-1 ℓ
) with c d > 0 is conveniently chosen (such a constant exists since d is the metric associated to a norm and X has finite dimension). Then (A-8)(ii) and (A-8)(iii) hold and, since for all y ∈ Y and x ∈ X, ψ y (x) has all its entries bounded from below by the positive entries of ω, W (ψ y (x), ψ y (x ′ )) is uniformly bounded over (x, x ′ , y) ∈ X × X × Y and (A-8)(iv) holds. This completes Step 1.

Step 2 We now show that Assumptions (B-2), (B-3) and (B-4) hold.

Clearly, (B-2) and (B-3) hold by definitions of ψ θ y (x) and g θ (x; y). It remains to show (B-4). Since Θ is compact, then Lemma 12], we note that this implies that, for all δ ∈ (ρ, 1), there exists C > 0 such that for all k ≥ 1 and all θ ∈ Θ,

ω ≤ min 1≤ℓ≤d ω ℓ , |ω| ≤ ω, b ≤ |b| ≤ b, |λ| max (A + bγ T ) ≤ ρ, A + bγ T ≤ L for some (ω, ω, b, b, ρ) ∈ (0, ∞) 4 × (0, 1) and L > 0. By [16,
A + bγ T k ≤ C δk . (36)
We set

X 1 = [ω, ∞) d ⊂ X so that (B-4)(i) holds. Moreover, for all (θ, x, y) ∈ Θ × X 1 × Y, g θ (x; y) ≤ (2πω) -1/2
. Thus, Condition (B-4)(ii) holds. Now let

x 1 ∈ X. Using (33), ( 36) and (35), we have, for all x ∈ X, y 1:n ∈ Y n and θ ∈ Θ,

ψ θ y 1:n (x 1 ) -ψ θ y 1:n (x) = |A n (x 1 -x)| ≤ C δn |x 1 -x| .
Using that the norm defining d is equivalent to the norm

|•|, we get (B-4)(iii) with ψ(x) = C′ |x 1 -x| ,
for some positive constant C′ . Hence (B-4)(iv) holds and since

x 1 -ψ θ y (x 1 ) ≤ (L + 1) |x 1 | + ω + y 2 b , we also get (B-4)(v) provided that φ(y) ≥ (L + 1) |x 1 | + ω + y 2 b . (37)
It is straightforward to show that, for all θ ∈ Θ, x ∈ X 1 , y ∈ R, and ℓ ∈ {1, . . . , d},

∂ ln g θ ∂x ℓ (x; y) ≤ 1 2 y 2 ω 2 + 1 ω .
Thus, by the mean value theorem, for all θ ∈ Θ, (x, x ′ ) ∈ X 1 × X 1 and y ∈ Y, ln g θ (x; y)ln g θ (x ′ ; y) ≤ 1 2

y 2 ω 2 + 1 ω |x -x ′ | .
We thus obtain (B-4)(v), (B-4)(vi) and (B-4)(vii) by setting C = 0,

H(u) = sup d(x,x ′ )≤u |x -x ′ | , and φ(y) = (L + 1) |x 1 | + ω + 1/(2ω) + y 2 ( b + ω 2 ) .
In addition, for all θ ∈ Θ and x ∈ X, we have

y 2 G θ (x, dy) = γ T x .
Hence, using (32) with the above definitions, we obtain (B-4)(viii) and the proof is concluded.

The Threshold INGARCH model

The threshold INGARCH(1, 1) in Example 3 is a specific case of Definition 1 where ν is the counting measure on Y = Z + ,

ψ θ y (x) = ω + ax + by , (38)
g θ (x; y) = e -(x∧τ ) (x ∧ τ ) y y! , (39)
with θ = (ω, a, b, τ ) in a compact subset Θ of (0, ∞) 4 and X = (0, ∞). In this model, if a < 1, we then have the ergodicity and consistency results as stated in Theorem 8 below.

Theorem 8. Suppose that all θ = (ω, a, b, τ ) in Θ satisfy a < 1. Then Assumptions (A-1), (A-2), (B-2), (B-3) and (B-4) hold with V being defined as the identity function on X and with any x 1 ∈ X.

Proof. As in the proofs of the two theorems above, for convenience, we divide the proof into two steps.

Step 1. We first prove Assumptions (A-1) and (A-2) by applying Theorem 4. We set V (x) = V (x) = x and thus we only need to check (A-3), (A-4), (A-5), (A-6), (A-7) and (A-8). Condition (A-3) holds with the usual metric on R.

We have for all θ ∈ Θ,

RV (x) = ω + ax + b(x ∧ τ ) ≤ aV (x) + (ω + bτ ),
which yields (A-4). The fact that the kernel R is weak Feller easily follows by observing that, as x → x ′ , P(x) converges weakly to P(x ′ ) and the map x → x ∧ τ is continuous, so (A-5) holds. The proof of (A-6) is similar to the NBIN-GARCH case of Theorem 6 and is thus omitted. Assumption (A-7) holds since we have for all (x, x ′ , y) ∈

X 2 × Y with x = x ′ , |ψ y (x) -ψ y (x ′ )| |x -x ′ | = a < 1 .
To prove (A-8), we apply Lemma 5 with C = X, S = {1} (so µ boils down to the Dirac measure on {1}). For all (x, y) ∈ X × Y, let j(x) = e -(x∧τ ) and h(x; y) = (x∧τ ) y y! . Then h indeed satisfies (F-1). Thus by Lemma 5, for all (x, x ′ ) ∈ X 2 and y ∈ Y, we get that

α(x, x ′ ) = e -(x∨x ′ )∧τ +(x∧x ′ )∧τ ∈ (0, 1] and φ(x, x ′ ) = x ∧ x ′ satisfy (A-8)(i).
Let W (x, x ′ ) = 1 for all (x, x ′ ) ∈ X 2 , which is a constant function. Thus (A-8)(ii) and (A-8)(iv) clearly hold. Moreover, (A-8)(iii) holds since for all (x, x ′ ) ∈ X 2 , we have that

1 -α(x, x ′ ) ≤ x ∨ x ′ -x ∧ x ′ = |x -x ′ | = W (x, x ′ )|x -x ′ | .
Therefore, (A-8) holds, which completes Step 1.

Step 2. We now prove (B-2), (B-3) and (B-4). By assumption on Θ, then there exists (ω, ω, b, b, τ , τ , α, ᾱ) ∈ (0, ∞) 6 × (0, 1) 2 such that

ω ≤ ω ≤ ω, ≤ b ≤ b, τ ≤ τ ≤ r, α ≤ a ≤ ᾱ .
Clearly, (B-2) and (B-3) hold by definitions of ψ θ y (x) and g θ (x; y). It remains to check (B-4) for a well-chosen closed subset X 1 and any x 1 ∈ X. Let X 1 = [ω, ∞) ⊂ X so that (B-4)(i) holds. By noting that for all (θ, x, y) ∈ Θ×X×Y, g θ (x; y) ≤ 1, we have (B-4)(ii). From ( 9) and (38), we have for all s ≤ t, y s:t ∈ Y t-s+1 , x ∈ X and θ ∈ Θ,

ψ θ y s:t (x) = ω 1 -a t-s+1 1 -a + a t-s+1 x + b t-s j=0 a j y t-j . (40)
Using (40), we have, for all θ ∈ Θ, x ∈ X and y 1:n ∈ Y n ,

ψ θ y 1:n (x 1 ) -ψ θ y 1:n (x) = a n |x 1 -x| ≤ ᾱn |x 1 -x| .
(x, x ′ ) ∈ X 2 1 = [ω, ∞) 2 , ln g θ (x; y) -ln g θ (x ′ ; y) = x ′ ∧ τ -x ∧ τ + y ln(x ∧ τ ) -ln(x ′ ∧ τ ) ≤ 1 + (ω ∧ τ ) -1 y |x -x ′ | .
Setting φ(y) = 1+ ω+(1+ ᾱ)x 1 + b ∨ (ω ∧ τ ) -1 y, H(x) = x and C = 0 then yield Conditions (B-4)(v), (B-4)(vi) and (B-4)(vii). Now (B-4)(viii) follows from ln

+ y G θ (x, dy) ≤ y G θ (x, dy) = x ∧ τ ≤ V (x) .
This concludes the proof.

5 Numerical experiments

Numerical procedure

In this part we provide a numerical method for computing the (conditional) MLE θx,n for the parameter θ = (ω, a, b, r) in the NBIN-GARCH(1, 1) model introduced in Example 1 and studied in Section 4.1. It is convenient to write θ = (ϑ, r) with ϑ = (ω, a, b) and then to write ψ ϑ y (x) and g r (x; y) instead of ψ θ y (x) and g θ (x; y) in ( 25) and (26), respectively. In contrast to the approach used in [START_REF] Zhu | A negative binomial integer-valued GARCH model[END_REF], we allow the component r to be any positive real number, rather than a discrete one and to be unknown as well. We thus maximize jointly with respect to the parameters ϑ and r the log-likelihood function L θ

x,n y 1:n = L (ϑ,r)

x,n y 1:n . In practice, one does not rely on a compact set Θ of parameters as in Theorem 6. Instead the maximization is performed over all parameters ω > 0, a > 0, b > 0, r > 0 such that the stability constraint a + br < 1 holds (taken from ( 27)). We use the constrained nonlinear optimization function auglag (Augmented Lagrangian Minimization Algorithm) from the package alabama (Augmented Lagrangian Adaptive Barrier Minimization Algorithm) in R. For this purpose we provide an initial parameter point and a numerical computation of the normalized log-likelihood function L θ x,n y 1:n and of its gradient. The initial point is obtained by applying a conditional least square (CLS) estimation based on an ARMA(1, 1) representation of the model, see [START_REF] Zhu | A negative binomial integer-valued GARCH model[END_REF]Section 3]. The computation of the log-likelihood and of its derivatives are derived as follows. For all x ∈ X, denoting u ϑ k = ψ ϑ y 1:k-1 (x) for all k ≥ 2 and u ϑ 1 = x, we have

L (ϑ,r) x,n y 1:n = n -1 n k=1 ln g r ψ ϑ y 1:k-1 (x); y k = n -1 ln g r (x, y 1 ) + n -1 n k=2 ln g r u ϑ k ; y k .
The computation of u ϑ k for all k ≥ 2 is done iteratively by observing that u ϑ k = ψ ϑ y k-1 (u ϑ k-1 ) and the computation of L (ϑ,r)

x,n y 1:n is deduced. The computation of the derivatives with respect to parameter θ = (ϑ, r) of the function L (ϑ,r)

x,n y 1:n are then obtained in two steps. First, for k ≥ 2, the derivative of u ϑ k with respect to ϑ are obtained iteratively by ∂u ϑ 1 /∂ϑ = 0 and

∂u ϑ k ∂ϑ = (1, u ϑ k-1 , Y k-1 ) + a ∂u ϑ k-1 ∂ϑ .
Then the derivatives of L (ϑ,r)

x,n y 1:n with respect to ϑ and r are given by

∂L (ϑ,r) x,n ∂ϑ = n -1 n k=1 ∂ln g r ∂x u ϑ k ; y k ∂u ϑ k ∂ϑ = n -1 n k=2 y k u ϑ k - y k + r 1 + u ϑ k ∂u ϑ k ∂ϑ and ∂L (ϑ,r) x,n ∂r = n -1 n k=1 ∂ln g r ∂r u ϑ k ; y k = n -1 n k=1 Γ 2 (r + y k ) -ln(1 + u ϑ k ) -Γ 2 (r) ,
respectively, where Γ 2 is the digamma function Γ 2 (r) = d dr ln •Γ(r), r > 0.

Simulation study

We consider two NBIN-GARCH(1, 1) models with parameters:

(M.1) θ ⋆ = (ω ⋆ , a ⋆ , b ⋆ , r ⋆ ) = (3, .2, .2, 2) and (M.2) θ ⋆ = (ω ⋆ , a ⋆ , b ⋆ , r ⋆ ) = (3, .35, .1, 1.5).
We simulated m = 200 data sets for each sample size n = 2 7 , 2 8 , 2 9 and 2 10 . In Figure 2, we display the obtained boxplots of the difference of the normalized log-likelihood function evaluated respectively at MLE and at the true value θ ⋆ . As predicted by the theory, this difference appears to converge to 0 as the number of observations n → ∞. For the NBIN-GARCH(1, 1) model, it can be shown that Θ ⋆ = {θ ⋆ }, which implies the convergence of the MLE to the true parameter. We can observe this behavior for each component of the MLE for the two models in Figure 3 and Figure 4. We also report the Monte Carlo mean along with the mean absolute deviation error (MADE): MADE = m -1 m j=1 | θj x,nθ ⋆ j | as an evaluation criterion for the estimated parameter in Table 1. which is in (0, 1) by (A-7). Then

d(X, X ′ ) d(x, x ′ ) = d(ψ Y (x), ψ Y (x ′ )) d(x, x ′ ) ≤ ρ .
Therefore, Condition (50) holds for all (x, x ′ ) ∈ X 2 with ρ as in (53).

We conclude this section with the postponed Proof of Lemma 14. Let (x, x ′ ) ∈ X 2 . We define Ḡ((x, x ′ ); •) as the distribution of (Y, Y ′ , U ) drawn as follows. We first draw a random variable Ȳ taking values in Y with density g(φ(x, x ′ ); •) with respect to ν. Then we define (Y, Y ′ , U ) by separating the two cases, α(x, x ′ ) = 1 and α(x, x ′ ) < 1.

• Suppose that α(x, x ′ ) = 1. Then from (A-8)(i), we have

G(x; •) = G(x ′ ; •) = G(φ(x, x ′ ); •) .
In this case, we set (Y, Y ′ , U ) = ( Ȳ , Ȳ , 1).

• Suppose now that α(x, x ′ ) < 1. Then, using (22), the functions (1α(x, x ′ )) -1 g(x; •)α(x, x ′ )g(φ(x, x ′ ); •)

and

(1α(x, x ′ )) -1 g(x ′ ; •)α(x, x ′ )g(φ(x, x ′ ); •) , are probability density functions with respect to ν and we let Λ and Λ ′ be two independent random variables taking values in Y drawn with these two density functions, respectively. In this case we draw U independently according to a Bernoulli variable with mean α(x, x ′ ) and set

(Y, Y ′ ) = ( Ȳ , Ȳ ) if U = 1 , (Λ, Λ ′ ) if U = 0 .
One can easily check that the so defined kernel Ḡ satisfies (42) and (43). Moreover, for all (x, x ′ ) ∈ X 2 , Ḡ((x, x ′ ); Y 2 × {1}) = P(U = 1) = α(x, x ′ ) , which is compatible with (46). The kernel R is defined by setting R((x, x ′ ); •) as the conditional distribution of (X, X ′ ) = (ψ Y (x), ψ Y (x ′ )) given that U = 1. To complete the proof of Lemma 14, observe that for any measurable f : X 2 → R + , we have, for all (x, x ′ ) ∈ X 2 such that α(x, x ′ ) > 0,

Note that for all m ≥ 1, ψ θ⋆ Y -m:0 (X -m ) = X 1 P θ⋆ -a.s., hence we get that

ψ θ⋆ Y -∞:0 = X 1 P θ⋆ -a.s. ( 63
)
To complete the proof of (B-1)(ii), we need to show that, under Pθ⋆ , y → g θ⋆ (ψ θ⋆ Y -∞:0 ; y) = g θ⋆ (X 1 ; y) is the conditional density of Y 1 given Y -∞:0 , that is, for any B ∈ Y, But since (63) implies that X 1 is σ(Y -∞:0 )-measurable, X 1 can be removed in the last conditioning, which concludes the proof (B-1)(ii). Finally, it remains to show the uniform convergence (41) in (B-5). By (B-3) and (57), we have, for all θ, θ ⋆ ∈ Θ, k ∈ Z + , Then Lemma 15 implies that, Pθ⋆ -a.s., A k (2) = O(η -k ) for any η ∈ (0, 1). The same property applies similarly to A k (3) by using (60) in place of (61). This yields (41) in the case where C > 0, which concludes the proof.

ψ θ Y -∞:k-1 = ψ θ

Figure 1 :

 1 Figure 1: Graphical representation of the observation-driven model.

  ; y) = Γ(y+r) y!Γ(r)x 1+x y

  This gives (B-4)(iii) and (B-4)(iv) by setting ̺ = ᾱ < 1 and ψ(x) = |x 1 -x|. Next we set φ, H and C to meet Conditions (B-4)(v) and (B-4)(vi) and (B-4)(vii). Let us write, for all θ ∈ Θ and y ∈ Y, x 1ψ θ y (x 1 ) ≤ ω + (1 + a)x 1 + by ≤ ω + (1 + ᾱ)x 1 + by and, for all

Figure 2 :Figure 3 :

 23 Figure2: Boxplots of the differences of log-likelihood functions evaluated at the estimated MLE and the true value for Models (M.1) and (M.2) with sample sizes n = 2 7 , 2 8 , 29 and n = 2 10 , respectively. The red "continuous" line indicates the position of zero.

Figure 4 :

 4 Figure 4: Same as Figure 3 but for Model (M.2).

R 1 =

 1 ((x, x ′ ); f ) = E f (ψ Y (x), ψ Y (x ′ )) | U = E f (ψ Ȳ (x), ψ Ȳ (x ′ )) = G(φ(x, x ′ ); f ) , where f (y) = f (ψ y (x), ψ y (x ′ )) for all y ∈ Y. so that, letting k → ∞, sup θ∈Θ d(ψ θ Y -∞:0 , x 1 ) ≤ m≥0 ̺ m φ (Y -m ) < ∞ , Pθ⋆ -a.s. (61)Let us now prove (B-1). Relation (56) directly yields (B-1)(i). Let us prove (B-1)(ii), hence consider the case θ = θ ⋆ . Using (58), we haved(ψ θ⋆ Y -m:0 (x 1 ), ψ θ⋆ Y -m:0 (X -m )) ≤ ̺ m+1 ψ(X -m ) P θ⋆ -a.s.Since { ψ(X -m )} m≥0 is stationary under P θ⋆ , it is bounded in probability, and since ̺ < 1, for all ǫ > 0, we have lim m→∞ P θ⋆ d ψ θ⋆ Y -m:0 (X -m ), ψ θ⋆ Y -m:0 (x) > ǫ = 0 .

½

  B (y)g θ⋆ (X 1 ; y) ν(dy) = P θ⋆ (Y 1 ∈ B | Y -∞:0 ) .Now, note that, by defintion of P θ⋆ , ½ B (y)g θ⋆ (X 1 ; y) ν(dy) = P θ⋆ (Y 1 ∈ B | X 1 ) = P θ⋆ (Y 1 ∈ B | X 1 , Y -∞:0 ) .

  Y 1:k-1 ψ θ Y -∞:0 , Pθ⋆ -a.s. (64)From (B-4)(iii) and (64), we getd(ψ θ Y 1:k-1 (x 1 ), ψ θ Y -∞:k-1 ) ≤ ̺ k-1 ψ ψ θ Y -∞:0 , Pθ⋆ -a.s.On the other hand (B-4)(iv) and (61) implysup θ∈Θ ψ ψ θ Y -∞:0 < ∞ , Pθ⋆ -a.s. ,(65)which, with the previous display, yields,sup θ∈Θ d(ψ θ Y 1:k-1 (x 1 ), ψ θ Y -∞:k-1 ) = O k→∞ ̺ k Pθ⋆ -a.s. (66)Since X 1 is closed and satisfies Condition (B-4)(i), we have that,ψ θ Y 1:k-1 (x 1 ) and ψ θ Y -∞:k-1 are in X 1 for all k ≥ 2. Thus Condition (B-4)(vi) gives that sup θ∈Θ ln g θ (ψ θ Y 1:k-1 (x 1 ); Y k ) g θ (ψ θ Y -∞:k-1 ; Y k ) ≤ A k (1) × A k (2) × A k (3) × A k (4) Pθ⋆ -a.s. ,whereA k (1) = sup θ∈Θ H d(ψ θ Y 1:k-1 (x 1 ), ψ θ Y -∞:k-1 ) A k (2) = sup θ∈Θ e C d(x 1 ,ψ θ Y -∞:k-1 ) A k (3) = sup θ∈Θ e C d(x 1 ,ψ θ Y 1:k-1 (x 1 )) A k (4) = φ(Y k ) .By (66) and (B-4)(vii), we haveA k (1) = O k→∞ ̺ k Pθ⋆ -a.s.With (59), this yields (41) in the case where C = 0. For C > 0, we further observe that, by (61) and (55), we have, for all θ ⋆ ∈ Θ and k ∈ Z + , Ẽθ⋆ ln + A k (2) ≤ Ẽθ⋆

Table 1 :

 1 Mean of estimates, MADEs (within parentheses) for the NBIN-GARCH(1, 1) models

	Sample size n

Acknowledgement

We are thankful to the editor-in-charge and anonymous referee for the insightful comments and the helpful suggestions that lead to improve this paper.

6 Postponed proofs 6.1 Convergence of the MLE Assumptions (A-1) and (A-2) are supposed to hold throughout this section. The proof of Theorem 3 relies on the approach introduced in [START_REF] Pfanzagl | On the measurability and consistency of minimum contrast estimates[END_REF], which was already used in [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF] for a restricted class of observation-driven models. Our main contribution here is to provide the handy conditions listed in Assumption (B-4). We first show that our conditions imply (B-1) and the following one.

(B-5) There exists x 1 ∈ X such that, for all θ, θ ⋆ ∈ Θ, p θ (Y 1 | Y -∞:0 ) defined as in ( 14) with x = x 1 is finite Pθ⋆ -a.s. Moreover, for all θ ⋆ ∈ Θ, we have

Indeed we have the following lemma.

Lemma 9. Assumptions (B-2), (B-3) and (B-4) imply (B-5) and (B-1).

Proof. See 6.3.

Now the proof of Theorem 3 directly follows from the following lemma.

Lemma 10. Assume that (B-2), (B-3) and (B-4)(i)-(ii) hold and that x 1 satisfies (B-5). Then Θ ⋆ defined by ( 15) is a non-empty closed subset of Θ and (16) holds.

Proof. By [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF]Theorem 33], to obtain [START_REF] Moulines | On recursive estimation for time varying autoregressive processes[END_REF], it is sufficient to show that, for all θ ⋆ ∈ Θ, the two following assertions hold.

(a) Ẽθ⋆ sup θ∈Θ ln

In (B-5), p θ (Y 1 | Y -∞:0 ) is defined Pθ⋆ -a.s. as the limit in ( 14) with x = x 1 . So, Pθ⋆ -a.s., by (B-4

) is bounded by the finite constant appearing in (B-4)(ii). Hence Condition (a) holds. Condition (b) then follows from (41). Since almost sure convergence implies the convergence in probability and Pθ⋆ is shift invariant, the random sequence

converges to zero in Pθ⋆ -probability. Then there exists a subsequence of (U m ) which converges Pθ⋆ -a.s. to zero. Hence, interpreting this convergence as a uniform (in θ) convergence of ln

) is continuous for all m Pθ⋆ -a.s. This is indeed the case by (B-2) and (B-3) and since g θ (x; y) is positive.

Ergodicity

For proving Theorem 4, we first recall a more general set of conditions derived in [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF], which are based on the following definition.

Definition 11. Let Ḡ be a probability kernel from X 2 to Y ⊗2 ⊗ P({0, 1}) satisfying the following marginal conditions, for all (x, x ′ ) ∈ X 2 and B ∈ Y,

and such that the following coupling condition holds

Define the following quantities successively.

• The trace measure of Ḡ((x, x ′ ); •) on the set {(y, y) :

• The measurable function α from X 2 to [0, 1] is defined by

• The kernel R is defined for all (x, x ′ ) ∈ X 2 and A ∈ X ⊗2 by

We can now introduce the so-called contracting condition which yields ergodicity.

(A-9) There exists a kernel Ḡ yielding α and R as in Definition 11, a measurable function W :

such that for all (x, x ′ ) ∈ X 2 and, for all n ≥ 1,

Under Conditions (A-3), (A-4), (A-5), (A-6) and (A-9) and by combining Theorem 6, Proposition 8 and Lemma 7 in [START_REF] Douc | Ergodicity of observationdriven time series models and consistency of the maximum likelihood estimator[END_REF], we immediately obtain the following result.

Theorem 12. Assume (A-3), (A-4), (A-5), (A-6) and (A-9). Then the Markov kernel K admits a unique invariant distribution π and π 1 ( V ) < ∞ for any V :

Assumptions (A-3), (A-4), (A-5) and (A-6) are quite usual and easy to check. The key point to obtain ergodicity is thus to construct Ḡ satisfying (A-9). For this, we can also rely on the following result which is quoted from [7, Lemma 9]. Lemma 13. Assume that there exists (ρ, β) ∈ (0, 1) × R such that for all

Then, (48) and (49) hold.

Now we can prove that our set of conditions is sufficient.

Proof of Theorem 4. We only need to show that (A-7) and (A-8) imply (A-9). We preface our proof by the following lemma. Definition 11 with the same α given in (46). Moreover, the kernel R defined by (47) satisfies, for all (x, x ′ ) ∈ X 2 such that α(x, x ′ ) > 0 and all measurable functions f :

Lemma 14. Assume (A-8)(i). Then one can define a kernel Ḡ as in

Let us conclude the proof of Theorem 4 before proving this lemma. By Lemma 14 and Lemma 13, it remains to check that (50) and (51) hold for all (x, x ′ ) ∈ X 2 . Observe that by definition of R, Condition (A-8)(iv) is equivalent to sup

so we can find β ∈ R such that (51) holds for all (x, x ′ ) ∈ X 2 . Now, let (x, x ′ ) ∈ X 2 and let (X, X ′ ) be distributed according to R((x, x ′ ); •) which is defined in (52). When x = x ′ , then d(X, X ′ ) = 0, implying that Condition (50) holds with any nonnegative ρ. For x = x ′ , let ρ be defined by ρ = sup

Proof of Lemma 9

Under (A-2), Assumptions (B-4)(viii) implies that for all θ ∈ Θ,

and if moreover C > 0,

For proving Lemma 9, we will also make use of [7, Lemma 34] which we quote here for convenience.

Lemma 15. Let {U n } n∈Z + be a stationary sequence of real-valued random variables on (Ω, F, P). Assume that E(ln

Then, for all η ∈ (0, 1), lim k→∞ η k U k = 0 , P-a.s.

Proof of Lemma 9. We first show that p θ (y|Y -∞:0 ) in ( 14) is finite for x = x 1 Pθ⋆ -a.s. By (B-2), this follows by writing

if, for all θ, θ ⋆ ∈ Θ, the limit

For all θ ∈ Θ, m ≥ 0, x ∈ X and y -m:0 ∈ Y m+1 , using (B-4)(iii), we have d(ψ θ y -m:0 (x 1 ), ψ θ y -m:0 (x)) ≤ ̺ m+1 ψ(x) .

Taking x = ψ θ y -m-1 (x 1 ) and using (B-4)(v), we obtain, for all y -m-1:0 ∈ Y m+2 , d(ψ θ y -m:0 (x 1 ), ψ θ y -m-1:0 (x 1 )) ≤ ̺ m+1 φ (y -m-1 ) .

Using (54) and Lemma 15, we have that ∀η ∈ (0, 1),

and thus ψ θ Y -m:0 (x 1 ) m≥0 is a Cauchy sequence Pθ⋆ -a.s. Its limit exists Pθ⋆ -a.s., since (X, d) is assumed to be complete, which defines the X-valued random variable ψ θ Y -∞:0 for all θ, θ ⋆ ∈ Θ when Y has distribution Pθ⋆ -a.s. Thus (57) holds and we further obtain that


domain: spi phys

text:
O Llopis
email: llopis@laas.fr

Z Abdallah

V Auroux

A Fernandez

High spectral purity laser characterization with a self-heterodyne frequency discriminator

Keywords: laser, frequency noise, phase noise, metrology, frequency discriminator I

The performance of a delay line discriminator dedicated to the characterization of high spectral purity lasers is investigated. The system noise floor is obtained using a symmetrical delay. The vibrations are one of the main problems with this approach at the low frequency offsets. However, the system is able to characterize commercially available narrow linewidth lasers, such as external cavity lasers and fiber lasers.

INTRODUCTION

The today availability of small size narrow linewidth lasers, such as semiconductor lasers optically coupled to high quality factor optical resonators (WGM resonators, fiber resonators…), is an opportunity to develop new embedded applications for which the laser phase noise level is mandatory. This is the case, as an example, of high precision interferometers and of microwave to terahertz signal generation using optics. Before being used in these systems, this type of laser has to be characterized in terms of phase or frequency noise. However, phase noise measurement in the optical domain is not as developed as it is for radio frequency (RF) sources, although some commercial systems are already available (ex : OEwaves). The reason is that it is much more difficult to design a high quality optical frequency synthesizer than a low phase noise RF synthesizer. Therefore, the easier way to measure this type of sources is still to set up the frequency discriminator technique.

In this approach, the laser signal is split in two paths, one path is delayed with a delay smaller than the coherence length of the laser (contrarily to linewidth measurement case) and the signals are combined on a photodiode [START_REF] Derickson | Fiber optic Test and Measurement[END_REF][START_REF] Ludvigsen | Laser linewidth measurements using self-homodyne detection with short delay[END_REF][START_REF] Camatel | Narrow linewidth CW laser phase noise characterization methods for coherent transmission system applications[END_REF][START_REF] Llopis | Phase noise measurement of a narrow linewidth CW laser using delay line approaches[END_REF]. In the selfheterodyne case, an acousto-optic modulator is added in one of the paths in order to shift the signal frequency around a few tens of MHz. In this case, the noise is analyzed using an RF phase noise test bench at this frequency [START_REF] Camatel | Narrow linewidth CW laser phase noise characterization methods for coherent transmission system applications[END_REF][START_REF] Llopis | Phase noise measurement of a narrow linewidth CW laser using delay line approaches[END_REF][START_REF] Kefelian | Ultralow frequency noise stabilization of a laser by locking to an optical fiber delay line[END_REF], which prevents the detection of some parasitic contributions of low frequency noise, such as the laser 1/f amplitude noise (but also increase the system cost and complexity).

The main difficulty in this approach is to evaluate the measurement noise floor. In microwave frequency discriminators, it is measured using a bypass of the delay line. However, in the optical case, part of the noise may come from nonlinear effects in the fiber [START_REF] Okusaga | Optical scattering induced noise in RF-photonic systems[END_REF], or from vibration sensitivity of the fiber spool [START_REF] Li | Low vibration sensitivity fiber spools for laser stabilization[END_REF][START_REF] Hati | Cancellation of vibration induced phase noise in optical fibers[END_REF] or of the whole system, and the noise floor has to be measured with the fiber maintained in the system.

In this paper, measurements of two types of high spectral purity lasers using a self-heterodyne system are presented. Then, the noise floor of the system is evaluated using a set of two identical fiber spools. Finally, we focus on the choice of the delay and on the importance of the reduction of vibration sensitivity.

II. MEASUREMENT SYSTEM

The measurement system is depicted in Figures 1 and3. The laser signal is split in two paths. On one path, an acoustooptic modulator shifts the laser frequency of 80 MHz. On the other path, the signal is delayed using a fiber spool, and then recombined to the first path. A polarization controller is added after the fiber spool in order to improve the output signal level. The beat frequency is recovered using a photodiode, amplified and, finally, feeds an RF phase noise measurement bench (in our case, an Agilent E5052B). Using this approach, the measured RF phase noise at 80 MHz is proportional to the laser frequency noise, at least at offset frequencies close to the carrier [START_REF] Camatel | Narrow linewidth CW laser phase noise characterization methods for coherent transmission system applications[END_REF][START_REF] Llopis | Phase noise measurement of a narrow linewidth CW laser using delay line approaches[END_REF]. Far from the carrier, the output spectrum is attenuated with a classical sin x/x shape. This is clearly visible in Figure 2, in which the phase noise of the 80 MHz output is depicted in case of an external cavity semiconductor laser (RIO Inc.) feeding the bench, and using a 2 km fiber spool in the bench. At 100 kHz, the measurement sensitivity locally drops down to zero because of the first zero of the sin x/x function.

Because of this response, we have limited our measurements on the 1 Hz -100 kHz offset range, and we have corrected the measured data using the reverse function of the sin x/x response of the bench and the time delay of the spools we are using. More precisely, the laser frequency noise can be computed from the measured 80 MHz RF phase noise [START_REF] Camatel | Narrow linewidth CW laser phase noise characterization methods for coherent transmission system applications[END_REF]

using equation 1. 20 / ) f ( L m m Hz / Hz in f m RF 10 ) f sin( 2 f S τ π = Δ (1)
were f m is the offset frequency, τ is the time delay between the two paths and L RF (f m ) is the measured single sideband phase noise, in dB c /Hz, of the RF signal.

Another important result of equation 1 is that it gives the proportionality coefficient between the RF phase noise and the laser close to carrier frequency noise, and that this coefficient is itself proportional to the delay τ. Thus, even if a large delay limits the measurement bandwidth because of the sin x/x effect, it enhances the measurement sensitivity and, as we will see, reduces the noise floor.

Equation 1 is used to compute the laser frequency noise data. However, when a noise floor is measured instead of a true laser frequency fluctuation, this noise floor is not affected by the sin x/x response, because it is either a noise superimposed on the system output, or the demodulation of the laser amplitude noise. Therefore, to convert the measured noise floor into frequency noise, we use a constant coefficient which is the limit of the above coefficient at low frequencies. This frequency discriminator has been used to measure two different high spectral purity lasers: a semiconductor laser stabilized on an external resonator, available in a small butterfly package (RIO), and a laboratory high quality fiber laser (Koheras). Both lasers are used at a relatively low output power: 7 mW (60 mA) for the RIO laser and 16 mW for the Koheras laser. The power dependence of the spectral purity is relatively weak for this type of lasers, and such an output power is largely sufficient to drive our bench (however, a small improvement of about 2 to 3 dB is observed on the RIO laser phase noise when biased at 120 mA instead of 60 mA).

Figure 4 depicts the result of the measurements performed on these two lasers, under these conditions. The measurement is performed in a Faraday's shielded room using the E5052 system and the frequency noise is computed using equation 1. In Figure 5, the same measurements are plotted in terms of single sideband phase noise, which shows the difference in dB between these two lasers, and also give an idea of the laser spectral shape, taking into account that when the phase noise is sufficiently low compared to 0 dB, it represents the lateral noise wings of the spectrum (we do not want to compute the linewidth, because this parameter is too much dependent on the integration time and has, in our opinion, no rigourous meaning for this type of high spectral purity optical sources). In a preceding paper [START_REF] Llopis | Phase noise measurement of a narrow linewidth CW laser using delay line approaches[END_REF], the phase noise of the Koheras laser had been a few dB overestimated in the 10 Hz-100 Hz range, because of a vibration noise sensitivity of the fiber spool we were using. Indeed, like in any type of noise measurement, it is essential to determine if the measured data are not too close from the measurement bench noise floor. The next paragraph is dedicated to the evaluation and the optimisation of this noise floor.

III. MEASUREMENT NOISE FLOOR

An obvious technique to measure a frequency discriminator noise floor is to cancel the delay τ by removing the fiber spool. In this case, the system sensitivity to the laser frequency noise is null, and only remains the other noise contributions. However, such an evaluation of the noise floor is quite optimistic. Indeed, it also cancels any noise contribution coming from the fiber spool. There are mainly three types of low frequency noise that can be generated in a fiber spool: noise related to optical nonlinear effects [START_REF] Okusaga | Optical scattering induced noise in RF-photonic systems[END_REF], such as Stimulated Brillouin Scattering (SBS) ; diffusion noise on impurities, such as Rayleigh scattering noise ; noise due to a vibration sensitivity of the spool [START_REF] Li | Low vibration sensitivity fiber spools for laser stabilization[END_REF][START_REF] Hati | Cancellation of vibration induced phase noise in optical fibers[END_REF].

Therefore, in order to get a more precise evaluation of the noise floor, two identical fiber spools have been used to cancel the delay τ, one in each path of the interferometer. In this case, the above described effects can still be observed, and are included in the real noise floor.

The plastic rod spool used in our first experiments [START_REF] Llopis | Phase noise measurement of a narrow linewidth CW laser using delay line approaches[END_REF] has been replaced by a packaged spool (Fig. 3). The fiber is free of constraint in this package, but it is however mechanically maintained thanks to the relative small size of the package and to an inside plastic protection. Two set of spools have been used, one of 1 km (τ = 5.1 μs) and one of 2 km (τ = 10.05 μs).

Figure 6 depicts the results of the noise floor measurements in the 2 km fiber spool case, with the battery biased RIO laser feeding the bench (7 mW optical input power). A clear contribution of vibration noise is observed between 10 Hz and 2 kHz when the system is set up on a normal table (green curve). When the system is put on a low vibration table (but still inside the Faraday shielded room), most of the vibration contribution is reduced and the noise floor measured with the two spools technique (red curve) is close to the noise floor measured with the spool removed from the bench (two harms with short delay ; black curve). These observations allow us to clarify the importance of the vibration noise contribution to the noise floor, and to demonstrate the absence of non-linear optical noise contribution. Indeed, in both cases (Koheras or RIO lasers), we are feeding the fiber spool with optical power levels much below the SBS threshold, which is in the range of 30 mW for a 2 km fiber spool [START_REF]Fiber-Optic Communication Systems[END_REF].

The same experiments have been performed with the 1 km spools. On Figure 7, it is clear that the vibration problem is still observed with the 1 km spool. Finally, with this configuration, the noise floor is increased of about 6 dB, which is directly related to the reduction of the delay τ of a factor 2, and thus of the same reduction in sensitivity.

Finally, the 2 km spool appears to be a good compromise between a good sensitivity, a high SBS threshold and a large offset frequency bandwidth, and this fiber length is selected for further experiments.

The last question is: what is limiting the noise floor when the fiber spool is removed. At high frequency offsets, it is clear that the signal to noise ratio may have an influence, and this signal to noise ratio is mainly determined by the laser RIN at the IF frequency (80 MHz). However, close to the carrier, we also found an influence of the oscillator which is driving the acousto-optic modulator. This oscillator features a -43.5 dBc/Hz phase noise level at 1 Hz, with a 30 dB/dec slope in this region. It determines the optical frequency discriminator noise floor between 1 Hz and 5 Hz offset, as shown in Figure 8. However, between 5 Hz and 2 kHz, the noise floor seems still to be determined by vibration sensitivity. Finally, between 2 kHz and 100 kHz, the noise floor changes with the laser type and power, which means that it should be mostly related to the laser RIN. The next step is now to improve the experimental set up to overcome these limitations. The problem of the acousto-optic modulator driver is easy to solve: this source will be replaced by a high spectral purity oscillator (a high quality OCXO) followed by a high power amplifier, and this should lead to a 30 dB improvement of the 1 Hz phase noise (case of a vibration free experiment). However, the noise between 10 Hz and 2 kHz will be more difficult to improve: it is directly related to the interferometer sensitivity to vibrations. A new set up will be designed, featuring higher compactness, and on which mechanical isolation will be easier to provide.

IV. CONCLUSION

The measurement of two different high spectral purity lasers has been performed with a delay line frequency discriminator. The problem of the measurement noise floor has been carefully studied and it has been found that it is mainly limited by the mechanical vibration sensitivity of the fiber spool, and also of the whole interferometer set-up. The system sensitivity and dynamic range is however sufficient to investigate on the two examples of lasers presented in this paper. However, for ultra-high spectral purity lasers characterization, the set-up must still be improved.

Fig 1 :

 1 Fig 1 : Self heterodyne frequency noise measurement system

Fig 2 :Fig 3 :Fig 4 :

 234 Fig 2 : RF phase noise data at 80 MHz measured for a semiconductor laser stabilised on an external resonator with a 2 km discriminator

Fig 5 :

 5 Fig 5 : Phase noise in dBc/Hz computed for the two lasers of Fig. 4

Fig 6 :

 6 Fig 6 : Noise floor measured for the 2 km delay bench with the RIO laser, using three different configurations : no spool (black) ; two spools and normal table to support the bench (green) ; two spools and low vibration table (red).

Fig 7 :

 7 Fig 7 : Noise floor measured for the 1 km delay system with the RIO laser, using two different configurations : no spool (black) and two spools and normal table to support the experiment (green).

Fig 8 :

 8 Fig 8 : Noise floor measured with two spools of 2 km on a low vibration table for two different lasers, and contribution of the 80 MHz oscillator.