Parameter Estimation of Linear Stochastic Differential Equations with Sparse Observations

Abstract

We consider parameter estimation for linear stochastic differential equations with independent experiments observed at infrequent and irregularly spaced follow-up times. The maximum likelihood method is used to obtain an asymptotically consistent estimator. A kernel-weighted score function is proposed for the parameter in drift terms. The strong consistency and the rate of convergence of the estimator are obtained. The numerical results show that the proposed estimator performs well with moderate sample sizes.

1. Introduction

To simulate the dynamic behavior of a complex system, linear stochastic differential equations (LSDEs) are frequently used. In many real-world applications, it is customary that the parameters that define the system must be estimated from the data. As an example, geometric Brownian motion (GBM) is one of the most popular stochastic processes and undoubtedly an effective instrument in modeling and predicting random changes in stock prices [1,2]. Both deterministic and stochastic components contribute to the pharmacokinetic and pharmacodynamic models: Although there are predictable trends in drug concentrations, it is not always possible to establish the precise concentration at any particular time [3].

In biometrics, a GBM model and an estimation procedure were developed for predicting the height growth of even-aged forest stands as part of a methodology for modeling growth in forest plantations [4].

Due to its growing use in a variety of domains, parameter estimation problems involving stochastic differential equations have received a lot of attention lately. Using the data, one should estimate the parameters that characterize the system. Several methods are proposed to evaluate the parameters, such as the least squares method [5,6,7,8,9], the maximum likelihood method [10,11,12,13,14], and the numerical approximation approach [15]. Several other methods, such as the generalized method of moments procedures [16,17], local linearization method [11,18], and MCMC methods [19], are also proposed.

Assume that n identical and independently distributed paths are observed. When a number of patients can be watched, for example, this situation can arise in pharmacokinetics. For each patient, a bolus of the medication is given, and the “path” of its diffusion through the body can be observed [3]. Such observations are typically sparse and only observed at infrequent and irregularly spaced follow-up times; the above methods are no longer applicable. In this case, we develop a computationally efficient method to deal with the observations with infrequent and irregularly spaced follow-up times. In this paper, we apply kernel methods to the parameter estimation of LSDEs. At the heart of the proposed approach is to “smooth” an individual’s contributions to the likelihood based on the distance of their observed time to the time of interest. The smoothing methods employed are where smoothing happens on an individual basis as compared to the population level, where all individuals are given the same weights. With a suitable choice of bandwidth, the consistency and asymptotic normality for the proposed estimator can be obtained. One can refer to [20,21] for statistical inference of diffusion processes.

In future research, we may consider parameter estimation for a more generalized drift term such as $f(X_t,\theta )$ and the nonparametric estimation for the drift function $f\left(X_t\right)$ where $f(·)$ is a measurable function based on identical and independently distributed observations.

Our paper is organized as follows. In Section 2, we propose an estimator for the drift parameter of the LSDE; we also obtain the consistency of $\mu$ and $\sigma$ and show their asymptotic normality where the convergence rate is ${\left(n{h}_n\right)}^{1/2}$, which is slower than $n^{1/2}$. In Section 3, numerical simulations are performed. Simulation findings show that the large sample approximations are suitable for usage in practice. Section 4 is the conclusion.

2. Models and Methods

2.1. Description of Models

Consider an LSDE model as follows

$$\left\{\begin{array}{c}\mathrm{d}S\left(t\right)=\mu S\left(t\right)\mathrm{d}t+\sigma S\left(t\right)\mathrm{d}W\left(t\right),\hfill \\ S\left(0\right)=s_0,\hfill \end{array}\right.$$

(1)

where $\mu$ is an unknown parameter, $\sigma$ is a constant, and $W\left(t\right)$ are independent standard Brownian motions. $W\left(t\right)$ is characterized by the following properties: (1) $W\left(0\right)=0$; (2) $W\left(t\right)$ has independent increments, which is, for every $t>0$, the future increments $W(t+u)-W\left(t\right)$, $u⩾0$, are independent of the past values $W\left(s\right),s⩽t$, and $W(t+u)-W\left(t\right)\sim \mathcal{N}(0,u)$; (3) $W\left(t\right)$ is continuous in t. Under the condition of Lipschitz and linear growth, the LSDE (1) has a unique strong solution $S\left(t\right)$,

$$S\left(t\right)=s_{0}exp\left\{\left(\mu -\frac{1}{2}{\sigma }^2\right)t+\sigma W_j\left(t\right)\right\}.$$

(2)

Let $\left\{S_i\left(t\right),s_{0i},{\mu }_i,{\sigma }_i;i=1,\cdots ,n\right\}$ be n dependent copies of $\{S\left(t\right),s_0,\mu ,\sigma \}$. As is known to all (see, e.g., [22]), for given t, $S_i\left(t\right),i=1,\cdots ,n$ follow a lognormal distribution

$$log{S}_i\left(t\right)\sim \mathcal{N}\left(log{s}_0+\mu t-\frac{1}{2}{\sigma }^{2}t,{\sigma }^{2}t\right).$$

(3)

We are aiming to use the observations to estimate $\mu$, where the observations consist of $\mathit{y}=\{y_i\left(t_{ik}\right);i=1,\cdots ,n;k=1,\cdots ,d_i\}$, and $d_i<\infty$. The probability of the ith subject at time point $t_{ik}$ is

$$p\left(y_i\left(t_{ik}\right);\mu \right)=\frac{1}{\sqrt{2\pi {\sigma }^2{t}_{ik}}}exp\left\{-\frac{(log{y}_i\left(t_{ik}\right)-\varphi {\left(t_{ik}\right)}^2}{2{\sigma }^2{t}_{ik}}\right\},$$

where $\varphi \left(t\right)=log{s}_0+\mu t-\frac{1}{2}{\sigma }^{2}t$.

If the observations are continuous, consider a time point $t^{*}$. We know that $y_i^{*}=y\left(t^{*}\right)$ are independent and identically distributed variables on the same probability space. Then we get the likelihood function

$$\begin{array}{cc}\hfill L_n({\mathit{y}}^{*};\mu )& =\prod _{i=1}^{n}p(y_i^{*};\mu )\hfill \\ \hfill & =\prod _{i=1}^n\frac{1}{\sqrt{2\pi {\sigma }_i^2{t}^{*}}}exp\left\{-\frac{(log{y}_i^{*}-\varphi {\left(t^{*}\right)}^2}{2{\sigma }_i^2{t}^{*}}\right\},\hfill \end{array}$$

where ${\mathit{y}}^{*}=\left(y_1^{*},y_2^{*},\cdots ,y_n^{*}\right)$. The log-likelihood function is

$$l_n({\mathit{y}}^{*};\mu )=\frac{1}{n}\sum _{i=1}^n\left[-\frac{{(log{y}_i^{*}-\varphi \left(t^{*}\right))}^2}{2{\left({\Sigma }_i\left(t^{*}\right)\right)}^2}-log\sqrt{2\pi }-log{\Sigma }_i\left(t^{*}\right)\right],$$

where ${\Sigma }_i\left(t^{*}\right)=\sqrt{{\sigma }_i^2{t}^{*}}$ and the score function is

$$\begin{array}{cc}\hfill B_n({\mathit{y}}^{*},\mu )=& \frac{\partial l_n({\mathit{y}}^{*};\mu )}{\partial \mu }.\hfill \end{array}$$

(4)

2.2. Kernel Estimation with Forward and Lagged Observation

The data $y_i\left(t\right),i=1,\cdots ,n$ are usually not observed continuously, and it is almost impossible for each individual to be observed at $t^{*}$. Hence $l_n({\mathit{y}}^{*},\theta )$ is not computable from the observations. We propose a method that formalizes the forwarding and lagging strategy, with kernel weighting enabling the use of all available forward and lagged observations. We “smooth” the observations’ contribution to the likelihood based on the distance of their observation time to the time of interest. If data continues to be collected on subjects for which observation has occurred, as in the case of the recurrent event, we use the kernel to impute missing values using both forward and backward-lagged observations. We construct a smoothed log-likelihood function by using kernel estimation

$$\begin{array}{cc}\hfill & l_n({\mathit{y}}^{*};\mu )=-\sum _{i=1}^n\int \left[\frac{{\left(K_{h_n}(s-t_k)log{y}_i\left(s\right)-{\varphi }_i\left(t^{*}\right)\right)}^2}{2{\Sigma }_i^2\left(t^{*}\right)}+log\sqrt{2\pi }+log{\Sigma }_i^2\left(t^{*}\right)\right]\mathrm{d}{N}_i\left(s\right),\hfill \end{array}$$

(5)

where ${\Sigma }_i^2\left(t^{*}\right)<\infty$ is the variance of $y_i^{*}$, $K_{h_n}\left(t\right)=K((t-t^{*})/h_n)/h_n$, $h_n$ is the bandwidth, and the kernel function $K\left(t\right)$ is a symmetric probability density with support $[-1,1]$ and mean 0 that bound the first derivative. In addition, $E\left[\mathrm{d}N\right(t\left)\right]=\lambda \left(t\right)\mathrm{d}t$, where $\lambda \left(t\right)$ is twice continuously differentiable and strictly positive for $t\in [0,T]$. The scoring function is given by

$$U_n\left(\mu \right)=\frac{1}{n}\sum _{i=1}^n\int \frac{K_{h_n}(t^{*}-t)log{y}_i\left(t\right)-(log{s}_0+(\mu -\frac{1}{2}{\Sigma }_i^2\left(t^{*}\right))t^{*})}{{\Sigma }_i^2\left(t^{*}\right)}\mathrm{d}{N}_i\left(t\right).$$

(6)

Assume that the following conditions hold:

(A.1)

${\Theta }_{{\mu }_0}$ is an open sets of $\mathbb{R}$, and ${\Theta }_{{\mu }_0}=\{\mu :|\mu -{\mu }_0|<\rho \}$ for some $\rho >0$ and ${\mu }_0$ is the true parameter.

(A.2)

${\lambda }^{*}\left(t\right)$ is twice continuously differentiable.

(A.3)

$K\left(z\right)$ is a symmetric density function satisfying ${\int }_{-\infty }^{\infty }K{\left(z\right)}^{2}dz<\infty .$ In addition, $h_n\to 0$, $n{h}_n\to \infty$, $n{h}_n^5\to 0$.

Condition (A.1) is a usual assumption for the proof of consistency, and condition (A.2) ensures the Taylor expansion of the score function to the second order. Our methods depend on a proper choice of bandwidth, which is shown in condition (A.3). The estimator ${\widehat{\mu }}_n$ is obtained based on solving Equation (5) with a kernel bandwidth selected to obtain the consistency.

Lemma 1.

Under conditions (A.1)–(A.3), we have

$$\mathrm{E}\left[{\left(n{h}_n\right)}^{\frac{1}{2}}{U}_n\left({\mu }_0\right)\right]=0,$$

as$n\to \infty$.

Proof of Lemma 1.

From the smoothed likelihood function (5), we have the smoothed scoring function

$$\begin{array}{cc}\hfill & {\left(n{h}_n\right)}^{\frac{1}{2}}{U}_n\left({\mu }_0\right)\hfill \\ \hfill =& {\left(n{h}_n\right)}^{\frac{1}{2}}\frac{1}{n}\sum _{i=1}^n\left[\int \frac{K_{h_n}(t^{*}-t)log{y}_i\left(t\right)-(log{s}_0+(\mu -\frac{1}{2}{\Sigma }_i^2\left(t^{*}\right))t^{*})}{{\Sigma }_i^2\left(t^{*}\right)}\mathrm{d}{N}_i\left(t\right)\right]\hfill \\ \hfill =& {\left(n{h}_n\right)}^{\frac{1}{2}}\frac{1}{n}\sum _{i=1}^n\left[\int \frac{\frac{1}{h_n}K\left(\frac{t^{*}-t}{h_n}\right)log{y}_i\left(t\right)-(log{s}_0+(\mu -\frac{1}{2}{\Sigma }_i^2\left(t^{*}\right))t^{*})}{{\Sigma }_i^2\left(t^{*}\right)}\mathrm{d}{N}_i\left(t\right)\right]\hfill \\ \hfill =& I.\hfill \end{array}$$

Let $z=\frac{t^{*}-t}{h_n}$. We have

$$\begin{array}{c}\hfill I=h_n^{\frac{1}{2}}{n}^{-\frac{1}{2}}\sum _{i=1}^n\left[\int \frac{K\left(z\right)log{y}_i(t^{*}-h_{n}z)-(log{s}_0+(\mu -\frac{1}{2}{\Sigma }_i^2\left(t^{*}\right))t^{*})}{{\Sigma }_i^2\left(t^{*}\right)}\mathrm{d}{N}_i\left(z\right)\right].\end{array}$$

Define $F(r,s)=\mathrm{E}[y(r-s)-(log{s}_0+(\mu -{\Sigma }_i^2\left(t^{*}\right)/2)r)]/{\Sigma }_i^2\left(t^{*}\right)$. Obviously, $F(t^{*},0)=0$. By taking expectations together with Taylor expansion, ${\int }_{-\infty }^{\infty }K\left(z\right)\mathrm{d}z=1$ and ${\int }_{-\infty }^{\infty }zK\left(z\right)\mathrm{d}z=0$, we have

$$\begin{array}{cc}& \mathrm{E}\left[I\right]\hfill \\ \hfill =& h_n^{\frac{1}{2}}{n}^{-\frac{1}{2}}\sum _{i=1}^n\left[\int \frac{K\left(z\right)log{y}_i(t^{*}-h_{n}z)-(log{s}_0+(\mu -\frac{{\Sigma }_i^2\left(t^{*}\right)}{2})t^{*})}{{\Sigma }_i^2\left(t^{*}\right)}\right.\hfill \\ & \left.\lambda (t^{*}-h_{n}z)\mathrm{d}\left(z\right)\right]\hfill \\ \hfill =& h_n^{\frac{1}{2}}{n}^{\frac{1}{2}}\frac{1}{\sigma }\left(-\frac{\partial F}{\partial s}{|}_{s=0}{h}_n\frac{\partial \lambda }{\partial t}{{|}_{t=t^{*}{h}_n}+\frac{1}{2}\frac{{\partial }^{2}F}{\partial s^2}|}_{s=0}{h}_n^2\right)\hfill \\ \hfill =& o\left(n^{\frac{1}{2}}{h}_n^{\frac{5}{2}}\right).\hfill \end{array}$$

From condition (A.3), we have $\mathrm{E}\left[I\right]=o\left(1\right)$. □

The following theorem shows the consistency of the proposed estimator ${\widehat{\mu }}_n$ obtained based on solving Equation (5).

Theorem 1.

Under conditions (A.1)–(A.3), ${\widehat{\mu }}_n$ admits the consistency as $n\to \infty$.

Proof of Theorem 1.

Solving Equation (5), we have

$${\widehat{\mu }}_n{t}^{*}=\frac{1}{n{h}_n}\sum _{i=1}^{n{h}_n}\frac{\int K_{h_n}(t-t^{*})logy\left(t\right)\mathrm{d}{N}_i\left(t\right)}{\int K_{h_n}(t-t^{*})\mathrm{d}{N}_i\left(t\right)}-log{s}_0+\frac{1}{2}{\sigma }^2{t}^{*}.$$

By properties of the kernel function $K(·)$, we have

$$\begin{array}{cc}\hfill |{\widehat{\mu }}_n{t}^{*}-\mu t^{*}|=& \left|\frac{1}{n{h}_n}\sum _{i=1}^{n{h}_n}\frac{\int K_{h_n}(t-t^{*})(logy(t)-log{s}_0-(\mu -\frac{{\sigma }^2}{2})t^{*})\mathrm{d}{N}_i(t)}{\int K_{h_n}(t-t^{*})\mathrm{d}{N}_i(t)}\right|\hfill \\ \hfill \le & \left|\frac{1}{n{h}_n}\sum _{i=1}^{n{h}_n}(logy(b_i)-log{s}_0-(\mu -\frac{{\sigma }^2}{2})t^{*})\right|,\hfill \end{array}$$

where $b_i\in [t^{*}-h_n,t^{*}+h_n]$ is some constant. By Equation (3), we have

$$logy(t)\sim \mathcal{N}(log{s}_0+(\mu -\frac{{\sigma }^2}{2})t,{\sigma }^{2}t).$$

Hence, we have

$$logy(b_i)-log{s}_0-(\mu -\frac{{\sigma }^2}{2})t^{*}\sim \mathcal{N}(a,{\sigma }^{2}t),$$

where $a=O\left(h_n\right)$ is some constant. By the Wiener–Khinchin law of large numbers, we have

$$|{\widehat{\mu }}_n{t}^{*}-\mu t^{*}|=O\left(h_n\right),$$

which goes to zero, as $n\to \infty$. □

The following theorem shows the asymptotic normality of ${\widehat{\mu }}_n$.

Theorem 2.

Assume conditions (A.1)–(A.3) hold, ${\widehat{\mu }}_n$ is consistent, and the asymptotic distribution of ${\widehat{\mu }}_n$ satisfies

$${\left(n{h}_n\right)}^{\frac{1}{2}}({\widehat{\mu }}_n-\mu )\sim N\left(0,{\left(C_n\left({\mu }_0\right)\right)}^2\Sigma \left({\mu }_0\right)\right),$$

as $n\to \infty$, where

$$C_n\left({\mu }_0\right)=-{\left({\int }_0^1\frac{\partial U_n}{\partial \mu }({\mu }_0+\lambda (\mu -{\mu }_0))\mathrm{d}\lambda \right)}^{-1},$$

and

$$\Gamma \left({\mu }_0\right)=\int {\left(K\left(z\right)\right)}^2\frac{1}{{\sigma }^2}{t}^{*}\lambda (t^{*}-h_{n}z)\mathrm{d}z.$$

Proof of Theorem 2.

Let $\left\{{\mu }_n\right\}$ be a strongly consistent sequence of $\mu$, i.e., ${\mu }_n\stackrel{a.s.}{\to }{\mu }_0$, and $\left\{{\mu }_n\right\}\in {\Theta }_{{\mu }_0}$. We can seek a solution ${\widehat{\mu }}_n$ of the log-likelihood function $l_n({\mathit{y}}^{*},{\widehat{\mu }}_n)$, and ${\widehat{\mu }}_n$ is a strongly consistent sequence. Note that

$$U_n(\mu )=\frac{1}{n}\sum _{i=1}^n\int \frac{K_{h_n}(t^{*}-t)y_i(t)-(log{s}_0+(\mu -\frac{1}{2}{\Sigma }_i{(t^{*})}^2)t^{*})}{{\Sigma }_i^2(t^{*})}\mathrm{d}{N}_i(t),$$

we denote $U_n\left(\mu \right)={\sum }_{i=1}^n\psi (y_i^{*},\mu )$. Expand $U_n$ as

$$\begin{array}{cc}\hfill U_n\left(\mu \right)=& U_n\left({\mu }_0\right)+{\int }_{{\mu }_0}^{\mu }\frac{\partial U_n\left(u\right)}{\partial \mu }\mathrm{d}u\hfill \\ \hfill =& U_n\left({\mu }_0\right)+{\int }_0^1\frac{\partial U_n}{\partial \mu }\left({\mu }_0+\lambda (\mu -{\mu }_0)\right)\mathrm{d}\lambda (\mu -{\mu }_0).\hfill \end{array}$$

Let $\mu ={\widehat{\mu }}_n$, we have $U_n\left({\widehat{\mu }}_n\right)=0$. Then

$${\widehat{\mu }}_n-\mu =-{\left({\int }_0^1\frac{\partial U_n}{\partial \mu }({\mu }_0+\lambda (\mu -{\mu }_0))\mathrm{d}\lambda \right)}^{-1}{U}_n\left({\mu }_0\right).$$

(7)

Multiply both sides of Equation (7) by ${\left(n{h}_n\right)}^{\frac{1}{2}}$, we denote

$$C_n\left({\mu }_0\right)=-{\left({\int }_0^1\frac{\partial U_n}{\partial \mu }({\mu }_0+\lambda (\mu -{\mu }_0))\mathrm{d}\lambda \right)}^{-1}.$$

Then

$${\left(n{h}_n\right)}^{\frac{1}{2}}({\widehat{\mu }}_n-\mu )=C_n\left({\mu }_0\right){\left(n{h}_n\right)}^{\frac{1}{2}}{U}_n\left({\mu }_0\right).$$

Hence we give the variance of ${\left(n{h}_n\right)}^{\frac{1}{2}}{U}_n\left({\mu }_0\right)$,

$${\left(n{h}_n\right)}^{\frac{1}{2}}{U}_n\left({\mu }_0\right)=\frac{{\left(n{h}_n\right)}^{\frac{1}{2}}}{n}\sum _{i=1}^n\psi (y_i^{*},{\mu }_0)=n^{\frac{1}{2}}\left(\frac{1}{n}\sum _{i=1}^n{h}_n^{\frac{1}{2}}\psi (y_i^{*},{\mu }_0)\right).$$

From Lemma 1, we have that $\mathrm{E}[h_n^{\frac{1}{2}}\psi (y_i^{*},{\mu }_0)]=0,i=1,\cdots ,n$. By central limit theorem, $\mathrm{var}\left[{\left(n{h}_n\right)}^{\frac{1}{2}}{U}_n\left({\mu }_0\right)\right]=\mathrm{var}\left[h_n^{\frac{1}{2}}\psi (y_i^{*},{\mu }_0)\right]$, and we denote ${\varphi }_i\left(t\right)=y_i\left(t\right)-(log{s}_0+({\mu }_0-\frac{1}{2}{\sigma }^2)t^{*})$. Then

$$\begin{array}{cc}& \mathrm{var}\left[h_n^{\frac{1}{2}}\psi (y_i^{*},{\mu }_0)\right]\hfill \\ \hfill =& \mathrm{E}\left[{(h_n^{\frac{1}{2}}\psi (y_i^{*},{\mu }_0))}^2\right]\hfill \\ \hfill =& h_n\mathrm{E}\left[{\left(\sum _{k=1}^n\frac{K_{h_n}(t^{*}-t_{ik})(y_i\left(t_{ik}\right)-(log{s}_0+({\mu }_0-\frac{1}{2}{\Sigma }_i^2\left(t^{*}\right))t^{*}))}{{\Sigma }_i^2\left(t^{*}\right)}\right)}^2\right]\hfill \\ \hfill =& h_n\frac{1}{{\sigma }^4}\left({\int }_{t_1\ne t_2}{K}_{h_n}(t^{*}-t_1)K_{h_n}(t^{*}-t_2)\mathrm{E}\left[{\varphi }_i\left(t_1\right){\varphi }_i\left(t_2\right)\right]\mathrm{E}\left[\mathrm{d}N\left(t_1\right)\mathrm{d}N\left(t_2\right)\right]\right.\hfill \\ & +\left.{\int }_{t_1=t_2}{\left(K_{h_n}(t^{*}-t_1)\right)}^2\lambda \left(t_1\right)\mathrm{E}\left[{\varphi }_i{\left(t_1\right)}^2\right]\mathrm{d}{t}_1\right).\hfill \end{array}$$

Assume that for $t_1\ne t_2$, $pr(\mathrm{d}N\left(t_1\right)=1|N\left(t_2\right)-N\left(t_2^{-}\right)=1)=g(t_1,t_2)\mathrm{d}{t}_1$, where $pr$ means the probability, $g(t_1,t_2)$ is continuous for $t_1\ne t_2$, and $g(t_1\pm ,t_2\pm )$ exists. Then

$$\begin{array}{cc}& \mathrm{E}\left[{\left(h_n^{\frac{1}{2}}\psi (y_i^{*},{\mu }_0)\right)}^2\right]\hfill \\ \hfill =& h_n\left(\frac{1}{{\sigma }^4}{\int }_{t_1\ne t_2}{K}_{h_n}(t^{*}-t_1)K_{h_n}(t^{*}-t_2)E\left[{\varphi }_i\left(t_1\right){\varphi }_i\left(t_2\right)\right]g(t_1,t_2)E\left[\mathrm{d}N\left(t_2\right)\right]\mathrm{d}{t}_1\right.\hfill \\ \hfill +& \left.\frac{1}{{\sigma }^4}{\int }_{t_1=t_2}{\left(K_{h_n}(t^{*}-t_1)\right)}^2\lambda \left(t_1\right)E\left[{\varphi }_i{\left(t_1\right)}^2\right]\mathrm{d}{t}_1\right)\hfill \\ \hfill =& h_n(I_1+I_2).\hfill \end{array}$$

Using a change of variables, we have $h_n{I}_1=O\left(h_n\right)$.

With notation $\mathrm{E}\left[{\left({\varphi }_i(t^{*}-h_{n}z)\right)}^2\right]={\sigma }^2{t}^{*}$, we have

$$\begin{array}{cc}\hfill h_n{I}_2=& h_n\frac{1}{{\sigma }^4}\int {\left(h_n\right)}^{-2}{\left(K\left(z\right)\right)}^2\mathrm{E}\left[{\left({\varphi }_i(t^{*}-h_{n}z)\right)}^2\right]\lambda (t^{*}-h_{n}z)h_n\mathrm{d}z\hfill \\ \hfill =& \int {\left(K\left(z\right)\right)}^2\frac{1}{{\sigma }^2}{t}^{*}\lambda (t^{*}-h_{n}z)\mathrm{d}z.\hfill \end{array}$$

From the Lyapunov central limit theorem, we have that ${\left(n{h}_n\right)}^{\frac{1}{2}}{U}_n\left({\mu }_0\right)$ converges to a continuous Gaussian process $Z\sim N\left(0,\Gamma \left({\mu }_0\right)\right)$. Hence, we have

$${\left(n{h}_n\right)}^{\frac{1}{2}}({\widehat{\mu }}_n-\mu )=C_n\left({\mu }_0\right){\left(n{h}_n\right)}^{\frac{1}{2}}{U}_n\left({\mu }_0\right)\sim N\left(0,{\left(C_n\left({\mu }_0\right)\right)}^2\Gamma \left({\mu }_0\right)\right).$$

□

Remark 1.

When there are several observations $d_i>1$, one can estimate the drift parameter μ by a standard maximum likelihood method. Let

$$z_{ik}=log\left(\frac{y_i\left(t_{ik+1}\right)}{y_i\left(t_{ik}\right)}\right),$$

for $k=1,\cdots ,d_i$ and $i=1,2,\cdots ,n$. Thus, conditional on the observation times,

$$z_{ik}\sim \mathcal{N}\left(\left(\mu -{\sigma }^2/2\right)\left(t_{ik+1}-t_{ik}\right),{\sigma }^2\left(t_{ik+1}-t_{ik}\right)\right)$$

and they are independent. For example, if we reparameterize μ as $\nu =\mu -{\sigma }^2/2$, the MLE for ν would be given by

$$\widehat{\nu }=\frac{{\sum }_{i=1}^n{\sum }_{j=1}^{d_i-1}{z}_j^{\left(i\right)}/{\Delta }_j^{\left(i\right)}}{{\sum }_{i=1}^n(d_i-1)},$$

(8)

where ${\Delta }_k^{\left(i\right)}=t_{ik+1}-t_{ik}$. Then μ is estimated by $\widehat{\mu }=\widehat{\nu }+{\sigma }^2/2$.

Remark 2.

When there is only one observation $d_i=1$, the estimator proposed in Equation (8) is not effective. Our estimator performs reasonably well in this extreme case, and we could give an explicit asymptotic variance for our estimator, which is ${\sigma }^2{t}^{*}$.

3. Simulation

In this section, utilizing both forward and backward-lagged observation, we examine the kernel estimator. We generate 1000 datasets, and each dataset consists of $n=100,400,900$ subjects with different bandwidths (BD). The process is generated through model (1); we set the initial condition $s_0=1$, $\mu =2$ and $\sigma =0.02$. Then the solution is

$$S\left(t\right)=s_{0}exp\{1.9998t+0.02w\left(t\right)\},$$

(9)

where $w\left(t\right)$ is a standard Brownian motion. The number of observation times for each subject is a Poisson distributed with an intensity rate of 5. The time points of each individual’s observation are generated from a uniform distribution, Unif$(0,1)$. The outcomes from other models’ parameter selections are not mentioned because they are essentially identical. All simulations were performed on a laptop running R 4.2.9 with 8G of RAM.

Based on Theorems 1 and 2, we obtain a kernel estimator with asymptotically negligible bias and employ bandwidths in the range $(n^{-1},n^{-\frac{1}{5}})$ when calculating ${\widehat{\mu }}_n$ using the smoothed likelihood score function (6). The kernel function we choose is the Epanechnikov kernel, which is $K\left(x\right)=\frac{3}{4}{(1-x^2)}_{+}$. The usage of additional kernel functions has little effect on the estimator’s empirical performance, according to additional simulations (not published).

The simulation results show that the estimates for the parameter in the model are accurate.

Table 1. Simulation results with different n and bandwidths.

n	BD	$\widehat{\mathit{\mu }}$	Bias ($\widehat{\mathit{\mu }}$)	RB ($\widehat{\mathit{\mu }}$)
100	$n^{-0.25}$	2.095	0.095	0.048
	$n^{-0.3}$	2.065	0.065	0.033
400	$n^{-0.25}$	2.052	0.052	0.026
	$n^{-0.3}$	2.030	0.030	0.015
900	$n^{-0.25}$	2.020	0.020	0.010
	$n^{-0.3}$	2.020	0.020	0.010

summarizes the main findings from over 1000 simulations. We note that the bias diminishes and is minor as the sample size grows. The performance improves the larger sample sizes and smaller bandwidths. The overall parameter estimates are evaluated by the bias and relative bias (RB), which are defined as

$$Bias\left(\widehat{\mu }\right)={\mu }_0-\widehat{\mu }, RB\left(\widehat{\mu }\right)=\frac{\left|Bias\left(\widehat{\mu }\right)\right|}{\left|{\mu }_0\right|},$$

where ${\theta }_0$ denotes the true parameter.

4. Conclusions

In this paper, we have presented kernel-weighting methods for the estimation of the LSDE model (1) in repeatable experiments when the observation time is a random variable, and the number of observations of each individual is uncertain or even sparse. This is a real improvement because the past literature usually supposed that observation intervals are equally spaced and could not deal with the sparse observations. We consider the maximum likelihood estimation of the drift parameter. This method has some assumptions, and we give the asymptotic normality of the proposed estimator. In numerical studies, we set the true parameter ${\mu }_0=2$, ${\sigma }_0=0.02$, and the initial condition $s_0=1$ for each individual with sparse observations (the frequency of observation follows a Poisson distribution with mean 5). Using the smoothed scoring function, we obtain the estimation of the drift parameter.