Regression estimation for continuous-time functional data processes with missing at random response

Abstract

In this paper, we are interested in nonparametric kernel estimation of a generalised regression function based on an incomplete sample $(X_t,Y_t,{\zeta }_t{)}_{t\in [0,T]}$ copies of a continuous-time stationary and ergodic process $(X,Y,\zeta )$. The predictor X is valued in some infinite-dimensional space, whereas the real-valued process Y is observed when the Bernoulli process $\zeta =1$ and missing whenever $\zeta =0$. Uniform almost sure consistency rate as well as the evaluation of the conditional bias and asymptotic mean square error are established. The asymptotic distribution of the estimator is provided with a discussion on its use in building asymptotic confidence intervals. To illustrate the performance of the proposed estimator, a first simulation is performed to compare the efficiency of discrete-time and continuous-time estimators. A second simulation is conducted to discuss the selection of the optimal sampling mesh in the continuous-time case. Then, a third simulation is considered to build asymptotic confidence intervals. An application to financial time series is used to study the performance of the proposed estimator in terms of point and interval prediction of the IBM asset price log-returns. Finally, a second application is introduced to discuss the usage of the initial estimator to impute missing household-level peak electricity demand.

Keywords:

• Asymptotic mean square error
• continuous-time ergodic processes
• confidence intervals
• functional data
• generalised regression
• missing at random

1. Introduction

Let $(\mathcal{E},d)$ be an infinite-dimensional space equipped with a semi-metric $d(\cdot ,\cdot )$. Consider $(X_t,Y_t,{\zeta }_t{)}_{t\in {\mathbb{R}}^{+}}$ a stationary and ergodic continuous-time process valued in the space $\mathcal{X}:=\mathcal{E}\times \mathbb{R}\times \{0,1\}$, such that each triplet $(X_t,Y_t,{\zeta }_t)$ has the same probability distribution as the random variable (r.v.) $(X,Y,\zeta )$ defined on probability space $(\mathrm{\Omega },\mathcal{F},\mathbb{P})$. Let S be a compact interval in $\mathbb{R}$ and ψ be a measurable function defined on the space $S\times \mathbb{R}$, $(y,Y)\mapsto \psi (y,Y)$, where y is a real variable such that $\mathbb{E}(|\psi (Y,y)|)<\mathrm{\infty }$. Consider the following regression model: (1) $\psi (Y,y)=m_{\psi }(X,y)+\epsilon ,$(1) where $m_{\psi }(X,y)$ is the conditional expectation of $\psi (Y,y)$ given the r.v. X. That is, for any $y\in S$ and a fixed $x\in \mathcal{E}$, $\mathbb{E}(\psi (Y,y)|X=x)=m_{\psi }(x,y)$. The error term ε is independent of X such $\mathbb{E}(\epsilon |X)=0$ almost surely (a.s.)

Usually, when no data are missing, a sample of stationary and ergodic process $(X_t,Y_t{)}_{0\le t\le T}$ is observed. Here, we allow the response variable $Y_t$ to be Missing At Random (MAR) any time t. To check whether an observation is complete or missing, a new variable ζ is introduced into the model as an indicator of the missing observations. Thus, for any $t\in [0,T]$, ${\zeta }_t=1$ if $Y_t$ is observed and 0 if $Y_t$ is missing. We suppose that the Bernoulli random variable ζ satisfies $\mathbb{P}(\zeta =1|X=x,Y=y)=\mathbb{P}(\zeta =1|X=x):=p(x)$. Here, $p(x)>0$ is the conditional probability of observing the response variable and is usually unknown. This assumption allows to conclude that ζ and Y are conditionally independent given X. Note that the above assumption says that the response variable does not provide additional information, on top of that given by the explanatory variable, to predict whether an individual will present a missing response.

In this paper, we are interested in the estimation of the regression function $m_{\psi }(\cdot ,y)$ based on the observed data $(X_t,Y_t,{\zeta }_t{)}_{0\le t\le T}$. Note that, for any $t\in [0,T]$, $t\mapsto \{X_t(\omega ),\omega \in \mathrm{\Omega }\}$ is an element in the space $\mathcal{E}$, which means that, for any fixed time $t=t_0$, $X_{t_0}$ is a curve. Specifically, if $\mathcal{E}={\mathcal{C}}_{[0,1]}$ is the space of square integrable functions defined on $[0,1]$, then the predictor $X_{t_0}:=\{X_{t_0}(s):s\in [0,1]\}$ describes a trajectory in the functional space $\mathcal{E}$ observed at the fixed time $t_0$.

In real life, there are several situations where the response variable might be missing at random. For instance, in survey sampling studies the non-response is an increasingly common problem, where the missing response reaches rates of $25\mathrm{\%}$–$30\mathrm{\%}$ or even higher (see, e.g. Sikov 2018). In such cases, the missing data become a real source of bias in survey sampling estimation. Another case where the response may be subject to the MAR phenomena is the household electricity consumption monitoring. Indeed, the real-time collection of intra-day electricity consumption is now possible after the deployment of smart meters at the household level. The transmission of the information from the smart meter towards the information system goes usually through WIFI or optical fibre networks which are significantly dependent on the weather conditions, among other factors. Therefore, a response variable such as the daily total electricity consumption might be subject to missing at random mechanism due to bad weather conditions (for more details , see Section 5.2). In financial market, despite the modern technology, which allows to collect data at a very fine time scale, financial data can still be missing. For instance, there are some regular holidays, such as Thanksgiving Day and Christmas, for which stock price data are missing. There are many other technical reasons (such as breakdown in devises recording data, computers' sudden shutdowns, …) that make stretches of data missing (see Section 5.1 for more details about this application). For further examples and details about missing at random data the reader is referred to Chapter 1 in Little and Rubin (2002).

Whenever $(X_i,Y_i,{\zeta }_i{)}_{1\le i\le n}$ is an independent and identically distributed (i.i.d) random sample several authors investigated nonparametric and semiparametric estimations of the regression function. In the framework where X is a finite dimensional, one can quote (Cheng 1994; Little and Rubin 2002; Nittner 2003; Tsiatis 2006; Liang et al. 2007; Efromovich 2011). See also, Ferraty et al. (2013) when the predictor is infinite-dimensional. However, less attention has been given to the case of dependent data including an infinite-dimensional covariate, except Ling et al. (2015), where a local constant estimation of the regression operator with discrete-time ergodic processes was considered.

In the case of continuous-time finite dimensional processes $(X,Y)\in {\mathbb{R}}^d\times {\mathbb{R}}^{d^{\mathrm{\prime }}}$ (d and $d^{\mathrm{\prime }}\ge 1$) satisfying a strong mixing condition, the estimation of the regression function based on completely observed data was considered by several authors, see for instance, the monograph by Bosq (1998) and the references therein.

Some of these results were extended by Didi and Louani (2014) and Bouzebda and Didi (2017) for stationary and ergodic processes. Chaouch and Laïb (2019) studied the asymptotic mean square error of the kernel regression estimator for MAR stationary and ergodic process and obtained an explicit upper bound of it.

It is worth noting that, even though a continuous-time functional processes framework is considered in this paper, in practice data are often collected according to some sampling scheme and the continuous-time process is discretised. Our results are then valid when considering a discrete-time ergodic stationary process $(X_{t_k},Y_{t_k},{\zeta }_{t_k}{)}_{1\le k\le n}$ sampled from a continuous-time process $\{(X_t,Y_t,{\zeta }_t){\}}_{0\le t\le T}$ with a regular sampling mesh $\delta =T/n.$ A discussion on the optimal choice of the sampling mesh in the continuous-time case will be then of great practical interest (see Section 3.5 and Simulation 2).

In the setting where $(X_t,Y_t{)}_{t\in [0,T]}\in \mathcal{E}\times \mathbb{R}$ is an α-mixing continuous-time process with Y is completely observed, Maillot (2008) established the convergence with rates of the regression operator, and a super-optimal mean square convergence rate was obtained in Chesneau and Maillot (2014).

This paper aims to complete and extend (Maillot 2008; Chesneau and Maillot 2014; Ling et al. 2015) work at several levels. First, we suppose that the continuous-time process satisfies an ergodic assumption rather than an α-mixing one. Therefore, the dependence structure considered here is more general and involves several processes which do not satisfy the mixing property. Indeed, our results are valid for α-mixing and non α-mixing as well as long memory and Bernoulli shift processes (for more details, see examples used to discuss Assumption A3 below). Moreover, our results are stated and proved without assuming neither mixing condition nor imposing a particular covariance structure on the process. This is due to the fact that the main technical tools used here are martingale difference devices and sequence of projections on appropriate σ-fields. Second, we complete and extend results established in Ling et al. (2015) for discrete-time functional data processes, to the continuous-time functional framework.

We estimate a general operator $m_{\psi }(x,y)$ which includes conditional mean, conditional distribution function and conditional quantiles. It is worth noting that such extension is not obvious since it requires an appropriate definition of σ-fields adapted to continuous-time context. Such adaptation is crucial when using martingale difference tools to establish asymptotic properties of the estimator. Third, the response variable considered here is affected by the MAR mechanism and therefore is not completely observed as in Maillot (2008). Moreover, in contrast to Maillot (2008) and Chesneau and Maillot (2014), we do not limit our study to the mean square convergence, but we provide a more exhaustive inference on the regression operator estimator including pointwise and uniform almost sure convergence rate, identification of the limiting distribution of our estimator, and provide method to build confidence intervals. Fourth, simulation study is also carried out to investigate the selection of the ‘optimal’ sampling mesh, which is one of the most important topics in nonparametric estimation with continuous-time processes.

The rest of this paper is organised as follows. In Section 2, we present the framework adapted to continuous-time ergodic processes and introduce assumptions needed for establishing asymptotic results. The main asymptotic properties of the estimator are discussed in Section 3. An illustration of the performance of the proposed estimator is discussed through simulated data in Section 4. Section 5 is devoted to an application of the proposed estimator to financial time series. Section 6 discusses the application of our theoretical results to continuous-time conditional quantiles estimation. Finally technical proofs are given in the Appendix.

2. Framework and assumptions

To define the framework of our study, we need to introduce some definitions. Let $X=(X_t{)}_{t\in [0,\mathrm{\infty })}$ be a continuous-time process defined on a probability space $(\mathrm{\Omega },\mathcal{F},\mathbb{P})$ and observed at any time $t\in [0,T]$. For more details about the definition of continuous-time ergodic processes, the reader is referred to Didi and Louani (2014). From now on, we consider $\{{\mathcal{F}}_t,t\ge 0\}$ the filtration defined on $(\mathrm{\Omega },\mathcal{F})$, that is $\{{\mathcal{F}}_t,t\ge 0\}$ is an increasing sequence of sub-σ-algebras of $\mathcal{F}$.

For a positive real number δ such that $n=\frac{T}{\delta }\in \mathbb{N}$ and $j\in \mathbb{N}\cap [1,n]$, consider the δ-partition $(T_j=\mathrm{j\delta }{)}_{1\le j\le n}$ of the interval $[0,T]$. Furthermore, for t>0 and $1\le j\le n$, we define the following σ-fields: (2) $\begin{array}{rl}{\mathcal{F}}_{t-\delta }& :=\sigma ((X_s,Y_s):0\le s<t-\delta ),{\mathcal{F}}_j:={\mathcal{F}}_{T_j}=\sigma ((X_s,Y_s),0\le s<T_j),\\ {\mathcal{S}}_{t,\delta }& :=\sigma ((X_s,Y_s),(X_r):0\le s<t,t\le r\le t+\delta ).\end{array}$(2) If t<0 we take ${\mathcal{F}}_t$ the trivial σ-field. Note that, for any $\delta >0$ and t>0, we have ${\mathcal{F}}_{t-\delta }\subset {\mathcal{S}}_{t-\delta ,\delta }\subset {\mathcal{S}}_{t,\delta }$. Moreover, for any $j\ge 2$, such that $T_{j-1}\le t\le T_j,$ we have ${\mathcal{F}}_{j-2}\subseteq {\mathcal{F}}_{t-\delta }\subset {\mathcal{S}}_{t,\delta }$.

Let $B(x,u)$ be a ball centred at $x\in \mathcal{E}$ with radius $u>0.$ Denote $D_t:=d(x,X_t)$ a nonnegative real-valued continuous-time process and let $F_x(u)=\mathbb{P}(D_t\le u):=\mathbb{P}(X_t\in B(x,u))$ and $F_x^{{\mathcal{F}}_{t-\delta }}=\mathbb{P}(X_t\in B(x,u)|{\mathcal{F}}_{t-\delta })$ be the distribution function and conditional distribution function of $(D_t{)}_{0\le t\le T}$ given the σ-field ${\mathcal{F}}_{t-\delta }$, respectively.

To define an estimator of the regression function adapted to the MAR, multiply Equation (1(1) $\psi (Y,y)=m_{\psi }(X,y)+\epsilon ,$(1) ) by ζ to get $\mathrm{\zeta \psi }(Y,y)=\zeta m_{\psi }(X,y)+\mathrm{\zeta \epsilon }.$ Taking conditional expectations with respect to X = x, one gets $\mathbb{E}(\mathrm{\zeta \psi }(Y,y)|X=x)=m_{\psi }(x,y)\mathbb{E}(\zeta |X=x).$ Thus we have $m_{\psi }(x,y)={\displaystyle \frac{\mathbb{E}(\mathrm{\zeta \psi }(Y,y)|X=x)}{\mathbb{E}(\zeta |X=x)}}.$ Given a random sample $(X_t,Y_t,{\zeta }_t{)}_{0\le t\le T}$ one can therefore define a kernel-type estimator of $m_{\psi }(x,y)$, say ${\hat{m}}_{\psi ,T}(x,y)$, adapted to the MAR response framework. Note that if there are missing observations in the response variable, a simple way to estimate $m_{\psi }(x,y)$ is to consider a kernel smoothing-type estimator which only considers observed data, in other words, those for which ${\zeta }_t=1$. Therefore, one gets (3) ${\hat{m}}_{\psi ,T}(x,y):=\{\begin{array}{ll}{\displaystyle \frac{{\int }_0^T{\zeta }_t\psi (Y_t,y){\mathrm{\Delta }}_t(x)\mathrm{d}t}{{\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t}},& \text{if}{\displaystyle {\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t\ne 0}\\ {\displaystyle \frac{1}{T}}{\displaystyle {\int }_0^T{\zeta }_t\psi (Y_t,y)\mathrm{d}t,}& \text{otherwise},\end{array}$(3) where ${\mathrm{\Delta }}_t(x)=K(\frac{D_t}{h_T})$, $K(\cdot )$ is a kernel density function, $h_T$ is the smoothing parameter tending to zero as T goes to infinity.

Remark 2.1

When the sample has missing observations in the response variable, two strategies can be followed to estimate $m_{\psi }(x,y)$. The first one, ${\hat{m}}_{\psi ,T}(x,y)$ given in Equation (3(3) ${\hat{m}}_{\psi ,T}(x,y):=\{\begin{array}{ll}{\displaystyle \frac{{\int }_0^T{\zeta }_t\psi (Y_t,y){\mathrm{\Delta }}_t(x)\mathrm{d}t}{{\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t}},& \text{if}{\displaystyle {\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t\ne 0}\\ {\displaystyle \frac{1}{T}}{\displaystyle {\int }_0^T{\zeta }_t\psi (Y_t,y)\mathrm{d}t,}& \text{otherwise},\end{array}$(3) ), called simplified estimator which only uses complete observations. The second approach consists in using the simplified estimator ${\hat{m}}_{\psi ,T}(x,y)$ to impute the missing values of the response variable $Y_t$ according to the following expression: $\tilde{\psi }(Y_t,y):={\zeta }_t\psi (Y_t,y)+(1-{\zeta }_t){\hat{m}}_{\psi ,T}(X_t,y)$, (see, e.g. Chu and Cheng 2003 or González-Manteiga and Pérez-González 2004). Consequently, an estimator, say ${\tilde{m}}_{\psi ,T}(x,y)$, based on imputed data may be defined as follows: ${\tilde{m}}_{\psi ,T}(x,y)=\frac{{\int }_0^T\tilde{\psi }(Y_t,y){\mathrm{\Delta }}_t(x)\mathrm{d}t}{{\int }_0^T{\mathrm{\Delta }}_t(x)\mathrm{d}t}.$

From now on, we set $Z_1(x):={\int }_0^{\delta }{\mathrm{\Delta }}_t(x)\mathrm{d}t$ and define the conditional bias as (4) $B_T(x,y):=\frac{{\overline{m}}_{\psi ,T,2}(x,y)}{{\overline{m}}_{\psi ,T,1}(x)}-m_{\psi }(x,y):=C_T(x,y)-m_{\psi }(x,y),$(4) where ${\overline{m}}_{\psi ,T,1}(x):={\overline{m}}_{\psi ,T,1}(x,1)$ and for $i=1,2,$ (5) ${\overline{m}}_{\psi ,T,i}(x,y^{i-1}):=\frac{1}{n\mathbb{E}(Z_1(x))}{\int }_0^T\mathbb{E}\{{\zeta }_t(\psi (Y_t,y){)}^{i-1}{\mathrm{\Delta }}_t(x)|{\mathcal{F}}_{t-\delta }\}\mathrm{d}t.$(5) Before introducing the assumptions under which we establish our asymptotic results, we add the following notations. Let $o_{a.s.}(u)$ denote a real random function ℓ such that $\ell (u)/u$ converges to zero almost surely (a.s.) as $u\to 0$ and denote ${\mathcal{O}}_{a.s.}(u)$ a real random function ℓ such that $\ell (u)/u$ is almost surely bounded.

1. (Assumptions on the kernel function). Let K be a nonnegative bounded kernel of class ${\mathcal{C}}^1$ over its support $[0,1]$ such that $K(1)>0$. The derivative $K^{\mathrm{\prime }}$ exists on $[0,1)$ and satisfies the condition $K^{\mathrm{\prime }}(v)<0$ for all $v\in [0,1)$ and $|{\int }_0^1(K^j{)}^{\mathrm{\prime }}(v)\mathrm{d}v|<\mathrm{\infty }$ for $j=1,2.$

2. (Assumptions related to the continuous-time functional ergodic processes)

   Let ${\alpha }_0$ be a nonnegative real number and $x\in \mathcal{E}$. Suppose, for any $0\le s<t\le T$ such that $t-s\le {\alpha }_0$, there exists a nonnegative continuous random function $f_{t,s}(x):=f_{X_t,s}(x)$ a.s. bounded by a deterministic function $b_{s,{\alpha }_0}(x)$¹.

   Moreover, let $g_{t,s,x}(\cdot )$ be a random function defined on $\mathbb{R}$, $f(x)$ is a deterministic nonnegative bounded function and $\varphi (\cdot )$ a nonnegative real function tending to zero (as its argument tends to 0), and assume that:

   1. $F_x(u):=\mathbb{P}(d(x,X_t)\le u)=\varphi (u)f(x)+o(\varphi (u))$ as $u\to 0.$

   2. For any $0\le s\le t$, $F_x^{{\mathcal{F}}_s}(u):={\mathbb{P}}^{{\mathcal{F}}_s}(d(x,X_t)\le u)=\mathbb{P}(d(x,X_t)\le u|{\mathcal{F}}_s)=\varphi (u)f_{t,s}(x)+g_{t,s,x}(u)$ with $g_{t,s,x}(u)=o_{a.s.}(\varphi (u))$ as $u\to 0,$ $g_{t,s,x}(u)/\varphi (u)$ a.s. bounded and $T^{-1}{\int }_0^T{g}_{t,t-\delta ,x}(u)\mathrm{d}t=o_{a.s.}(\varphi (u))$ as $T\to \mathrm{\infty }$ and $u\to 0.$

   3. For any $x\in \mathcal{E}$: $\underset{T\to \mathrm{\infty }}{lim}\frac{1}{T}{\int }_0^T{f}_{t,t-\delta }(x)\mathrm{d}t=f(x)$, a.s.

   4. There exists a nondecreasing bounded function ${\tau }_0:[0,1]\to [0,1]$ such that, uniformly in $u\in [0,1]$, $\frac{\varphi (\mathrm{hu})}{\varphi (h)}={\tau }_0(u)+o(1)\text{as}h\downarrow 0\text{and}{\int }_0^1(K(v){)}^{\prime }{\tau }_0(v)\mathrm{d}v<\mathrm{\infty }.$

   5. $T^{-1}{\int }_0^T{b}_{t,{\alpha }_0}(x)\mathrm{d}t\to D_{{\alpha }_0}(x)<\mathrm{\infty }$ as $T\to \mathrm{\infty }.$

3. (Local smoothness and continuity conditions)

   Suppose for any $(y,t)\in S\times [0,T]$ and r>0 such that $t\le r\le t+\delta$:

   1. $\mathbb{E}(\psi (Y_r,y)|{\mathcal{S}}_{t,\delta })=\mathbb{E}(\psi (Y_r,y)|X_r)=m_{\psi }(X_r,y)$ a.s.

   2. $\mathbb{E}({\zeta }_r\psi (Y_r,y)|{\mathcal{S}}_{t,\delta })=\mathbb{E}({\zeta }_r\psi (Y_r,y)|X_r)$ a.s. and for $\kappa \ge 2$,

      $\mathbb{E}({\zeta }_r|\psi (Y_r,y){|}^{\kappa }|{\mathcal{S}}_{t,\delta })=\mathbb{E}({\zeta }_r|\psi (Y_r,y){|}^{\kappa }|X_r)$ a.s.

   3. $\mathrm{\exists }\beta >0$ and a constant c>0 such that, for any $(x^{\prime },x^{″})\in {\mathcal{E}}^2$, $|m_{\psi }(x^{\prime },y)-m_{\psi }(x^{″},y)|\le c{d}^{\beta }(x^{\prime },x^{″}).$

   4. For any ${\kappa }_1\ge 2$, $\mathbb{E}(|\psi (Y_r,y){|}^{{\kappa }_1}|{\mathcal{S}}_{t,\delta })=\mathbb{E}(|\psi (Y_r,y){|}^{{\kappa }_1}|X_r)$ a.s.

      The functions $W_{{\kappa }_1}(x,y):=\mathbb{E}(|\psi (Y,y){|}^{{\kappa }_1}|X=x)$ and ${\overline{W}}_{{\kappa }_1}(x,y):=\mathbb{E}(|\psi (Y,y)-m_{\psi }(x,y){|}^{{\kappa }_1}|X=x)$ are continuous in the neighbourhood of x and $\underset{x\in \mathcal{C},y\in S}{sup}|W_{{\kappa }_1}(x,y)|<\mathrm{\infty }$ a.s.

   5. $\mathbb{E}({\zeta }_r|{\mathcal{S}}_{t,\delta })=\mathbb{E}({\zeta }_r|X_r)=p(X_r)$ a.s. and for any $x\in \mathcal{E}$,

      $\underset{\{z:d(x,z)\le u\}}{sup}|p(z)-p(x)|=o(1)$ a.s. as $u\to 0.$

For j = 1, 2, define the following moments, which are independent of $x\in \mathcal{E}$ (6) $M_j=K^{(j)}(1)-{\int }_0^1{\left(K^j\right)}^{\mathrm{\prime }}(u){\tau }_0(u)\mathrm{d}u,$(6) where $K^{(j)}(\cdot )$ denotes the jth derivative of the kernel $K(\cdot )$ and $(K^j{)}^{\mathrm{\prime }}$ the first derivative of K raised to the power j.

2.1. Comments on the assumptions

Condition (A1) is related to the choice of the kernel K, which is very usual in nonparametric functional estimation. Note that Parzen symmetric kernel is not adequate in this context since the random process $D_t=d(x,X_t)$ is positive, therefore we consider K with support $[0,1]$. This is a natural generalisation of the assumption usually made on the kernel in the multivariate case where K is supposed to be a spherically symmetric density function. The assumptions $K(1)>0$ and $K^{\mathrm{\prime }}<0$ guarantee that $M_1>0$ for all limit functions ${\tau }_0.$ In the case of non-smooth processes, ${\tau }_0$ may be equal to the Dirac δ-function at 1, the condition $K(1)>0$ is needed to define the moments $M_j$ which are, in this case, determined by the value $K(1)$.

Conditions (A2)(i) –(ii) reflect the ergodicity property assumed on the continuous-time functional process. It plays an important role in studying the asymptotic properties of the estimator. The functions $f_{t,s}$ and f play the same role as the conditional and unconditional densities in finite dimensional case, whereas $\varphi (u)$ characterises the impact of the radius u on the small ball probability as u goes to 0. Several examples to satisfy these conditions are given in Laïb and Louani (2010) for discrete-time functional data process. Some other examples satisfying this condition are also given in Didi and Louani (2014) when observations $(X_t,Y_t{)}_{0\le t\le T}$ are sampled from an ergodic continuous-time process taking values in ${\mathbb{R}}^d\times \mathbb{R}$ space.

Condition (A2)-(iii) involves the ergodic nature of the process where the random function $f_{t,t-\delta }$ belongs to the space of continuous functions. Note that approximating the integral ${\int }_0^T{f}_{t,t-\delta }(x)\mathrm{d}t$ by its Riemann's sum: $T^{-1}{\int }_0^T{f}_{t,t-\delta }(x)\mathrm{d}t\simeq n^{-1}\sum _{j=1}^n{f}_{\mathrm{j\delta },(j-1)\delta }(x)$ allows to easily prove that the sequence $(f_{\mathrm{j\delta },(j-1)\delta }(x){)}_{j\ge 1}$ is stationary and ergodic (see Didi and Louani 2014). (A2)-(iv) is a usual condition when dealing with functional data, whereas (A2)-(v) is a consequence of ergodic assumption.

(A.3)(ii) is a Hölder-type assumption that requires a certain smoothness of the regression operator $m_{\psi }(\cdot ,y).$ Such assumption is commonly used in nonparametric estimation. (A.3)(iii) is a smoothness condition on the κth centred conditional moments of $\psi (Y,y).$ (A3)(iv) assumes the continuity of the conditional probability of observing a missing response. Finally, note that the moments $(M_j{)}_{j=1,2}$ are linked to the small probability function through ${\tau }_0$. One can refer to Ferraty et al. (2007) for a discussion on the choice of ${\tau }_0$, the Kernel K and the positively of $(M_j{)}_{j=1,2}$.

Discussion on the assumptions (A3)(i) –(i${}^{\mathrm{\prime }}$). These hypotheses are Markov-type conditions that characterise the conditional moments of $\psi (Y,y)$. They are satisfied for a general class of processes including the α-mixing and non α-mixing as well as long memory and the Bernoulli shift processes. As pointed in Doukhan and Louhichi (1999), the main attraction of Bernoulli shift processes is that they provide examples of processes that are weakly dependent, but not mixing. According to the discussion made in the introduction we consider below some examples in both context (continuous and discretised processes) where the predictor X is a stationary and ergodic Markovian process that might be α-mixing or not and satisfies the conditions (A.3)(i) –($i^{\mathrm{\prime }}$).

First of all let us recall the following definitions.

Definition 2.1

see Doukhan and Louhichi (1999)

Let $({ϵ}_i{)}_{i\in \mathbb{Z}}$ be a sequence of independent real-valued r.v.s and F be a measurable function defined on ${\mathbb{R}}^{\mathbb{Z}}.$ A Bernoulli shift is a sequence $(U_i{)}_{i\in \mathbb{Z}}$ defined by $U_i=F({ϵ}_{i-j},j\in \mathbb{Z}).$

Definition 2.2

see Doukhan (2018, p. 60)

A centred second-order stationary process $(X_n)$ is called long-range dependent (LRD) if $\sum _{k=0}^{\mathrm{\infty }}{r}_k^2<\mathrm{\infty }$ and $\sum _{k=0}^{\mathrm{\infty }}|r_k|=\mathrm{\infty }$, where $r_k=\mathrm{Cov}(X_k,X_0).$

Definition 2.3

A process $(B_t^H{)}_{t\in \mathbb{R}}$ is called a fractional Brownian motion (fBm) with Hurst exponent $H\in (0,1]$; if it is almost surely continuous, centred Gaussian process with covariance ${\mathrm{\Gamma }}_H(s,t)=\mathrm{Cov}(B_t^H,B_s^H)=\frac{1}{2}(|t{|}^{2H}+|s{|}^{2H}(|t-s{|}^{2H}),\mathrm{\forall }s,t\in \mathbb{R}.$

Definition 2.4

see Lemma 4.2 in Maslowski and Pospíšil (2008)

A strictly stationary centred Gaussian process $(Y_t{)}_{t\ge 0}$ is ergodic if $\underset{t\to \mathrm{\infty }}{lim}R(t):=\underset{t\to \mathrm{\infty }}{lim}\mathbb{E}(Y(0)Y(t))=0.$

Example 2.5

Continuous-time long memory processes

Let $\lambda ,\sigma >0$ and consider the Langevin equation with fBM noise $B_t^H$ and initial condition $Z_0$: (7) $Z_t=Z_0-\lambda {\int }_0^t{Z}_s\mathrm{d}s+\sigma B_t^H,t\ge 0.$(7) Then, for each $H\in (0,1)$, the following Gaussian stationary Markovian fractional Ornstein–Uhlenbeck process $(X_t^H{)}_{t\ge 0}$ defined as $X_t^H:=\sigma {\int }_{-\mathrm{\infty }}^t{e}^{-\lambda (t-u)}\mathrm{d}{B}_u^H,t\ge 0,$ is the unique (a.s.) solution of (7(7) $Z_t=Z_0-\lambda {\int }_0^t{Z}_s\mathrm{d}s+\sigma B_t^H,t\ge 0.$(7) ) with initial condition $Z_0=X_0^H$ (for more details, see Section 2, p. 5, in Cheridito et al. 2003).

Note that for $H\in (\frac{1}{2},1)$, the auto-covariance function of $(X_t^H{)}_{t\ge 0}$ is similar to that of the increments of $(B_t^H{)}_{t\ge 0}$. Therefore, $X_t^H$ is ergodic (by Definition 2.3) and exhibits long-range dependence as detailed in Theorem 2.3 and the discussion in the end of page 8 in Cheridito et al. (2003).

Now, to check the condition (A.3)(i), consider the model: $\psi (Y_t,y)=m_{\psi }(X_t^H,y)+{ϵ}_t,\text{where }{ϵ}_t\text{'s are centred, i.i.d. and independent of }{X}_t^H.$ Let ${\mathcal{S}}_{t,\delta }$ be the σ-field generated by $\sigma ((X_s^H,{ϵ}_s),(X_r^H):0\le s<t,t\le r\le t+\delta )$.

It follows that, for any $r\ge t$, $\mathbb{E}[\psi (Y_r,y)|{\mathcal{S}}_{t,\delta }]=\mathbb{E}[m_{\psi }(X_r^H,y)+{ϵ}_r|{\mathcal{S}}_{t,\delta }]$.

Since $(X_r^H,{ϵ}_r)$ are Markovian then $\mathbb{E}[\psi (Y_r,y)|{\mathcal{S}}_{t,\delta }]=m(X_r^H)$ almost surely. Thus condition (A.3)(i) is satisfied.

Example 2.6

Discrete-time processes

As discussed above, in real life we do not observe the process continuously at any time $t\in [0,T]$. We rather observe a discretised version of it based on some sampling scheme. The following examples are used to show that Assumption (A3)(i) is also satisfied for discretised processes as well.

(i) Long-memory discrete-time processes. Let $({ϵ}_t{)}_{t\in \mathbb{Z}}$ be a white noise process with variance ${\sigma }^2$, and let I and B be the identity operator and the backshift operator, respectively. Giraitis and Leipus (1995) have proved (see Theorem 1 p. 55) that the k-factor Gegenbauer process $\prod _{i\le i\le k}(I-2{\nu }_{i}B+B^2{)}^{d_i}{X}_t={ϵ}_t,$ where $0<d_i<1/2$ if $|{\nu }_i|<1$ or $0<d_i<1/4$ if $|{\nu }_i|=1$, for $i=1,\dots ,k$, is long memory, stationary, causal and invertible and has the moving average representation. That is $X_t=\sum _{j\ge 0}{\psi }_j(d,\nu ){ϵ}_{t-j}$ with $\sum _{j=0}^{\mathrm{\infty }}{\psi }_j^2(d,\nu )<\mathrm{\infty }.$

On the other hand, Guégan and Ladoucette (2001) have shown that, if $({ϵ}_t{)}_{t\in \mathbb{Z}}$ is a Gaussian process, then the above process is not strong mixing whereas the moving average representation of $(X_t)$ confirms that it is a stationary Gaussian and ergodic process.

(ii) The stationary solution of the linear Markov AR(1) process: $X_i=\frac{1}{2}{X}_{i-1}+{ϵ}_i$, where $({ϵ}_i)$ are independent symmetric Bernoulli random variables taking values $-1$ and 1, is not α-mixing (see Andrews 1984). However, $(X_i)$ is a Markovian stationary and ergodic process.

(iii) Let $(u_i)$ be an i.i.d. sequence uniformly distributed on $\{1,\dots ,9\}$, and set $X_t:=\sum _{i=0}^{\mathrm{\infty }}{10}^{-i-1}{u}_{t-i}$, where the sequence $u_t,u_{t-1},\dots ,$ represents the decimals of $X_t$. The process $X=(X_t)$ is stationary and admits the following AR(1) representation: $X_t=\frac{1}{10}{X}_{t-1}+\frac{1}{10}{u}_t=\frac{1}{10}{X}_{t-1}+\frac{1}{2}+{ϵ}_t$ where ${ϵ}_t=\frac{1}{10}{u}_t-\frac{1}{2}$ is a strong white noise. This process is not α-mixing (see Francq and Zakoïan 2010, Example A.3, p. 349), but it is ergodic.

To check the hypothesis (A3)(i) for Examples (i), (ii) and (iii), consider the regression model $\psi (Y_i,y)=m_{\psi }(X_i,y)+{\sigma }_{\psi }(X_i){\eta }_i,$ where $({\eta }_i)$ is a white noise process independent of $(X_i)$ and define the σ-field: ${\mathcal{G}}_i=\sigma ((X_1,{\eta }_1,{\zeta }_1),\dots ,(X_i,{\eta }_i,{\zeta }_i),X_{i+1})$. It is then easy to see that condition (A3)(i) is fulfilled. The discrete-time processes in examples (i)–(iii) are still valid for the regression model developed in Section 3.5 under the context of sampling schemes.

3. Main results

In this section, we investigate several asymptotic properties of the continuous-time generalised regression estimator. Some particular cases, related to specific choices of the function $\psi (\cdot ,y)$, including the conditional cumulative distribution function and the conditional quantiles will also be discussed.

3.1. Almost sure consistency rates

3.1.1. Pointwise consistency

The following theorem establishes an almost sure pointwise consistency rate of ${\hat{m}}_{\psi ,T}(x,y).$

Theorem 3.1

Pointwise consistency

Assume that (A1)–(A3) hold true and the following conditions are satisfied: (8) $\underset{T\to \mathrm{\infty }}{lim}\mathrm{T\varphi }(h_T)=\mathrm{\infty }\text{and}\underset{T\to \mathrm{\infty }}{lim}\frac{\log T}{\mathrm{T\varphi }(h_T)}=0.$(8) Then, for T sufficiently large, we have (9) ${\hat{m}}_{\psi ,T}(x,y)-m_{\psi }(x,y)={\mathcal{O}}_{a.s.}(h_T^{\beta })+{\mathcal{O}}_{a.s.}\left(\sqrt{\frac{\log T}{\mathrm{T\varphi }(h_T)}}\right).$(9)

The proof of Theorem 3.1 is detailed in the supplementary material in Chaouch and Laïb (2023).

Remark 3.1

Theorem 3.1 generalises Theorem 1 of Laïb and Louani (2011) established in the context of discrete-time stationary and ergodic processes, and Theorem 3.4 of Ferraty et al. (2005) stated under a strong mixing assumption with completely observed response where the support of y is reduced to one point. Moreover, the function $\varphi (h_T)$ can decrease to zero at an exponential rate, whenever $h_T$ goes to zero, therefore $h_T$ should be chosen to decrease to zero at a logarithmic rate.

3.1.2. Uniform consistency

To establish the uniform consistency with rate of the regression operator, we need some additional definitions and assumptions that allow to express the uniform convergence rate as a function of the entropy number. Let $\mathcal{C}$ and S be compact sets in $\mathcal{E}$ and $\mathbb{R}$, respectively. Consider, for any $ϵ>0$, the ϵ-covering number of the compact set $\mathcal{C}$, say $N_{ϵ}:=\mathcal{N}(ϵ,\mathcal{C},d)$, defined by $\begin{array}{rl}{N}_{ϵ}& :=min\{n:\text{there exist}{c}_1,\dots ,c_n\in \mathcal{C}\text{such that}\mathrm{\forall x}\in \mathcal{C}\text{we can find }\\ & 1\le i\le n\text{for which}d(x,c_i)<ϵ\}.\end{array}$ The number $N_{ϵ}$ measures how full is the class $\mathcal{C}$. The finite set of points $c_1,c_2,\dots ,c_{N_{ϵ}}$ is called an ϵ-net of $\mathcal{C}$ if $\mathcal{C}\subset {\cup }_{k=1}^{N_{ϵ}}B(c_k,ϵ)$, where $B(c_k,ϵ)$ is the ball, centred at $c_k$ and of radius ϵ, with respect to the topology induced by the semi-metric $d(\cdot ,\cdot )$. The quantity ${\phi }_{\mathcal{C}}(ϵ)=\log (N_{ϵ})$ is called the Kolmogorov's ϵ-entropy of the set $\mathcal{C}$ that may be seen as a tool to measure the complexity of the subset $\mathcal{C}$, in the sense that high entropy means that a large amount of information is needed to describe an element of $\mathcal{C}$ with an accuracy ϵ. Several examples of ${\phi }_{\mathcal{C}}(ϵ)$ covering special cases of functional processes are given in Ferraty et al. (2010) and Laïb and Louani (2011).

(U0)	Assume that (A2) holds uniformly in the following sense: (A2)(i) and (A2)(ii) hold true with the remaining term $o(\varphi (u))$ is uniform in x. For any $x\in \mathcal{C}$, $\underset{T\to \mathrm{\infty }}{lim}\underset{x\in \mathcal{C}}{sup}|\frac{1}{T}{\int }_0^T{f}_{t,t-\delta }(x)\mathrm{d}t-f(x)|=0a.s.$ $T^{-1}{\int }_0^T{b}_{t,{\alpha }_0}(x)\mathrm{d}t\to D_{{\alpha }_0}(x)$ as $T\to \mathrm{\infty }$ with $0<\underset{x\in \mathcal{C}}{sup}{D}_{{\alpha }_0}(x)<\mathrm{\infty }$. $b_0<\underset{x\in \mathcal{C}}{inf}f(x)\le \underset{x\in \mathcal{C}}{sup}f(x)<\mathrm{\infty }$ for some nonnegative real number $b_0$. $\underset{x\in \mathcal{C}}{inf}p(x)>b_1$ for some nonnegative real number $b_1$.
(U1)	The kernel function K satisfies the following conditions: K is a Hölder function of order 1 with a constant $a_K$. There exist two constants $a_2$ and $a_3$ such that $0<a_2\le K(x)\le a_3<\mathrm{\infty }$, for any $x\in \mathcal{C}.$
(U2)	For $1\le \ell \le 2$, the sequence of random variables $({\psi }^{\ell }(Y_t,y){)}_t$ is ergodic and $\mathbb{E}(|{\psi }^{\ell }(Y_0,y)|)<\mathrm{\infty }.$
(U3)	There exist $c_{\psi }>0$ and nonnegative real number γ such that for any $y\in S$ $\underset{y^{\mathrm{\prime }}\in [y-u,y+u]\cap S}{sup}|\psi (Y_t,y)-\psi (Y_t,y^{\mathrm{\prime }})|\le c_{\psi }{u}^{\gamma }.$
(U4)	Let $T_n$ be the integer part of T and suppose for T large enough $\frac{(\log T{)}^2}{\mathrm{T\varphi }(h_T)}<{\phi }_C({ϵ}_n)<\frac{\mathrm{T\varphi }(h_T)}{\log T}\mathrm{with}{ϵ}_n=\frac{\log T_n}{T_n}.$

Conditions in (U0) are standard in this context to get uniform consistency rate. Condition (U1) is usually used when we deal with nonparametric estimation for functional data, (U2) requires the existence of the moments up to order 2 of $\psi (Y,y)$. (U3) is a regularity condition upon the function $\psi (\cdot ,y)$ which is necessary to obtain the uniform consistency result over the compact S. (U4) allows to cover the subset $\mathcal{C}$ with a finite number of balls and to express the convergence rate in terms of the Kolmogorov's entropy of this subset. Similar condition has been used in Ferraty et al. (2010), where the authors have pointed out that, for a radius not too large, one requires the quantity ${\phi }_C(\log T_n/T_n)$ to be not too small and not too large. This condition seems to satisfy this exigence, since it implies that $\frac{{\phi }_C(\log T_n/T_n)}{\mathrm{T\varphi }(h_T)}$ goes to 0 for sufficiently large T. Examples given in Ferraty et al. (2010) and Laïb and Louani (2011) satisfy (U4).

Theorem 3.2 states uniform consistency rate of the Kernel regression estimator. It generalises Theorem 2 in Ferraty et al. (2010) in the i.i.d. case and that in Laïb and Louani (2011) established in the context of discrete-time stationary and ergodic processes with completely observed response.

Theorem 3.2

Uniform consistency

Assume (A1), (U0)–(U4), (A3) hold true. Moreover, suppose conditions in (8(8) $\underset{T\to \mathrm{\infty }}{lim}\mathrm{T\varphi }(h_T)=\mathrm{\infty }\text{and}\underset{T\to \mathrm{\infty }}{lim}\frac{\log T}{\mathrm{T\varphi }(h_T)}=0.$(8) ) are satisfied and (10) $\sum _{n\ge 1}{n}^{\gamma }\exp \{(1-\eta ){\phi }_{\mathcal{C}}\left(\frac{\log n}{n}\right)\}<\mathrm{\infty }\mathrm{forsome}\eta >0\mathrm{where}\gamma \mathrm{isasin}(U3).$(10) Then we have (11) $\underset{y\in S}{sup}\underset{x\in \mathcal{C}}{sup}|{\hat{m}}_{\psi ,T}(x,y)-m_{\psi }(x,y)|={\mathcal{O}}_{a.s.}\left(h_T^{\beta }\right)+{\mathcal{O}}_{a.s.}\left(\sqrt{\frac{{\phi }_{\mathcal{C}}({ϵ}_T)}{\mathrm{T\varphi }(h_T)}}\right)\mathrm{as}T\to +\mathrm{\infty }.$(11)

The proof of Theorem 3.2 is detailed in the supplementary material in Chaouch and Laïb (2023).

3.2. Asymptotic conditional bias and risk evaluation

Before evaluating the conditional bias, let us introduce some additional notations. Consider for i = 1, 2, the following assumption:

(BC1) Recall that $D_t(x)=d(X_t,x)$ and, for any $t\ge 0$, denote $\begin{array}{rl}& \mathbb{E}[m_{\psi }(X_t,y)-m_{\psi }(x,y)|D_t(x),{\mathcal{F}}_{t-\delta }]\\ & =\mathbb{E}[m_{\psi }(X_t,y)-m_{\psi }(x,y)|D_t(x)]=:{\mathrm{\Psi }}_y(D_t(x)).\end{array}$ Assume that the function ${\mathrm{\Psi }}_y$ is differentiable at 0 and satisfies ${\mathrm{\Psi }}_y(0)=0$ and ${\mathrm{\Psi }}_y^{\mathrm{\prime }}(0)\ne 0$ for any $y\in \mathbb{R}$. This condition was introduced in Ferraty et al. (2007) and used by Laïb and Louani (2010) to evaluate the conditional bias. The introduction of ${\mathrm{\Psi }}_y(\cdot )$ allows to make an integration with respect to the real random variable $D_t(x)$ rather than the couple of random variables $(D_t(x),X_t)$, where $X_t$ being functional continuous random variable.

The following proposition gives asymptotic expression of the conditional bias term, which generalises Proposition 1 in Laïb and Louani (2010) for discrete-time estimator to our setting. Its proof is similar to the one in the discrete-time framework and therefore is omitted.

Proposition 3.3

Conditional Bias

Under assumptions (A1)–(A3), (BC1) and conditions in (8(8) $\underset{T\to \mathrm{\infty }}{lim}\mathrm{T\varphi }(h_T)=\mathrm{\infty }\text{and}\underset{T\to \mathrm{\infty }}{lim}\frac{\log T}{\mathrm{T\varphi }(h_T)}=0.$(8) ), we have $B_T(x,y)=\frac{h_T{\mathrm{\Psi }}_y^{\mathrm{\prime }}(0)}{M_1}[K(1)-{\int }_0^1(\mathrm{sK}(s){)}^{\mathrm{\prime }}{\tau }_0(s)\mathrm{d}s+o_{a.s.}(1)]+{\mathcal{O}}_{a.s.}\left(h_T\sqrt{\frac{\log T}{\mathrm{T\varphi }(h_T)}}\right).$

The next result gives an explicit expression of the asymptotic quadratic risk associated to the estimator ${\hat{m}}_{\psi ,T}(x,y)$.

Theorem 3.4

Quadratic risk

Suppose that Assumptions (A1)–(A3) hold true. Then, whenever $p(x)>0$ and $f(x)>0$, we have, for a fixed $(x,y)\in \mathcal{E}\times \mathbb{R}$, that $\begin{array}{rl}& \mathrm{MSE}(x,y):=\mathbb{E}\left[{({\hat{m}}_{\psi ,T}(x,y)-m_{\psi }(x,y))}^2\right]\\ & =A_1{h}_T^2[A_1+\mathcal{O}\left(\sqrt{\frac{\log T}{\mathrm{T\varphi }(h_T)}}\right)]+\frac{A_2(x,y)}{\mathrm{T\varphi }(h_T)},\end{array}$ where $\begin{array}{rl}{A}_1& =\frac{{\mathrm{\Psi }}_y^{\mathrm{\prime }}(0)}{M_1}[K(1)-{\int }_0^1(\mathrm{sK}(s){)}^{\mathrm{\prime }}{\tau }_0(s)\mathrm{d}s+o(1)]\mathrm{and}\\ A_2(x,y)& =\frac{4(W_2(x,y)+(m_{\psi }(x,y){)}^2)M_2}{p(x)M_1^{2}f(x)}.\end{array}$

Remark 3.2

1. Note that, for sufficiently large T, the expression of MSE becomes $A_1^2{h}_T^2+\frac{A_2(x,y)}{\mathrm{T\varphi }(h_T)}$. This result generalises the one in Chaouch and Laïb (2019) established in the framework of real-valued continuous-time processes. Note, however that, for finite-dimensional continuous-time processes with MAR response, the bias term obtained in Chaouch and Laïb (2019) is of order $h_T^2$ which is smaller than $h_T$ given in Proposition 3.3. The increase in the bias term is because of the infinite dimensional characteristic of the functional space.

2. The mean squared error can be used as a theoretical guidance to select the ‘optimal’ bandwidth by minimising the quantity $A_1^2{h}_T^2+\frac{A_2(x,y)}{\mathrm{T\varphi }(h_T)}$ with respect to $h_T$. However, $A_1$ and $A_2(x,y)$ depend on some unknown quantities which should be replaced by their empirical consistent estimators, namely ${\mathrm{\Psi }}_{y,T}^{\mathrm{\prime }}(0),(M_{j,T}{)}_{j=1,2},{\tau }_{0,T},p_T,W_{2,T},$ and $f_T$. Note that ${\mathrm{\Psi }}_y^{\mathrm{\prime }}(0)$ may be viewed as real regression function with response variable $m_{\psi }(X,y)-m_{\psi }(x,y)$ and predictor $d(X,x)$. It may be then estimated by a kernel regression estimate ${\mathrm{\Psi }}_{y,T}^{\mathrm{\prime }}(0)$ by replacing $m_{\psi }(x,y)$ by its estimator ${\hat{m}}_{\psi ,T}(X,y)$.

3.3. Asymptotic normality

The following theorem establishes the asymptotic distribution of the estimator.

Theorem 3.5

Assume that conditions (A1)–(A3) are fulfilled. Suppose that, for β as defined in (A3)(ii), the following conditions hold true: (12) $\underset{T\to \mathrm{\infty }}{lim}\mathrm{T\varphi }(h_T)=\mathrm{\infty },h_T^{\beta }\sqrt{\mathrm{T\varphi }(h_T)}=o(1)\mathrm{and}{h}_T^{\beta }\log T^{1/2}=o(1)\mathrm{as}T\to \mathrm{\infty }.$(12) Then, for any $(x,y)\in \mathcal{E}\times S$ such that $f(x)>0,$ we have $\sqrt{\mathrm{T\varphi }(h_T)}({\hat{m}}_{\psi ,T}(x,y)-m_{\psi }(x,y))\stackrel{d}{\to }N(0,{\sigma }^2(x,y)),$ where (13) ${\sigma }^2(x,y)\le \frac{1}{f(x)}\frac{M_2}{M_1^{2}p(x)}{\overline{W}}_2(x,y)=:\frac{1}{f(x)}\tilde{V}(x,y)\mathrm{as}T⟶\mathrm{\infty }$(13) and ${\overline{W}}_2(x,y):=\mathbb{E}[{(\psi (Y,y)-m_{\psi }(x,y))}^2|X=x].$

Note that the statement (13(13) ${\sigma }^2(x,y)\le \frac{1}{f(x)}\frac{M_2}{M_1^{2}p(x)}{\overline{W}}_2(x,y)=:\frac{1}{f(x)}\tilde{V}(x,y)\mathrm{as}T⟶\mathrm{\infty }$(13) ) gives only an upper bound of the asymptotic variance ${\sigma }^2(x,y)$. The following proposition gives an estimate of $\tilde{V}(x,y)$ that will be needed to construct confidence intervals for the unknown operator $m_{\psi }(x,y)$.

Proposition 3.6

Suppose conditions of Theorem 3.5 hold and ${\sigma }^2(x,y)>0$, then (14) ${\hat{V}}_T(x,y):=\frac{\sqrt{M_{T,2}}}{M_{1,T}}\sqrt{\frac{{\overline{W}}_{2,T}(x,y)}{T{F}_{x,T}(h_T)p_T(x)}},$(14) is a consistent estimator for $\tilde{V}(x,y)$. The quantities $M_{1,T}$, $M_{2,T}$, ${\overline{W}}_{2,T}$, $p_T(x)$ and $F_{x,T}$ are empirical versions of $M_1$, $M_2$, ${\overline{W}}_2$, $p(x)$ and $F_x,$ respectively.

$M_{1,T}$ and $M_{2,T}$ are calculated by replacing ${\tau }_0$, given in (A2)(iv), by its empirical version ${\tau }_{0,T}(u)=\frac{F_{x,T}(u{h}_T)}{F_{x,T}(h_T)}\mathrm{where}{F}_{x,T}(u)=\frac{1}{T}{\int }_0^{T}1{\mathrm{l}}_{\{d(x,X_t)\le u\}}\mathrm{d}t.$ On the other hand ${\overline{W}}_{2,T}(x,y)$ and $p_T(x)$ are given by ${\overline{W}}_{2,T}(x,y)=\frac{{\int }_0^T{\zeta }_t(\psi (Y_t,y){)}^2{\mathrm{\Delta }}_t(x)\mathrm{d}t}{{\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t}-({\hat{m}}_{\psi ,T}(x,y){)}^2\text{and}{p}_T(x)={\displaystyle \frac{{\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t}{{\int }_0^T{\mathrm{\Delta }}_t(x)\mathrm{d}t}}.$

3.4. Continuous-time confidence intervals

Using the non-decreasing property of the cumulative standard Gaussian distribution function, the estimator ${\hat{V}}_T(x,y)$, defined in (14(14) ${\hat{V}}_T(x,y):=\frac{\sqrt{M_{T,2}}}{M_{1,T}}\sqrt{\frac{{\overline{W}}_{2,T}(x,y)}{T{F}_{x,T}(h_T)p_T(x)}},$(14) ), with the help of Proposition 3.6 and Theorem 3.5, the following corollary provides estimated confidence intervals for $m_{\psi }(x,y)$ at any x fixed.

Corollary 3.7

Assume conditions of Theorem 3.5 are fulfilled and the conditions in (12(12) $\underset{T\to \mathrm{\infty }}{lim}\mathrm{T\varphi }(h_T)=\mathrm{\infty },h_T^{\beta }\sqrt{\mathrm{T\varphi }(h_T)}=o(1)\mathrm{and}{h}_T^{\beta }\log T^{1/2}=o(1)\mathrm{as}T\to \mathrm{\infty }.$(12) ) are replaced by (15) $\underset{T\to \mathrm{\infty }}{lim}T{F}_{x,T}(h_T)=\mathrm{\infty }\mathrm{and}\underset{T\to \mathrm{\infty }}{lim}{h}_T^{\beta }\sqrt{T{F}_{x,T}(h_T)}=0.$(15) Then, for any $0<\alpha <1$, the $(1-\alpha )$ confidence intervals for $m_{\psi }(x,y)$ are (16) $I_{1-\alpha }(x,y)=[{\hat{m}}_{\psi ,T}(x,y)-c_{\alpha /2}\sqrt{{\hat{V}}_T(x,y)};{\hat{m}}_{\psi ,T}(x,y)+c_{1-\alpha /2}\sqrt{{\hat{V}}_T(x,y)}],$(16) where $c_{\alpha }$ is the αth quantile of the standard normal distribution.

These intervals are similar to those given in Remark 2 in Laïb and Louani (2010) for discrete-time ergodic context with complete data, and those obtained in Ling et al. (2015) for discrete-time stationary ergodic data with missing at random response.

3.5. Sampling schemes and computation of the confidence intervals

In the previous section, the process was supposed to be observable over $[0,T]$. However, in practice the data are often collected according to a sampling scheme since it is difficult to observe a path continuously at any time t over the interval $[0,T]$. Hereafter, we briefly discuss the effect of a sampling scheme on the construction of confidence intervals for the regression function $m_{\psi }(x,y)$. Assume that the data are sampled, either regularly, irregularly or even randomly, from an underlying continuous-time process at instants $(t_k{)}_{k=1,\dots ,n}$. For a sake of simplicity, we consider here the case where the instants $(t_k)$ are irregularly spaced, that is $\underset{1\le k\le n}{inf}|t_{k+1}-t_k|=\delta >0.$ Now, for $k\in \{1,\dots n\}$, we define the following increasing families of σ-algebra: ${\mathcal{F}}_k:={\mathcal{F}}_{t_k}=\sigma ((X_{t_1},Y_{t_1}),\dots ,(X_{t_k},Y_{t_k}))$ and ${\mathcal{G}}_k:={\mathcal{G}}_{t_k}=\sigma ((X_{t_1},Y_{t_1}),\dots ,(X_{t_k},Y_{t_k});X_{t_{k+1}}).$ The purpose then consists in estimating $m_{\psi }(x,y)$ given the discrete-time ergodic stationary process $(X_{t_k},Y_{t_k},{\zeta }_{t_k}{)}_{1\le k\le n}$ sampled from the underlying continuous-time process $\{(X_t,Y_t,{\zeta }_t){\}}_{0\le t\le T}$. In case of a regular sampling scheme, that is $T=\mathrm{n\delta }$, the regression function is $m_{\psi }(x,y):=\frac{\mathbb{E}({\zeta }_{t_k}\psi (Y_{t_k},y)|X_{t_k}=x)}{\mathbb{E}(\zeta |X_{t_k}=x)},1\le k\le n$ and its estimator ${\hat{m}}_{\psi ,T}(x,y)$ defined in (3(3) ${\hat{m}}_{\psi ,T}(x,y):=\{\begin{array}{ll}{\displaystyle \frac{{\int }_0^T{\zeta }_t\psi (Y_t,y){\mathrm{\Delta }}_t(x)\mathrm{d}t}{{\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t}},& \text{if}{\displaystyle {\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t\ne 0}\\ {\displaystyle \frac{1}{T}}{\displaystyle {\int }_0^T{\zeta }_t\psi (Y_t,y)\mathrm{d}t,}& \text{otherwise},\end{array}$(3) ) becomes (17) ${\hat{m}}_{\psi ,n}(x,y)=\frac{\sum _{k=1}^n{\zeta }_{t_k}\psi (Y_{t_k},y)K\left(\frac{d(x,X_{t_k})}{h_n}\right)}{\sum _{k=1}^n{\zeta }_{t_k}K\left(\frac{d(x,X_{t_k})}{h_n}\right)},t_k=\mathrm{k\delta }(1\le k\le n).$(17) Note that Theorem 3.5 holds for the estimate ${\hat{m}}_{\psi ,n}(x,y)$ when replacing T by $\mathrm{n\delta }$. The limiting law is a Gaussian random variable with mean zero and variance function ${\sigma }^2(x,y)=\frac{1}{f(x)}\frac{M_2}{M_1^{2}p(x)}{\overline{W}}_2(x,y).$ Making use of Corollary 3.7 and considering similar steps as in Laïb and Louani (2010), it follows that, for any $0<\alpha <1$, the $(1-\alpha )$ asymptotic confidence intervals of $m_{\psi }(x,y)$ are (18) ${\hat{m}}_{\psi ,n}(x,y)\pm c_{1-\alpha /2}\frac{\sqrt{M_{n,2}}}{M_{n,1}}\sqrt{\frac{{\overline{W}}_{n,2}(x,y)}{n{F}_{x,n}(x)p_n(x)}},\mathrm{as}n\to \mathrm{\infty },$(18) where $c_{1-\alpha /2}$ is the quantile of standard normal distribution.

4. Simulation study

This section aims to discuss numerically some aspects related to continuous-time processes that might affect the quality of estimation of the operator $m_{\psi }(x,y)$. Here we consider $\psi (Y_t,y)=Y_t,$ therefore $m_{\psi }(x,y)=m(x)=\mathbb{E}(Y_t|X_t=x)$, where $m(x)$ is the conditional expectation of $Y_t$ given $X_t=x.$ The first simulation aims to compare the quality of estimation of $m(x)$ based on the continuous-time and discrete-time processes. In the second simulation, we discuss the choice of the ‘optimal’ sampling mesh δ in the case of continuous-time processes and assess its sensitivity to the missing at random mechanism. Finally, the third simulation discusses the effect of the MAR rate on the coverage rate and length of the estimated confidence intervals.

4.1. Simulation 1: continuous-time versus discrete-time estimators

In this first simulation, we try to compare the estimation of the regression operator when discrete- and continuous-time processes are considered. We want to know whether considering a continuous-time processes may improve the quality of the predictions or not. We suppose that the functional space $\mathcal{E}=L^2([-1,1])$ endowed with its natural norm. The generation of continuous-time processes $(\{X_t(s):s\in [-1,1]\},Y_t{)}_{t\in [0,T]}$ is obtained by considering the following steps:

1. First, we simulate an Ornstein–Uhlenbeck (OU) process $(Z_t{)}_{t\ge 0}$ solution of the following stochastic differential equation: (19) $\mathrm{d}{Z}_t=2(5-Z_t)\mathrm{d}t+7\mathrm{d}{\text{W}}_t,$(19) where ${\text{W}}_t$ denotes a Wiener process. Here, we take $\mathrm{d}t=\mathrm{0.005.}$

2. Let $\mathrm{\Gamma }(\cdot )$ be the operator mapping $\mathbb{R}$ into $L^2([-1,1])$ defined, for any $z\in \mathbb{R}$, as follows: $\mathrm{\Gamma }(z):=(1+\lfloor z\rfloor -z){\text{P}}_{\text{num}(\lfloor z\rfloor )}+(z-\lfloor z\rfloor ){\text{P}}_{\text{num}(\lfloor z+1\rfloor )},$ where $P_j$ is the Legendre polynomials of degree j and $\text{num}(z):=1+2z\text{sign}(z)-\text{sign}(z)(1+\text{sign}(z))/2$ and $\lfloor \cdot \rfloor$ denotes the floor function.

3. We consider that curves are sampled at 400 equispaced values in $[-1,1]$ and defined, for any $t\in [0,T]$, as $X_t(s)=\mathrm{\Gamma }(Z_t)(s),s\in [-1,1].$

4. To generate the real-valued process $(Y_t{)}_{t\in [0,T]}$, the following nonlinear functional regression model is considered: (20) $Y_t=m(X_t)+{\epsilon }_t,$(20) where $m(x):={\int }_{-1}^1{x}^2(s)\mathrm{d}s,$ ${\epsilon }_t=U_t-U_{t-1}$ where $U_t$ is a Wiener process independent of $X_t$.

Observe that the OU process $\{Z_t:t\in [0,T]\}$ is a real-valued continuous-time process (since dt tends to zero). The operator $\mathrm{\Gamma }(\cdot )$ has a role to transform each observation in the process $Z_t$ into a curve through the Legendre polynomials. In such way, the functional variable X is generated continuously as is the process $(Z_t)$. Moreover, note that steps 1, 2 and 3 are devoted to simulate the continuous-time functional process $\{X_t(s):s\in [-1,1]{\}}_{t\in [0,T]}$, whereas in step 4 the real-valued continuous-time process $(Y_t{)}_{t\in [0,T]}$ is generated. A sample of 20 simulated curves is displayed in Figure 1(left) and an example of the real-valued process $(Y_t)$ is given in Figure 1(right).

Figure 1. Left: A sample of 20 simulated curves $\{X_t(s):s\in [-1,1]\}$. Right: A realisation of the process $(Y_t{)}_{t\in [0,200]}$.

Now, our purpose is to compare, in terms of estimation accuracy, the continuous-time estimator with the discrete-time one for different values of T = 50, 200, 1000 and several missing at random rates. It is worth noting that the continuous-time process $(X_t,Y_t)$ is observed at every instant $t=\delta ,2\delta ,\dots ,\mathrm{n\delta }$, where $\delta =0.005$ and $n=T/\delta$. However, the discrete-time process is observed only at the instants $t=1,2,\dots ,n.$

As in Ferraty et al. (2013) and Ling et al. (2015), we consider that the missing at random mechanism is led by the following probability distribution: (21) $p(x)=\mathbb{P}({\zeta }_t=1|X_t=x)=\text{expit}({\int }_{-1}^1{x}^2(s)\mathrm{d}s),$(21) where $\text{expit}(u)=e^u/(1+e^u)$, for $u\in \mathbb{R}.$ Now, we specify the tuning parameters on which depend our estimator given in (3(3) ${\hat{m}}_{\psi ,T}(x,y):=\{\begin{array}{ll}{\displaystyle \frac{{\int }_0^T{\zeta }_t\psi (Y_t,y){\mathrm{\Delta }}_t(x)\mathrm{d}t}{{\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t}},& \text{if}{\displaystyle {\int }_0^T{\zeta }_t{\mathrm{\Delta }}_t(x)\mathrm{d}t\ne 0}\\ {\displaystyle \frac{1}{T}}{\displaystyle {\int }_0^T{\zeta }_t\psi (Y_t,y)\mathrm{d}t,}& \text{otherwise},\end{array}$(3) ). We choose the quadratic kernel defined as $K(u)=\frac{3}{4}(1-u^2)1{\mathrm{l}}_{(0,1)}(u)$ and because curves are smooth enough we choose as semi-metric the $L_2$-norm of the second derivatives of the curves, that is for $t_1\ne t_2$, (22) $\mathrm{d}(X_{t_1},X_{t_2})={({\int }_{-1}^1[X_{t_1}^{(2)}(s)-X_{t_2}^{(2)}(s){]}^2\mathrm{d}s)}^{1/2}.$(22) We used the local cross-validation method on the κ-nearest neighbours introduced in Ferraty and Vieu (2006) page 116 to select the optimal bandwidth for both discrete- and continuous-time regression estimators. The accuracy of the discrete- and continuous-time regression estimators is evaluated over M = 500 replications. The accuracy is measured, at each replication $j=1,\dots ,M$, by using the squared errors ${\text{SE}}_T^j:=({\hat{m}}_T^j(x)-m(x){)}^2$ and ${\text{SE}}_n^j:=({\hat{m}}_n^j(x)-m(x){)}^2$ for the continuous-time and discrete-time estimators, respectively. Observe that the discrete-time estimator of the regression operator is defined as ${\hat{m}}_n(x):={\displaystyle \frac{\sum _{t=1}^n{\zeta }_t{Y}_t{\mathrm{\Delta }}_t(x)}{\sum _{t=1}^n{\zeta }_t{\mathrm{\Delta }}_t(x)}}.$

To get a better idea about the variability of the errors, Table 1 summarises the distribution of the squared errors (multiplied by ${10}^{-2}$) $({\mathrm{SE}}_T^j{)}_{j=1,\dots ,M}$ and $({\mathrm{SE}}_n^j{)}_{j=1,\dots ,M}$. It shows that continuous-time regression estimator is more accurate than the discrete-time one. Moreover, when T increases the squared errors decrease faster when working with the continuous-time process.

Table 1. Summary statistics of $({\text{SE}}^j{)}_{j=1,\dots ,500}$ for discrete- and continuous-time estimators of the regression function.

		Continuous ($\times {10}^{-2}$)			Discrete ($\times {10}^{-2}$)
MAR rate	T	50	200	1000	50	200	1000
$p=20\mathrm{\%}$	$Q_{25\mathrm{\%}}$	0.56	0.54	0.18	1.6	1.3	0.7
	Median	2.4	2.58	1.25	7.5	3.7	2.1
	Mean	5.39	4.62	4	15.11	11.8	9.1
	$Q_{75\mathrm{\%}}$	6.21	7.13	4.7	16.06	10.9	4.7
$p=40\mathrm{\%}$	$Q_{25\mathrm{\%}}$	0.49	0.6	0.3	1.7	1	0.5
	Median	4.77	2.7	2.2	4.9	6.1	2
	Mean	6.93	8.8	4.6	12.4	9.9	8.3
	$Q_{75\mathrm{\%}}$	9.39	10.6	4.5	16.05	11.3	11.1

4.2. Simulation 2: optimal sampling mesh selection

The purpose of this simulation is to investigate another aspect related to continuous-time processes. The selection of the ‘optimal’ sampling mesh is one of the most important topics in continuous-time processes.

First of all, we generate a continuous-time functional data process according to the following equation: $X_t(s)=Z_t(1-\sin (s-\pi /3)),s\in [0,\pi /3]\text{and}t\in [0,T],$ where $Z_t$ is an OU process solution of the stochastic differential equation (19(19) $\mathrm{d}{Z}_t=2(5-Z_t)\mathrm{d}t+7\mathrm{d}{\text{W}}_t,$(19) ) and practically observed at the instants $t=\delta ,2\delta ,\dots ,\mathrm{n\delta }$ with n = 200 fixed. Here, we take different values of sampling mesh δ, calculate the corresponding empirical version of the Mean Integrated Square Error ($\mathrm{MISE}(\delta )$) and identify the optimal mesh, say ${\delta }^{\star }$, that minimises $\mathrm{MISE}(\delta )$. Note that each curve observed at the instant t is discretised at 100 equidistant points over the interval $[0,\pi /3].$ The response variable is obtained following the hereafter nonlinear functional regression model (20(20) $Y_t=m(X_t)+{\epsilon }_t,$(20) ), where the operator $m(\cdot )$ is defined as $m(X_t)={({\int }_0^{\pi /3}{X}_t^{\mathrm{\prime }}(s)\mathrm{d}s)}^2\mathrm{and}{\epsilon }_t\sim N(0,0.075).$ Moreover, the missing at random mechanism in this simulation is also supposed to be the same as described in the first simulation as per Equation (21(21) $p(x)=\mathbb{P}({\zeta }_t=1|X_t=x)=\text{expit}({\int }_{-1}^1{x}^2(s)\mathrm{d}s),$(21) ). For the tuning parameters used to build the estimator, we considered the quadratic kernel and given the shape of the true regression operator, which depends on the first derivative of the functional predictor, the Euclidean distance between the first-order derivatives of the curves is adopted as a semi-metric. Finally the bandwidth is selected according to the local cross-validation method based on the κ-nearest neighbours as detailed in Ferraty and Vieu (2006, p. 116).

For each value of sampling mesh δ, the regression operator $m(\cdot )$ is estimated over a grid of 50 different fixed curves and the whole procedure is repeated over M = 500 replications. Finally, the empirical MISE is calculated, for each sampling mesh δ, according to the following equation: $\mathrm{MISE}(\delta ):={\displaystyle \frac{1}{M}}\sum _{k=1}^M{\displaystyle \frac{1}{50}}\sum _{j=1}^{50}{(m(x_j)-{\hat{m}}_{n,\delta }^k(x_j))}^2.$ Observe that ${\hat{m}}_{n,\delta }^k(\cdot )$ is the estimator of $m(\cdot )$, obtained at the kth iteration, depends on the sampling mesh δ, so is the MISE.

Figure 2 displays the values of $\mathrm{MISE}(\delta )$ obtained for different values of sampling mesh δ and missing at random rate of 10%, 50% and 0% (complete data), respectively. One can observe that higher is the missing at random rate, higher are the errors in estimating the regression operator.

Figure 2. The $\mathrm{MISE}(\delta )$ obtained for different values of sampling mesh δ and several missing at random rates.

Table 2 reports the optimal sampling mesh ${\delta }^{\star }$, which minimises $\mathrm{MISE}(\delta )$, for different missing at random rates. It also provides some summary statistics to have an idea about the distribution of $\mathrm{MISE}(\delta )$ for several values of δ. One can observe, from Table 2, that higher is the missing at random rate longer we need to observe the underlying process to collect the n = 200 observations to be able to reasonably estimate the regression operator. Indeed, when the data is complete the optimal time interval $T^{\star }=n{\delta }^{\star }=200\times 0.3=60$. However, when $\mathrm{MAR}=10\mathrm{\%}$ (resp. 50%) the optimal time interval is equal to $T^{\star }=200\times 0.36=72$ (resp. $T^{\star }=200\times 0.38=76$). Consequently, it can be concluded that when the missing at random mechanism is heavily affecting the response variable, we need to collect data over a longer period of time. This allows to get sufficient information about the dynamic of the underlying continuous-time process and therefore get a better estimate of the regression operator.

Table 2. The optimal sampling mesh (${\delta }^{\star }$) obtained for different MAR rates and some summary statistics of the MISE(δ).

	${\delta }^{\star }$	MISE(${\delta }^{\star }$)	$\overline{\mathrm{MISE}(\delta )}$	$Q_{25\mathrm{\%}}$	$Q_{50\mathrm{\%}}$	$Q_{75\mathrm{\%}}$
MAR=0%	0.30	0.0476	0.0562	0.0510	0.0524	0.0550
MAR=10%	0.36	0.0555	0.0643	0.0587	0.0612	0.0635
MAR=50%	0.38	0.0635	0.0767	0.0714	0.0750	0.0770

4.3. Simulation 3: asymptotic confidence intervals

In this section, we are interested in evaluating the coverage rate, as well as the length, of the asymptotic confidence intervals given in (16(16) $I_{1-\alpha }(x,y)=[{\hat{m}}_{\psi ,T}(x,y)-c_{\alpha /2}\sqrt{{\hat{V}}_T(x,y)};{\hat{m}}_{\psi ,T}(x,y)+c_{1-\alpha /2}\sqrt{{\hat{V}}_T(x,y)}],$(16) ). The effect of the sampling mesh on the coverage rate will also be discussed numerically. Since this paper aims to extend results in Delsol (2009) about confidence intervals to continuous-time functional data, we consider the same simulation framework.

Let $X_t(s)=\cos (Z_t+\pi (2s-1)),s\in [0,1]\mathrm{and}t\in [0,T],$ where $Z_t$ is an OU process solution of the stochastic differential equation (19(19) $\mathrm{d}{Z}_t=2(5-Z_t)\mathrm{d}t+7\mathrm{d}{\text{W}}_t,$(19) ) observed at the instants $t=\delta ,2\delta ,\dots ,\mathrm{n\delta }$ with n = 100, 200 fixed. Here, for comparison purpose, we consider two sampling mesh $\delta =0.1,\mathrm{0.7.}$ The regression operator is defined as $m(x)=\frac{1}{2\pi }{\int }_{1/2}^{3/4}(x^{\mathrm{\prime }}(s){)}^2\mathrm{d}s$, while the errors $\{{\epsilon }_t\}$ are independent centred normal random variable with variance $0.1\times s_n^2$, where $s_n^2$ is the empirical variance of $\{m(X_1),\dots ,m(X_n)\}.$ Because the regression operator is defined as a function of the derivative of the functional random variable, the appropriate semimetric to be used in such case is based on the first derivative of the curve (see (22(22) $\mathrm{d}(X_{t_1},X_{t_2})={({\int }_{-1}^1[X_{t_1}^{(2)}(s)-X_{t_2}^{(2)}(s){]}^2\mathrm{d}s)}^{1/2}.$(22) )). Moreover, the quadratic kernel is used to perform this simulation. The optimal bandwidth is selected based on local cross-validation method on the κ-nearest neighbours. The missing at random rate is simulated according the conditional probability distribution given in (21(21) $p(x)=\mathbb{P}({\zeta }_t=1|X_t=x)=\text{expit}({\int }_{-1}^1{x}^2(s)\mathrm{d}s),$(21) ).

For a fixed $\alpha \in (0,1)$, the asymptotic $(1-\alpha )$-confidence intervals for $m(x)$ with $x\in \mathrm{\Xi }$ are computed and compared for several values of sample size n and sampling mesh δ. Here $\mathrm{\Xi }:=\{x_1,\dots ,x_{n_{\mathrm{\Xi }}}\}$ is a grid of $n_{\mathrm{\Xi }}=50$ independently simulated curves where the regression operator is estimated. For every fixed curve $x\in \mathrm{\Xi },$ a number of M = 500 replications is considered to approximate the coverage rate. In this simulation $1-\alpha =0.95,0.90$ were considered.

As expected, Table 3 shows that the average coverage rate varies with the sample size n, the sampling mesh and the MAR Rate. Higher are the sample size and the sampling mesh and smaller is the MAR rate closer will be the average coverage rate to $1-\alpha .$ Moreover, one can also observe that the asymptotic confidence intervals length decreases when the sample size increases and the MAR rate decreases.

Table 3. Average coverage over the grid Ξ and average confidence interval length appears in brackets.

		0.95		0.90
$1-\alpha$	MAR rate	n = 100	n = 200	n = 100	n = 200
$\delta =0.1$	$\mathrm{MAR}=0\mathrm{\%}$	0.762 (0.148)	0.778 (0.104)	0.743 (0.126)	0.737 (0.088)
	$\mathrm{MAR}=25\mathrm{\%}$	0.750 (0.149)	0.767 (0.104)	0.725 (0.125)	0.731 (0.088)
	$\mathrm{MAR}=45\mathrm{\%}$	0.722 (0.152)	0.750 (0.106)	0.692 (0.127)	0.697 (0.089)
$\delta =0.7$	$\mathrm{MAR}=0\mathrm{\%}$	0.817 (0.179)	0.851 (0.148)	0.850 (0.144)	0.81 (0.124)
	$\mathrm{MAR}=25\mathrm{\%}$	0.811 (0.178)	0.841 (0.149)	0.831 (0.144)	0.798 (0.123)
	$\mathrm{MAR}=45\mathrm{\%}$	0.766 (0.179)	0.791 (0.148)	0.771 (0.145)	0.738 (0.124)

Figure 3 (resp. Figure 4) displays an example of asymptotic confidence intervals obtained for the 50 curves in the testing sample when n = 100, the MAR rate = 0%, 25%, 45%, $\delta =0.7$ and $1-\alpha =0.95$ (resp. $1-\alpha =0.9$). One can observe that the coverage rate decreases with an increase in the MAR rate. Similar results are also obtained when $\delta =0.1$.

Figure 3. Asymptotic $95\mathrm{\%}$ confidence intervals when $\delta =0.7$ and n = 100. Red dot represents the true regression function value at any fixed curve $x\in \mathrm{\Xi }$ and the black dot is its estimation. The vertical lines represent the confidence intervals.

Figure 4. Asymptotic $90\mathrm{\%}$ confidence intervals when $\delta =0.7$ and n = 100. Red dot represents the true regression function value at any fixed curve $x\in \mathrm{\Xi }$ and the black dot is its estimation. The vertical lines represent the confidence intervals.

5. Applications to real data

5.1. Application 1: prediction of financial asset returns

In Financial market, despite the modern technology, which allows to collect data at a very fine time scale, financial data can still be missing. For instance, there are some regular holidays, such as Thanksgiving Day and Christmas, for which stock price data are missing. There are many other technical reasons (such as breakdown in devises recording data, computers' sudden shutdowns, …) that make stretches of data missing.

This section aims to assess the performance of the estimator proposed in this paper on missing at random financial functional time series. The International Business Machine cooperation (IBM) asset price is considered as the response variable and the Standard & Poor's 500 (SP500) stock market index as predictor. While the IBM asset price is observed at a daily frequency from 24 March 2016 to 28 September 2016, the SP500 is observed every minute during the same period. Note that the daily trading activity lasts for 7 hours excluding the weekend. Since in this paper, we are interested in stationary processes, a first-order differentiation of the IBM daily asset price and the SP500 stock market index was considered to make the original time series stationary.

Our sample here can be denoted as follows: $(X_d,Y_d{)}_{d=1,\dots 129}$, where the sample size n = 129 is the total number of trading days from 24 March 2016 to 28 September 2016 after first differentiation of the original time series, $Y_d=\mathrm{\Delta }{\mathrm{IBM}}_d$, $X_d=\{X_d(s):=\mathrm{\Delta }{\mathrm{SP}500}_d(s):1\le s\le 420\}.$ Originally, the data are completely observed. Therefore, to validate our estimator, we artificially create missing observations. We assume here that the missing data are generated according to the conditional probability distribution given in (21(21) $p(x)=\mathbb{P}({\zeta }_t=1|X_t=x)=\text{expit}({\int }_{-1}^1{x}^2(s)\mathrm{d}s),$(21) ). We split the original sample into training and testing subsets. Our purpose then is to predict the IBM asset price in the testing subset using the regression operator. Three MAR rates $0\mathrm{\%}$ (complete data), $20\mathrm{\%}$ and $45\mathrm{\%}$ were considered to test the performance of the estimator in terms of prediction. Similarly, as in the simulation section, we considered here the quadratic kernel and the bandwidth was selected using the cross-validation method on the κ-nearest neighbours. For the semi-metric, because the curves are not smooth (as it can be seen in Figure 5, right panel) we use the PCA-semi-metric, say $d_4^{\mathrm{PCA}}(\cdot ,\cdot )$, based on the projection on the four eigenfunctions, $v_1(\cdot ),\dots ,v_4(\cdot )$, associated to the four largest eigenvalues of the empirical covariance operator of the functional predictor X: $d_4^{\mathrm{PCA}}(X_{t_1},X_{t_2})=\sqrt{\sum _{k=1}^4{({\int }_0^{420}(X_{t_1}(s)-X_{t_2}(s))v_k(s)\mathrm{d}s)}^2}.$ As criteria to measure the accuracy of the estimator in predicting 30 observations in the testing subset, we considered the Absolute Error: ${\mathrm{AE}}_t:=|Y_t-m_n(X_t)|$ for $t=1,\dots ,30$. Figure 6 displays the distribution of the absolute errors obtained when MAR rate is $0\mathrm{\%},20\mathrm{\%}$ and $45\mathrm{\%}$, respectively. One can clearly observe the effect of the MAR rate on the quality of the prediction. Higher is the MAR rate lower is the quality of prediction.

Figure 5. Left: First-order differentiated IBM asset price. Right: First-order differentiated SP500 intraday (minute frequency) stock market index curves.

Figure 6. The $(\mathrm{AE}{)}_{t=1,\dots ,30}$ obtained for different values of missing at random rates.

Moreover, we build a 95% prediction interval for the IBM asset price in the testing subset. Figure 7 shows that the coverage rate is sensitive to the percentage MAR data in the training subset.

Figure 7. A 95% prediction interval of the IBM asset price in the testing subset. Red dot represents the true values of IBM asset price and the black dot is their prediction using the regression operator. The vertical lines represent the prediction intervals.

5.2. Application 2: Daily peak electricity demand imputation

By accurately predicting household peak load, utility companies can better balance the overall electricity demand and supply. This information helps in optimising power generation and distribution, ensuring a stable and reliable electricity grid. It also allows to plan for peak demand periods and avoid potential blackouts or overloading of the grid. Moreover, predicting peak loads empowers consumers with information about their electricity consumption patterns. Thus, with this knowledge, households can make informed decisions to manage their energy usage more effectively, reduce electricity bills and contribute to energy conservation efforts. Furthermore, peak load predictions enable demand response programs, where utility companies offer incentives for consumers to adjust their energy consumption during peak periods, thereby reducing strain on the grid.

For these reasons, electricity companies deployed smart meters to replace the mechanical one. This new generation of smart meters allows to record the electricity demand of any household at very fine time scale and send it to the information system. The transmission of the information from the smart meter towards the information system goes usually through WIFI or optical fibre networks which are significantly dependent on the weather conditions, among several other factors. Therefore, the calculation of the daily peak electricity demand might be subject to missing at random mechanism due to bad weather conditions.

Figure 8. Daily peak electricity demand of a household containing 10% missing data.

Figure 9. Intraday temperature curves. Coloured curves are for the missing peak load days.

Figure 10. Daily peak electricity demand process. Red dots represent values of imputed missing data.

Figure 11. 95% confidence intervals for the imputed values.

Figure 8 displays the daily peak load $\{Y_{t_k}{\}}_{k=1,\dots ,1009}$ obtained from a household smart meter from 24 September 1996 to 29 June 1999 (leading to a total of n = 1009 days). The original data contains 10% of missing observations. Here, we assume that the intraday (3-hour frequency) temperature curve $\{X_{t_k}(s):s=3,6,\dots ,24\}$ explains the missingness mechanism in the daily peak demand. Figure 9 displays the intraday, 3-hour frequency, temperature curves. Our purpose in this application is to impute the missing data in the peak demand process using the initial estimator of the regression operator ${\hat{m}}_n(x):={\hat{m}}_{\psi ,n}(x,y)$ defined in (17(17) ${\hat{m}}_{\psi ,n}(x,y)=\frac{\sum _{k=1}^n{\zeta }_{t_k}\psi (Y_{t_k},y)K\left(\frac{d(x,X_{t_k})}{h_n}\right)}{\sum _{k=1}^n{\zeta }_{t_k}K\left(\frac{d(x,X_{t_k})}{h_n}\right)},t_k=\mathrm{k\delta }(1\le k\le n).$(17) ) with $\psi (Y_{t_k},y)=Y_{t_k}.$ Figure 10 displays the imputed peak electricity demand process obtained according to the following formula: ${\tilde{Y}}_{t_k}={\delta }_{t_k}{Y}_{t_k}+(1-{\delta }_{t_k}){\hat{m}}_n(X_{t_k}).$ If $Y_{t_k}$ is observed (that is ${\delta }_{t_k}=1$), then ${\tilde{Y}}_{t_k}=Y_{t_k}$, otherwise $Y_{t_k}$ is missing (i.e. ${\delta }_{t_k}=0$) and will be imputed by ${\hat{m}}_n(X_{t_k}).$ The red dots in Figure 10 represent the imputed values of the missing observations in the peak electricity demand process. Figure 11 shows the 95% confidence intervals around the missing values of the peak load.

6. Discussion of a special case: conditional quantiles

Let $x\in \mathcal{E}$ be fixed and $y\in \mathbb{R}$, then if $\psi (Y,y)=1{\mathrm{l}}_{\{]-\mathrm{\infty },y]\}}(Y)$ the operator $m_{\psi }(x,y)$ is the conditional cumulative distribution function (df) of Y given X = x, namely $F(y|x)=\mathbb{P}(Y\le y|X=x)$ which may be estimated by ${\hat{F}}_T(y|x):={\hat{m}}_{\psi ,T}(x,y).$ For a given $\alpha \in (0,1),$ the $\alpha \mathrm{th}$ -order conditional quantile of the distribution of Y given X = x is defined as $q_{\alpha }(x)=inf\{y\in \mathbb{R}:F(y|x)\ge \alpha \}.$

Notice that, whenever $F(\cdot |x)$ is strictly increasing and continuous in a neighbourhood of $q_{\alpha }(x)$, the function $F(\cdot |x)$ has a unique quantile of order α at a point $q_{\alpha }(x)$, that is $F(q_{\alpha }(x)|x)=\alpha .$ In such case $q_{\alpha }(x)=F^{-1}(\alpha |x)=inf\{y\in \mathbb{R}:F(y|x)\ge \alpha \},$ which may be estimated uniquely by ${\hat{q}}_{T,\alpha }(x)={\hat{F}}_T^{-1}(\alpha |x)$. Conditional quantiles have been widely studied in the literature when the predictor X is of finite dimension, see for instance, Gannoun et al. (2003) and Ferraty et al. (2005) for dependent functional data.

(a) Almost sure pointwise and uniform convergence

Under the same conditions of Theorem 3.1, the statement (9(9) ${\hat{m}}_{\psi ,T}(x,y)-m_{\psi }(x,y)={\mathcal{O}}_{a.s.}(h_T^{\beta })+{\mathcal{O}}_{a.s.}\left(\sqrt{\frac{\log T}{\mathrm{T\varphi }(h_T)}}\right).$(9) ) still holds for the estimator of the cumulative conditional distribution function ${\hat{F}}_T(y|x)$. That is ${\hat{F}}_T(\alpha |x)$ converges, almost surely, towards $F(y|x)$ with a rate $\mathcal{O}(h_T^{\beta })+\mathcal{O}(\sqrt{\log T/(\mathrm{T\varphi }(h_T))}).$

Consequently, since $F(q_{\alpha }(x)|x)=\alpha ={\hat{F}}_T({\hat{q}}_{T,\alpha }(x)|x)$ and ${\hat{F}}_T(\cdot |x)$ is continuous and strictly increasing, then we have $\mathrm{\forall }ϵ>0,\mathrm{\exists \eta }(ϵ)>0,\mathrm{\forall y},|{\hat{F}}_T(y|x)-{\hat{F}}_T(q_{\alpha }(x)|x)|\le \eta (ϵ)\Rightarrow |y-q_{\alpha }(x)|\le ϵ$ which implies that, $\mathrm{\forall }ϵ>0,\mathrm{\exists \eta }(ϵ)>0,$ (23) $\begin{array}{rl}\mathbb{P}(|{\hat{q}}_{T,\alpha }(x)-q_{\alpha }(x)|\ge \eta (ϵ))& \le \mathbb{P}(|{\hat{F}}_T({\hat{q}}_{T,\alpha }(x)|x)-{\hat{F}}_T(q_{\alpha }(x)|x)|\ge \eta (ϵ))\\ & =\mathbb{P}(|F(q_{\alpha }(x)|x)-{\hat{F}}_T(q_{\alpha }(x)|x)\ge \eta (ϵ)).\end{array}$(23) Therefore, the statement (9(9) ${\hat{m}}_{\psi ,T}(x,y)-m_{\psi }(x,y)={\mathcal{O}}_{a.s.}(h_T^{\beta })+{\mathcal{O}}_{a.s.}\left(\sqrt{\frac{\log T}{\mathrm{T\varphi }(h_T)}}\right).$(9) ) still holds for the conditional quantile estimator ${\hat{q}}_{T,\alpha }(x)$ whenever conditions of Theorem 3.1 are satisfied. Ferraty et al. (2005) derived similar pointwise convergence rate by inverting the estimator of the conditional cumulative distribution function. Their result has been obtained under mixing condition and additional assumptions on the joint distribution, and the Lipschitz condition on $F(y|x)$ and its derivatives with respect to y.

Regarding the almost sure uniform convergence, observe that under conditions of Theorem 3.2, the statement (11(11) $\underset{y\in S}{sup}\underset{x\in \mathcal{C}}{sup}|{\hat{m}}_{\psi ,T}(x,y)-m_{\psi }(x,y)|={\mathcal{O}}_{a.s.}\left(h_T^{\beta }\right)+{\mathcal{O}}_{a.s.}\left(\sqrt{\frac{{\phi }_{\mathcal{C}}({ϵ}_T)}{\mathrm{T\varphi }(h_T)}}\right)\mathrm{as}T\to +\mathrm{\infty }.$(11) ) still holds true for the $\underset{y\in S}{sup}\underset{x\in \mathcal{C}}{sup}|{\hat{F}}_T(y(x)-F_T(y(x)|$, when $\psi (Y,y)$ is replaced by $1{\mathrm{l}}_{\{]-\mathrm{\infty },y]\}}(Y)$. Moreover, assume that, for fixed $x_0\in \mathcal{C}$, $F(y|x_0)$ is differentiable at $q_{\alpha }(x_0)$ with $\frac{\mathrm{\partial }}{\mathrm{\partial }y}F(y|x_0){|}_{y=q_{\alpha }(x_0)}:=g(q_{\alpha }(x_0)|x_0)>\nu >0$, where ν is a real number, and $g(\cdot |x)$ is uniformly continuous for all $x\in \mathcal{C}$. Knowing that ${\hat{F}}_T({\hat{q}}_{T,\alpha }(x)|x)=F(q_{\alpha }(x)|x)=\alpha$ and making use of a Taylor's expansion of the function $F({\hat{q}}_{T,\alpha }(x)|x)$ around $q_{\alpha }(x)$, we can write (24) $F({\hat{q}}_{T,\alpha }(x)|x)-F(q_{\alpha }(x)|x)=({\hat{q}}_{T,\alpha }(x)-q_{\alpha }(x))g(q_{T,\alpha }^{\ast }(x)|x)$(24) where $q_{T,\alpha }^{\ast }(x)$ lies between $q_{\alpha }(x)$ and ${\hat{q}}_{T,\alpha }(x)$. It follows then from (24(24) $F({\hat{q}}_{T,\alpha }(x)|x)-F(q_{\alpha }(x)|x)=({\hat{q}}_{T,\alpha }(x)-q_{\alpha }(x))g(q_{T,\alpha }^{\ast }(x)|x)$(24) ) that the inequality (23(23) $\begin{array}{rl}\mathbb{P}(|{\hat{q}}_{T,\alpha }(x)-q_{\alpha }(x)|\ge \eta (ϵ))& \le \mathbb{P}(|{\hat{F}}_T({\hat{q}}_{T,\alpha }(x)|x)-{\hat{F}}_T(q_{\alpha }(x)|x)|\ge \eta (ϵ))\\ & =\mathbb{P}(|F(q_{\alpha }(x)|x)-{\hat{F}}_T(q_{\alpha }(x)|x)\ge \eta (ϵ)).\end{array}$(23) ) still holds true uniformly in x and y. Moreover, the fact that ${\hat{q}}_{T,\alpha }(x)$ converges a.s. towards $q_{\alpha }(x)$ as T goes to infinity, combined with the uniformly continuity of $g(\cdot |x)$, allow to write that (25) $\underset{x\in \mathcal{C}}{sup}|{\hat{q}}_{T,\alpha }(x)-q_{\alpha }(x)|\underset{x\in \mathcal{C}}{sup}|g(q_{\alpha }(x)|x)|=O_{a.s.}(\underset{y\in S}{sup}\underset{x\in \mathcal{C}}{sup}|{\hat{F}}_T(y|x)-F(y|x)|).$(25) Since $g(q_{\alpha }(x)|x)$ is uniformly bounded from below, we can then claim that the estimator ${\hat{q}}_{T,\alpha }(x)$ converges uniformly towards $q_{\alpha }(x)$ with the same convergence rate given in (11(11) $\underset{y\in S}{sup}\underset{x\in \mathcal{C}}{sup}|{\hat{m}}_{\psi ,T}(x,y)-m_{\psi }(x,y)|={\mathcal{O}}_{a.s.}\left(h_T^{\beta }\right)+{\mathcal{O}}_{a.s.}\left(\sqrt{\frac{{\phi }_{\mathcal{C}}({ϵ}_T)}{\mathrm{T\varphi }(h_T)}}\right)\mathrm{as}T\to +\mathrm{\infty }.$(11) ), as T goes to infinity.

(b) Continuous-time confidence intervals

Confidence intervals for the conditional quantiles $q_{\alpha }(x)$ may be obtained according to the following steps. First, consider a Taylor's expansion of ${\hat{F}}_T(\cdot |x)$ around $q_{\alpha }(x)$ and making use of the fact that ${\hat{q}}_{T,\alpha }(x)$ converges a.s. towards $q_{\alpha }(x)$ as T goes to infinity, one gets (26) ${\hat{q}}_{T,\alpha }(x)-q_{\alpha }(x)=-\frac{1}{{\hat{g}}_T(q_{\alpha }(x)|x)}({\hat{F}}_T(q_{\alpha }(x)|x)-F(q_{\alpha }(x)|x)),$(26) where ${\hat{g}}_T(\cdot |x)$ is a consistent estimator of $g(\cdot |x)$. Then, replacing $\psi (Y,y)$ by the indicator function, we get under conditions of Corollary 3.7, the following $(1-\alpha )$ confidence intervals for $q_{\alpha }(x):$ (27) ${\hat{q}}_{T,\alpha }(x)\pm c_{1-\alpha /2}\frac{\sqrt{M_{T,2}}}{M_{T,1}{\hat{g}}_T({\hat{q}}_{T,\alpha }(x)|x)}\sqrt{\frac{\alpha (1-\alpha )}{T{F}_{x,T}(h_T)p_T(x)}},\mathrm{as}T\to \mathrm{\infty }.$(27)

Acknowledgments

Open Access funding provided by the Qatar National Library.

We thank the editor, associate editor and the two referees for their valuable and constructive comments which helped improve the manuscript substantially.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Notes

1 For any $0\le s<t\le T$ such that $t-s\le {\alpha }_0$, there exists a non-negative continuous random function $f_{t,s}(\omega ,x)$ such that $\underset{s\ge 0,\omega \in \mathrm{\Omega }}{sup}|f_{t,s}(\omega ,x)|\le b_{s,{\alpha }_0}(x)p.s.$ where $b_{s,{\alpha }_0}(x)$ is a deterministic function.

Appendix. Proofs of main results

In this section and for sake of simplification will denote ${\hat{m}}_T(x,y)$ and $m(x,y)$ for ${\hat{m}}_{\psi ,T}(x,y)$ and $m_{\psi }(x,y)$, respectively, and ${\psi }_y(Y_t)$ for $\psi (Y_t,y)$. Consider now the following quantities: (A1) $\begin{array}{rl}{Q}_T(x,y)& :=({\hat{m}}_{T,2}(x,y)-{\overline{m}}_{T,2}(x,y))-m(x,y)({\hat{m}}_{T,1}(x)-{\overline{m}}_{T,1}(x)),\end{array}$(A1) (A2) $\begin{array}{rl}{R}_T(x,y)& :=-B_T(x,y)({\hat{m}}_{T,1}(x)-{\overline{m}}_{T,1}(x)).\end{array}$(A2) We have then (A3) ${\hat{m}}_T(x,y)-m(x,y)={\hat{m}}_T(x,y)-C_T(x,y)+B_T(x,y)=B_T(x,y)+\frac{Q_T(x,y)+R_T(x,y)}{{\hat{m}}_{T,1}(x)}.$(A3) We start first by stating some technical lemmas that will be used later.

Lemma A.1

Assume that assumptions (A1)–(A2) are satisfied, then we have for any $j\ge 1$ and $\ell \ge 1$

1. $\frac{1}{\varphi (h_T)}\mathbb{E}({\mathrm{\Delta }}_t^{\ell }(x)|{\mathcal{F}}_{j-2})=M_{\ell }{f}_{t,T_{j-2}}(x)+{\mathcal{O}}_{a.s.}(\frac{g_{t,T_{j-2},x}(h_T)}{\varphi (h_T)}),$

2. $\frac{1}{\varphi (h_T)}\mathbb{E}({\mathrm{\Delta }}_t^{\ell }(x))=M_{\ell }f(x)+o(1)$.

Proof.

The proof is similar to the proof of Lemma 1 of Laïb and Louani (2010).

Lemma A.2

Let $(Z_n{)}_{n\ge 1}$ be a sequence of real martingale differences with respect to the sequence of σ-fields $({\mathcal{F}}_n=\sigma (Z_1,\dots Z_n){)}_{n\ge 1}$ where $\sigma (Z_1,\dots ,Z_n)$ is the sigma-field generated by the random variables $Z_1,\dots ,Z_n$. Set $S_n=\sum _{i=1}^n{Z}_i.$ For any $p\ge 2$ and any $n\ge 1$, assume that there exist some nonnegative constants C and $d_n$ such that $\mathbb{E}(Z_n^p|{\mathcal{F}}_{n-1})\le C^{p-2}p!d_n^2$ almost surely. Then, for any $ϵ>0,$ we have $\mathbb{P}(|S_n|>ϵ)\le 2\exp \{-\frac{{ϵ}^2}{2(D_n+\mathrm{Cϵ})}\},$ where $D_n=\sum _{i=1}^n{d}_i^2.$

Proof.

See Theorem 8.2.2 of de la Peña and Giné (1999).

Proof

Proof of Theorem 3.4

From (A3(A3) ${\hat{m}}_T(x,y)-m(x,y)={\hat{m}}_T(x,y)-C_T(x,y)+B_T(x,y)=B_T(x,y)+\frac{Q_T(x,y)+R_T(x,y)}{{\hat{m}}_{T,1}(x)}.$(A3) ) and Lemma 1.2 (in Chaouch and Laïb 2023), we have for T large enough (A4) $\begin{array}{rl}\mathrm{MSE}(x,y)& =\mathbb{E}{({\hat{m}}_T(x,y)-m(x,y))}^2\simeq \mathbb{E}{(B_T(x,y)+\frac{Q_T(x,y)+R_T(x,y)}{p(x)})}^2\\ & \simeq \mathbb{E}(B_T^2(x,y))+\frac{1}{p^2(x)}[\mathbb{E}(Q_T^2(x,y))+\mathbb{E}(R_T^2(x,y))],\end{array}$(A4) where the products $2\mathbb{E}[B_T(x,y)(Q_T(x,y)+R_T(x,y))]$ and $2\mathbb{E}[Q_T(x,y)\times R_T(x,y)]$ have been ignored because by the Cauchz–Schwarz inequality $\begin{array}{rl}\mathbb{E}[B_T(x,y)(Q_T(x)+R_T(x,y))]& \le \mathbb{E}{(B_T^2(x,y))}^{1/2}\times \mathbb{E}{((Q_T(x,y)+R_T(x,y))}^{1/2}\\ & \le max\{\mathbb{E}(B_T^2(x,y)),\mathbb{E}[(Q_T(x,y)+R_T(x,y){)}^2]\}.\end{array}$ We have the same inequality for the second product. The proof of Theorem 3.4 results from Proposition 3.3 and Lemma A.3 below, which gives an upper bound of the expectation of $Q_T^2(x,y)$ and $R_T^2(x,y)$, respectively.

Lemma A.3

Assume that (A1)–(A3) hold true, then we have (A5) $\mathbb{E}(Q_T^2(x,y))\simeq \frac{4p(x)(W_2(x,y)+(m(x,y){)}^2)M_2}{\mathrm{T\varphi }(h_T)M_1^{2}f(x)}.$(A5)

Proof.

Ignoring the product term as above, one may write $\begin{array}{rl}\mathbb{E}[Q_T^2(x,y)]& \simeq \mathbb{E}{[{\hat{m}}_{T,2}(x,y)-{\overline{m}}_{T,2}(x,y)]}^2\\ & +m^2(x,y)\mathbb{E}{[{\hat{m}}_{T,1}(x)-{\overline{m}}_{T,1}(x)]}^2:=I_{T,1}+m^2(x,y)I_{T,2}.\end{array}$ The terms $I_{T,1}$ and $I_{T,2}$ can be handled similarly. Let us now evaluate the first one. Since $(T_j=\mathrm{j\delta }{)}_{0\ge j\ge n}$ is a δ-partition of $[0,T]$, we have (A6) $\begin{array}{rl}{\hat{m}}_{T,2}(x,y)-{\overline{m}}_{T,2}(x,y)& =\frac{1}{n\mathbb{E}(Z_1(x))}\sum _{j=1}^n{\int }_{T_{j-1}}^{T_j}[{\zeta }_t{\psi }_y(Y_t){\mathrm{\Delta }}_t(x)-\mathbb{E}\{{\zeta }_t{\psi }_y(Y_t){\mathrm{\Delta }}_t(x)|{\mathcal{F}}_{t-\delta }|\}]\mathrm{d}t\\ & :=\frac{1}{n\mathbb{E}(Z_1(x))}\sum _{j=1}^n{\mathcal{L}}_{T,j}(x,y).\end{array}$(A6) Since $({\mathcal{L}}_{T,j}(x,y){)}_{j\ge 1}$ is a sequence of martingale differences with respect to the family $({\mathcal{F}}_{j-1}{)}_{j\ge 1}$, then

$\mathbb{E}({\mathcal{L}}_{T,j}(x,y)\times {\mathcal{L}}_{T,k}(x,y))=0$ for every $j,k\in \{1,\dots n\}$ such that $j\ne k.$ Therefore (by ignoring the product term), we have (A7) $I_{T,1}=\mathbb{E}{[{\hat{m}}_{T,2}(x,y)-{\overline{m}}_{T,2}(x,y)]}^2\simeq \frac{1}{n^2(\mathbb{E}(Z_1(x)){)}^2}\sum _{j=1}^n\mathbb{E}({\mathcal{L}}_{T,j}(x,y){)}^2.$(A7) Using Jensen inequality and a double conditioning with respect to ${\mathcal{S}}_{t-\delta ,\delta }$ combined with (A3)(iii) –(iv), $I_{T,1}$ may bounded as $\begin{array}{rl}{I}_{T,1}& \le \frac{4}{n^2(\mathbb{E}(Z_1(x)){)}^2}\sum _{j=1}^n{\int }_{T_{j-1}}^{T_j}\mathbb{E}[{\mathrm{\Delta }}_t^2(x)p(X_t)W_2(X_t,y)]\mathrm{d}t.\\ & =\frac{4(p(x)+o(1))(W_2(x,y)+o(1))[M_{2}f(x)+o(1)]}{\mathrm{T\varphi }(h_T)[M_{1}f(x)+o(1){]}^2}.\end{array}$ Similarly, we have $I_{T,2}\le \frac{4(p(x)+o(1))[M_{2}f(x)+o(1)]}{\mathrm{T\varphi }(h_T)[M_{1}f(x)+o(1){]}^2}.$ Therefore $\mathbb{E}(Q_T^2(x,y))\simeq I_{T,1}+m^2(x,y)I_{T,2}=\frac{4p(x)(W_2(x,y)+(m(x,y){)}^2)M_2}{\mathrm{T\varphi }(h_T)M_1^2{f}^2(x)}.$ Moreover, using the decomposition (A2(A2) $\begin{array}{rl}{R}_T(x,y)& :=-B_T(x,y)({\hat{m}}_{T,1}(x)-{\overline{m}}_{T,1}(x)).\end{array}$(A2) ), Theorem 3.3 and Lemma 1.1 (in Chaouch and Laïb 2023) one can see that $\mathbb{E}(Q_T^2(x,y))$ is negligible with respect to $\mathbb{E}(Q_T^2(x,y)).$ This completes the proof.

Proof

Proof of Theorem 3.5

The proof of Theorem 3.5 is based essentially on Lemma A.4 established below, which gives the normality asymptotic of the principal term $Q_T(x,y)$ in (A3(A3) ${\hat{m}}_T(x,y)-m(x,y)={\hat{m}}_T(x,y)-C_T(x,y)+B_T(x,y)=B_T(x,y)+\frac{Q_T(x,y)+R_T(x,y)}{{\hat{m}}_{T,1}(x)}.$(A3) ). Indeed, we have from (A3(A3) ${\hat{m}}_T(x,y)-m(x,y)={\hat{m}}_T(x,y)-C_T(x,y)+B_T(x,y)=B_T(x,y)+\frac{Q_T(x,y)+R_T(x,y)}{{\hat{m}}_{T,1}(x)}.$(A3) ) that (A8) $\sqrt{\mathrm{T\varphi }(h_T)}({\hat{m}}_T(x,y)-m(x,y))=\sqrt{\mathrm{T\varphi }(h_T)}{B}_T(x,y)+\frac{\sqrt{\mathrm{T\varphi }(h_T)}{Q}_T(x,y)+\sqrt{\mathrm{T\varphi }(h_T)}{R}_T(x,y)}{{\hat{m}}_{T,1}(x)}.$(A8) Under (A1)–(A3), Lemma 1.2 (in Chaouch and Laïb 2023) implies that ${\hat{m}}_{T,1}(x)$ converges, almost surely, to $p(x)$ as $T\to \mathrm{\infty }.$ Moreover, using Lemma 1.3 (in Chaouch and Laïb 2023), we get under (A3)(i) –(ii) combined with conditions (12(12) $\underset{T\to \mathrm{\infty }}{lim}\mathrm{T\varphi }(h_T)=\mathrm{\infty },h_T^{\beta }\sqrt{\mathrm{T\varphi }(h_T)}=o(1)\mathrm{and}{h}_T^{\beta }\log T^{1/2}=o(1)\mathrm{as}T\to \mathrm{\infty }.$(12) ) that $\sqrt{\mathrm{T\varphi }(h_T)}{B}_T(x,y)={\mathcal{O}}_{a.s.}(h_T^{\beta }\sqrt{\mathrm{T\varphi }(h_T)})=o_{a.s}(1),$ and $\sqrt{\mathrm{T\varphi }(h_T)}{R}_T(x,y)={\mathcal{O}}_{a.s.}\left(\sqrt{\mathrm{T\varphi }(h_T)}{h}_T^{\beta }{\left(\frac{\log T}{\mathrm{T\varphi }(h_T)}\right)}^{1/2}\right)={\mathcal{O}}_{a.s.}(h_T^{\beta }\log T^{1/2})=o_{a.s.}(1).$ The proof may be then achieved by Lemma A.4 and Slutsky's Theorem.

Lemma A.4

Under conditions (A1)–(A3), we have $\begin{array}{rl}& \sqrt{\mathrm{T\varphi }(h_T)}({\hat{m}}_T(x,y)-m(x,y))\stackrel{d}{\to }N(0,{\tilde{\sigma }}^2(x,y))\mathrm{as}T⟶\mathrm{\infty }\\ \mathrm{where}& {\tilde{\sigma }}^2(x,y)\le \frac{M_2}{M_1^{2}f(x)p(x)}{\overline{W}}_2(x,y).\end{array}$

Proof

Proof of Lemma A.4

We have (A9) $\begin{array}{rl}\sqrt{\mathrm{T\varphi }(h_T)}{Q}_T(x,y)& =\sum _{i=1}^n{\xi }_{T,i}(x,y),\mathrm{with}{\xi }_{T,i}(x,y)={\eta }_{T,i}(x,y)-\mathbb{E}[{\eta }_{T,i}(x,y)|{\mathcal{F}}_{t-\delta }]\mathrm{and}\\ {\eta }_{T,i}(x,y)& =\frac{1}{\mathbb{E}{Z}_1}\sqrt{\frac{{\varphi }_T(h)}{n}}{\int }_{T_{i-1}}^{T_i}{\zeta }_t{\mathrm{\Delta }}_t(x)[{\psi }_y(Y_t)-m(x,y)]\mathrm{d}t.\end{array}$(A9) Since for any $i\ge 1$ and $t\in [T_{i-1},T_i]$, ${\mathcal{F}}_{i-2}\subset {\mathcal{F}}_{t-\delta }\subset {\mathcal{F}}_{i-1}$, then $({\xi }_{T,i}(x,y){)}_{i\ge 1}$ is ${\mathcal{F}}_{i-1}$-measurable, $\mathbb{E}(|{\xi }_{T,i}|)<\mathrm{\infty }$ provided $\mathbb{E}({\zeta }_t^2)<\mathrm{\infty }$ and $\mathbb{E}(X_t^2)<\mathrm{\infty }$. Moreover, we have for any $1\le i\le n$, $\mathbb{E}({\xi }_{T,i}|{\mathcal{F}}_{i-2})=\mathbb{E}\{\mathbb{E}[{\eta }_i|{\mathcal{F}}_{t-\delta }]|{\mathcal{F}}_{i-2}\}-\mathbb{E}\{\mathbb{E}[{\eta }_i|{\mathcal{F}}_{t-\delta }]|{\mathcal{F}}_{i-2}\}=0$ a.s. Hence $({\xi }_{T,i}(x,y){)}_{i\ge 1}$ is a sequence of martingale differences with respect to the σ-fields $({\mathcal{F}}_{i-1}{)}_{i\ge 1}$. To prove the asymptotic normality, it suffices to check both following conditions (see Corollary 3.1, p. 56, Hall and Heyde 1980):

(a) $\sum _{i=1}^n\mathbb{E}[{\xi }_{T,i}^2(x,y)|{\mathcal{F}}_{i-2}]\stackrel{\mathbb{P}}{\to }{\tilde{\sigma }}^2(x,y)$ and (b) $n\mathbb{E}[{\xi }_{T,i}^2(x,y)1{\mathrm{l}}_{\{|{\xi }_{T,i}(x,y)|>ϵ\}}]=o(1)$ holds for any $ϵ>0$.

Proof of (a) Observe now that $|\sum _{i=1}^n\mathbb{E}[{\eta }_{T,i}^2(x,y)|{\mathcal{F}}_{i-2}]-\sum _{i=1}^n\mathbb{E}[{\xi }_{T,i}^2(x,y)|{\mathcal{F}}_{i-2}]|\le |\sum _{i=1}^n{\left(\mathbb{E}[{\eta }_{T,i}(x,y)|{\mathcal{F}}_{i-2}]\right)}^2|.$ Using (A1), (A3)(i),(i${}^{\mathrm{\prime }}$), (ii) and (iv) with Lemma A.1, and a double conditioning with respect to the σ-field ${\mathcal{S}}_{t-\delta ,\delta }$ and the fact that $n\mathbb{E}{Z}_1(t)=O(\mathrm{T\varphi }(h))$, we have $\begin{array}{rl}& \left|\mathbb{E}\left({\eta }_{T,i}|{\mathcal{F}}_{i-2}\right)\right|\\ & =\frac{1}{\mathbb{E}{Z}_1}\sqrt{\frac{{\varphi }_T(h)}{n}}|{\int }_{T_{i-1}}^{T_i}\mathbb{E}(p(X_t){\mathrm{\Delta }}_t(x)[m(X_t,y)-m(x,y)]|{\mathcal{F}}_{i-2})\mathrm{d}t|\\ & \le \sqrt{n{\varphi }_T(h)}\frac{p(x)}{n\mathbb{E}{Z}_1}\underset{u\in \mathcal{B}(x,h)}{sup}|m(u,y)-m(x,y)|\underset{u\in \mathcal{B}(x,h)}{sup}|p(u)-p(x)||{\int }_{T_{i-1}}^{T_i}\mathbb{E}({\mathrm{\Delta }}_t(x)|{\mathcal{F}}_{i-2})\mathrm{d}t|\\ & =O\left(\sqrt{n{\varphi }_T(h)}{h}^{\beta }\right)O(\varphi (h_T){\int }_{T_{i-1}}^{T_i}{f}_{t,T_{i-2}}(x)\mathrm{d}t+o(1))\frac{1}{n\mathbb{E}{Z}_1}\\ & =O\left(\sqrt{n{\varphi }_T(h)}{h}^{\beta }\right)O(\frac{1}{T}{\int }_{T_{i-1}}^{T_i}{b}_{t,{\alpha }_0}(x)\mathrm{d}t)).\end{array}$ It follows by (A2)-(iii) and the Cauchy–Schwarz inequality that $\sum _{i=1}^n{\left(\mathbb{E}[{\eta }_{T,i}(x,y)|{\mathcal{F}}_{i-2}]\right)}^2=O(h^{2\beta }\varphi (h))=o(1).$ Thus we have only to show that $\sum _{i=1}^n\mathbb{E}[{\eta }_{T,i}^2(x,y)|{\mathcal{F}}_{i-2}]\stackrel{\mathbb{P}}{\to }{\sigma }^2(x,y)$. Using again the Cauchy–Schwarz inequality, one may write (A10) $\begin{array}{rl}\sum _{i=1}^n\mathbb{E}[{\eta }_{T,i}^2(x,y)|{\mathcal{F}}_{i-2}]& =\frac{1}{(\mathbb{E}{Z}_1{)}^2}\frac{{\varphi }_T(h)}{n}\sum _{i=1}^n\mathbb{E}\left[{({\int }_{T_{i-1}}^{T_i}{\zeta }_t{\mathrm{\Delta }}_t(x)[{\psi }_y(Y_t)-m(x,y)])}^2|{\mathcal{F}}_{i-2}\right]\\ & \le \frac{\delta }{(\mathbb{E}{Z}_1{)}^2}\frac{{\varphi }_T(h)}{n}\sum _{i=1}^n\mathbb{E}\left[{\int }_{T_{i-1}}^{T_i}{\zeta }_t^2{\mathrm{\Delta }}_t^2(x)[{\psi }_y(Y_t)-m(x,y){]}^2\mathrm{d}t|{\mathcal{F}}_{i-2}\right]\\ & =\frac{\delta }{(\mathbb{E}{Z}_1{)}^2}\frac{{\varphi }_T(h)}{n}\sum _{i=1}^n\mathbb{E}\left[{\int }_{T_{i-1}}^{T_i}{\zeta }_t^2{\mathrm{\Delta }}_t^2(x)[{\psi }_y(Y_t)-m(X_t,y){]}^2\mathrm{d}t|{\mathcal{F}}_{i-2}\right]\\ & +\frac{\delta }{(\mathbb{E}{Z}_1{)}^2}\frac{{\varphi }_T(h)}{n}\sum _{i=1}^n\mathbb{E}\left[{\int }_{T_{i-1}}^{T_i}{\zeta }_t{\mathrm{\Delta }}_t^2(x)[m(X_t,y)-m(x,y){]}^2\mathrm{d}t|{\mathcal{F}}_{i-2}\right]\\ & =:A_n+C_n.\end{array}$(A10) Now, let us evaluate the term $A_n$. Conditioning three times with respect to ${\mathcal{F}}_{t-\delta }$ and ${\mathcal{S}}_{t-\delta }$, and making use of Conditions (A3)(i${}^{\mathrm{\prime }}$), (iii), (iv) and the fact that $T=\mathrm{n\delta }$, to get from Lemma A.1 that (A11) $\begin{array}{rl}{A}_n& =\frac{\delta }{(\mathbb{E}{Z}_1{)}^2}\frac{{\varphi }_T(h)}{n}\sum _{i=1}^n\mathbb{E}[{\int }_{T_{i-1}}^{T_i}p(X_t){\mathrm{\Delta }}_t^2(x){\overline{W}}_2(X_t,y)\mathrm{d}t|{\mathcal{F}}_{i-2}]\\ & \le \frac{\delta }{(\mathbb{E}{Z}_1{)}^2}\frac{{\varphi }_T(h)}{n}(p(x)+o(1))|({\overline{W}}_2(x,y)+o(1))\sum _{i=1}^n\mathbb{E}[{\int }_{T_{i-1}}^{T_i}{\mathrm{\Delta }}_t^2(x)\mathrm{d}t|{\mathcal{F}}_{i-2}]\\ & \le (\delta +o(1))p(x){\overline{W}}_2(x,y)\frac{{\varphi }_T^2(h)}{(\mathbb{E}{Z}_1{)}^2}\frac{1}{n}\sum _{i=1}^n{\int }_{T_{i-1}}^{T_i}\mathbb{E}[\frac{1}{{\varphi }_T(h)}{\mathrm{\Delta }}_t^2(x)|{\mathcal{F}}_{i-2}]\mathrm{d}t\\ & \le \delta (\delta +o(1))p(x){\overline{W}}_2(x)\frac{{\varphi }_T^2(h)}{(\mathbb{E}(Z_1{)}^2}{M}_2\{\frac{1}{T}\sum _{i=1}^n{\int }_{T_{i-1}}^{T_i}{f}_{t,T_{i-2}}(x)\mathrm{d}t\\ & +O_{a.s.}\left[\frac{1}{T}\sum _{i=1}^n{\int }_{T_{i-1}}^{T_i}\frac{g_{t,T_{i-2},x}(h_T)}{{\varphi }_T(h)}\mathrm{d}t\right]\}\end{array}$(A11) The Riemann's sum combined with condition (A2)(iii) gives that $\frac{1}{T}\sum _{i=1}^n{\int }_{T_{i-1}}^{T_i}{f}_{t,T_{i-2}}(x)\mathrm{d}t\le \frac{1}{T}{\int }_0^T{f}_{t,T_{t-\delta }}(x)\mathrm{d}t⟶f(x)\mathrm{a}.\mathrm{s}.\mathrm{as}T⟶\mathrm{\infty }.$ Moreover, by (A2)(ii) one gets $\frac{g_{t,T_{i-2},x}(h_T)}{{\varphi }_T(h)}=o(1)\mathrm{as}T⟶\mathrm{\infty }$. Therefore, we have (A12) $\begin{array}{rl}{A}_n& \le \delta (\delta +o(1))p(x){\overline{W}}_2(x,y)\frac{{\varphi }_T^2(h)}{(\delta {\varphi }_T(h_T)M_{1}f(x)+o(1){)}^2}{M}_2[f(x)+o(1)]\\ & =\frac{M_2}{M_1^{2}f(x)}p(x){\overline{W}}_2(x,y):={\tilde{\sigma }}^2(x,y)\mathrm{as}T⟶\mathrm{\infty }.\end{array}$(A12) On the other hand, by the same arguments as above combined with the fact that $\underset{u\in \mathcal{B}(x,h)}{sup}|m(x)-m(u)|=O(h^{2\beta }{)}_{a.s.}$, we get $C_n=o_{a.s.}(1)$.

Proof of part (b) Using successively Hölder, Markov, Jensen and Minkowski inequalities combined with conditions (A3)(iii), (A3)(iv) and Lemma A.1, we get, for any $ϵ>0$ and fixed real numbers p>1 and q>1 such that $1/p+1/q=1$, $n\mathbb{E}[{\xi }_{T,i}^2(x,y)1{\mathrm{l}}_{\{|{\xi }_{T,i}(x,y)|>ϵ\}}]\le 4n(ϵ/2{)}^{-2q/p}\mathbb{E}[|{\eta }_{T,i}{|}^{2q}]=O\left({(\mathrm{T\varphi }(h_T))}^{-\gamma /2}\right)=o_{a.s}(1)$ by taking $2q=2+\gamma$ ($0<\gamma <1$), since $\mathrm{T\varphi }(h_T)$ towards to infinity as T goes to infinity.

Proof

Proof of Corollary 3.7

Observe that (A13) $\begin{array}{rl}& \sqrt{\frac{T{F}_{x,T}(h_T)}{{\tilde{V}}_T^2(x,y)}}({\hat{m}}_T(x,y)-m(x,y))\\ & =\sqrt{\frac{F_{x,T}(h_T)}{\varphi (h_T)f(x)}}\sqrt{\frac{{\sigma }^2(x,y)f(x)}{{\tilde{V}}_T^2(x,y)}}\sqrt{\frac{\mathrm{T\varphi }(h_T)}{{\sigma }^2(x,y)}}({\hat{m}}_T(x,y)-m(x,y)).\end{array}$(A13) It follows from the consistency of $F_{x,T}(h_T)$ and (A2)(i) that $\frac{F_{x,T}(h_T)}{\varphi (h_T)f(x)}$ goes to 1 a.s. as T goes to infinity. By Theorem 3.5, the quantity $\sqrt{\frac{\mathrm{T\varphi }(h_T)}{\sigma (x,y)}}({\hat{m}}_T(x,y)-m(x,y))$ converges to $N(0,1)$ as $T\to \mathrm{\infty }$. Then using the non-decreasing property of the cumulative standard normal distribution function Ψ, we get, for a given risk $0<\alpha <1$, the $(1-\alpha )$- pseudo-confidence interval (A14) $\sqrt{\frac{\mathrm{T\varphi }(h_T)}{{\sigma }^2(x,y)}}|{\hat{m}}_T(x,y)-m(x,y)|\le {\mathrm{\Psi }}^{-1}(1-\frac{\alpha }{2}).$(A14) Considering now the statement (13(13) ${\sigma }^2(x,y)\le \frac{1}{f(x)}\frac{M_2}{M_1^{2}p(x)}{\overline{W}}_2(x,y)=:\frac{1}{f(x)}\tilde{V}(x,y)\mathrm{as}T⟶\mathrm{\infty }$(13) ) combined with Proposition 3.6, it holds that (A15) $\underset{T\to \mathrm{\infty }}{lim}\frac{{\sigma }^2(x,y)f(x)}{{\tilde{V}}_n^2(x,y)}\le \underset{T\to \mathrm{\infty }}{lim}\frac{f(x)V^2(x,y)}{{\tilde{V}}_n^2(x,y)}=\underset{T\to \mathrm{\infty }}{lim}\frac{{\tilde{V}}^2(x,y)}{{\tilde{V}}_n^2(x,y)}=1\mathrm{a}.\mathrm{s}.,$(A15) since ${\tilde{V}}_n^2(x,y)$ is a consistent estimator of ${\tilde{V}}^2(x,y)$. The proofs follows then from the statements (A13(A13) $\begin{array}{rl}& \sqrt{\frac{T{F}_{x,T}(h_T)}{{\tilde{V}}_T^2(x,y)}}({\hat{m}}_T(x,y)-m(x,y))\\ & =\sqrt{\frac{F_{x,T}(h_T)}{\varphi (h_T)f(x)}}\sqrt{\frac{{\sigma }^2(x,y)f(x)}{{\tilde{V}}_T^2(x,y)}}\sqrt{\frac{\mathrm{T\varphi }(h_T)}{{\sigma }^2(x,y)}}({\hat{m}}_T(x,y)-m(x,y)).\end{array}$(A13) ), (A14(A14) $\sqrt{\frac{\mathrm{T\varphi }(h_T)}{{\sigma }^2(x,y)}}|{\hat{m}}_T(x,y)-m(x,y)|\le {\mathrm{\Psi }}^{-1}(1-\frac{\alpha }{2}).$(A14) ) and (A15(A15) $\underset{T\to \mathrm{\infty }}{lim}\frac{{\sigma }^2(x,y)f(x)}{{\tilde{V}}_n^2(x,y)}\le \underset{T\to \mathrm{\infty }}{lim}\frac{f(x)V^2(x,y)}{{\tilde{V}}_n^2(x,y)}=\underset{T\to \mathrm{\infty }}{lim}\frac{{\tilde{V}}^2(x,y)}{{\tilde{V}}_n^2(x,y)}=1\mathrm{a}.\mathrm{s}.,$(A15) ).