Empirical likelihood method for longitudinal data generated from unequally-spaced Lèvy processes

Abstract

By introducing the notion of “empirical likelihood function of observing sums”, unequally-spaced time series data and longitudinal data generated from Lévy processes can be analyzed. Characteristic function is further incorporated to handle the situations where the density function of the increments is difficult to obtain. In the situations where both characteristic function and the density function are available, it is shown through the simulation examples that the proposed empirical maximum likelihood method does not suffer from significant information loss comparing to the maximum likelihood estimation method. The performances of the proposed method is tested for both equally-spaced and unequally-spaced observations.

Appendix

1.1 Integration via adaptive Gauss Hermite (AGH) quadrature method

In this appendix, adaptive Gauss Hermite (AGH) quadrature method (see Liu and Pierce 1994) is described for the evaluation of the integral in (18). Denote by $\exp [-G( {b})]$ the integrand. This means that

$$\begin{aligned} G(b)=-\sum ^q_{\ell =1}\mathcal {L}(\varvec{\theta }_{i\ell }(b);Y_{i\ell })-\log f(b). \end{aligned}$$

The gradient and Hessian of $\mathcal {L}(\varvec{\theta }_{i\ell }(b);Y_{i\ell })$ will be given later on. Choose $b_0 = \arg max_{b} \left[ -G( b) \right]$. Let $y=(b- b_0)\sqrt{ G{''}(b_0 )/2 }$ and $\sigma =\sqrt{ 1/G{''}(b_0 )}$. Then, $b =b_0+ \sqrt{2} \sigma \cdot y$ and $db = \sqrt{2} \sigma dy$. Consider

$$\begin{aligned}&\int exp \left( -G(b) \right) db \nonumber \\&\quad = exp (-G(b_0)) \int \{ {exp \left( - G(b)+G(b_0) + \frac{1}{2} G^{''}(b_0)(b -b_0)^2 \right) } \nonumber \\&\qquad \cdot {exp \left( -\frac{1}{2} G^{''}(b_0)(b -b_0)^2 \right) \}} db \nonumber \\&\quad = exp (-G(b_0)) \sqrt{2} \sigma \cdot \int {exp \left( - G(b_0 + \sqrt{2}\sigma y)+G(b_0) + y^2 \right) }\cdot {exp ( -y^2 )}d{y} \nonumber \\&\quad \simeq \sqrt{2} \sigma \cdot \sum _{s=1}^r \omega _j exp ( y_s^2 ) \cdot exp \left( - G(b_0 + \sqrt{2}\sigma y_s) -y_s^2 \right) . \end{aligned}$$

(A.1)

On the right-hand-side of (A.1), the values of $w_s$ and $y_s$ for $j=1,2,\dots ,r$ can be obtained easily using R package fastGHQuad. In the computation, $r=4$ nodes are used.

Below, the gradient and Hessian of $\mathcal {L}(\varvec{\theta }_{i\ell }(b);Y_{i\ell })$ are given. For brevity, we suppress the notation of i and $\ell$ and write $\mathcal {H}_j(b)$ for $\mathcal {H}(U,\frac{\Delta Y_{ij\ell }}{\Delta T_{ij}},\varvec{\theta }_{i\ell }(b))\,.$ The notation $\partial$ refers to the first order derivative with respect to a component of b and $\partial ^2$ refers to the second order derivative with respect to one or two components of b. Then, $\mathcal {L}(\varvec{\theta }_{i\ell }(b);Y_{i\ell })$ can be rewritten as

$$\begin{aligned} \mathcal {L}(\varvec{\theta }(b);Y) = \sum _{j=1}^n {\frac{\Delta T_j}{T}} \cdot log \left[ 1+\mathcal {H}^T_j(b) \cdot \eta (b) \right] . \end{aligned}$$

(A.2)

Here, $\eta (b)$ is obtained from Eq. (15) and thus fulfills

$$\begin{aligned} \sum _{i=1}^n {\frac{\Delta T_j}{T}} \cdot \left[ \frac{\mathcal {H}_j (b)}{ 1+\mathcal {H}_j^T(b) \cdot \eta (b)} \right] =0. \end{aligned}$$

(A.3)

Then, the first order derivative is

$$\begin{aligned} \partial \mathcal {L}(\varvec{\theta }(b);Y) = \left[ \sum _{j=1}^n {\frac{\Delta T_j}{T}} \cdot \frac{\partial \mathcal {H}_j(b)}{ 1+\mathcal {H}_j^T(b) \cdot \eta (b)} \right] \cdot \eta (b) \end{aligned}$$

The second order derivative is

$$\begin{aligned} \partial ^2 \mathcal {L}(\varvec{\theta }(b);Y)&=\left( \sum _{j=1}^n {\frac{\Delta T_j}{T}} \cdot \frac{\partial \mathcal {H}^T_j (b)}{ 1+\mathcal {H}_j^T(b) \cdot \eta (b)} \right) \cdot \partial \eta (b)\\&\quad + \left( \sum _{j=1}^n {\frac{\Delta T_j}{T}} \cdot \frac{\partial ^2 \mathcal {H}_j^T (b)}{ 1+\mathcal {H}_j^T(b) \cdot \eta (b)} \right) \cdot \eta (b) \\&\quad - \eta ^T(b) \cdot \left( \sum _{j=1}^n {\frac{\Delta T_j}{T}} \cdot \frac{ \partial \mathcal {H}_j(b)\partial ^T \mathcal {H}_j(b)}{ (1+\mathcal {H}_j^T(b) \cdot \eta (b))^2} \right) \cdot \eta (b) \\&\quad - [\partial \eta (b)]^T \cdot \left( \sum _{j=1}^n {\frac{\Delta T_j}{T}} \cdot \frac{\mathcal {H}_j(b) \partial ^T \mathcal {H}_j (b))}{ (1+\mathcal {H}_j^T(b) \cdot \eta (b))^2} \right) \cdot \eta (b) \end{aligned}$$

Differentiating (A.3) with respect to b yields

$$\begin{aligned}&\left( \sum _{j=1}^n {\frac{\Delta T_i}{T}} \cdot \frac{\mathcal {H}_j (b)\cdot \mathcal {H}_j^T (b)}{ (1+\mathcal {H}_j^T(b) \cdot \eta (b))^2}\right) \partial \eta (b)\\&\quad = \left\{ \sum _{j=1}^n {\frac{\Delta T_j}{T}} \cdot \left[ \frac{\partial \mathcal {H}_j (b)}{ 1+\mathcal {H}_j^T(b) \cdot \eta (b)} \right] - \sum _{j=1}^n {\frac{\Delta T_j}{T}} \cdot \frac{\mathcal {H}_j(b) \cdot \partial \mathcal {H}_j (b)^T}{ (1+\mathcal {H}_j^T(b) \cdot \eta (b))^2} \cdot \eta (b)\right\} . \end{aligned}$$

$\partial \eta (b)$ can therefore be solved from the above system of linear equations.