Quantile estimation in fractional Levy Ornstein-Uhlenbeck processes

Abstract

First we study estimation of the drift parameter in the fractional Ornstein-Uhlenbeck process whose marginal distribution is Student $t$ -distribution. We obtain Spearman’s correlation based estimator, quantile estimator and Brownian excursion based estimator of the drift parameter. Then we study method of moments estimator and quantile estimator in fractional inverse Gaussian and fractional gamma Ornstein-Uhlenbeck processes.

Keywords

Fractional Brownian motion fractional student process fractional Levy Ornstein-Uhlenbeck process heavy tail long-memory Spearman’s correlation quantile estimator Brownian excursion fractional inverse Gaussian process fractional gamma process

1. Introduction and preliminaries

Levy processes are processes with stationary independent increments. Levy Ornstein-Uhlenbeck (LOU) process generalizes the Ornstein-Uhlenbeck process to include jumps, see Jacod and Shiryayev (1987). The Levy Ornstein-Uhlenbeck (LOU) process, is an extension of Ornstein-Uhlenbeck process with Levy process driving term. Levy driven processes of Ornstein-Uhlenbeck type have been extensively studied over the last decade and widely used in finance, see Barndorff-Neilsen and Shephard (2001). In finance, it is useful as a generalization of Vasicek model, as one-factor short-term interest rate model which could take into account the jump of the interest rate. It also generalizes stochastic volatility model where the volatility has jumps.

Jump processes are of two types: Finite activity processes and infinite activity processes. Finite activity processes have finite number of jumps in a finite time interval, e.g., a Poisson process and infinite activity processes have infinite number of jumps in a finite time interval, e.g., gamma process, inverse Gaussian process and tempered stable process.

It is well known that the suitably parametrized autoregresive (AR) process with Gaussian error has the continuous limit the Vasicek model. Wolfe (1982) studied continuous analogue of the stochastic difference equation of AR type with Levy type innovations whose limit is a Levy driven OU Process. Gourieroux and Jasiak (2006) studied autoregressive gamma (ARG) process and showed that its continuous time limit is the Cox-Ingersoll-Ross (CIR) model. Thus the stationary ARG process is a discretized version of the CIR process. Gourieroux and Jasiak (2006) studied pseudo-maximum likelihood estimation in autoregressive gamma (ARG) process. This process can also be used for application in series of squared returns and intertrade durations for high-frequency data, i.e., it is a stochastic duration model. ARG model also fits a series of volumes per trade, which is an alternative proxy for liquidity. This is different from gamma autoregressive process (GAR) process studied in Sim (1990) and Gaver and Lewis (1980) where just the noise of the linear autoregressive process is Gamma distributed. For intertrade durations, the most popular model is autoregressive conditional duration (ACD) model introduced by Engle and Russell (1998).

For the exponential AR (1) model, the ratio estimator was introduced in Davis and McCormick (1989). This estimator is motivated by the extreme value theory for the correlation parameter of an AR (1) process whose innovation distribution is positive. In the case of exponential AR (1) process, it coincides with the maximum likelihood estimator. See Neilsen and Shephard (2003).

The weak consistency of the ratio estimator in the LOU process was studied in Jongbloed et al. (2006). The strong consistency and asymptotic Weibullness was studied in Brockwell, Davis and Yang (2007) in the case of Gamma innovations.

Recently long memory processes and processes with jumps have received attention in finance, engineering and physics. The simplest continuous time long memory process is the fractional Brownian motion discovered by Kolmogorov (1940) and later on studied by Levy (1948) and Mandelbrot and van Ness (1968). Hence fractional Brownian motion can be named as the Kolmogorov process.

A normalized fractional Brownian motion $\{W_{t}^{H},t\geqslant 0\}$ with Hurst parameter $H\in(0,1)$ is a centered Gaussian process with continuous sample paths whose covariance kernel is given by

$\displaystyle E(W_{t}^{H}W_{s}^{H})=\frac{1}{2}(s^{2H}+t^{2H}-|t-s|^{2H}),s,t% \geqslant 0.$

The process is self similar (scale invariant) and it can be represented as a stochastic integral with respect to standard Brownian motion. For $H=\frac{1}{2}$ , the process is a standard Brownian motion. For $H\neq\frac{1}{2}$ , the fBm is not a semimartingale and not a Markov process, but a Dirichlet process. The increments of the fBm are negatively correlated for $H<\frac{1}{2}$ and positively correlated for $H<\frac{1}{2}$ and in this case they display long-range dependence. The parameter $H$ which is also called the self similarity parameter, measures the intensity of the long range dependence. The ARIMA (p, d, q) with autoregressive part of order $p$ , moving average part of order $q$ and fractional difference parameter $d\in(0,0.5)$ process converge in Donsker sense to fBm.

Symmetric Fractional Student Diffusion

The fractional student diffusion is given by

$\displaystyle dX_{t}=-\alpha X_{t}dt+\sqrt{\frac{2\alpha}{\nu-1}(\nu+X_{t}^{2}% )}dW_{t}^{H}.$

The invariant density is given by Student t distribution with $\nu$ degrees of freedom:

$\displaystyle f(x)=\frac{\Gamma((\nu+1)/2)}{\sqrt{\nu\pi}\Gamma(\nu/2)}\left(1% +\frac{x^{2}}{\nu}\right)^{-(\nu+1)/2},x\in{\mathbb{R}},\nu>0=\frac{1}{\sqrt{% \nu}{\cal B}(\nu/2,1/2))[1+x^{2}/\nu]^{(\nu+1)/2}]}$

where

$\displaystyle{\cal B}(\nu/2,1/2))=\frac{\Gamma(1/2)\Gamma(\nu/2)}{\Gamma((\nu+% 1)/2)}$

is the beta function.

As $\nu\to\infty$ , this SDE becomes the classical OU process

$\displaystyle dX_{t}=-\alpha X_{t}dt+\sqrt{2\alpha}dW_{t}.$

Skewed Fractional Student Diffusion

$\displaystyle dX_{t}=-\alpha X_{t}dt+\sqrt{\frac{2\alpha}{\nu-1}[X_{t}^{2}+2% \gamma\nu^{1/2}X_{t}+(1+\gamma^{2})\nu]}dW_{t}^{H}.$

Here $\gamma$ is the skewness parameter. $\gamma=0$ gives the fractional symmetric Student diffusion.

Next we focus on non-Gaussian Levy processes. A fractional Poisson process $\{N_{H}(t),t>0\}$ with Hurst parameter $H\in(1/2,1)$ is defined as

$\displaystyle N_{H}(t)=\frac{1}{\Gamma(H-\frac{1}{2})}\int_{0}^{t}u^{\frac{1}{% 2}-H}\left(\int_{u}^{t}\tau^{H-\frac{1}{2}}(\tau-u)^{H-\frac{3}{2}}d\tau\right% )dq(u)$

where $q(u)=\frac{N(u)}{\sqrt{\lambda}}-\sqrt{\lambda}u$ and $N(u)$ is a homogeneous Poisson process with intensity $\lambda>0$ . See Araya et al. (2019). The covariance of $N^{H}$ is given by

$\displaystyle E(N_{H}(t)N_{H}(s))=\frac{V_{H}^{2}}{2}(t^{2H}+s^{2H}-|t-s|^{2H})$

where $V_{H}^{2}:=-\frac{\Gamma(2-2H)\cos(\pi H)}{(2H-1)\pi H}$ . The process is self-similar in the wide sense, has wide sense stationary increments, has fat-tailed non-Gaussian distribution, and exhibits long range dependence. The process converges to fractional Brownian motion in distribution. The process is self-similar in the asymptotic sense.

Strict sense, wide sense and asymptotic sense self-similarity are equivalent for fractional Brownian motion. Stock returns and interest rates are far from being self-similar in strict sense.

Hawkes processes are an efficient generalization of the Poisson processes to model a sequence of arrivals over time of some types of events, that present self-exciting feature, in the sense that each arrival increases the rate of future arrivals for some period of time. This class of counting processes allows one to capture self-exciting phenomena in a more accurate way compared to inhomogeneous Poisson processes or Cox processes. In finance, they are accurate to model for example credit risk contagion, order book or microstructure noises’s feature of financial markets.

A Hawkes process is a counting process $A_{t}$ with stochastic intensity $\lambda_{t}$ given by

$\displaystyle\lambda_{t}=\mu+\int_{0}^{t}\Phi(t-s)dA_{s}$

where $\mu>0$ and $\Phi:{\mathbb{R}}\to{\mathbb{R}}_{+}$ are two parameters. The parameter $\mu>0$ is called the background intensity and the function $\Phi$ is called the excitation function. When $\Phi=0$ , this a homogeneous Poisson process.

A fractional Hawkes process $\{A_{H}(t),t>0\}$ with Hurst parameter $H\in(1/2,1)$ is defined as

$\displaystyle A_{H}(t)=\frac{1}{\Gamma(H-\frac{1}{2})}\int_{0}^{t}u^{\frac{1}{% 2}-H}\left(\int_{u}^{t}\tau^{H-\frac{1}{2}}(\tau-u)^{H-\frac{3}{2}}d\tau\right% )dR(u)$

where $R(u)=\frac{A(u)}{\sqrt{\lambda}_{t}}-\sqrt{\lambda}_{t}u$ and $A(u)$ is a Hawkes process with stochastic intensity $\lambda_{t}$ .

As a generalization of these processes, we have the following process. Fractional Levy Process (FLP) is defined as

$\displaystyle M_{H,t}=\frac{1}{\Gamma(H+\frac{1}{2})}\int_{\mathbb{R}}[(t-s)^{% H-1/2}_{+}-(-s)^{H-1/2}_{+}]dL_{s},t\in{\mathbb{R}}$

where $\{L_{t},t\in{\mathbb{R}}\}$ is a Levy process on ${\mathbb{R}}$ with $E(L_{1})=0$ , $E(L_{1}^{2})<\infty$ .

Here are some properties of the fractional Levy process:

1) the covariance of the process is given by

$\displaystyle\text{cov}(M_{H,t},M_{H,s})=\frac{E(L_{1}^{2})}{2\Gamma(2H+1)\sin% (\pi H)}[|t|^{2H}+|s|^{2H}-|t-s|^{2H}].$

2) $M_{H}$ is not a martingale. For a large class of Levy processes, $M_{H}$ is neither a semimartingale. 3) $M_{H}$ is Hölder continuous of any order $\beta$ less than $H-\frac{1}{2}$ . 4) $M_{H}$ has stationary increments. 5) $M_{H}$ is symmetric. 6) $L$ is self-similar, but $M_{H}$ is not self-similar. 7) $M_{H}$ has infinite total variation on compacts.

Thus FLP is a generalization and a natural counterpart of FBM. Fractional stable motion is a special case of FLP. First we discuss estimation in partially observed models and then we discuss estimation in directly observed model in finite dimensional set up. In finance, the log-volatility process can be modeled as a fractionally integrated moving average (FIMA) process which is defined as

$\displaystyle Y_{H}(t)=\int_{-\infty}^{t}g_{H}(t-u)dM_{u},t\in{\mathbb{R}},g_{% H}(t)=\frac{1}{\Gamma(H-\frac{1}{2})}\int_{0}^{t}g(t-s)s^{H-\frac{3}{2}}ds,t% \in{\mathbb{R}}$

which is the Riemann-Liouville fractional integral of order $H$ and the kernel $g$ is the kernel of a short memory moving average process. The log-volatility process will have slow (hyperbolic rate) decay of the auto-correlation function (acf).

The process $Y_{H}(t)$ can be written as

$\displaystyle Y_{H}(t)=\int_{-\infty}^{t}g(t-u)dM_{H,u},t\in{\mathbb{R}}.$

We assume the following conditions on the kernel $g:{\mathbb{R}}\rightarrow{\mathbb{R}}$ , namely 1) $g(t)=0$ for all $t<0$ (causality), 2) $|g(t)|\leqslant Ce^{-ct}$ for some constants $C>0$ and $c>0$ (short memory).

The FIMA process is stationary and is infinite divisible. It has long memory and jumps which agree empirically with stochastic volatility models. The asset return can be modeled as a COGARCH process

$\displaystyle dX(t)=\sqrt{e^{Y_{H}(t)}}dL_{t}$

where $(L_{t},t\in{\mathbb{R}}$ is another Levy process and the initial value $Y_{H}(0)$ is independent of $L$ .

Consider the kernel

$\displaystyle g(t-s)=\sigma e^{-\theta(t-s)}I_{(0,\infty)}(t-s),\theta>0$

then

$\displaystyle g_{H}(t)=\frac{\sigma}{\Gamma(H-\frac{1}{2})}\int_{0}^{\infty}e^% {\theta(t-s)}I_{(0,\infty)}(t-s)s^{H-\frac{3}{2}}ds,t\in{\mathbb{R}}.$

Note that

$\displaystyle U_{t}^{H,\theta,\sigma}=\int_{\mathbb{R}}g_{H}(t-u)dM_{u},t\in{% \mathbb{R}}$

is the fractional Levy Ornstein-Uhlenbeck (FLOU) process satisfying the fractional Langevin equation

$\displaystyle dU_{t}=-\theta U_{t}dt+\sigma dM_{H,t},t\in{\mathbb{R}}.$

The process has long memory. Levy driven processes of Ornstein-Uhlenbeck type have been extensively studied over the last few years and widely used in finance, see Barndorff-Neilsen and Shephard (2001). FLOU process generalizes FOU process to include jumps. Maximum quasi-likelihood estimation in fractional Levy stochastic volatility model was studied in Bishwal (2011b). Berry-Esseen inequalities for the discretely observed Ornstein-Uhlenbeck-Gamma process was studied in Bishwal (2011c). Minimum contrast estimation in fractional Ornstein-Uhlenbeck process based on both continuous and discrete observations was studied in Bishwal (2011d). Bishwal (2011a) studied sufficiency problem in Vasicek model.

Now we focus on the fundamental semimartingale behind the O-U model. Define

$\displaystyle\kappa_{H}:=2H\Gamma(3/2-H)\Gamma(H+1/2),$ $\displaystyle k_{H}(t,s):=\kappa_{H}^{-1}(s(t-s))^{\frac{1}{2}-H},\eta_{H}:=% \frac{2H\Gamma(3-2H)\Gamma(H+\frac{1}{2})}{\Gamma(3/2-H)},$ $\displaystyle v_{t}\equiv v_{t}^{H}:=\eta_{H}^{-1}t^{2-2H},{\cal M}_{t}^{H}:=% \int_{0}^{t}k_{H}(t,s)dM_{s}^{H}.$

For using Girsanov theorem for Brownian motion, since a Radon-Nikodym derivative process is always a martingale, a central problem is how to construct an appropriate martingale which generates the same filtration, up to sets of measure zero, as the non-semimartingale called the fundamental martingale.

Extending Norros et al. (1999) it can be shown that ${\cal M}_{t}^{H}$ is a martingale, called the fundamental martingale whose quadratic variation $\langle{\cal M}^{H}\rangle_{t}$ is $v_{t}^{H}$ . Moreover, the natural filtration of the martingale ${\cal M}^{H}$ coincides with the natural filtration of the FLP $M^{H}$ since

$\displaystyle M_{t}^{H}:=\int_{0}^{t}K(t,s)d{\cal M}_{s}^{H}$

holds for $H\in(1/2,1)$ where

$\displaystyle K_{H}(t,s):=H(2H-1)\int_{s}^{t}r^{H-\frac{1}{2}}(r-s)^{H-\frac{3% }{2}}dr,0\leqslant s\leqslant t$

and for $H=1/2$ , the convention $K_{1/2}\equiv 1$ is used.

Define

$\displaystyle Q(t):=\frac{d}{dv_{t}}\int_{0}^{t}k_{H}(t,s)X(s)ds.$

It is easy to see that

$\displaystyle Q(t)=\frac{\eta_{H}}{2(2-2H)}\left\{t^{2H-1}Z_{i}(t)+\int_{0}^{t% }r^{2H-1}dZ(s)\right\}.$

Define the process $Z=(Z(t),t\in[0,T])$ by

$\displaystyle\widetilde{Z}(t):=\int_{0}^{t}k_{H}(t,s)dX(s).$

Extending Kleptsyna and Le Breton (2002), we have:

(i) $\widetilde{Z}$ is the fundamental semimartingale associated with the process $X$ .

(ii) $\widetilde{Z}$ is a $({\cal F}_{t})$ -semimartingale with the decomposition

$\displaystyle\widetilde{Z}(t)=\mu(\theta)\int_{0}^{t}X(s)dv_{s}+{\cal M}_{t}^{% H}.$

(iii) $X$ admits the representation

$\displaystyle X(t)=\int_{0}^{t}K_{H}(t,s)d\widetilde{Z}(s).$

(iv) The natural filtration $({\cal\widetilde{Z}}(t))$ of $\widetilde{Z}$ and $({\cal X}(t))$ of $X$ coincide.

We describe on our observations now. Note that for equally spaced data (homoscedastic case) $v_{t_{k}}-v_{t_{k-1}}=\eta_{H}^{-1}\linebreak\left(\frac{T}{n}\right)^{2-2H}[k% ^{2-2H}-(k-1)^{2-2H}],k=1,2,\ldots,n$ . For $H=0.5$ , $v_{t_{k}}-v_{t_{k-1}}=\eta_{H}^{-1}\left(\frac{T}{n}\right)^{2-2H}[k^{2-2H}-(k% -1)^{2-2H}]=\frac{T}{n},k=1,2,\ldots,n$ . We have

$\displaystyle Q(t)=\frac{d}{dv_{t}}\int_{0}^{t}k_{H}(t,s)X(s)ds=\kappa_{H}^{-1% }\frac{d}{dv_{t}}\int_{0}^{t}s^{1/2-H}(t-s)^{1/2-H}X(s)ds=\kappa_{H}^{-1}\eta_% {H}t^{2H-1}\frac{d}{dt}\int_{0}^{t}s^{1/2-H}(t-s)^{1/2-H}X(s)ds=\kappa_{H}^{-1% }\eta_{H}t^{2H-1}\int_{0}^{t}\frac{d}{dt}s^{1/2-H}(t-s)^{1/2-H}X(s)ds=\kappa_{% H}^{-1}\eta_{H}t^{2H-1}\int_{0}^{t}s^{1/2-H}(t-s)^{-1/2-H}X(s)ds.$

The process $Q$ depends continuously on $X$ and therefore, the discrete observations of $X$ does not allow one to obtain the discrete observations of $Q$ . The process $Q$ can be approximated by

$\displaystyle\widetilde{Q}(n)=\kappa_{H}^{-1}\eta_{H}n^{2H-1}\sum_{j=0}^{n-1}j% ^{1/2-H}(n-j)^{-1/2-H}X(j).$

It is easy to show that $\widetilde{Q}(n)\rightarrow Q(t)$ almost surely as $n\rightarrow\infty$ , see Tudor and Viens (2007).

Define a new partition $0\leqslant r_{1}<r_{2}<r_{3}<\ldots<r_{m_{k}}=t_{k},k=1,2,\ldots,n$ . Define

$\displaystyle\widetilde{Q}(t_{k})=\kappa_{H}^{-1}\eta_{H}t_{k}^{2H-1}\sum_{j=1% }^{m_{k}}r_{j}^{1/2-H}(r_{m_{k}}-r_{j})^{-1/2-H}{X(r_{j})}(r_{j}-r_{j-1}),$

$k=1,2,\ldots,n$ .

It is easy to show that $\widetilde{Q}(t_{k})\rightarrow Q(t)$ almost surely as $m_{k}\rightarrow\infty$ for each $k=1,2,\ldots,n$ .

We use this approximate observation in the calculation of our estimators. Thus our observations are

$\displaystyle X(t)\approx\int_{0}^{t}K_{H}(t,s)d\widetilde{Z}(s)\quad\text{% where}\quad\widetilde{Z}(t)=\theta\int_{0}^{t}\widetilde{Q}(s)dv_{s}+{\cal M}_% {t}^{H}.$

2. Quantile estimation in fractional student Ornstein-Uhlenbeck process

Consider the fractional Ornstein-Uhlenbeck process satisfying the SDE

$\displaystyle dX_{t}=-\alpha X_{t}dt+dL^{H}_{t}$

where $(L_{t}^{H})_{t\geqslant 0}$ is a fractional Levy process with driving term $(L_{t})_{t\geqslant 0}$ being a Levy process with cumulant function

$\displaystyle\ln E(e^{\textit{iuL}_{1}})=iu\mu-\lambda\frac{K_{\nu/2-1}(% \lambda|u|)}{K_{\nu/2}(\lambda|u|)},u\in{\mathbb{R}}\backslash\{0\}$

where $K$ is the modified Bessel function of the third kind defined as

$\displaystyle K_{\nu}(x)=\frac{1}{2}\int_{0}^{\infty}y^{\nu-1}e^{-\frac{x}{2}(% y+y^{-1})}dy\equiv\frac{1}{2}\int_{-\infty}^{\infty}e^{\nu u-x\cosh x}du,x>0$

where $\nu>0$ is arbitrary. Let $\mu\in{\mathbb{R}}$ and $\lambda>0$ be arbitrary. The weakly stationary solution satisfies

$\displaystyle X_{t}=e^{-\alpha t}X_{0}+e^{-\alpha t}\int_{0}^{t}e^{\alpha s}dL% _{s}$

and $X_{t}$ has marginal t-distribution $T(\nu,\lambda,\mu)$ with pdf

$\displaystyle f(x)=\frac{1}{\lambda{\cal B}(\nu/2,1/2)[1+((x-\mu)/\lambda)^{2}% ]^{(\nu+1)/2}},$

${\cal B}(\cdot,\cdot)$ being the beta function. If $\nu>1$ , then $E(X_{t})=\mu$ , and if $\nu>2$ , then $\text{corr}(X_{t+\tau},X_{t})=e^{-\alpha|\tau|}$ .

The fSOU process $(X_{t})$ has power function decay of the distribution tails. For a strictly increasing continuous mapping $h:{\mathbb{R}}\to{\mathbb{R}}$ consider the process $Y_{t}=h(X_{t})$ . It is obvious that $\{Y_{t}\}_{t\geqslant 0}$ is a stationary Markov process. Further, if $h$ is twice continuously differentiable then $\{Y_{t}\}_{t\geqslant 0}$ is a diffusion process.

2.1 Spearman’s estimator

The Spearman’s rank correlation coefficient $\rho_{S}$ for a pair of random variables $U$ and $V$ is

$\displaystyle\rho_{S}(U,V)=\rho(F_{U}(U),F_{V}(V))$

where $F_{U}$ and $F_{V}$ are (marginal) distribution functions of $U$ and $V$ respectively and $\rho$ is the linear correlation. Spearman’s rho is the linear correlation of the probability transformed random variables.

For any bivariate random vector $(Q_{1},Q_{2})$ following meta-Gaussian distribution

$\displaystyle\rho_{S}(Q_{1},Q_{2})=\frac{6}{\pi}\arcsin\frac{\rho}{2}.$

Since Spearman’s rho is invariant under strictly increasing transformation of data and $(Q_{0},Q_{k})$ is Gaussian, we have

$\displaystyle\rho_{S}(Y_{0},Y_{k})=\rho_{S}(Q_{0},Q_{k})=\frac{6}{\pi}\arcsin% \frac{\rho(Q_{0},Q_{k})}{2}=\frac{6}{\pi}\arcsin\frac{e^{-\alpha k}}{2}.$

Hence by inverting and plugging in gives the plug-in estimator of $\alpha$ based on Spearman’s rank correlation coefficient as

$\displaystyle\widehat{\alpha}_{k}:=-\frac{1}{k}\ln\left(2\sin\left(\frac{\pi}{% 6}\widehat{\rho}_{S}(Y_{0},Y_{k})\right)\right),k=1,2,\ldots$

where

$\displaystyle\widehat{\rho}_{S}(Y_{0},Y_{k}):=\frac{12}{n(n^{2}-1)}\sum_{k=1}^% {n}\left(\textit{rank}(Y_{0})-\frac{n+1}{2}\right)\left(\textit{rank}(Y_{k})-% \frac{n+1}{2}\right)$

2.2 Quantile estimator

Let

$\displaystyle g(x):=-\ln[\sin(2\pi(x-1/4)],x\in(1/4,1/2)$

and

$\displaystyle\widehat{p}_{k}^{0}:=\frac{1}{n-k}\sum_{j=0}^{n-k}I(Q_{j}>0,Q_{j+% k}>0)$

be an empirical estimate of $p_{k}:=P(Q_{0}>0,Q_{k}>0)$ . Then define

$\displaystyle\widehat{\alpha}_{k}:=\frac{1}{k}g(\widehat{p}_{k}^{0})=-\frac{1}% {k}\ln[\sin(2\pi(\widehat{p}_{k}^{0}-1/4)],k=1,2,\ldots.$

Borovkov and Decrouez (2013) showed that $\widehat{\alpha}_{k}$ is strongly consistent and asymptotically normal as $k\to\infty$ . Since $Q_{t}$ is unobservable, $\widehat{p}_{k}^{0}$ needs to be replaced with an estimator of $p_{k}$ from $Y_{t}$ .

Let $\widehat{m}_{n}(Y)$ is the sample median of the sample $Y_{0},Y_{1},\ldots,Y_{n}$ and $m$ is the median of $F$ . We have $\widehat{m}_{n}(Y)\to m$ a.s. as $n\to\infty$ .

Since $\{Q_{j}>0\}=\{Y_{j}>m\}$ , we have

$\displaystyle\widehat{p}_{k}^{0}:=\frac{1}{n-k}\sum_{j=0}^{n-k}I(Y_{j}>m,Y_{j+% k}>m).$

Since we do not know $m$ , we use the sample median and obtain

$\displaystyle\widehat{p}_{k}:=\frac{1}{n-k}\sum_{j=0}^{n-k}I(Y_{j}>\widehat{m}% _{n}(Y),Y_{j+k}>\widehat{m}_{n}(Y)).$

Therefore

$\displaystyle\check{\alpha}_{k}:=\frac{1}{k}g(\widehat{p}_{k}),k=1,2,\ldots$

will also be strongly consistent estimator of $\alpha$ . Further as $n\to\infty$

$\displaystyle\sqrt{n}(\widehat{p}_{k}^{0}-p_{k})\to{\cal N}(0,v_{k}^{2}),$

where

$\displaystyle v_{k}^{2}=p_{k}(1-p_{k})+2\sum_{l=1}^{\infty}[P(Q_{k+l}>0,Q_{l}>% 0,Q_{k}>0,Q_{0}>0)-p_{k}^{2}])$

and

$\displaystyle\sqrt{n}(\widehat{\alpha}_{k}^{0}-\alpha_{k})/k\to{\cal N}(0,v_{k% }^{2}(g^{\prime}(p_{k}))^{2}k^{-2}).$

follows from the continuity of $g$ .

2.3 Brownian excursion based estimator

The process $(Q_{t},t\geqslant 0)$ is unobservable, but the transformed process $Y_{t}=h(Q_{t})$ is observable. Even if, $Y_{t}$ is partially observed, that is, if we know the number of times the process’ trajectory crossed the band between two given levels $A$ and $B$ , during a given time interval, then one can estimate $\alpha$ , provided we also know the values of the distribution function $F(y)$ of the stationary distribution of $Y_{t}$ at the points $A$ and $B$ .

Let

$\displaystyle\tau_{a,b}:=\inf\{t>0:Q_{t}=a|Q_{0}=b\},a,b\in{\mathbb{R}}.$

$\tau_{a,b}$ is the time it takes the fOU process $(Q_{t},t\geqslant 0)$ to transit from $b$ to $a$ .

The excursion time $\tau_{(a),b}$ from $a$ to $b$ is the time it takes the fOU process starting from $b$ to reach the point $a$ and return back to $b$ . Clearly,

$\displaystyle\tau_{(a),b}=^{D}\tau_{a,b}+\tau_{b,a}$

where $\tau_{a,b}$ and $\tau_{b,a}$ are assumed to be independent.

Set $a:=h^{-1}(A),b:=h^{-1}(B)$ . Since $Q_{t}\sim{\cal N}(0,1/(2\alpha))$ in the stationary regime, one has $\Phi(\sqrt{2\alpha}a)=P(Q_{t}\leqslant a)=P(h(Q_{t})\leqslant A)=F(A)$ and likewise $\Phi(\sqrt{2\alpha}b)=P(Q_{t}\leqslant b)=P(h(Q_{t})\leqslant B)=F(B)$ . Therefore,

$\displaystyle a=\frac{\Phi^{-1}(F(A))}{\sqrt{2\alpha}}=:\frac{a^{*}}{\sqrt{% \alpha}},b=\frac{\Phi^{-1}(F(B))}{\sqrt{2\alpha}}=:\frac{b^{*}}{\sqrt{\alpha}}$

The standard scaling argument shows that the process $Q_{t}^{*}:=\alpha^{1/2}Q_{t/\alpha}$ satisfies the OU equation with $\alpha=1$ and another Brownian motion, namely $W_{t}^{*}=\alpha^{1/2}W_{t/\alpha}$ . Therefore

$\displaystyle\tau_{(a),b}=^{D}\alpha^{-1}\tau^{*}_{(a^{*}),b^{*}}$

where the asterisk indicates that the excursion time on the right side is for the process $X_{t}^{*}$ , the distribution of $\tau^{*}_{(a^{*}),b^{*}}$ being independent of $\alpha$ and having a mean $m_{a^{*},b^{*}}$ .

Now strong Markov property implies that subsequent excursion from $b$ to $a$ in $\{Q_{t}\}$ (or, equivalently, from $B$ to $A$ in $\{Y_{t}\}$ ) form an i.i.d. sequence. Therefore, if $N_{A,B}(T)$ denotes the number of times $\{Y_{t}\}$ crosses the band between the levels $A$ and $B$ in a given direction during the time interval $[0,T]$ , the integral renewal theorem asserts that, as $T\to\infty$ ,

$\displaystyle\frac{N_{A,B}(T)}{T}\to^{a.s.}\frac{1}{E(\tau_{(a),b})}=\frac{% \alpha}{m_{a^{*},b^{*}}}$

which gives another consistent estimator of $\alpha$ given by

$\displaystyle\widetilde{\alpha}_{T}:=\frac{m_{a^{*},b^{*}}N_{A,B}(T)}{T}.$

3. Quantile estimation in fractional inverse gaussian-Ornstein-Uhlenbeck stochastic volatility model

We consider the SDEs

$\displaystyle dY_{t}=(\mu+\beta X_{t})dt+\sqrt{X_{t}}dW_{t}+\rho dZ_{\theta t}$ $\displaystyle dX_{t}=-\theta X_{t}dt+dZ_{\theta t}^{H}$

where $Z_{t}$ is a Levy process and $Z_{t}^{H}$ is a fractional Levy process independent of $X_{0}$ with ${\cal L}(X_{0})=IG(\delta,\gamma)$ . We suppose that $\delta$ and $\gamma$ are known. Here $\theta>0$ and $\rho<1$ .

When the process $Z$ is inverse-Gaussian, the model is fIGOU process. In fIGOU model, calculation of conditional cummulants of the integrated volatility conditioned on the initial value is enough to be able to compute European style options.

Note that the cumulative process or the integrated process $I_{t}=\int_{0}^{t}X_{u}du$ has long range dependence or long memory, see Barndorff-Neilsen and Shephard (1998).

The cumulant functions of IGOU process are given by

$\displaystyle k(u)=\log E[e^{-uZ(1)}]=-u\delta\gamma^{-1}(1+2u\gamma^{-2})^{-1% /2},$ $\displaystyle k^{\prime}(u)=\log E(e^{-uX_{t}})=\delta\gamma-\delta\gamma(1+2u% \gamma^{-2})^{1/2},u\in{\mathbb{R}}.$

The process $Z$ is the sum of two independent Levy processes $Z=Z^{(1)}+Z^{(2)}$ where ${\cal L}(X_{0})=IG(\delta/2,\gamma)$ and $Z^{(2)}$ is a compound Poisson process given by $Z^{(2)}=\gamma^{-2}\sum_{j=1}^{N_{t}}u_{j}$ with $N$ being a Poisson process with intensity $\delta\gamma/2$ and $u_{j}$ is a sequence of independent and identically $\chi_{1}^{2}$ -distributed random variables independent of $N$ , see Barndorff-Neilsen and Shephard (1998).

The processes $Z$ and $X$ have infinitely many jumps in any finite time interval, hence they are infinite activity processes.

Invariant distribution is Generalized Inverse Gaussian (GIG).

$\displaystyle{\cal L}(X_{0})={\cal L}(X_{t})={\cal GIG}(\lambda,\delta,\sqrt{% \alpha^{2}-\beta^{2}}).$

Mixture distribution is $Q_{\theta}$ Generalized Hyperbolic (GH).

$\displaystyle{\cal L}(Y_{t})={\cal GH}\left(\theta,\frac{\alpha}{\sqrt{t}},% \frac{\beta}{\sqrt{t}},\sqrt{t}\delta,\mu t\right).$ $\displaystyle{\cal L}(Y_{t+1}-Y_{t})={\cal GH}(\lambda,\alpha,\beta,\delta,\mu).$

First we study the LAD estimators in the stable OU case for small $\Delta$ . Consider the model

$\displaystyle dX_{t}=-\theta X_{t}dt+dZ_{t}$

where $Z_{t}$ is a Levy process independent of $X_{0}$ . The Least Absolute Deviation (LAD) estimator is robust to “outlying data”. The LAD estimation has a long history and is one of popular estimation procedure robust to outlers. The LAD estimation is based on the Laplacian $L^{1}$ -loss while the LSE is on the Laussian $L^{2}$ -loss. We refer to Portnoy and Koenker (1997), Knight (1998) and Koenker (2005) as well as the references therein for a detailed account and historical background on LAD estimation. For time series literature, see Davis and Dunsmuir (1997) and Davis, Knight and Liu (1992).

The least absolute deviation (LAD) estimator is defined as the minimizer $\hat{\theta}_{n}$ of the contrast function $\theta\to\sum_{i=1}^{n}|X_{t_{i}}-e^{-\theta\Delta}X_{t_{i-1}}|$ . For fixed $T$ , asymptotic normality of $\hat{\theta}_{n}$ is achieved at the rate $n^{1/\beta-1/2}$ where $\beta$ stands for the activity index of the driving Levy process, also known as the Blumenthal-Getoor (1961) activity index defined as

$\displaystyle\beta=\inf\left\{r>0:\int_{|z|\leqslant 1}|z|^{r}\nu(dz)<\infty\right\}$

which is the degree of small-jump fluctuations.

Now consider another least absolute deviation (LAD) estimator which is defined as the minimizer $\tilde{\theta}_{n}$ of the contrast function $\theta\to\sum_{i=1}^{n}|X_{t_{i}}+\theta X_{t_{i-1}}\Delta|$ .

Under infill and large time sampling design, that is when $\Delta\to 0$ and $n\Delta\to\infty$ , asymptotic normality of $\tilde{\theta}_{n}$ is achieved at the rate $\sqrt{n}\Delta^{1-1/\beta}$ where $\beta$ stands for the activity index of the driving Levy process. Note that $\sqrt{n}\Delta^{1-1/\beta}=T^{1-1/\beta}n^{(2-\beta)/(2\beta)}$ . This implies that the rate of convergence is determined by the most active part of the driving Levy process, the presence of a driving Wiener part leads to $\sqrt{n\Delta}$ , which is familiar in the context of asymptotically efficient estimation of diffusions with compound Poisson jumps, while a pure-jump driving Levy process leads to a faster one. As a result, when $Z$ is a pure jump Levy process, we have a faster rate of convergence that the familiar rate $\sqrt{n\Delta}$ . It is interesting to note that rate of convergence is faster by only changing the type of loss from $L^{2}$ to $L^{1}$ .

Using self-weighted LAD (SLAD) contrast function $\theta\to\sum_{i=1}^{n}w(X_{t_{i-1}})|X_{t_{i}}+\theta X_{t_{i-1}}|$ , for a bounded continuous weight function $w$ , the rate of convergence can be improved to the conventional rate $\sqrt{n}$ , see Masuda (2010) which extended autoregressive process with infinite variance studied in Ling (2005). An example of weight function is $w(x)=\exp(-|x|)$ . The unweighted weight function corresponds to the case $w\equiv 1$ .

In the Wiener case, the unweighted SLAD estimator leans to the asymptotic variance $\pi\theta_{0}$ where as the asymptotic variance of the exact MLE is $2\theta_{0}$ . Hence the asymptotic efficiency of the SLAD estimator relative to the MLE is $2/\pi$ . This is same as asymptotic relative efficiency for the same the sample median over sample mean in estimating the mean of i.i.d. normal samples. For the asymptotic normality of SLAD estimator one does not need the rapidly increasing experimental design $n\Delta^{2}\to 0$ which is quite inevitable while adopting contrast function based on Euler-type approximation. For SLAD estimator, the weaker condition $n\Delta^{3}\to 0$ is sufficient. The SLAD estimator converges more rapidly than the LSE as soon as $\Delta=o(n^{-1/2}(\log n)^{1/(2-\beta)})$ . The sampling design condition for $h$ is $\Delta=n^{-\tau}$ where $\tau\in(0,1)$ .

An interval estimator of $\tau$ is given by

$\displaystyle\frac{\beta}{2(2\beta-1)}<\tau<\frac{\beta}{2\beta-1},\text{if }% \beta>1,\text{and }\frac{\beta}{2(2\beta-1)}<\tau<1,\text{if }\beta\in(2/3,1).$

An example of infill interval is $h=n^{-3/5}$ .

Next we consider the SDEs

$\displaystyle d\bar{Y}_{t}=(\mu+\beta X_{t})dt+\sqrt{X_{t}}dW_{t}+\rho dZ_{% \theta t},$ $\displaystyle dX_{t}=-\theta X_{t}dt+dZ_{\theta t}^{H}$

where $Z_{t}$ is a Levy process independent of $X_{0}$ and ${\cal L}(X_{0})=IG(\delta,\gamma)$ . We suppose that $\delta$ and $\gamma$ are known and we are interested in estimating $\vartheta_{0}=(\theta,\rho)$ , where $\theta>0$ and $\rho<0$ .

Note that if $Q_{j\Delta},1\leqslant j\leqslant n$ were observed, then

$\displaystyle\theta_{n}:=-\frac{1}{\Delta}\ln\left(\min_{1\leqslant j\leqslant n% }\frac{Q_{j\Delta}}{Q_{(j-1)\Delta}}\right)$

is a weakly consistent estimator of $\theta$ as $n\rightarrow\infty$ similar to the non-fractional case Jongbloed et al. (2005). In our case, $Q_{j\Delta},1\leqslant j\leqslant n$ are unobserved.

We can estimate the state process $Q$ by $\hat{Q}$ given by the recursion

$\displaystyle\hat{Q}_{t_{i}}=e^{-\theta\Delta}\hat{Q}_{t_{i-1}}+Z_{t_{i}}-Z_{t% _{i-1}}.$

Then we estimate $\theta$ by

$\displaystyle\check{\theta}_{n}:=-\frac{1}{\Delta}\ln\left(\min_{1\leqslant j% \leqslant n}\frac{\hat{Q}_{j\Delta}}{\hat{Q}_{(j-1)\Delta}}\right).$

Consider

$\displaystyle dY_{t}=(\mu+\beta Q_{t})dt+\sqrt{Q_{t}}dW_{t}+\rho dZ_{\theta t}.$

In order to construct the estimating functions, we use the first and second cumulants which are given respectively by

$\displaystyle\kappa_{y_{1}}^{(1)}=\theta\rho\Delta\kappa_{IG}^{(1)},\kappa_{y_% {1}}^{(2)}=\Delta\kappa_{IG}^{(1)}+2\theta\rho^{2}\Delta\kappa_{IG}^{(2)}.$

Inverting these cumulants and replacing the cumulants by their sample quantities, we obtain the explicit the moment estimators of $\rho$ and $\lambda$ .

The moment estimators of $\rho$ and $\theta$ are given by

$\displaystyle\hat{\rho}_{n}:=\frac{\gamma(\gamma s_{y}^{2}-\Delta\delta)}{2% \bar{y}},\hat{\theta}_{n}:=\frac{\gamma\bar{y}}{\Delta\delta\hat{\rho}_{n}}$

where

$\displaystyle s_{y}^{2}:=\frac{1}{n}\sum_{j=1}^{n}(y_{j}-\bar{y})^{2}=\frac{1}% {n}\sum_{j=1}^{n}y_{j}^{2}-(\bar{y})^{2},$ $\displaystyle\bar{y}:=\frac{1}{n}\sum_{j=1}^{n}y_{j},y_{j}:=Y_{j\Delta}-Y_{(j-% 1)\Delta}$

Let $\vartheta=(\rho,\theta).$ and $\hat{\vartheta}_{n}=(\hat{\rho}_{n},\hat{\theta}_{n})$ . We have the following properties of the estimators following Masuda (2005):

Proposition 3.1. For fixed $\Delta>0$ as $n\to\infty$ ,

$\displaystyle(\text{a})\hat{\vartheta}_{n}\rightarrow\vartheta_{0}\ \ \text{a.% s.}\ \ \text{as }n\rightarrow\infty.(\text{b})\sqrt{n}(\hat{\vartheta}_{n}-% \vartheta_{0})\rightarrow^{\cal D}{\cal N}_{2}(0,(2\theta\rho^{2}\Delta^{2}% \delta^{2}\gamma^{-4})^{-2}D(\vartheta_{0}))\text{ as }n\rightarrow\infty.$

where $D(\vartheta_{0})$ is the limiting covariance matrix.

Remark: In the fIG-OU stochastic volatility model, for the case of $\rho=-1,\gamma=1,\delta=1$ the moment estimator is given by $\hat{\theta}_{n\Delta}:=-\frac{Y_{n\Delta}-Y_{0}}{{n\Delta}}$ . Thus the parameter $\theta$ can be estimated by just the two terminal observations.

The moment estimators are sensitive to outliers since they are based on mean and standard deviation of the sample data. In order to incorporate outliers and model misspecifications, we consider robust estimators. The robust estimators of $\rho$ and $\lambda$ are given by

$\displaystyle\tilde{\rho}_{n}:=\frac{\gamma(\gamma a_{y}-\Delta\delta)}{2% \tilde{y}},\tilde{\theta}_{n}:=\frac{\gamma\tilde{y}}{\Delta\delta\tilde{\rho}% _{n}}\ \ \text{where}\ \ a_{y}:=\frac{1}{n}\sum_{j=1}^{n}|y_{j}-\tilde{y}|$

is the sample mean absolute deviation from median,

$\displaystyle\tilde{y}:=\text{median of }\{y_{j},1\leqslant j\leqslant n\}$

which is defined as

$\displaystyle\tilde{y}=\left\{\begin{array}[]{r@{\quad:\quad}l}\frac{y_{k}+y_{% k+1}}{2}\quad:\quad&n=2k\\ y_{k+1}\quad:\quad&n=2k+1\end{array}\right.$

Let $\vartheta=(\rho,\theta).$ and $\tilde{\vartheta}_{n}=(\tilde{\rho}_{n},\tilde{\theta}_{n}).$ We have the following properties of the estimators. By using the standard theory of order statistics (see Theorem 5.9 and 5.21 in Van der Vaart (2000)) and mixing property of the process, along with Glivenko-Cantelli argument and Delta method, we obtain:

Proposition 3.2. For fixed $\Delta>0$ as $n\to\infty$ ,

$\displaystyle(\text{a})\tilde{\vartheta}_{n}\rightarrow\vartheta_{0}\ \ \text{% a.s.}\ \ \text{as }n\rightarrow\infty.(\text{b})\sqrt{n}(\tilde{\vartheta}_{n}% -\vartheta_{0})\rightarrow^{\cal D}{\cal N}_{2}(0,\frac{\pi}{2}(2\theta\rho^{2% }\Delta^{2}\delta^{2}\gamma^{-4})^{-2}D(\vartheta_{0}))\text{ as }n\rightarrow\infty.$

where $D(\vartheta_{0})$ is the limiting covariance matrix.

4. Method of moments estimation in fractional-Gamma-Ornstein-Uhlenbeck stochastic volatility model

We generalize Ornstein-Uhlenbeck process to include non-normal innovations. First we study the asymptotic behavior of the ratio estimator of the drift parameter in Gamma-Ornstein-Uhlenbeck (GOU) volatility process based on observations of the asset price process. This model captures the stylized facts as it preserves jumps in the volatility process. We study the behavior of the moment estimators.

Bishwal (2011c) studied estimation for the discretely observed Ornstein-Uhlenbeck-Gamma (OUG) process. We note that this is different from Gamma-Ornstein-Uhlenbeck (GOU) process we study here. OUG process is an Ornstein-Uhlenbeck process with an additive Brownian noise and gamma noise. GOU process is a pure jump process. OUG process may be compared with the GAR process.

Recall that the autoregressive gamma (ARG) process can be written as

$\displaystyle X_{t}=\sum_{j=1}^{N_{t}}U_{j}+\epsilon_{t}$

where $N_{t}$ is a Poisson process, $U_{j}$ are independent and identically distributed Gamma distributed random variables with shape parameter 1 and scale parameter $c_{t}$ . independent of $N$ and $\epsilon_{t}$ are gamma distributed random variables with shape parameter $\nu$ and scale parameter $c_{t}$ .

We propose the discretized version of the Heston model as:

$\displaystyle y_{i}=\mu+\sqrt{x_{i}}\xi_{i},x_{i}=\sum_{j=1}^{N_{t}}U_{i,j}+% \epsilon_{i}$

where $y_{i}$ is the return, $y_{i}$ is the volatility, $\xi_{i}$ are normally distributed, $\epsilon_{i}$ are gamma distributed, $(\xi_{i},\epsilon_{i})$ has the correlation $\rho$ and $U_{j}$ are independent $\chi_{1}^{2}$ distributed.

We consider the special cases when $X$ is an fGOU process. $Z$ has the gamma density with parameters $\delta$ and $\gamma$ . The stochastic volatility model has the form

$\displaystyle d\bar{Y}_{t}=(\mu+\beta X_{t})dt+\sqrt{X_{t}}dW_{t}+\rho dZ_{t},$ $\displaystyle dX_{t}=-\theta X_{t}dt+dZ_{t}^{H}$

where $Z_{t}^{H}$ is a fractional gamma process independent of $X_{0}$ and $Z_{t}$ is a gamma process.

Then we consider the fractional Gamma-OU process. Let $({\Omega,{\cal F},\{{{\cal F}_{t}}\}_{t\geqslant 0},P})$ be a stochastic basis on which is defined the Ornstein-Uhlenbeck process ${X_{t}}$ satisfying the Itô stochastic differential equation

$\displaystyle dX_{t}=-\theta X_{t}dt+dZ_{t}^{H},t\geqslant 0,$

where $\{Z_{t}^{H}\}$ is a fractional Gamma process with the filtration $\{{\cal F}_{t}\}_{t\geqslant 0}$ and $\theta>0$ is the unknown parameter to be estimated on the basis of continuous observation of the process $\{Y_{t}\}$ on the time interval $[0,T]$ . The solution of the above SDE is given by

$\displaystyle X_{t}=\int_{-\infty}^{t}e^{-\theta(t-s)}dZ_{s}^{H}.$

This process is stationary. In fact, it can be shown that $X_{t_{i}}$ is a stationary discrete time AR (1) process with autoregression coefficient $\phi\in(0,1)$ with the following representation

$\displaystyle X_{t_{i}}=\phi X_{t_{i-1}}+\epsilon_{t_{i-1}}$

where

$\displaystyle\phi=e^{-\theta\Delta},\epsilon_{t_{i-1}}=\int_{t_{i-1}}^{t_{i}}e% ^{-\theta(t_{i}-u)}dZ_{u}^{H}.$

Consider

$\displaystyle dY_{t}=(\mu+\beta Q_{t})dt+\sqrt{Q_{t}}dW_{t}+\rho dZ_{t}.$

The invariant distribution is Gamma:

$\displaystyle{\cal L}(Q_{0})={\cal L}(Q_{t})={\cal G}(\theta,\delta,\sqrt{% \alpha^{2}-\beta^{2}}).$

Mixture distribution is Variance-Gamma denoted by $V G$ .

$\displaystyle{\cal L}(Y_{t})={\cal VG}\left(\theta,\frac{\alpha}{\sqrt{t}},% \frac{\beta}{\sqrt{t}},\sqrt{t}\delta,\mu t\right).$ $\displaystyle{\cal L}(Y_{t+1}-Y_{t})={\cal VG}(\theta,\alpha,\beta,\delta,\mu).$

Let $q(x,\theta)$ is a density function of $Y_{t}$ . The contrast function is defined as

$\displaystyle u(x,\theta)=-\log q(x,\theta).$

Minimum contrast estimator (MCE) is defined as

$\displaystyle\hat{\theta}_{T}:=\arg\inf_{\theta}u(x,\theta).$

Asymptotic properties of MCE was studied in Bishwal (2006, 2010, 2011e). The density function of the GH distribution is given by

$\displaystyle d_{{\cal GH}(\theta,\alpha,\beta,\delta,\mu)}(x)=\frac{(\alpha^{% 2}-\beta^{2})^{\frac{\theta}{2}}}{\sqrt{2\pi}\alpha^{\theta-\frac{1}{2}}\delta% ^{\theta}K_{\theta}(\delta\sqrt{\alpha})}(\delta^{2}+(x-\mu)^{2})^{(\theta-% \frac{1}{2})/2}e^{\beta(x-\mu)}K_{\theta-\frac{1}{2}}(\alpha\sqrt{\delta^{2}+(% x-\mu)^{2}}).$

Here $K_{\nu}$ denotes the modified Bessel function of third kind with index $\nu$ . When $\theta=-\frac{1}{2}$ , GIG becomes $IG(\delta,\gamma)$ .

The moment estimators are defined as

$\displaystyle\hat{\theta}_{n}:=\frac{\frac{1}{n^{2}}\left[\sum_{i=1}^{n}(Y_{i% \Delta}-Y_{(i-1)\Delta})\right]^{2}}{\frac{1}{n^{2}}\sum_{i=1}^{n}(Y_{i\Delta}% -Y_{(i-1)\Delta})^{2}-\frac{\Delta}{n}\left[\sum_{i=1}^{n}(Y_{i\Delta}-Y_{(i-1% )\Delta})\right]}\frac{2a^{3}(a+1)}{b^{4}\Delta},$ $\displaystyle\hat{\rho}_{n}:=\frac{\frac{1}{n^{2}}\sum_{i=1}^{n}(Y_{i\Delta}-Y% _{(i-1)\Delta})^{2}-\frac{\Delta}{n}\left[\sum_{i=1}^{n}(Y_{i\Delta}-Y_{(i-1)% \Delta})\right]}{\frac{1}{n^{2}}\left[\sum_{i=1}^{n}(Y_{i\Delta}-Y_{(i-1)% \Delta})\right]}\frac{b^{3}\Delta}{2a^{2}(a+1)}.$

For the exponential AR (1) model, the ratio estimator of $\theta$ is defined as

$\displaystyle\hat{\theta}_{n}:=-\frac{1}{\Delta}\ln\left[\min_{1\leqslant i% \leqslant n}\frac{Q_{i\Delta}}{Q_{(i-1)\Delta}}\right].$

This estimator is motivated by the extreme value theory for the correlation parameter of an AR (1) process whose innovation distribution is positive. See Davis and McCormick (1989). In the case of exponential AR (1) process, it coincides with the maximum likelihood estimator. See Neilsen and Shephard (2003).

The weak consistency of the ratio estimator in the LOU process was studied in Jongbloed et al. (2005). The strong consistency and asymptotic Weibullness was studied in Brockwell, Davis and Yang (2007) in the case of Gamma innovations.

Using the techniques of AR (1) type model with exponential innovations (Davis and McCormick (1989)), we also obtain the following estimator of the drift $\theta$ which are defined as

$\displaystyle\check{\theta}_{n}:=-\frac{1}{\Delta}\ln\left(\min_{1\leqslant i% \leqslant n}\frac{\hat{Q}_{i\Delta}}{\hat{Q}_{(i-1)\Delta}}\right)$

and

$\displaystyle\tilde{\theta}_{n}:=-\frac{1}{\Delta}\ln\left[\frac{\sum_{i=1}^{n% }\hat{Q}_{i\Delta}}{\sum_{i=1}^{n}\hat{Q}_{(i-1)\Delta}}\right]$

where $\hat{Q}$ is the estimator of $Q$ based on observations of $Y$ which could be obtained, for example, by Kitagawa algorithm. The limit distribution of the first estimator would be Weibull which can be useful for extreme value theory in finance. Bishwal (2011f) studied extreme value theory in finance.

Let $\vartheta=(\rho,\theta)$ and $\hat{\vartheta}_{n}=(\hat{\rho}_{n},\hat{\theta}_{n})$ . By using Theorem 2.2 in Masuda (2005) (see also Theorem 4.1 Van der Vaart (2000)), we obtain the strong consistency and asymptotic normality of the method of moments (MM) estimators:

Proposition 4.1. For fixed $\Delta>0$ as $n\to\infty$ ,

$\displaystyle(\text{a})\hat{\vartheta}_{n}\rightarrow\vartheta_{0}\ \ \text{a.% s.}\ \ \text{as }n\rightarrow\infty.(\text{b})\sqrt{n}(\hat{\vartheta}_{n}-% \vartheta_{0})\rightarrow^{\cal D}{\cal N}_{2}(0,(I^{-1}(\vartheta_{0}))\text{% as }n\rightarrow\infty.$

where $I(\vartheta_{0})$ is the Fisher information matrix.

5. Numerical example

Consider the example satisfying the SDEs

$\displaystyle dY_{t}=\sqrt{X_{t}}dW_{t}+\rho dZ_{\theta t},$ $\displaystyle dX_{t}=-\theta X_{t}dt+dZ_{\theta t}^{H}$

To generate the trajectories of $Y$ , one can apply the Rosinski’s method, see Rosinski (2000). By virtue of this method, we do not need any discretization scheme like the Euler-Maruyama scheme, see Glasserman (2004). For simulation, suppose that $\delta=1$ and $\gamma=1$ are known, and the true values are $\theta=3$ and $\rho=-1$ and $H=0.5$ . Hence the simulated model is given by

$\displaystyle dY_{t}=\sqrt{X_{t}}dW_{t}-dZ_{3t},$ $\displaystyle dX_{t}=-3X_{t}dt+dZ_{3t}.$

The estimator is given by

$\displaystyle\hat{\theta}_{n}=-\frac{\bar{y}}{\Delta}=-\frac{1}{n\Delta}\sum_{% j=1}^{n}y_{j}=-\frac{1}{n\Delta}\sum_{j=1}^{n}(Y_{j\Delta}-Y_{(j-1)\Delta})=-% \frac{1}{n\Delta}(Y_{n\Delta}-Y_{0}).$

For $\Delta=1$ , the estimator is just the negative of the sample mean. The mean and standard deviation of the estimator based on large independent trajectories of $Y$ can be calculated. The estimator becomes better when $n$ and $\Delta$ increase.

The similarity of this estimator can be compared with the classical Brownian model. In a vanishing interval setup when $\Delta\to 0$ such that $n\Delta=T$ , we have $\hat{\theta}_{T}=-\frac{W_{T}}{T}$ which is the MLE of $\theta=0$ . By the strong law of large numbers for Brownian motion, $\hat{\theta}_{T}\to 0$ almost surely as $T\to\infty$ . The example should help the use of Monte Carlo techniques to fit empirical log-return data.

Remarks. For the Spearman’s rank correlation based estimator, with $H=0.5$ , the estimator values are different for different values of $k$ . The root mean square errors of the estimators $\hat{p}_{k}$ are smaller than those for $\hat{p}_{k}^{0}$ . As $\alpha$ increases, the value of the estimator ratio decreases to some limiting value. The limiting value corresponds to the case of independent observations $Y_{j},j=0,1,\ldots,n$ . The quality of the estimators $\hat{p}_{k}^{0}$ and $\hat{p}_{k}$ is determined by how often the difference $Y_{j}-a$ changes its sign for $a=0$ and $a=\hat{m}_{n}(Y)$ , respectively, and in the later case this occurs more often. For the robust estimator, standard deviation of the unweighted LAD estimator is much smaller than the LSE. Hence LAD estimator is a better estimator.

References

Araya

Bahamonde

Torres

, & Viens

(2019). Donsker type theorem for fractional Poisson process. Statist. Probab. Letters, 150, 1-8.

Barndorff-Neilsen

O.E.

, & Shephard

(1998). Processes of normal inverse Gaussian type. Finance and Stochastics, 2, 41-68.

Barndorff-Nielsen

O.E.

, & Shephard

(2001). Non-Gaussian Ornstein-Uhlenbeck-based models and some of their uses in financial economics (with discussion). Journal of the Royal Statistical Society, Series B, 63, 167-241.

Bishwal

J.P.N.

(2006). Rates of weak convergence of the approximate minimum contrast estimators for the discretely observed Ornstein-Uhlenbeck process. Statistics and Probability Letters, 76(13), 1397-1409.

Bishwal

J.P.N.

(2008). Parameter Estimation in Stochastic Differential Equations, Lecture Notes in Mathematics, 1923, Springer-Verlag.

Bishwal

J.P.N.

(2010). Uniform rate of weak convergence for the minimum contrast estimator in the Ornstein-Uhlenbeck process. Methodology and Computing in Applied Probability, 12(3), 323-334.

Bishwal

J.P.N.

(2011a). Sufficiency and Rao-Blackwellization of Vasicek model. Theory of Stochastic Processes, 17(1), 12-15.

Bishwal

J.P.N.

(2011b). Maximum quasi-likelihood estimation in fractional Levy stochastic volatility model. Journal of Mathematical Finance, 1(3), 12-15.

Bishwal

J.P.N.

(2011c). Berry-Esseen inequalities for the discretely observed Ornstein-Uhlenbeck-Gamma process. Markov Processes and Related Fields, 17(1), 119-150.

10.

Bishwal

J.P.N.

(2011d). Minimum contrast estimation in fractional Ornstein-Uhlenbeck process: Continuous and discrete sampling. Fractional Calculus and Applied Analysis, 14(3), 375-410.

11.

Bishwal

J.P.N.

(2011e). Some new estimators of integrated volatility. American Open Journal of Statistics, 1(2), 74-80.

12.

Bishwal

J.P.N.

(2011f). Financial extremes: A short review. Advances and Applications in Statistics, 25(1), 1-14.

13.

Bishwal

J.P.N.

(2022a). Quasi-likelihood Estimation in fractional Levy SPDEs from Poisson sampling. European Journal of Mathematical Analysis, 2(2), 1-16.

14.

Bishwal

J.P.N.

(2022b). Parameter Estimation in Stochastic Volatility Models, Springer Nature, Cham.

15.

Bishwal

J.P.N.

(2023). Parameter estimation for SPDEs driven by cylindrical stable processes. Eur. J. Math. Anal, 3(4), 1-25. doi: 10.28924/ada/ma.3.4.

16.

Blumenthal

R.M.

, & Getoor

R.K.

(1961). Sample functions of stochastic processes with stationary independent increments. J. Math. Mech, 10, 493-516.

17.

Brockwell

P.J.

Davis

R.A.

, & Yang

(2007). Estimation of non-negative Levy driven Ornstein-Uhlenbeck processes. Journal of Applied Probability, 44, 977-989.

18.

Borovkov

, & Decrouez

(2013). Ornstein-Uhlenbeck processes with heavy distribution tails. Theory of Probability and its Applications, 57(3), 396-418.

19.

Davis

R.A.

, & McCormick

W.P.

(1989). Estimation for first order autoregressive processes with positive or bounded innovations, Stochastic Process. Appl, 31, 237-250.

20.

Davis

R.A.

Knight

, & Liu

(1992). M-estimation for autoregressions with infinite variance. Stochastic Processes and their Applications, 40, 145-180.

21.

Davis

R.A.

, & Dunsmuir

W.T.M.

(1997). Least absolute deviation estimation for regression with ARMA errors. Journal of Theoretical Probability, 10, 481-497.

22.

Davis

R.A.

, & McCormick

W.P.

(1989). Estimation for first order autoregressive processes with positive or bounded innovations. Stochastic Process. Appl, 31, 237-250.

23.

Engle

, & Russell

(1998). The autoregressive conditional duration model: A new model for irregularly spaced data. Econometrica, 66, 1127.1162.

24.

Gaver

D.P.

, & Lewis

P.A.W.

(1980). First order autoregressive gamma sequences and point processes. Adv. Appl. Probab, 12, 727-745.

25.

Glasserman

(2004). Monte Carlo Methods in Financial Engineering, Springer, New York.

26.

Gourieroux

, & Jasiak

(2006). Autoregressive gamma process. Journal of Forecasting, 25(2), 129-152.

27.

Grigelionis

(2012). Student’s t-distribution and Related Stochastic Processes, Springer.

28.

Jacod

, & Shiryev

A.N.

(1987). Limit Theorems for Stochastic processes, Springer-Verlag, Berlin.

29.

Jongbloed

van der Meulen

F.H.

, & van der Vaart

A.W.

(2005). Nonparametric inference for Levy driven Ornstein-Uhlenbeck processes. Bernoulli, 11(5), 759-791.

30.

Kleptsyna

M.L.

, & Le Breton

(2002). Statistical analysis of the fractional Ornstein-Uhlenbeck type process. Statistical Inference for Stochastic Processes, 5(3), 229-248.

31.

Knight

(1998). Limit distributions for L1 regression estimators under general conditions. Annals of Statistics, 26, 755-770.

32.

Koenker

(2005). Quantile Regression, Cambridge University Press, Cambridge.

33.

Kolmogorov

A.N.

(1940). Wiener skewline and other interesting curves in Hilbert space. Doklady Akad. Nauk, 26, 115-118.

34.

Levy

(1948). Processus stochastiques et movement Brownien, Paris.

35.

Ling

(2005). Self-weighted least absolute deviation estimation for infinite variance autoregressive models. Journal of Royal Statistical Society, Series B, Statistiacl Methodology, 67, 381-393.

36.

Leonenko

N.N.

Papic

Sikorskii

, & Suvak

(2020). Approximation of heavy-tailed fractional Pearson diffusions in Skorohod topology. Journal of Mathematical Analysis and Applications, 486, 1-22.

37.

Leonenko

N.N.

, & Suvak

(2010). Statistical inference for student diffusion process. Stoch. Anal. Appl, 28, 972-1002.

38.

Mandelbrot

, & Van Ness

J.W.

(1968). Fractional Brownian motions, fractional noises and applications. SIAM Review, 10, 422-437.

39.

Masuda

(2005). Classical method of moments for partially and discretely observed ergodic models. Statistical Inference for Stochastic Processes, 8, 25-50.

40.

Masuda

(2010). Asymptotic mixed normality in the self-weighted LAD estimation of Levy OU process. Electronic Journal of Statistics, 4, 525-565.

41.

Neilsen

, & Shephard

(2003). Likelihood analysis of a first order autoregressive model with exponential innovations. J. Time Series Anal, 24(3), 337-344.

42.

Norros

Valkeila

, & Virtamo

(1999). An elementary approach to a Girsanov formula and other analytical results on fractional Brownian motion. Bernoulli, 5, 571-587.

43.

Portnoy

, & Koenker

(1997). The Gaussian hare and Laplacian tortoise: Computability of squared-error versus absolute-error estimators. Statistical Science, 12, 279-300.

44.

Rosinski

(2000). Stochastic series representations of inverse Gaussian and other exponentially tempered stable processes, Research Report 42, Center for Mathematics, Physics and Stochastics, University of Aarhus, Aarhus.

45.

Tudor

C.A.

, & Viens

(2007). Statistical aspects of the fractional stochastic calculus. Annals of Statistics, 35(3), 1183-1212.

46.

Van der Vaart

A.W.

(2000). Asymptotic Statistics, Cambridge University Press, Cambridge.

47.

Wolfe

A.J.

(1982). On a continuous analogue of the stochastic difference equation Xn=ρ⁢Xn-1+Bn. Stochastic Process Appl, 12, 301-312.