Modeling and Estimation of Multivariate Discrete and Continuous Time Stationary Processes

1. Introduction

Stationary stochastic processes provide a significant instrument for modeling numerous temporal phenomena related to different fields of science. In particular, due to the evidence of long dependence structures in the real financial data, stationary processes possessing long-memory have been widely applied in mathematical finance.

When discrete time is considered, stationary data is typically modeled by applying autoregressive–moving-average (ARMA) processes or their extensions. One focal reason for popularity of ARMA processes is that for every stationary process with a vanishing autocovariance γ(·) and for every n ∈ ℕ there exists an ARMA process X such that γ_X(k) = γ(k) for |k| ≤ n. For a comprehensive overview of ARMA processes we mention [1–3]. The immense ARMA family include for example SARIMA (seasonal autoregressive integrated moving-average) models, where a seasonal ARMA process is obtained by differencing the original data. We mention also different GARCH (generalized autoregressive conditional heteroskedasticity) models originating from [4] and [5] that are commonly used in financial modeling taking account of the time-dependent volatility. ARMA processes, their extensions and estimation in these models have been concerned e.g., in [6–16], to name but a few. Moreover, in [17] we showed that all univariate strictly stationary processes indexed by the integers are characterized by the AR(1) (autoregressive model of order 1) equation

$\begin{array}{l}{X}_t-\varphi X_{t-1}=Z_t,t\in \mathbb{Z}.\end{array}$

However, contrary to the classical AR(1), the noise Z belonging to a certain class of stationary processes is not necessarily white. Established on the characterization, we proposed an estimation method for ϕ in the case of a square integrable stationary process. This method has several advantages over conventional ones such as maximum likelihood and least squares fitting of ARMA models. Furthermore, in [18] we applied our method in estimation of a generalization of the ARCH model involving a covariate process that can be interpreted as the liquidity of an asset. We would like to emphasize that the proposed method is applicable in estimation of essentially any square integrable one-dimensional stationary process. Hence, it also covers stationary solutions, which are often of central interest, of models of the vast ARMA family.

In the case of continuous time, the Ornstein-Uhlenbeck process X given by the Langevin equation

$\begin{array}{ll}d{X}_t=-\theta X_{t}dt+d{B}_t,t\in \mathbb{Z},& (1)\end{array}$

where θ > 0 and B is a two-sided Brownian motion, can be seen as the analog of the discrete time AR(1) process. By posing a suitable initial condition, (1) yields a stationary solution. The foregoing can be generalized, for example, by replacing Brownian motion with other stationary increment processes satisfying certain integrability conditions. One popular option as the replacement is fractional Brownian motion recovering the fractional Ornstein-Uhlenbeck process introduced in [19]. This kind of generalized Ornstein-Uhlenbeck processes are applied e.g., in mathematical finance to describe mean-reverting systems under the influence of shocks, and they are a highly active topic of research. Equations of type (1) with varying driving forces, and estimation in such models have been concerned e.g., in [20–32], to mention but a few. Furthermore, in [33] we showed that a generalized version of (1) characterizes all multivariate strictly stationary processes with continuous paths. Consequently, we proposed an estimation method for the parameter matrix of (1). The method is based on continuous time algebraic matrix Riccati equations (CAREs) and it is applicable in estimation of essentially any square integrable multivariate stationary process.

Algebraic Riccati equations, occurring naturally e.g., in optimal control and filtering theory, is an intensively studied topic in the literature on its own. In many applications, real-valued CAREs often take the symmetric form

$\begin{array}{ll}{B}^{\top }A+AB-ACA+D=0,& (2)\end{array}$

where C and D are symmetric, and symmetric solutions A are to be found. For a general approach to algebraic Riccati equations the reader may consult for example [34]. The existence and uniqueness of a solution to (2) is a well-studied topic, especially when C and D are positive semidefinite [see e.g., [35–37]].

The main contribution of this paper is to extend the characterization and the consequent estimation method of discrete stationary processes of [17] to multivariate settings. Serving as an instrument for the characterization, we also define a multivariate discrete Lamperti transform giving a one-to-one correspondence between stationary and H-self-similar processes with a positive definite matrix H. The characterization leads to covariance based symmetric CAREs for the related parameter matrix providing us with a novel estimation method of multivariate discrete stationary processes. In outline, this results in the following correspondence between the continuous time case. When the concepts of noises are extended from the conventional ones, the reputed analogous AR(1) and Langevin equations characterize discrete and continuous time stationary processes, respectively. The characterizations provide us with models of stationary processes yielding symmetric CAREs for the corresponding parameter matrices. Furthermore, these equations can be applied in a similar manner in estimation in both cases. For readers' interest and to highlight the obtained connection between discrete and continuous time, we also give the key results of [33].

The rest of the paper is organized as follows. In section 2.1, we first give an AR(1) type of characterization covering all multivariate strictly stationary processes indexed by the integers. Consequently, under the assumption of square integrability, we obtain a set of symmetric CAREs for the model parameter matrix serving as a basis for estimation. Finally, we state theorems for consistency and asymptotic distribution of the parameter (matrix) estimator. In section 2.2, we present the main results of [33], while at the same time comparing them to the results obtained in discrete time. For the reader's convenience, all technical proofs are postponed to section 3.

2. Main Results

We begin with some preliminaries and a short notational introduction. The processes we consider in this paper are n-dimensional, real-valued and indexed by I ∈ {ℤ, ℝ}. For such a process Y we write Y = (Y_t)_t∈I, where the ith component of the random vector Y_t is denoted by $Y_t^{(i)}$. Equality of the distributions of two random vectors Y_t and Z_t is denoted by $Y_t\stackrel{\mathrm{\text{law}}}{=}{Z}_t$. Similarly, equality of two processes Y and Z in the sense of finite dimensional distributions is denoted by $Y={(Y_t)}_{t\in I}\stackrel{\mathrm{\text{law}}}{=}{(Z_t)}_{t\in I}=Z$. Throughout the paper, we investigate strictly stationary processes meaning that ${(X_{t+s})}_{t\in I}\stackrel{\mathrm{\text{law}}}{=}{(X_t)}_{t\in I}$ for every s ∈ I. Consequently, we omit the word “strictly” and simply say that X is stationary. By writing A ≥ 0 or A > 0 we mean that the matrix A is positive semidefinite or positive definite, respectively. We denote an eigendecomposition of a symmetric matrix by A = QΛQ^⊤, where Λ = diag(λ_i). Furthermore, the L² vector norm and the corresponding induced matrix norm is denoted by ‖·‖.

In sections 2.1 and 2.2, we introduce models for discrete and continuous time stationary processes, respectively. Consequently, this leads to symmetric CAREs for the model parameter matrix of the form

$\begin{array}{ll}{B}^{\top }A+AB-ACA+D=0,& (3)\end{array}$

where C, D ≥ 0, and we are solving the equation for a positive definite A. There exists a vast amount of literature on existence and uniqueness of a solution [see e.g., [35] or [36]] in the described setting. In particular, if C, D > 0, then there exists a unique positive semidefinite solution to (3). Furthermore, there exists several numerical methods for finding the positive semidefinite solution of (3) [see e.g., [38, 39] or monograph [40]]. Hence, provided that the solution to (3) is unique, a prospective estimation method of A can be based on the equation.

2.1. Discrete Time

In this subsection, we extend the characterization of discrete (I = ℤ) stationary processes of [17] to multivariate settings. Consequently, we derive quadratic equations for the corresponding model parameter matrix providing us with a natural way to define an estimator for the parameter. Finally, we pose theorems for consistency and asymptotic distribution. A strong analog with the continuous time case I = ℝ covered in [33] is obtained. We start by providing some definitions.

DEFINITION 2.1. Let G = (G_t)_t∈ℤ be an n-dimensional stationary increment process. We define a stationary process ΔG = (Δ_tG)_t∈ℤ by

$\begin{array}{l}{\Delta }_{t}G=G_t-G_{t-1}.\end{array}$

Next, we define a class of stationary increment processes having sub-exponentially deviating sample paths. These processes serve as the noise in the subsequent AR(1) type of characterization of stationary processes.

DEFINITION 2.2. Let H > 0 be a positive definite n × n matrix, and let G = (G_t)_t∈ℤ be an n-dimensional stochastic process with stationary increments and G₀ = 0. If

$\begin{array}{ll}{\displaystyle \underset{l\to -\infty }{\text{lim}}}{\displaystyle \sum _{k=l}^0}{e}^{kH}{\Delta }_{k}G& (4)\end{array}$

exists in probability and defines an almost surely finite random variable, we denote $G\in {{G}}_H$.

REMARK 2.3. Lemma 3.1 shows that existence of a logarithmic moment is sufficient for $G\in {{G}}_H$ for all H > 0. Particularly, this is the case if G is square integrable. On the other hand, an example of an one-dimensional stationary increment process G with G₀ = 0, but $G\notin {{G}}_H$ for any H > 0 was provided in [41]. Moreover, Definition 2.2 could also be stated without the assumption of positive definiteness. However, in the one-dimensional case with H ≤ 0, the convergence of (4) would imply G ≡ 0. As expected, we encounter a similar kind of dimensional degeneracy when considering e.g., symmetric matrices with non-positive eigenvalues. Hence, the assumption H > 0 may be regarded as natural. See also Remark 3.3.

The next theorem characterizes all multivariate stationary processes, including processes possessing long-memory.

THEOREM 2.4. Let H > 0 be a positive definite n × n matrix, and let X = (X_t)_t∈ℤ be an n-dimensional stochastic process. Then X is stationary if and only if ${\text{lim}}_{t\to -\infty }{e}^{tH}{X}_t\stackrel{\mathbb{P}}{=}0$ and

$\begin{array}{ll}{\Delta }_{t}X=(e^{-H}-I)X_{t-1}+{\Delta }_{t}G& (5)\end{array}$

for $G\in {{G}}_H$ and t ∈ ℤ. Moreover, the process $G\in {{G}}_H$ is unique.

COROLLARY 2.5. Let H > 0 be a positive definite n × n matrix, and let X be stationary. Then X admits an AR(1) type of representation

$\begin{array}{ll}{X}_t-\Phi X_{t-1}={\Delta }_{t}G,& (6)\end{array}$

where Φ = e^−H and $G\in {{G}}_H$.

By using (6) and the expression (17) from the proof of Theorem 2.4, it is straightforward to show that ΔG is centered and square integrable if and only if X is centered and square integrable, respectively. In what follows, we assume these two attributes and write $\gamma (t)=𝔼(X_t{X}_0^{\top })$ and $r(t)=𝔼\left[({\Delta }_{t}G){({\Delta }_{0}G)}^{\top }\right]$. Furthermore, since

$\begin{array}{l}{G}_t={\displaystyle \sum _{k=1}^t}{\Delta }_{k}G,t\ge 1,\end{array}$

in this case also G is centered and square integrable, and we denote $v(t)=\text{cov}(G_t)=𝔼(G_t{G}_t^{\top })$. We would like to point out that centeredness can be assumed without loss of generality (see Remark 2.11).

Under the discussed assumptions, we obtain an expression for γ(t) in terms of the noise process.

REMARK 2.6. The autocovariance function γ(t) is given by

$\begin{array}{l}\gamma (t)=e^{-tH}{\displaystyle \sum _{k=-\infty }^t}{\displaystyle \sum _{j=-\infty }^0}{e}^{kH}r(k-j)e^{jH}.\end{array}$

Furthermore, if G has independent components, we obtain

$\begin{array}{l}\gamma (t)=\frac{e^{-tH}}{2}{\displaystyle \sum _{k=-\infty }^t{\displaystyle \sum _{j=-\infty }^0{e}^{kH}}}(v(k-j+1)+v(k-j-1)\\ -2v(k-j))e^{jH}.\end{array}$

The following lemma writes the quadratic equations for the model parameter Φ = e^−H presented in [17] in our multivariate setting.

LEMMA 2.7. Let H > 0 be a positive definite n × n matrix, and let X be stationary of the form (6). Then

$\begin{array}{ll}r(t)=\Phi \gamma (t)\Phi -\gamma (t+1)\Phi -\Phi \gamma (t-1)+\gamma (t)& (7)\end{array}$

for every t ∈ ℤ.

REMARK 2.8. In the proof of the lemma, we utilize the increment process ΔG similarly as in [17] yielding quadratic equations for Φ in terms of r and γ. One could consider an estimation method for Φ based on the above equations. However, for a general stationary X, (7) is a symmetric CARE only if t = 0, and even in this case, existence of a unique positive semidefinite solution is not guaranteed.

By adapting the approach of [33], we obtain a set of symmetric CAREs on which we construct an estimator for the model parameter Φ = e^−H. For this, we define the following matrix coefficients.

DEFINITION 2.9. We set

$\begin{array}{l}{B}_t={\displaystyle \sum _{k=1}^t}\gamma (k-1)-\gamma {(k)}^{\top }\\ C_t={\displaystyle \sum _{k=1}^t}{\displaystyle \sum _{j=1}^t}\gamma (k-j)\\ D_t=v(t)-2\gamma (0)+\gamma (t)+\gamma {(t)}^{\top }\end{array}$

for every t ∈ ℕ.

THEOREM 2.10. Let H > 0 be a positive definite n × n matrix, and set Θ = I − e^−H. Let X = (X_t)_t∈ℤ be stationary of the form (5). Then the CARE

$\begin{array}{ll}{B}_t^{\top }\Theta +\Theta B_t-\Theta C_t\Theta +D_t=0& (8)\end{array}$

is satisfied for every t ∈ ℕ.

REMARK 2.11. Equations (7) and (8) are covariance based. Consequently, they hold also when X and G in Theorem 2.4 are not centered.

REMARK 2.12. It is worth to emphasize the interplay between the parameter matrix H and the noise G in the context of a fixed stationary X. The noise in Theorem 2.4 and Corollary 2.5 is unique only after H > 0 is fixed. On the other hand, if we fix some covariances of the noise process such that the equations of Theorem 2.10 (or Lemma 2.7) yield a positive definite solution, then we may set H according to this solution. Furthermore, the corresponding noise is now given by Corollary 2.5 and its covariance satisfies the set prerequisites. In practice, the parameter matrix H is estimated by assuming some information on the noise, and by estimating the autocovariance of the observed stationary process X. The estimation will discussed in more detail at the end of this subsection. See also the following examples illustrating how to obtain familiar noise terms in the case of two stationary ARMA type of processes.

We give a couple of examples on how some basic multivariate processes of ARMA type can be presented in the form (6), and how to derive the corresponding noise G together with its covariance function v.

EXAMPLE 2.13. Let X be an n-dimensional stationary AR(1) type of process given by

$\begin{array}{l}{X}_t-\varphi X_{t-1}={ϵ}_t,\end{array}$

with 0 < ϕ = QΛQ^⊤, ‖ϕ‖ < 1 and ϵ ~ iid(0, Σ). Then, we may set $H=-Q\text{diag}(log{\lambda }_i)Q^{\top }$ giving Φ = ϕ. Now ΔG = ϵ and $G_t=\sum _{k=1}^t{ϵ}_k$. Furthermore, $v(t)=\sum _{k=1}^t\text{cov}({ϵ}_k)=t\Sigma$ for t ≥ 1.

EXAMPLE 2.14. Let X be an n-dimensional stationary ARMA(1, q) type of process given by

$\begin{array}{l}{X}_t-\varphi X_{t-1}={ϵ}_t+{\theta }_1{ϵ}_{t-1}+\dots +{\theta }_q{ϵ}_{t-q},\end{array}$

with 0 < ϕ, ‖ϕ‖ < 1 and ϵ ~ iid(0, Σ). Similarly as above, we may set Φ = ϕ and now ΔG equals to the MA(q) process on the right. Consequently, for t ≥ 1,

$\begin{array}{l}{G}_t={\displaystyle \sum _{k=1}^t}{ϵ}_k+{\theta }_1{ϵ}_{k-1}+\dots +{\theta }_q{ϵ}_{k-q}\end{array}$

and

$\begin{array}{l}v(t)={\displaystyle \sum _{i,j=0}^q}\text{max}(0,t-|i-j|){\theta }_i\Sigma {\theta }_j^{\top },\end{array}$

where θ₀ = I.

In [17], we proposed an estimation method of one-dimensional stationary processes based on equations (7). In particular, we showed that the method is applicable except in some special class of stationary processes. In [42], we provided a comprehensive analysis of the class, and proved that it consists of highly degenerate processes. On the other hand, due to the strong dependence structure, the failure of different estimation methods is expected. Fundamentally, a stationary process X belongs to the class if there exists two values H and $\tilde{H}$ such that the corresponding processes ΔG and $\Delta \tilde{G}$ in (6) have identical autocovariance functions. Next, we state a lemma showing that these degenerate processes have a special characteristic also under the new set of Equations (8).

LEMMA 2.15. Let X be a one-dimensional stationary process and let H > 0 be fixed. Set Φ = e^−H and

$\begin{array}{l}{X}_t-\Phi X_{t-1}={\Delta }_{t}G,t\in \mathbb{Z}.\end{array}$

If the equation

$\begin{array}{ll}\gamma (t){\Phi }^2-(\gamma (t+1)-\gamma (t-1))\Phi +\gamma (t)-r(t)=0& (9)\end{array}$

yields the same two solutions $\Phi ,\tilde{\Phi }>0$ for every t ∈ ℤ, then also the equation

$\begin{array}{ll}{C}_t{\Theta }^2-2B_t\Theta -D_t=0& (10)\end{array}$

yields the same two solutions Θ = 1 − Φ and $\tilde{\Theta }=1-\tilde{\Phi }$ for every t ∈ ℕ.

The following remark illuminates the connection between the coefficient matrices of Definition 2.9 and uniqueness of the solution to (8), which we discussed in the beginning of section 2.

REMARK 2.16. By Lemma 3.8, the matrix Θ is positive definite. Since

$\begin{array}{l}{C}_t=𝔼\left[{\displaystyle \sum _{k=1}^t}{X}_{k-1}{\left({\displaystyle \sum _{k=1}^t}{X}_{k-1}\right)}^{\top }\right]=\text{cov}\left({\displaystyle \sum _{k=1}^t}{X}_{k-1}\right),\end{array}$

the matrix C_t is positive semidefinite. Furthermore, if the smallest eigenvalue of v(t) grows enough in time, D_t becomes positive definite [see [33]]. This is the case e.g., when the noise G has independent components with growing variances.

For estimation, it is desirable that (8) admits a unique positive semidefinite solution. In this case, the solution is the correct parameter matrix Θ by the construction. Moreover, it guarantees that a convergent numerical scheme has the desired limit Θ. In the sequal, we simply assume that t is chosen in such a way that C_t, D_t > 0 ensuring the uniqueness of the solution, and we omit the subindex t from Equation (8). We have justified the assumption of positive definiteness in detail in continuous time [see section 2.1 and Remark 2.11 in [33]]. Furthermore, we assume that v(t) is known and the stationary process X is observed up to the time T > t, and the coefficient matrices B, C, and D are estimated from these observations by replacing the autocovariances γ(·) with some estimators ${\widehat{\gamma }}_T(·)$. The coefficient estimators are denoted by ${\widehat{B}}_T,{\widehat{C}}_T$, and ${\widehat{D}}_T$, and we set

$\begin{array}{l}{\Delta }_{T}B={\widehat{B}}_T-B,{\Delta }_{T}C={\widehat{C}}_T-C,{\Delta }_{T}D={\widehat{D}}_T-D.\end{array}$

Next, we define an estimator ${\widehat{\Theta }}_T$ for the matrix Θ = I − Φ = I − e^−H. The proofs of the related asymptotic results allow a certain amount of flexibility in the definition. Thus, we give a definition that probably is the most convenient from the practical point of view. Consistency and the rate of convergence of ${\widehat{\Theta }}_T$ are inherited from autocovariance estimators ${\widehat{\gamma }}_T(·)$ of the observed stationary process. In addition, the limiting distribution is obtained as a linear function of the limiting distribution of the autocovariance estimators.

DEFINITION 2.17. The estimator ${\widehat{\Theta }}_T$ is defined as the unique positive semidefinite solution to the perturbed CARE

$\begin{array}{l}{\widehat{B}}_T^{\top }{\widehat{\Theta }}_T+{\widehat{\Theta }}_T{\widehat{B}}_T-{\widehat{\Theta }}_T{\widehat{C}}_T{\widehat{\Theta }}_T+{\widehat{D}}_T=0\end{array}$

whenever ${\widehat{C}}_T,{\widehat{D}}_T>0$. Otherwise, we set ${\widehat{\Theta }}_T=0$.

THEOREM 2.18. Let C, D > 0. Assume that

$\begin{array}{l}{\displaystyle \underset{s\in \left\{0,1,\dots ,t\right\}}{\text{max}}}\Vert {\widehat{\gamma }}_T(s)-\gamma (s)\Vert \stackrel{\mathbb{P}}{\to }0.\end{array}$

Then

$\begin{array}{l}\Vert {\widehat{\Theta }}_T-\Theta \Vert \stackrel{\mathbb{P}}{\to }0,\end{array}$

where ${\widehat{\Theta }}_T$ is given by Definition 2.17.

THEOREM 2.19. Let l(T) be a rate function. If

$\begin{array}{l}l(T)\left[\begin{array}{c}\text{vec}({\widehat{\gamma }}_T(0)-\gamma (0))\\ \text{vec}({\widehat{\gamma }}_T(1)-\gamma (1))\\ ⋮\\ \text{vec}({\widehat{\gamma }}_T(t)-\gamma (t))\end{array}\right]\stackrel{law}{\to }Z,\end{array}$

where Z is a (t + 1)n²-dimensional random vector, then:

(1) Let $\tilde{Z}$ be the permutation of elements of Z corresponding to the order of elements of

$\begin{array}{l}\left[\begin{array}{c}\text{vec}\left({({\widehat{\gamma }}_T(0)-\gamma (0))}^{\top }\right)\\ \text{vec}\left({({\widehat{\gamma }}_T(1)-\gamma (1))}^{\top }\right)\\ ⋮\\ \text{vec}\left({({\widehat{\gamma }}_T(t)-\gamma (t))}^{\top }\right).\end{array}\right]\end{array}$

Define a linear mapping $L_1:{ℝ}^{(t+1)n^2}\to {ℝ}^{3n^2}$ by

$\begin{array}{l}{L}_1(Z)=\left[\begin{array}{c}{\displaystyle {\sum }_{k=0}^{t-1}}(t-k)\left[\begin{array}{c}{Z}^{(k{n}^2+1)}\\ ⋮\\ Z^{((k+1)n^2)}\end{array}\right]+{\displaystyle {\sum }_{k=1}^{t-1}}(t-k)\left[\begin{array}{c}{\tilde{Z}}^{(k{n}^2+1)}\\ ⋮\\ {\tilde{Z}}^{((k+1)n^2)}\end{array}\right]\\ {\displaystyle {\sum }_{k=1}^t}\left[\begin{array}{c}{Z}^{((k-1)n^2+1)}-{\tilde{Z}}^{(k{n}^2+1)}\\ ⋮\\ Z^{(k{n}^2)}-{\tilde{Z}}^{((k+1)n^2)}\end{array}\right]\\ \left[\begin{array}{c}{Z}^{(t{n}^2+1)}+{\tilde{Z}}^{(t{n}^2+1)}\\ ⋮\\ Z^{((t+1)n^2)}+{\tilde{Z}}^{((t+1)n^2)}\end{array}\right]-2\left[\begin{array}{c}{Z}^{(1)}\\ ⋮\\ Z^{(n^2)}\end{array}\right]\end{array}\right],\end{array}$

where $\sum _1^0$ is an empty sum. Then

$\begin{array}{l}l(T)\text{vec}({\Delta }_{T}C,{\Delta }_{T}B,{\Delta }_{T}D)\stackrel{law}{\to }{L}_1(Z).\end{array}$

(2) If D, C > 0 and ${\widehat{\Theta }}_T$ is given by Definition 2.17, then

$\begin{array}{l}l(T)\text{vec}({\widehat{\Theta }}_T-\Theta )\stackrel{law}{\to }{L}_2(L_1(Z)),\end{array}$

where $L_2:{ℝ}^{3n^2}\to {ℝ}^{n^2}$ is a linear mapping expressible in terms of Θ, t and r.

2.2. Continuous Time

We have collected the main results (Theorems 2.21, 2.28, 2.32 and 2.33) of [33] concerning continuous time multivariate stationary processes into this subsection. In addition, in order to complete the analog between discrete and continuous time, we derive quadratic equations for the model parameter (Proposition 2.24) that correspond to the equations of Lemma 2.7. Throughout the subsection, we assume that the considered processes have continuous paths almost surely and hence, the related stochastic integrals can be interpreted as pathwise Riemann-Stieltjes integrals. Again, we start by defining the class ${{G}}_H$ of stationary increment processes for H > 0.

DEFINITION 2.20. Let H > 0 be a positive definite n × n matrix, and let G = (G_t)_t∈ℝ be an n-dimensional stochastic process with stationary increments and G₀ = 0. If

$\begin{array}{l}{\displaystyle \underset{s\to -\infty }{\text{lim}}}{\displaystyle {\int }_s^0}{e}^{Hu}d{G}_u\end{array}$

exists in probability and defines an almost surely finite random variable, we denote $G\in {{G}}_H$.

As in discrete time, it can be shown that existence of some logarithmic moments ensure that $G\in {{G}}_H$ for all H > 0. In particular, square integrability of G suffices, which is the case in our second moment based estimation method.

The next theorem is the continuous time counterpart of Theorem 2.4 showing that all stationary processes are characterized by the Langevin equation, whereas in discrete time, the characterization was given by an AR(1) type of equation.

THEOREM 2.21. Let H > 0 be a positive definite n × n matrix, and let X = (X_t)_t∈ℝ be an n-dimensional stochastic process. Then X is stationary if and only if

$\begin{array}{l}{X}_0={\displaystyle {\int }_{-\infty }^0}{e}^{Hu}d{G}_u\end{array}$

and

$\begin{array}{ll}d{X}_t=-H{X}_{t}dt+d{G}_t,& (11)\end{array}$

for $G\in {{G}}_H$ and t ∈ ℝ. Moreover, the process $G\in {{G}}_H$ is unique.

COROLLARY 2.22. From Theorem 2.21 it follows that X is the unique stationary solution

$\begin{array}{ll}{X}_t=e^{-Ht}{\displaystyle {\int }_{-\infty }^t}{e}^{Hu}d{G}_u& (12)\end{array}$

to (11).

In order to apply Theorem 2.21 in estimation, we pose the assumption

$\begin{array}{l}{\displaystyle \underset{s\in \left[0,1\right]}{\text{sup}}}𝔼(\Vert G_s{\Vert }^2)<\infty .\end{array}$

This guarantees that $G\in {{G}}_H$ for all H > 0, and square integrability of X and G. On the other hand, if X is square integrable, then G is also. In addition and without loss of generality, we assume that the processes are centered. Again, we write $\gamma (t)=𝔼(X_t{X}_0^{\top })$ and $v(t)=𝔼(G_t{G}_t^{\top })$. Now, the autocovariance function of the following stationary process is well-defined.

DEFINITION 2.23. Let G = (G_t)_t∈ℝ be a centered square integrable stationary increment process and let δ > 0. We define a stationary process ${\Delta }^{\delta }G={({\Delta }_t^{\delta }G)}_{t\in ℝ}$ by

$\begin{array}{l}{\Delta }_t^{\delta }G=G_t-G_{t-\delta }\end{array}$

and the corresponding autocovariance function r_δ by

$\begin{array}{l}{r}_{\delta }(t)=𝔼\left[({\Delta }_t^{\delta }G){({\Delta }_0^{\delta }G)}^{\top }\right].\end{array}$

As in discrete time (Lemma 2.7), we obtain quadratic equations for the model parameter H in terms of r_δ and γ. The equations could potentially be used to construct an estimator for H and they might be also of independent interest.

PROPOSITION 2.24. Let H > 0 be a positive definite n × n matrix, and let X be of the form (12). Then

$\begin{array}{l}{r}_{\delta }(t)=2\gamma (t)-\gamma (t+\delta )-\gamma (t-\delta )+\left({\displaystyle {\int }_t^{t+\delta }\gamma }(s)ds-{\displaystyle {\int }_{t-\delta }^t\gamma }(s)ds\right)H\\ +H\left({\displaystyle {\int }_{t-\delta }^t\gamma }(s)ds-{\displaystyle {\int }_t^{t+\delta }\gamma }(s)ds\right)+H({\displaystyle {\int }_{t-\delta }^t(s-t+\delta )}\gamma (s)ds\\ +{\displaystyle {\int }_t^{t+\delta }(t-s+\delta )}\gamma (s)ds)H\end{array}$

for every t ∈ ℝ.

REMARK 2.25. The advantage of the equations above is that we have to consider γ(s) only for s ∈ [t − δ, t + δ], but as in the discrete case, for a general stationary X we obtain a symmetric CARE only when t = 0. In addition, similarly as above, we could set in discrete time ${\Delta }_t^{k}G:\; =G_t-G_{t-k},k\in \mathbb{N}$. However, this would lead to more complicated equations in Lemma 2.7.

A significant difference compared to the discrete time Equations (7) occurs in the univariate case. Namely, the first order term with respect to H vanishes.

COROLLARY 2.26. The univariate case yields

$\begin{array}{l}{r}_{\delta }(t)=2\gamma (t)-\gamma (t+\delta )-\gamma (t-\delta )+H^2({\displaystyle {\int }_{t-\delta }^t(s-t+\delta )}\gamma (s)ds\\ +{\displaystyle {\int }_t^{t+\delta }(t-s+\delta )}\gamma (s)ds)\end{array}$

for every t ∈ ℝ.

One could base a univariate estimation method on the above equations without the concern of existence of a unique positive solution. However, since we wish to treat also multivariate settings, we present the most central results of [33] that are obtained from Theorem 2.21 by considering the noise G directly. First, we define matrix coefficients corresponding to Definition 2.9. Consequently, we write symmetric CAREs for the parameter H that are similar to the CAREs (8) for the discrete time parameter Θ.

DEFINITION 2.27. We set

$\begin{array}{l}{B}_t={\displaystyle {\int }_0^t}\gamma (s)-\gamma {(s)}^{\top }ds\\ C_t={\displaystyle {\int }_0^t}{\displaystyle {\int }_0^t}\gamma (s-u)duds\\ D_t=v(t)-2\gamma (0)+\gamma (t)+\gamma {(t)}^{\top }\end{array}$

for every t ≥ 0.

THEOREM 2.28. Let H > 0 be a positive definite n × n matrix, and let X = (X_t)_t∈ℝ be stationary of the form (12). Then the CARE

$\begin{array}{ll}{B}_t^{\top }H+H{B}_t-H{C}_{t}H+D_t=0& (13)\end{array}$

is satisfied for every t ≥ 0.

REMARK 2.29. As in discrete time, the Equations (13) and (14) are covariance based and hence, they hold also for a non-centered stationary X.

REMARK 2.30. Contrary to the discrete time Equations (8), the first order term vanishes in the univariate setting as in (13).

Again, we assume that t is chosen in such a way that C_t, D_t > 0 ensuring the existence of a unique positive semidefinite solution. We have discussed this assumption in detail in [33]. We define an estimator $\widehat{H}$_T for the model parameter matrix H identically to the discrete time by replacing the autocovariances γ(·) in the matrix coefficients with their estimators ${\widehat{\gamma }}_T(·)$. The below given definition differs slightly from the definition in [33], but the same asymptotic results still apply.

DEFINITION 2.31. The estimator $\widehat{H}$_T is defined as the unique positive semidefinite solution to the perturbed CARE

$\begin{array}{l}{\widehat{B}}_T^{\top }{\widehat{H}}_T+{\widehat{H}}_T{\widehat{B}}_T-{\widehat{H}}_T{\widehat{C}}_T{\widehat{H}}_T+{\widehat{D}}_T=0\end{array}$

whenever ${\widehat{C}}_T,{\widehat{D}}_T>0$. Otherwise, we set $\widehat{H}$_T = 0.

As in discrete time, asymptotic properties of $\widehat{H}$_T are inherited from the autocovariance estimators. However, due to the continuous time setting, instead of pointwise convergence, we have to consider functional form of convergence of ${\widehat{\gamma }}_T(·)$. In [33], we have provided sufficient conditions in the case of Gaussian noise G with independent components, under which the assumptions of the following theorems are satisfied. In particular, the results are valid for fractional Brownian motion that is widely applied in the field of mathematical finance.

THEOREM 2.32. Let C, D > 0. Assume that

$\begin{array}{l}{\displaystyle \underset{s\in \left[0,t\right]}{\text{sup}}}\Vert {\widehat{\gamma }}_T(s)-\gamma (s)\Vert \stackrel{\mathbb{P}}{\to }0.\end{array}$

Then

$\begin{array}{l}\Vert \widehat{H_T}-H\Vert \stackrel{\mathbb{P}}{\to }0,\end{array}$

where $\widehat{H}$_T is given by Definition 2.31.

THEOREM 2.33. Let Y = (Y_s)_{s∈[0, t]} be an n²-dimensional stochastic process with continuous paths almost surely and let l(T) be a rate function. If

$\begin{array}{l}l(T)\text{vec}({\widehat{\gamma }}_T(s)-\gamma (s))\stackrel{law}{\to }{Y}_s\end{array}$

in the uniform topology of continuous functions, then:

(1) Let $\tilde{Y}$_s be the permutation of elements of Y_s that corresponds to the order of elements of $\text{vec}\left({(\widehat{\gamma }(s)-\gamma (s))}^{\top }\right)$. Then

$\begin{array}{l}l(T)\text{vec}({\Delta }_{T}C,{\Delta }_{T}B,{\Delta }_{T}D)\stackrel{law}{\to }\left[\begin{array}{c}{\int }_0^t(t-s)(Y_s+{\tilde{Y}}_s)ds\\ {\int }_0^t\left(Y_s-{\tilde{Y}}_s\right)ds\\ Y_t+{\tilde{Y}}_t-2Y_0\end{array}\right]\\ =:L_1(Y).\end{array}$

(2) If C, D > 0 and $\widehat{H}$_T is given by Definition 2.31, then

$\begin{array}{l}l(T)\text{vec}({\widehat{H}}_T-H)\stackrel{law}{\to }{L}_2(L_1(Y)),\end{array}$

where $L_2:{ℝ}^{3n^2}\to {ℝ}^{n^2}$ is a linear mapping expressible in terms of H, t and the covariance function of G.

3. Proofs

In the following, we denote the smallest eigenvalue of H > 0 by λ_min. Consequently $\Vert e^{kH}\Vert =\Vert Q\text{diag}(e^{{\lambda }_{i}k})Q^{\top }\Vert =e^{{\lambda }_{min}k}$ for a negative k.

3.1. Discrete Time

The proof of the next lemma follows the lines of the proof of Theorem 2.2. in [41] that concerns the one-dimensional continuous time case. However, in our setting, we obtain a weaker sufficient condition for $G\in {{G}}_H$ for all H > 0.

LEMMA 3.1. Let G = (G_t)_t∈ℤ be an n-dimensional stationary increment process with G₀ = 0. Assume that

$\begin{array}{l}𝔼\left(|\text{log}\Vert G_1\Vert {𝟙}_{\left\{\Vert G_1\Vert >1\right\}}{|}^{1+\delta }\right)<\infty \end{array}$

for some δ > 0. Then $G\in {{G}}_H$ for all positive definite n × n matrices H.

PROOF. Let H > 0. We apply the Borel-Cantelli lemma together with Markov's inequality to show that $\Vert e^{kH}{\Delta }_{k}G\Vert \to 0$ almost surely as k → −∞. Let ϵ > 0 be fixed below.

$\begin{array}{l}\mathbb{P}\left(\Vert e^{kH}{\Delta }_{k}G\Vert >ϵ\right)\le \mathbb{P}\left(e^{{\lambda }_{min}k}\Vert {\Delta }_{k}G\Vert >ϵ\right)\\ =\mathbb{P}\left(e^{{\lambda }_{min}k}\Vert G_1\Vert >ϵ\right)=\mathbb{P}\left(\Vert G_1\Vert >\frac{ϵ}{e^{{\lambda }_{min}k}}\right)\\ =\mathbb{P}\left(\text{log}\Vert G_1\Vert >\text{log }ϵ-{\lambda }_{min}k\right),\end{array}$

since ΔG is stationary and G₀ = 0. Furthermore, $k(\frac{\text{log }ϵ}{k}-{\lambda }_{min})\ge -Ck$ for some C > 0 and k ≤ k_ϵ. Thus, for k ≤ k_ϵ,

$\begin{array}{l}\mathbb{P}\left(\Vert e^{kH}{\Delta }_{k}G\Vert >ϵ\right)\le \mathbb{P}\left(\text{log}\Vert G_1\Vert \ge -Ck\right)\\ =\mathbb{P}\left(\text{log}\Vert G_1\Vert {𝟙}_{\left\{\Vert G_1\Vert >1\right\}}\ge -Ck\right)\\ \le \frac{𝔼|\text{log}\Vert G_1\Vert {𝟙}_{\left\{\Vert G_1\Vert >1\right\}}{|}^{1+\delta }}{{(-Ck)}^{1+\delta }}\le c\frac{1}{{(-k)}^{1+\delta }}\end{array}$

giving the wanted result. We conclude the proof by noting that

$\begin{array}{l}{\displaystyle \sum _{k=-\infty }^0}\Vert e^{kH}{\Delta }_{k}G\Vert \le {\displaystyle \sum _{k=-\infty }^0}\Vert e^{\frac{1}{2}kH}\Vert \Vert e^{\frac{1}{2}kH}{\Delta }_{k}G\Vert \le \\ {\displaystyle \underset{k}{\text{sup}}}\Vert e^{\frac{1}{2}kH}{\Delta }_{k}G\Vert {\displaystyle \sum _{k=-\infty }^0}{e}^{\frac{1}{2}{\lambda }_{min}k}<\infty \end{array}$

almost surely.

Next, we extend the concept of self-similarity to discrete time multivariate processes. While the following two definitions are natural, to the best of our knowledge they are not widely acknowledged in the literature.

DEFINITION 3.2. Let H > 0 be a positive definite n × n matrix, and let $Y={(Y_{e^t})}_{t\in \mathbb{Z}}$ be an n-dimensional stochastic process. Then Y is H-self-similar if

$\begin{array}{l}{(Y_{e^{t+s}})}_{t\in \mathbb{Z}}\stackrel{law}{=}{(e^{sH}{Y}_{e^t})}_{t\in \mathbb{Z}}\end{array}$

for every s ∈ ℤ.

REMARK 3.3. Similarly as in Definition 2.2, positive definite matrices serve as natural counterparts of the conventional positive scalar valued exponents of self-similarity. In addition, ${\text{lim}}_{l\to -\infty }{Y}_{e^l}$ does not converge (see the proof of the auxiliary Lemma 3.7) in the case of a general H-self-similar Y with non-positive definite H. It seems intuitive to expect the term to be convergent. However, posing such requirement leads to degeneracy of Y e.g., when H is symmetric with negative eigenvalues [cf. Remark 3.2 in [33]]. Furthermore, this indicates a possibility of dimensional reduction. For details on continuous time self-similar processes, we refer to [43]. See also the discussion on non-positive exponents of self-similarity in [44].

The following transform and the corresponding theorem giving one-to-one correspondence between self-similar and stationary processes were originally introduced by Lamperti in the univariate continuous time setting [45].

DEFINITION 3.4. Let H > 0 be a positive definite n × n matrix, and let X = (X_t)_t∈ℤ and $Y={(Y_{e^t})}_{t\in \mathbb{Z}}$ be n-dimensional stochastic processes. We define

$\begin{array}{l}{({{L}}_{H}X)}_{e^t}=e^{tH}{X}_t\end{array}$

and

$\begin{array}{l}{({{L}}_H^{-1}Y)}_t=e^{-tH}{Y}_{e^t}.\end{array}$

THEOREM 3.5. The operator ${{L}}_H$ together with its inverse ${{L}}_H^{-1}$ define a bijection between n-dimensional stationary processes and n-dimensional H-self-similar processes.

PROOF. First, let X be stationary and set $Z_{e^t}={({{L}}_{H}X)}_{e^t}$. Then

$\begin{array}{l}\left[\begin{array}{c}{Z}_{e^{t_1+s}}\\ Z_{e^{t_2+s}}\\ ⋮\\ Z_{e^{t_m+s}}\end{array}\right]=\left[\begin{array}{c}{e}^{(t_1+s)H}{X}_{t_1+s}\\ e^{(t_2+s)H}{X}_{t_2+s}\\ ⋮\\ e^{(t_m+s)H}{X}_{t_m+s}\end{array}\right]\stackrel{\text{law}}{=}\left[\begin{array}{c}{e}^{(t_1+s)H}{X}_{t_1}\\ e^{(t_2+s)H}{X}_{t_2}\\ ⋮\\ e^{(t_m+s)H}{X}_{t_m}\end{array}\right]=\left[\begin{array}{c}{e}^{sH}{Z}_{e^{t_1}}\\ e^{sH}{Z}_{e^{t_1}}\\ ⋮\\ e^{sH}{Z}_{e^{t_m}}\end{array}\right]\end{array}$

for every m ∈ ℕ, t ∈ ℤ^m and s ∈ ℤ. Hence, Z is H-self-similar.

Now, let Y be H-self-similar and set $Z_t={({{L}}_H^{-1}Y)}_t$. Then

$\begin{array}{l}\left[\begin{array}{c}{Z}_{t_1+s}\\ Z_{t_2+s}\\ ⋮\\ Z_{t_m+s}\end{array}\right]=\left[\begin{array}{c}{e}^{-(t_1+s)H}{Y}_{e^{t_1+s}}\\ e^{-(t_2+s)H}{Y}_{e^{t_2+s}}\\ ⋮\\ e^{-(t_m+s)H}{Y}_{e^{t_m+s}}\end{array}\right]\stackrel{\text{law}}{=}\left[\begin{array}{c}{e}^{-t_{1}H}{Y}_{e^{t_1}}\\ e^{-t_{2}H}{Y}_{e^{t_2}}\\ ⋮\\ e^{-t_{m}H}{Y}_{e^{t_m}}\end{array}\right]=\left[\begin{array}{c}{Z}_{t_1}\\ Z_{t_2}\\ ⋮\\ Z_{t_m}\end{array}\right]\end{array}$

for every m ∈ ℕ, t ∈ ℤ^m and s ∈ ℤ. Hence, Z is stationary completing the proof.

REMARK 3.6. The form of the conventional continuous time Lamperti transform is not directly applicable in the discrete time setting due to the scaling of time. Hence, in order to stay within the given discrete parameter sets, we use exponential clocks with self-similar processes in Definitions 3.2 and 3.4.

Before the proof of Theorem 2.4 we state an auxiliary lemma.

LEMMA 3.7. Let H > 0 be a positive definite n × n matrix, and let ${(Y_{e^t})}_{t\in \mathbb{Z}}$ be an n-dimensional H-self-similar process. We define a process G = (G_t)_t∈ℤ by

$\begin{array}{l}{G}_t=\{\begin{array}{cc}{\displaystyle {\sum }_{k=1}^t}{e}^{-kH}{\Delta }_k{Y}_{e^k},\hfill & t\ge 1\hfill \\ 0,\hfill & t=0\hfill \\ -{\displaystyle {\sum }_{k=t+1}^0}{e}^{-kH}{\Delta }_k{Y}_{e^k},\hfill & t\le -1.\hfill \end{array}\hfill \end{array}$

Then $G\in {{G}}_H$.

PROOF. It is straightforward to verify that ${\Delta }_{t}G=e^{-tH}{\Delta }_t{Y}_{e^t}$ for every t ∈ ℤ. In addition

$\begin{array}{l}{\displaystyle \underset{l\to -\infty }{\text{lim}}}{\displaystyle \sum _{k=l}^0}{e}^{kH}{\Delta }_{k}G={\displaystyle \underset{l\to -\infty }{\text{lim}}}{\displaystyle \sum _{k=l}^0}{\Delta }_k{Y}_{e^k}=Y_1-{\displaystyle \underset{l\to -\infty }{\text{lim}}}{Y}_{e^{l-1}},\end{array}$

where by self-similarity of Y

$\begin{array}{l}\mathbb{P}(\Vert Y_{e^{l-1}}\Vert \ge ϵ)=\mathbb{P}(\Vert e^{H(l-1)}{Y}_1\Vert \ge ϵ)\le \mathbb{P}(e^{{\lambda }_{min}(l-1)}\Vert Y_1\Vert \ge ϵ)\\ \to 0.\end{array}$

Hence, we set

$\begin{array}{l}{\displaystyle \underset{l\to -\infty }{\text{lim}}}{\displaystyle \sum _{k=l}^0}{e}^{kH}{\Delta }_{k}G=Y_1.\end{array}$

PROOF OF THEOREM 2.4. Assume that ${\text{lim}}_{t\to -\infty }{e}^{tH}{X}_t\stackrel{\mathbb{P}}{=}0$ and (5) holds for $G\in {{G}}_H$. Then, by using (5) repeatedly

$\begin{array}{l}{X}_t=e^{-H}{X}_{t-1}+{\Delta }_{t}G=e^{-(n+1)H}{X}_{t-n-1}+{\displaystyle \sum _{j=0}^n}{e}^{-jH}{\Delta }_{t-j}G\\ =e^{-(n+1)H}{X}_{t-n-1}+e^{-tH}{\displaystyle \sum _{k=t-n}^t}{e}^{kH}{\Delta }_{k}G=e^{-tH}\\ \left(e^{(t-n-1)H}{X}_{t-n-1}+{\displaystyle \sum _{k=t-n}^t}{e}^{kH}{\Delta }_{k}G\right)\end{array}$

for every n ∈ ℕ. Since, as n → ∞, the limit of the sum above is well-defined, and ${\text{lim}}_{n\to \infty }{e}^{(t-n-1)H}{X}_{t-n-1}\stackrel{\mathbb{P}}{=}0$, we obtain that

$\begin{array}{ll}{X}_t=e^{-tH}{\displaystyle \sum _{k=-\infty }^t}{e}^{kH}{\Delta }_{k}G.& (14)\end{array}$

Let m ∈ ℕ, t ∈ ℤ^m, and s ∈ ℤ. Then, by stationary increments of G, we have

$\begin{array}{l}\left[\begin{array}{c}{e}^{-(t_1+s)H}{\displaystyle {\sum }_{j=-M}^{t_1}}{e}^{(j+s)H}{\Delta }_{j+s}G\\ ⋮\\ e^{-(t_m+s)H}{\displaystyle {\sum }_{j=-M}^{t_m}}{e}^{(j+s)H}{\Delta }_{j+s}G\end{array}\right]\stackrel{\text{law}}{=}\left[\begin{array}{c}{e}^{-t_{1}H}{\displaystyle {\sum }_{j=-M}^{t_1}}{e}^{jH}{\Delta }_{j}G\\ ⋮\\ e^{-t_{m}H}{\displaystyle {\sum }_{j=-M}^{t_m}}{e}^{jH}{\Delta }_{j}G\end{array}\right]\end{array}$

for every −M < min{t_i}. Since the random vectors above converge in probability as M → ∞, we obtain that

$\begin{array}{l}\left[\begin{array}{c}{X}_{t_1+s}\\ ⋮\\ X_{t_m+s}\end{array}\right]=\left[\begin{array}{c}{e}^{-(t_1+s)H}{\displaystyle {\sum }_{j=-\infty }^{t_1}}{e}^{(j+s)H}{\Delta }_{j+s}G\\ ⋮\\ e^{-(t_m+s)H}{\displaystyle {\sum }_{j=-\infty }^{t_m}}{e}^{(j+s)H}{\Delta }_{j+s}G\end{array}\right]\\ \stackrel{\text{law}}{=}\left[\begin{array}{c}{e}^{-t_{1}H}{\displaystyle {\sum }_{j=-\infty }^{t_1}}{e}^{jH}{\Delta }_{j}G\\ ⋮\\ e^{-t_{m}H}{\displaystyle {\sum }_{j=-\infty }^{t_m}}{e}^{jH}{\Delta }_{j}G\end{array}\right]=\left[\begin{array}{c}{X}_{t_1}\\ ⋮\\ X_{t_m}\end{array}\right]\end{array}$

and hence, X is stationary.

Next, assume that X is stationary. Then, by Theorem 3.5 there exists a H-self-similar Y such that

$\begin{array}{l}{\Delta }_{t}X=e^{-tH}{Y}_{e^t}-e^{-(t-1)H}{Y}_{e^{t-1}}=(e^{-H}-I)X_{t-1}+e^{-tH}{\Delta }_t{Y}_{e^t}.\end{array}$

Defining G as in Lemma 3.7 completes the proof of the other direction.

To prove uniqueness, we use (15). Assume that, for $G,\tilde{G}\in {{G}}_H$,

$\begin{array}{l}{e}^{tH}{X}_t={\displaystyle \sum _{k=-\infty }^t}{e}^{kH}{\Delta }_{k}G={\displaystyle \sum _{k=-\infty }^t}{e}^{kH}{\Delta }_k\tilde{G}\end{array}$

for every t ∈ ℤ. Then

$\begin{array}{l}{e}^{tH}{X}_t-e^{(t-1)H}{X}_{t-1}=e^{tH}{\Delta }_{t}G=e^{tH}{\Delta }_t\tilde{G}.\end{array}$

Since e^tH is invertible and both processes start from zero, we conclude that $G=\tilde{G}$.

PROOF OF LEMMA 2.7. We have that

$\begin{array}{l}{\Delta }_{t}G{({\Delta }_{0}G)}^{\top }=(X_t-\Phi X_{t-1})(X_0^{\top }-X_{-1}^{\top }\Phi ).\end{array}$

Taking expectations yields

$\begin{array}{l}r(t)=\Phi \gamma (t)\Phi -\gamma (t+1)\Phi -\Phi \gamma (t-1)+\gamma (t).\end{array}$

PROOF OF THEOREM 2.10. Let

$\begin{array}{l}{\Delta }_{t}X=-\Theta X_{t-1}+{\Delta }_{t}G.\end{array}$

Then for t ∈ ℕ we have

$\begin{array}{l}{G}_t={\displaystyle \sum _{k=1}^t}{\Delta }_{k}G={\displaystyle \sum _{k=1}^t}{\Delta }_{k}X+\Theta {\displaystyle \sum _{k=1}^t}{X}_{k-1}=X_t-X_0+\Theta {\displaystyle \sum _{k=1}^t}{X}_{k-1}.\end{array}$

Hence

$\begin{array}{l}\text{ cov}(G_t)=\text{cov}(X_t-X_0)+𝔼\left[(X_t-X_0){\displaystyle \sum _{k=1}^t}{X}_{k-1}^{\top }\right]\Theta \\ +\Theta 𝔼\left[{\displaystyle \sum _{k=1}^t}{X}_{k-1}{(X_t-X_0)}^{\top }\right]+\Theta 𝔼\left[{\displaystyle \sum _{k=1}^t}{X}_{k-1}{\left({\displaystyle \sum _{k=1}^t}{X}_{k-1}\right)}^{\top }\right]\Theta \end{array}$

giving (8) since

$\begin{array}{l}𝔼\left[(X_t-X_0){\displaystyle \sum _{k=1}^t}{X}_{k-1}^{\top }\right]={\displaystyle \sum _{k=1}^t}\gamma (t-k+1)-\gamma (-k+1)\\ ={\displaystyle \sum _{k=1}^t}\gamma (k)-\gamma {(k-1)}^{\top }\end{array}$

and

$\begin{array}{l}\text{cov}(X_t-X_0)=𝔼\left[(X_t-X_0){(X_t-X_0)}^{\top }\right]=2\gamma (0)-\gamma (t)\\ -\gamma (-t)=2\gamma (0)-\gamma (t)-\gamma {(t)}^{\top }\end{array}$

PROOF OF LEMMA 2.15. Assume that $\tilde{\Phi }=e^{-\tilde{H}}$ satisfies (9) for every t ∈ ℤ and set

$\begin{array}{l}{X}_t-\tilde{\Phi }{X}_{t-1}={\Delta }_t\tilde{G},t\in \mathbb{Z}.\end{array}$

Consequently

$\begin{array}{l}\tilde{r}(t)=\gamma (t){\tilde{\Phi }}^2-(\gamma (t+1)-\gamma (t-1))\tilde{\Phi }+\gamma (t)=r(t),t\in \mathbb{Z},\end{array}$

where $\tilde{r}(t)$ is the autocovariance function of ${({\Delta }_t\tilde{G})}_{t\in \mathbb{Z}}$. Now, since $G_0={\tilde{G}}_0=0$, we obtain that

$\begin{array}{l}\text{var}(G_t)=\text{var}\left({\displaystyle \sum _{k=1}^t}{\Delta }_{k}G\right)={\displaystyle \sum _{k,j=1}^t}\text{cov}({\Delta }_{k}G,{\Delta }_{j}G)={\displaystyle \sum _{k,j=1}^t}r(k-j)\\ =\text{var}\left({\displaystyle \sum _{k=1}^t}{\Delta }_k\tilde{G}\right)=\text{var}({\tilde{G}}_t)\end{array}$

for all t ∈ ℕ. Hence, both Θ and $\tilde{\Theta }$ are solutions to (10).

LEMMA 3.8. The matrix Θ = I − e^−H is positive definite.

PROOF. Let a be a real vector of length n, and let H = QΛQ^⊤ be an eigendecomposition of H. Then

$\begin{array}{l}{a}^{\top }(I-e^{-H})a=\Vert a{\Vert }^2-a^{\top }{e}^{-H}a,\end{array}$

where

$\begin{array}{l}|a^{\top }{e}^{-H}a|\le \Vert a{\Vert }^2{e}^{-{\lambda }_{min}}<\Vert a{\Vert }^2\end{array}$

completing the proof.

In order to show that ${\widehat{\Theta }}_T$ is consistent, we simply need to find suitable bounds for Δ_TB, Δ_TC and Δ_TD in terms of the autocovariance estimators. After that, the same strategy as in [33] can be applied.

LEMMA 3.9. Set

$\begin{array}{l}{M}_{t,T}={\displaystyle \underset{s\in \left\{0,1,\dots ,t\right\}}{\text{max}}}\Vert {\widehat{\gamma }}_T(s)-\gamma (s)\Vert .\end{array}$

Then the coefficients of the perturbed CARE satisfy

$\begin{array}{l}\Vert {\Delta }_{T}D\Vert \le 4M_{t,T}\\ \Vert {\Delta }_{T}C\Vert \le t^2{M}_{t,T}\\ \Vert {\Delta }_{T}B\Vert \le 2t{M}_{t,T},\end{array}$

PROOF. First, we recall first that

$\begin{array}{l}\Vert {\widehat{\gamma }}_T(-s)-\gamma (-s)\Vert =\Vert {\widehat{\gamma }}_T{(s)}^{\top }-\gamma {(s)}^{\top }\Vert =\Vert {\widehat{\gamma }}_T(s)-\gamma (s)\Vert .\end{array}$

Now, since v(t) is known,

$\begin{array}{l}\Vert {\Delta }_{T}D\Vert \le 2\Vert {\widehat{\gamma }}_T(0)-\gamma (0)\Vert +\Vert {\widehat{\gamma }}_T(t)-\gamma (t)\Vert +\Vert {\widehat{\gamma }}_T{(t)}^{\top }\\ -\gamma {(t)}^{\top }\Vert \le 4M_{t,T}.\end{array}$

Moreover

$\begin{array}{l}\Vert {\Delta }_{T}C\Vert \le {\displaystyle \sum _{k=1}^t}{\displaystyle \sum _{j=1}^t}\Vert {\widehat{\gamma }}_T(k-j)-\gamma (k-j)\Vert \le t^2{M}_{t,T}.\end{array}$

Finally

$\begin{array}{l}\Vert {\Delta }_{T}B\Vert \le {\displaystyle \sum _{k=1}^t}\Vert {\widehat{\gamma }}_T(k-1)-\gamma (k-1)\Vert +\Vert \gamma {(k)}^{\top }-{\widehat{\gamma }}_T{(k)}^{\top }\Vert \\ \le 2t{M}_{t,T}.\end{array}$

PROOF OF THEOREM 2.18. The result follows by replacing ${\text{sup}}_{s\in \left[0,t\right]}\Vert {\widehat{\gamma }}_T(s)-\gamma (s)\Vert$ with M_t,T in Corollary 3.14 and in the proof of Theorem 2.9 of [33]. The details are left to the reader.

PROOF OF THEOREM 2.19. For the first part of the theorem, we notice that

$\begin{array}{l}{C}_t={\displaystyle \sum _{k=1}^t}{\displaystyle \sum _{j=1}^t}\gamma (k-j)={\displaystyle \sum _{k=1}^t}{\displaystyle \sum _{l=k-t}^{k-1}}\gamma (l)={\displaystyle \sum _{l=1-t}^{-1}}{\displaystyle \sum _{k=1}^{t+l}}\gamma (l)\\ +{\displaystyle \sum _{l=0}^{t-1}}{\displaystyle \sum _{k=l+1}^t}\gamma (l)={\displaystyle \sum _{1-t}^{-1}}(t+l)\gamma (l)+{\displaystyle \sum _{l=0}^{t-1}}(t-l)\gamma (l)\\ ={\displaystyle \sum _{l=0}^{t-1}}(t-l)\gamma (l)+{\displaystyle \sum _{l=1}^{t-1}}(t-l)\gamma {(l)}^{\top },\end{array}$

where $\sum _0^{-1}$ and $\sum _1^0$ are interpreted as empty sums. Now we have that

$\begin{array}{l}{\Delta }_{T}C={\displaystyle \sum _{k=0}^{t-1}}(t-k)({\widehat{\gamma }}_T(k)-\gamma (k))\\ +{\displaystyle \sum _{k=1}^{t-1}}(t-k)\left({\widehat{\gamma }}_T{(k)}^{\top }-\gamma {(k)}^{\top }\right)\\ {\Delta }_{T}B={\displaystyle \sum _{k=1}^t}{\widehat{\gamma }}_T(k-1)-\gamma (k-1)-{\widehat{\gamma }}_T{(k)}^{\top }+\gamma {(k)}^{\top }\\ {\Delta }_{T}D=2(\gamma (0)-{\widehat{\gamma }}_T(0))+{\widehat{\gamma }}_T(t)-\gamma (t)+{\widehat{\gamma }}_T{(t)}^{\top }-\gamma {(t)}^{\top }\end{array}$

and furthermore

$\begin{array}{l}l(T)\text{vec}({\Delta }_{T}C,{\Delta }_{T}B,{\Delta }_{T}D)=l(T)\\ \left[\begin{array}{c}{\displaystyle {\sum }_{k=0}^{t-1}(t-k)}\text{vec}(\widehat{\gamma }(k)-\gamma (k))+{\displaystyle {\sum }_{k=1}^{t-1}(t-k)}\text{vec}\left({(\widehat{\gamma }(k)-\gamma (k))}^{\top }\right)\\ {\displaystyle {\sum }_{k=1}^t\text{vec}}(\widehat{\gamma }(k-1)-\gamma (k-1))-\text{vec}\left({(\widehat{\gamma }(k)-\gamma (k))}^{\top }\right)\\ -2\text{vec}(\widehat{\gamma }(0)-\gamma (0))+\text{vec}(\widehat{\gamma }(t)-\gamma (t))+\text{vec}\left({(\widehat{\gamma }(t)-\gamma (t))}^{\top }\right)\end{array}\right]\\ =L_1\left[l(T)\left[\begin{array}{c}\text{vec}({\widehat{\gamma }}_T(0)-\gamma (0))\\ \text{vec}({\widehat{\gamma }}_T(1)-\gamma (1))\\ ⋮\\ \text{vec}({\widehat{\gamma }}_T(t)-\gamma (t))\end{array}\right]\right]\stackrel{\text{law}}{\to }{L}_1\text{(Z)}\end{array}$

by the continuous mapping theorem. For the second part of the theorem, the proof of the continuous time case of [33] can be applied just by replacing ${\text{sup}}_{s\in \left[0,t\right]}\Vert {\widehat{\gamma }}_T(s)-\gamma (s)\Vert$ with M_t,T in the definition of the set A_T.

3.2. Continuous Time

We only provide the proof of Proposition 2.24, while the other proofs can be found from [33]

PROOF OF PROPOSITION 2.24. Integrating (11) from 0 to t gives

$\begin{array}{l}{G}_t=X_t-X_0+H{\displaystyle {\int }_0^t}{X}_{s}ds.\end{array}$

Hence

$\begin{array}{l}{\Delta }_t^{\delta }G{({\Delta }_0^{\delta }G)}^{\top }=\left(X_t-X_{t-\delta }+H{\displaystyle {\int }_{t-\delta }^t}{X}_{s}ds\right)\\ \left(X_0^{\top }-X_{-\delta }^{\top }+{\displaystyle {\int }_{-\delta }^0}{X}_s^{\top }dsH\right).\end{array}$

Taking expectations yields

$\begin{array}{l}{r}_{\delta }(t)=2\gamma (t)-\gamma (t+\delta )-\gamma (t-\delta )+{\displaystyle {\int }_{-\delta }^0}\gamma (t-s)-\gamma (t-\delta -s)dsH\\ +H{\displaystyle {\int }_{t-\delta }^t}\gamma (s)-\gamma (s+\delta )ds+H{\displaystyle {\int }_{t-\delta }^t}{\displaystyle {\int }_{-\delta }^0}\gamma (s-u)dudsH,\end{array}$

where the first order terms can be treated with a simple change of variables. For the second order term we obtain that

$\begin{array}{l}{\displaystyle {\int }_{t-\delta }^t}{\displaystyle {\int }_{-\delta }^0}\gamma (s-u)duds={\displaystyle {\int }_{t-\delta }^t}{\displaystyle {\int }_s^{s+\delta }}\gamma (x)dxds\\ ={\displaystyle {\int }_{t-\delta }^t}{\displaystyle {\int }_{t-\delta }^x}\gamma (x)dsdx+{\displaystyle {\int }_t^{t+\delta }}{\displaystyle {\int }_{x-\delta }^t}\gamma (x)dsdx\\ ={\displaystyle {\int }_{t-\delta }^t}(x-t+\delta )\gamma (x)dx+{\displaystyle {\int }_t^{t+\delta }}(t-x+\delta )\gamma (x)dx.\end{array}$

Author Contributions

The author confirms being the sole contributor of this work and has approved it for publication.

Conflict of Interest

The author declares that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

I would like to thank Lauri Viitasaari for his comments and suggestions.

References

1. Brockwell PJ, Davis RA. Time Series: Theory and Methods. New York, NY: Springer Science & Business Media (1991). doi: 10.1007/978-1-4419-0320-4

CrossRef Full Text | Google Scholar

2. Hamilton JD. Time Series Analysis. Vol. 2. Princeton: Princeton University Press (1994).

3. Neusser K. Time Series Econometrics. Springer (2016). doi: 10.1007/978-3-319-32862-1

CrossRef Full Text | Google Scholar

4. Engle RF. Autoregressive conditional heteroskedasticity with estimates of the variance of United Kingdom inflation. Econometr J Econometr Soc. (1982) 50:987–1007. doi: 10.2307/1912773

CrossRef Full Text | Google Scholar

5. Bollerslev T. Generalized autoregressive conditional heteroskedasticity. J. Econometr. (1986) 31:307–27. doi: 10.1016/0304-4076(86)90063-1

CrossRef Full Text | Google Scholar

6. Hannan EJ. The asymptotic theory of linear time-series models. J Appl Probabil. (1973) 10:130–45. doi: 10.1017/S0021900200042145

CrossRef Full Text | Google Scholar

7. Tiao GC, Tsay RS. Consistency properties of least squares estimates of autoregressive parameters in ARMA models. Ann Stat. (1983) 11:856–71. doi: 10.1214/aos/1176346252

CrossRef Full Text | Google Scholar

8. Mikosch T, Gadrich T, Kluppelberg C, Adler RJ. Parameter estimation for ARMA models with infinite variance innovations. Ann Stat. (1995) 23:305–26. doi: 10.1214/aos/1176324469

CrossRef Full Text | Google Scholar

9. Mauricio JA. Exact maximum likelihood estimation of stationary vector ARMA models. J Am Stat Assoc. (1995) 90:282–91. doi: 10.1080/01621459.1995.10476511

CrossRef Full Text | Google Scholar

10. Ling S, McAleer M. Asymptotic theory for a vector ARMA-GARCH model. Econometr Theory. (2003) 19:280–310. doi: 10.1017/S0266466603192092

CrossRef Full Text | Google Scholar

11. Francq C, Zakoian JM. Maximum likelihood estimation of pure GARCH and ARMA-GARCH processes. Bernoulli. (2004) 10:605–37. doi: 10.3150/bj/1093265632

CrossRef Full Text | Google Scholar

12. Bollerslev T. Glossary to ARCH (GARCH). CREATES Res Pap. (2008) 49:1–44. doi: 10.2139/ssrn.1263250

CrossRef Full Text | Google Scholar

13. Han H, Kristensen D. Asymptotic theory for the QMLE in GARCH-X models with stationary and nonstationary covariates. J Bus Econ Stat. (2014) 32:416–29. doi: 10.1080/07350015.2014.897954

CrossRef Full Text | Google Scholar

14. Baillie RT, Bollerslev T, Mikkelsen HO. Fractionally integrated generalized autoregressive conditional heteroskedasticity. J Econometr. (1996) 74:3–10. doi: 10.1016/S0304-4076(95)01749-6

CrossRef Full Text | Google Scholar

15. Ling S, Li WK. On fractionally integrated autoregressive moving-average time series models with conditional heteroscedasticity. J Am Stat Assoc. (1997) 92:1184–94. doi: 10.1080/01621459.1997.10474076

CrossRef Full Text | Google Scholar

16. Zhang MY, Russell JR, Tsay RS. A nonlinear autoregressive conditional duration model with applications to financial transaction data. J Econometr. (2001) 104:179–207. doi: 10.1016/S0304-4076(01)00063-X

CrossRef Full Text | Google Scholar

17. Voutilainen M, Viitasaari L, Ilmonen P. On model fitting and estimation of strictly stationary processes. Modern Stochast Theory Appl. (2017) 4:381–406. doi: 10.15559/17-VMSTA91

CrossRef Full Text | Google Scholar

18. Voutilainen M, Viitasaari L, Ilmonen P, Torres S, Tudor C. On the ARCH model with stationary liquidity. (2020). Available online at: https://link.springer.com/content/pdf/10.1007/s00184-020-00779-x.pdf

Google Scholar

19. Cheridito P, Kawaguchi H, Maejima M. Fractional Ornstein-Uhlenbeck processes. Electron J Probabil. (2003) 8:1–44. doi: 10.1214/EJP.v8-125

CrossRef Full Text | Google Scholar

20. Hu Y, Nualart D. Parameter estimation for fractional Ornstein-Uhlenbeck processes. Stat Probabil Lett. (2010) 80:1030–8. doi: 10.1016/j.spl.2010.02.018

CrossRef Full Text | Google Scholar

21. Kleptsyna ML, Le Breton A. Statistical analysis of the fractional Ornstein-Uhlenbeck type process. Stat Inference Stochast Process. (2002) 5:229–48. doi: 10.1023/A:1021220818545

CrossRef Full Text | Google Scholar

22. Azmoodeh E, Viitasaari L. Parameter estimation based on discrete observations of fractional Ornstein-Uhlenbeck process of the second kind. Stat Inference Stochast Process. (2015) 18:205–27. doi: 10.1007/s11203-014-9111-8

CrossRef Full Text | Google Scholar

23. Bajja S, Es-Sebaiy K, Viitasaari L. Least squares estimator of fractional Ornstein-Uhlenbeck processes with periodic mean. J Korean Stat Soc. (2017) 46:608–22. doi: 10.1016/j.jkss.2017.06.002

CrossRef Full Text | Google Scholar

24. Balde MF, Es-Sebaiy K, Tudor CA. Ergodicity and drift parameter estimation for infinite-dimensional fractional Ornstein-Uhlenbeck process of the second kind. Appl Math Opt. (2020) 81:785–814. doi: 10.1007/s00245-018-9519-4

CrossRef Full Text | Google Scholar

25. Brouste A, Iacus SM. Parameter estimation for the discretely observed fractional Ornstein-Uhlenbeck process and the Yuima R package. Comput Stat. (2013) 28:1529–47. doi: 10.1007/s00180-012-0365-6

CrossRef Full Text | Google Scholar

26. Douissi S, Es-Sebaiy K, Tudor CA. Hermite Ornstein-Uhlenbeck processes mixed with a Gamma distribution. Publicationes Mathematicae Debrecen. (2020) 96:1–22. doi: 10.5486/PMD.2020.8443

CrossRef Full Text | Google Scholar

27. Hu Y, Nualart D, Zhou H. Parameter estimation for fractional Ornstein-Uhlenbeck processes of general Hurst parameter. Stat Inference Stochast Process. (2019) 22:111–42. doi: 10.1007/s11203-017-9168-2

CrossRef Full Text | Google Scholar

28. Nourdin I, Tran TD. Statistical inference for Vasicek-type model driven by Hermite processes. Stochast Process Appl. (2019) 129:3774–91. doi: 10.1016/j.spa.2018.10.005

CrossRef Full Text | Google Scholar

29. Sottinen T, Viitasaari L. Parameter estimation for the Langevin equation with stationary-increment Gaussian noise. Stat Inference Stochast Process. (2018) 21:569–601. doi: 10.1007/s11203-017-9156-6

CrossRef Full Text | Google Scholar

30. Tanaka K. Maximum likelihood estimation for the non-ergodic fractional Ornstein-Uhlenbeck process. Stat Inference Stochast Process. (2015) 18:315–32. doi: 10.1007/s11203-014-9110-9

CrossRef Full Text | Google Scholar

31. Applebaum D, et al. Infinite dimensional Ornstein-Uhlenbeck processes driven by Lévy processes. Probabil Surv. (2015) 12:33–54. doi: 10.1214/14-PS249

CrossRef Full Text | Google Scholar

32. Magdziarz M. Fractional Ornstein-Uhlenbeck processes. Joseph effect in models with infinite variance. Phys A Stat Mech Appl. (2008) 387:123–33. doi: 10.1016/j.physa.2007.08.016

CrossRef Full Text | Google Scholar

33. Voutilainen M, Viitasaari L, Ilmonen P, Torres S, Tudor C. Vector-valued generalised Ornstein-Uhlenbeck processes. arXiv[Preprint].arXiv:1909.02376. (2019).

Google Scholar

34. Lancaster P, Rodman L. Algebraic Riccati Equations. New York, NY: Oxford University Press (1995).

Google Scholar

35. Kucera V. A contribution to matrix quadratic equations. IEEE Trans Automat Control. (1972) 17:344–7. doi: 10.1109/TAC.1972.1099983

CrossRef Full Text | Google Scholar

36. Wonham WM. On a matrix Riccati equation of stochastic control. SIAM J Control. (1968) 6:681–97. doi: 10.1137/0306044

CrossRef Full Text | Google Scholar

37. Sun Jg. Perturbation theory for algebraic Riccati equations. SIAM J Matrix Anal Appl. (1998) 19:39–65. doi: 10.1137/S0895479895291303

CrossRef Full Text | Google Scholar

38. Byers R. Solving the algebraic Riccati equation with the matrix sign function. Linear Algebra Appl. (1987) 85:267–79. doi: 10.1016/0024-3795(87)90222-9

CrossRef Full Text | Google Scholar

39. Laub A. A Schur method for solving algebraic Riccati equations. IEEE Trans Automat Control. (1979) 24:913–21. doi: 10.1109/TAC.1979.1102178

CrossRef Full Text | Google Scholar

40. Bini DA, Iannazzo B, Meini B. Numerical Solution of Algebraic Riccati Equations. Vol. 9. Philadelphia: SIAM (2012). doi: 10.1137/1.9781611972092

CrossRef Full Text | Google Scholar

41. Viitasaari L. Representation of stationary and stationary increment processes via Langevin equation and self-similar processes. Stat Probabil Lett. (2016) 115:45–53. doi: 10.1016/j.spl.2016.03.020

CrossRef Full Text | Google Scholar

42. Voutilainen M, Viitasaari L, Ilmonen P. Note on AR(1)-characterisation of stationary processes and model fitting. Modern Stochast Theory Appl. (2019) 6:195–207. doi: 10.15559/19-VMSTA132

CrossRef Full Text | Google Scholar

43. Embrechts P, Maejima M. Selfsimilar Processes. Princeton: Princeton University Press (2002). doi: 10.1515/9781400825103

CrossRef Full Text | Google Scholar

44. Pipiras V, Taqqu MS. Long-Range Dependence and Self-Similarity. Vol. 45. Cambridge: Cambridge University Press (2017). doi: 10.1017/CBO9781139600347

CrossRef Full Text | Google Scholar

45. Lamperti J. Semi-stable stochastic processes. Trans Am Math Soc. (1962) 104:62–78. doi: 10.1090/S0002-9947-1962-0138128-7

CrossRef Full Text | Google Scholar