Nonparametric Functional Least Absolute Relative Error Regression: Application to Econophysics

1. Introduction

Examining how a functional predictor co-varies with a scalar response variable is a fundamental issue in functional statistics. This relationship is commonly evaluated through regression models, which is estimated using the least squares loss function. In this paper, we adopt another smoothing approach based on the Least Absolute Relative Error (LARE). Undoubtedly, the LARE criterion improves the robustness and the accuracy of the standard regression approaches.

The general framework of this paper is the nonparametric prediction in functional statistics, which has attracted considerable interest due to numerous valuable contributions in recent years. The foundational work in this field was introduced by [1], who established the almost complete pointwise consistency of the kernel-type regression estimator under the assumption of independent and identically distributed (i.i.d.) observations. [2] proved the convergence of this estimator in the $L^p$ norm. [3] derived the exact asymptotic expression for the corresponding $L^p$ errors. The asymptotic normality of the same estimator in the strong mixing case was later demonstrated by [4]. In [5], the authors established a uniform version of the almost complete convergence rate in the i.i.d. framework. In addition to these classical kernel-based approaches, recent research has explored M-estimation techniques and/or local linear smoothing methods. Significant results in this field of nonparametric functional prediction were presented by [6,7,8,9]., which provide comprehensive references on this topics.

Alternatively, in this paper, we focus on relative error regression (RE-regression), a methodology in which the the relative error is used. The latter constitutes a performance measure in many practical applications. Despite its utility, the literature on nonparametric RE-regression remains limited. The most existing studies adopting parametric approaches. In this context, [9] were among the first to investigate the use of relative squared error in estimation. From a practical point of view, RE-regression has been applied in diverse fields, including medicine [10] and finance [11]. In [12], authors studied estimation via RE-regression in multiplicative regression models, while [13] extended this framework to nonparametric estimation, establishing the convergence of the local linear estimator derived from relative error loss.

In last decade, RE-regression has been explored extensively for dependent data. For instance, [14] considered quasi-associated time series, and [15] studied spatial processes. The functional nonparametric RE-regression was first introduced by [16], who established strong consistency and derived the asymptotic distribution of the RE-regression estimator. [17] extended the method to incomplete functional data, examining kernel smoothing techniques for truncated observations. The quasi-associated functional time series case was investigated by[18]. For a comprehensive overview of the latest progress on this topic, see [19,20].

While all the cited studies estimate the regression operator using the mean least squares relative error loss function, in the present contribution, we employ the absolute relative error loss function. This approach is combined with kernel-based local weighting techniques to develop a robust estimator for the nonparametric functional regression model. This alternative estimator offers substantial gains in robustness and accuracy compared to the earlier approach. Such an improvement is especially important in functional statistics, where data are typically noisy and high-dimensional. Specifically, in functional statistics, extreme values and rapid fluctuations are common, which can disproportionately influence least squares-based estimators. Unlike the Least Squares Relative Error (LSRE), the Absolute Relative Error (ARE) treats deviations proportionally, making it more robust to outliers. This robustness ensures a more reliable estimation of the underlying functional relationship, making ARE particularly suitable for high-frequency applications in many applied areas such as economics, finance, signal processing, and environmental monitoring. From a theoretical point of view, the manipulation of the LARE-regression model is more challenging than LSRE-regression. This difficulty is because the LSRE-regression can be explicitly expressed as the ratio of the first and second inverted conditional moments, whereas the LARE-regression lacks such an explicit form. Its asymptotic properties must be derived through its Bahadur representation. Despite this additional complexity, we have established the almost complete convergence of the proposed estimator with precision of the convergence rate. The latter is stated under standard assumptions in nonparametric functional statistics. The assumed conditions explore the functional nature of the model and the structure of the data. An additional advantage of the LARE-regression model is its broader scope of applicability compared to the LSRE-regression. Indeed, the LSRE-regression is limited to strictly positive response variables; in contrast, the LARE regression overcomes this restriction allowing it to extend its potential applicability. Furthermre, the practical usefulness of the proposed estimator is demonstrated through data-driven selection procedures for tuning parameters, as well as through simulation studies and real data applications, which confirm its robustness and superior performance in practice.

The paper is organized as follows. In Section 2, we introduce the proposed estimation algorithms. Section 3 presents the required assumptions and the main asymptotic results. In Section 4, we discuss methods for selecting the smoothing parameter and the metric using the random matrix theory from econophysics. The performance of the constructed estimator on simulated and real data is evaluated in Section 5. Finally, Section 6 is devoted to summarize some concluding remarks. The technical proofs are provided in the appendix’s Appendix A.

2. Functional LARE Regression and Its Estimation

As outlined in the introduction, our main objective is to evaluate the relationship between an exogenous functional variable T and a real-valued endogenous variable S. More precisely, the variables under consideration belong to the space $\mathcal{X}\times \mathrm{I}\mathrm{R}$. The set $\mathcal{X}$ is a functional space endowed with an appropriate semi-norm N. Henceforth. we fix a point $\mathcal{T}$ in $\mathcal{X}$ and consider a neighborhood ${\mathcal{N}}_{\mathcal{T}}$ of $\mathcal{T}$ within this functional space. In our earlier contributions ([16,18]), we modeled the underlying relationship between T and S using the LSRE regression defined by

$$\begin{array}{c}\hfill REG\left(\mathcal{T}\right)=arg\underset{s}{min}\mathrm{I}\mathrm{E}\left({\left({\displaystyle \frac{S-s}{S}}\right)}^2|T=\mathcal{T}\right).\end{array}$$

The main feature of this loss function is its ability to provide an explicit definition of the predictor. Indeed, by simple differentiation with respect to $\theta$, we can show that

$$REG\left(\mathcal{T}\right))={\displaystyle \frac{\mathrm{I}\mathrm{E}\left[S^{-1}|T=\mathcal{T}\right]}{\mathrm{I}\mathrm{E}\left[S^{-2}|T=\mathcal{T}\right]}}.$$

(1)

Clearly, such an explicit definition of the model simplifies the analysis and derivation of its mathematical properties. Nevertheless, this approach is not without limitations. While the least squares relative error method reduces the effect of large outliers values, it can give excessive weight to large deviations which compromised its robustness against the very small values or highly dispersed observations. To address this limitation, we propose to proceed with the least absolute relative error (LARE) approach, which provides a more balanced treatment of both small and large deviations in the data. Formally, the LARE regression is defined as the estimator that minimizes

$$\begin{array}{c}\hfill RER\left(\mathcal{T}\right)=arg\underset{s}{min}\mathrm{I}\mathrm{E}\left(\left|{\displaystyle \frac{S-s}{S}}\right||T=\mathcal{T}\right).\end{array}$$

(2)

Evidently, this alternative loss function enhances the robustness of LSRE regression by treating all deviations proportionally in a linear manner, rather than using the traditional quadratic form. The LARE loss function evaluates the differences between predicted and true values more equitably, providing a more balanced, resilient, and reliable assessment of model performance. This advantage is particularly important in the presence of large deviations or extremely small observations, where the LSRE approach fails to handle. Overall, the LARE rule is more robust and effective, especially in datasets with irregular or heavy-tailed distributions.

For the estimation setup, we consider n independent and identically distributed (i.i.d.) pairs $(T_1,S_1),\dots ,(T_n,S_n)$, sampled from the distribution of $(T,S)$. The kernel estimator of $RER\left(\mathcal{T}\right)$ is then given by

$${\displaystyle \widehat{RER}\left(\mathcal{T}\right)=arg\underset{s}{min}\frac{{\sum }_{i=1}^{n}M\left(m^{-1}N(\mathcal{T},T_i)\right)\left|\frac{S_i-s}{S_i}\right|}{{\sum }_{i=1}^{n}M\left(m^{-1}N(\mathcal{T},T_i)\right)},}$$

where M is a kernel function and $m:=m_n$ (to simplify the notation) is a sequence of positive real numbers that goes to zero as n goes to infinity.

The primary objective of this paper is to establish the asymptotic properties of the estimator $\widehat{RER}\left(\mathcal{T}\right)$ of $RER\left(\mathcal{T}\right)$, where the explanatory variable T takes values in a semi-metric space $\mathcal{X}$. To the best of our knowledge, this work is the first to employ this approach to construct an estimator of the regression operator, even in multivariate settings. Although the finite-dimensional case can be viewed as a special case of this framework, the infinite-dimensional case is especially interesting, both for its theoretical challenges and the variety of practical applications it supports. For further motivation and detailed discussions on this topic, we refer the reader to [21,22,23] and the references therein.

3. The Strong Consistency of the Kernel Estimator

Now, to investigate the almost complete convergence, we denote by C or $C^{\prime }$ strictly positive constants and put

$$M_i=M\left(m^{-1}N(\mathcal{T},T_i)\right),i=1,\dots ,n.$$

We define the closed ball of radius r centered at $\mathcal{T}$ as

$$B(\mathcal{T},r)=\left\{\mathcal{Z}\in \mathcal{X}:N(\mathcal{Z},\mathcal{T})<r\right\},$$

and set

$$U_1(s,\mathcal{T})={\mathbb{E}\left[\right|S|}^{-1}1{\mathrm{I}}_{S\le s}\mid T=\mathcal{T}],.U_2(s,\mathcal{T})={\mathbb{E}\left[\right|S|}^{-1}1{\mathrm{I}}_{S\ge s}\mid T=\mathcal{T}].$$

The following assumptions are imposed:

(HM1)

The small-ball probability satisfies $\mathbb{P}(X\in B(\mathcal{T},r))=:{\mathcal{G}}_{\mathcal{T}}\left(r\right)>0$ for all $r>0$, with

$$\underset{r\to 0}{lim}{\mathcal{G}}_{\mathcal{T}}\left(r\right)=0.$$

(HM2)

For all ${\mathcal{T}}_1,{\mathcal{T}}_2\in {\mathcal{V}}_{\mathcal{T}}$, the functions $U_{\gamma }$ are of class $C^1$ (with respect to s ) and satisfy

$$|U_{\gamma }(s,{\mathcal{T}}_1)-U_{\gamma }(s,{\mathcal{T}}_2)|\le C{d}^{b_{\gamma }}({\mathcal{T}}_1,{\mathcal{T}}_2),b_{\gamma }>0.$$

(HM3)

The kernel M is measurable, supported on $(0,1)$, and bounded: $0<C_2\le M(·)\le C_3<\infty$.

(HM4)

The small-ball probability satisfies

$${\displaystyle \frac{n{\mathcal{G}}_{\mathcal{T}}\left(m\right)}{logn}}\to \infty \mathrm{as}n\to \infty .$$

(HM5)

The response variable is bounded by inverse moments:

$$\forall p\ge 2,\mathbb{E}[S^{-p}\mid T=\mathcal{T}]<C<\infty .$$

These assumptions are generally mild and consistent with those commonly adopted in the literature on functional mode estimation. The functional component is typically described by two main conditions. The first, (HM1), concerns the small-ball probability function, which captures the concentration of the probability measure of the functional regressor. The second, (HM2), addresses the time-series structure of the functional data, quantifying the local dependence among functional observations. Both of these conditions are standard in functional time series analysis (see [24] for detailed discussions and examples of functional data satisfying these assumptions). The nonparametric component is characterized by (Hy3), which is commonly used in nonparametric functional statistics and has appeared in several previous studies. Finally, assumptions (HM3)–(HM5) impose technical conditions that facilitate the proofs and ensure the validity of the asymptotic results.

Theorem 1.

Given assumptions (HM1)–(HM5), we obtain

$$|\widehat{RER}\left(\mathcal{T}\right)-RER\left(\mathcal{T}\right)\left)\right|=O\left(m^{b_1}\right)+O\left(m^{b_2}\right)+O_{a.co.}\left(\sqrt{{\displaystyle \frac{logn}{n{\mathcal{G}}_{\mathcal{T}}\left(m\right)}}}\right).$$

(3)

The proof of this Theorem employs the Bahadur representation of the $\widehat{RER}$ estimator. This representation is established through the following auxiliary results.

Lemma 1.

Let $\left\{{\mathcal{U}}_n\right\}$ be a sequence of decreasing real-valued random functions and $\left\{{\mathcal{D}}_n\right\}$ a real random sequence satisfying

$$\begin{array}{cc}\hfill {\mathcal{D}}_n& =o_{a.co.}\left(1\right),\hfill \\ \hfill \underset{\left|d\right|\le A}{sup}\left|{\mathcal{U}}_n\left(d\right)+\beta d-{\mathcal{D}}_n\right|& =o_{a.co.}\left(1\right),\hfill \end{array}$$

 for some constants $\beta ,A>0$. Then, for any real random sequence $d_n$ such that ${\mathcal{U}}_n\left(d_n\right)=o_{a.co.}\left(1\right)$, we have

$$\sum _{n=1}^{\infty }\mathbb{P}\left(|d_n|\ge A\right)<\infty .$$

4. Smoothing Parameter Selection and Random Matrix Theory for Metric Approximation

The practical performance of the estimator $\tilde{R}$ depends crucially on the choice of parameters involved in its construction. In particular, the bandwidth parameter m plays a central role in the implementation of this regression method. Various strategies have been proposed in the nonparametric regression literature to address this issue. In this study, we adopt two widely used approaches from classical regression: the cross-validation rule and the bootstrap algorithm.

4.1. Leave-One-Out Cross-Validation Principle

In regression analysis, the leave-one-out cross-validation (LOOCV) criterion is based on minimizing the mean squared error. This criterion is highly effective for bandwidth selection because it directly optimizes for predictive accuracy. The LOOCV algorithm was successfully extended to the domain of functional statistics by [25] who demonstrated its particular efficacy and suitability for functional data. The main advantage of this approach is its flexibility; for instance, it can be readily adapted to our model for functional LARE regression as

$$m^{opt}=arg\underset{m\in \mathcal{M}}{min}\sum _{i=1}^n\left|{\displaystyle \frac{(S_i-{\widetilde{RER}}^{-i}\left(T_i\right))}{S_i}}\right|,$$

(4)

where ${\widetilde{RER}}^{-i}\left(T_i\right)$ and ${\widehat{RER}}^{-i}\left(T_i\right)$ denote the leave-one-out estimators of $RER$, respectively, computed without the i-th observation $(T_i,S_i)$. The efficiency of these estimators depends on the choice of the subset $\mathcal{M}$, which determines the optimization of (4).

Concerning the subset $\mathcal{M}$, we point out that there are two primary strategies for choosing it. The local approach and the global one. In the local case, $\mathcal{M}$ is determined by the number of neighborhoods around the location point. In the global case, $\mathcal{M}$ is derived from the quantiles of the distance distribution between the functional regressors.

4.2. Bootstrap Approach

The bootstrap procedure is an alternative popular method to the cross-validation rule. Its principal feature is its data-driven nature, as it empirically approximates the sampling distribution of the estimator without external assumptions. This automatism analysis makes it particularly valuable in complex modeling scenarios, where asymptotic theory may be intractable or slow to converge. Furthermore, the bootstrap provides a direct and computationally feasible means to minimize a targeted loss function. This allows for a more adaptive for many loss functions. For the LARE- regression, we state the following algorithm :

Step 1:

Choose an arbitrary bandwidth $m_0$ and compute ${\widehat{RER}}_{m_0}\left(\mathcal{T}\right)$.

Step 2:

Compute the residual $\widehat{ϵ}=S-{\widehat{RER}}_{m_0}\left(\mathcal{T}\right)$.

Step 3:

Generate a bootstrap residual ${ϵ}^{*}$ using the distribution. $G^{*}=((\sqrt{5}+1)/2\sqrt{5})d_{\widehat{ϵ}(1-\sqrt{5})/2}-((1-\sqrt{5})/2\sqrt{5})d_{\widehat{ϵ}(\sqrt{5}+1)/2}$, (d is dirac measure ).

Step 4:

Construct a bootstrap sample ${(S_i^{*},T_i^{*})}_i=({ϵ}^{*}-{\widehat{RER}}_{m_0}\left(T_i\right)$

Step 5:

Compute the estimators ${\widehat{RER}}_m\left(\mathcal{T}\right)$ using the sample ${(S_i^{*},T_i^{*})}_i$.

Step 6:

Repeat the steps 3-6 $N_B$ times and we put ${\widehat{RER}}_m^r\left(\mathcal{T}\right)$ the estimator at the replication r.

Step 7:

Choose the optimal bandwidth m according to the following criterion

$$m_{op{t}^{Boo}}=arg\underset{m}{min}\sum _{r=1}^{N_B}\sum _{i=1}^n|{\widehat{RER}}_{m_0}\left(T_i\right)-{\widehat{RER}}_m\left(T_i\right)|.$$

4.3. Random Matrix Metric

The principal task of functional data analysis is to move from comparing two curves directly to comparing them within the context of a larger, noisy system. In this context, Random Matrix Theory (RMT) is used to define a stable distance (without noise ) between any two curves by eliminating the noise of the entire system. The main idea is based on the use of RMT to clean the covariance matrix estimated from the data. The eigenvalues and eigenvectors of the noisy covariance matrix are a mixture of true signal and random noise. RMT provides a null model to identify which parts of the spectrum are likely to be signal and which are noise. To do that, we follow the following algorithm:

• Step 1:
  Create the sample covariance matrix

  $$\mathbf{C}={\displaystyle \frac{1}{n}}{\mathbf{T}}^{tr}\mathbf{T}$$
  Here $\mathbf{T}$ is assumed to have been column-centered. This matrix $\mathbf{C}$ gives the observed covariances between all pairs of sampling points.
• Step 2:
  The Marčenko-Pastur law describes the eigenvalue distribution for a random matrix. We compare the eigenvalue spectrum of our empirical $\mathbf{C}$ to the spectrum predicted by the Marčenko-Pastur distribution for a random matrix with the same dimensions and variance.
  • Perform the eigen-decomposition: $\mathbf{C}=\mathbf{V}\mathrm{\Lambda }{\mathbf{V}}^{tr}$, where $\mathrm{\Lambda }$ is the diagonal matrix of eigenvalues ${\lambda }_i$.
  • Identify the eigenvalues that are outside the support of the random distribution. These are considered as signal eigenvalues. The eigenvalues within the random bounds are considered noise.
• Step 3:
  Create a cleaned covariance matrix ${\mathbf{C}}_{\mathrm{clean}}$ by retaining only the signal components. A common method is to replace the noise eigenvalues with a constant value:

  $${\mathrm{\Lambda }}_{\mathrm{clean}}=\mathrm{diag}\left(max({\lambda }_1,{\overline{\lambda }}_{\mathrm{noise}}),max({\lambda }_2,{\overline{\lambda }}_{\mathrm{noise}}),\dots \right)$$

  $${\mathbf{C}}_{\mathrm{clean}}=\mathbf{V}{\mathrm{\Lambda }}_{\mathrm{clean}}{\mathbf{V}}^{tr}$$

  where ${\overline{\lambda }}_{\mathrm{noise}}$ is a representative value for the noise eigenvalues .
• Step 4:
  With the cleaned covariance matrix, we can define a stable metric. The Mahalanobis distance is a natural choice. For two curves represented as vectors $\overrightarrow{\mathcal{T}}$ and $\overrightarrow{{\mathcal{T}}^{\prime }}$, the RMT-metric is:

  $$d_{\mathrm{RMT}}(\overrightarrow{\mathcal{T}},{\mathcal{T}}^{\prime })=\sqrt{{(\mathcal{T}-{\mathcal{T}}^{\prime })}^{tr}{\mathbf{C}}_{\mathrm{clean}}^{-1}(\mathcal{T}-{\mathcal{T}}^{\prime })}$$
  This metric is more stable and reliable than one calculated from the noisy covariance matrix ${\mathbf{C}}^{-1}$, which is highly distorted by spurious correlations.
  In conclusion, RMT acts as a filter to denoise the covariance structure of the dataset. This cleaned structure increases the accurate metric for comparing any two individual curves within that dataset.

5. Data-Driven Analysis

5.1. A Simulated Data Case

This section is devoted to demonstrate the practical implementation of the proposed estimator through two illustrative examples. First, we assess the feasibility of both selection methods (LOOCV and Bootstrap selectors). Second, we illustrate the superior efficiency of the LARE regression over the LSRE approach. For this purpose, we generate two types of functional variables:

$$\left\{\begin{array}{cc}\mathrm{Smooth}\mathrm{functional}\mathrm{regessor}& T_i\left(t\right)=2b_{i}cos\left(\pi t\right)+3b_i{sin}^2\left(2\pi t\right)+4b_i{t}^2\mathrm{for}\mathrm{all}t\in (-1/2,1/2)\\ \mathrm{Rough}\mathrm{functional}\mathrm{regessor}& T_i\left(t\right)=2b_{i}cos\left(\pi t\right)+3b_i{sin}^2\left(\pi t\right)+4b_{i}t+{\eta }_i\left(t\right)\mathrm{for}\mathrm{all}t\in (-1/2,1/2)\end{array}\right.$$

where $b_i$ is generated from a $N(0,1)$ distribution and ${\eta }_i\left(t\right)$ is drown from $N(t,1+t)$ distribution. The curves are then discretized over a grid of points given by 100 equispaced measurements in $(0,1)$. A visual comparison of the two functional samples is provided in (cf. Figure 1).

We model the scalar response Y using the following non-parametric regression

$$S_i=3{\int }_{-0.5}^{0.5}{\displaystyle \frac{dt}{1+T_i^2\left(t\right)}}+{ϵ}_i\mathrm{for}i=1,\dots ,150,$$

(5)

where the errors $\left({ϵ}_i\right)$ are assumed to be independent of $\left(T_i\right)$. For the simulation study, we consider four error types of white noise generated from a normal mixture distribution, as follows:

$$\left\{\begin{array}{cc}\mathrm{Standard}\mathrm{normal}\mathrm{distribution} \left(\mathrm{S.N.D.}\right)& N(0,1)\\ \mathrm{Skewed}\mathrm{unimodal}\mathrm{distribution} \left(\mathrm{S.U.D.}\right)& {\displaystyle \frac{1}{5}}N(0,1)+{\displaystyle \frac{1}{5}}N({\displaystyle \frac{1}{2}},{\displaystyle \frac{4}{9}})+{\displaystyle \frac{3}{5}}N({\displaystyle \frac{13}{12}},{\displaystyle \frac{25}{81}})\\ \mathrm{Skewed}\mathrm{bimodal}\mathrm{distribution} \left(\mathrm{C.B.D.}\right)& {\displaystyle \frac{3}{4}}N(0,1)+{\displaystyle \frac{1}{4}}N({\displaystyle \frac{3}{2}},{\displaystyle \frac{1}{9}})\\ \mathrm{Claw}\mathrm{distribution} (\mathrm{C}.\mathrm{D}.)& {\displaystyle \frac{1}{2}}N(0,1)+\sum _{i=1}^4{\displaystyle \frac{1}{10}}N({\displaystyle \frac{i}{2}}-1,{\displaystyle \frac{1}{100}})\end{array}\right.$$

Further details on normal mixture distributions can be found in [26]. Now, in order to compare the LOOCV to the Bootstrap selector, we optimize both rules over $\mathcal{M}$, the set of m such that the ball centered at $\mathcal{T}$ with radius m contains exactly k neighbors of $\mathcal{T}$. Specifically, the number k is selected from the subset $\{5,15,25,\dots ,75\}$. Observe that the Bootstrap selector requires a pilot bandwidth $m_0$ which is chosen using the R- routine h.default from the R-package fda.usc. The two selection procedures performance are tested by comparing their averaged absolute errors ASE, where

$${\displaystyle \mathrm{ASE}=\frac{1}{150}\sum _{i=1}^{150}\left|S_i-\widehat{RER}\left(T_i\right)\right|.}$$

The model performance also depends critically on using a metric N that is appropriate for the data type. The latter is related to the smoothness of the curves. Here, the shape of the curves (cf. Figure 1) confirms that the metric defined by the derivatives is adequate for the smooth curves,

$$N(\mathcal{T},{\mathcal{T}}^{\prime })={\left({\int }_0^1{({\mathcal{T}}^{\left(i\right)}\left(t\right)-{\mathcal{T}}^{\prime \left(i\right)}\left(t\right))}^{2}dt\right)}^{1/2},$$

where ${\mathcal{T}}^{\left(i\right)}$ denotes the ith derivative of the curve $\mathcal{T}$. While for the rough curves we employ the PCA semi metric associated to the q-greatest eigenvalues.

While the simulation study investigated numerous values of i, q, and n, we report here only the results for $i=2$, $q=1$, and $n=150$ for brevity. A quadratic kernel on $(0,1)$ was used for this empirical analysis . Table 1 summarizes the obtained ASE-results for both selectors.

We observe that the LOOCV strategy yields superior results compared to the Bootstrap method. This advantage is pronounced in the case of Skewed Bimodal Distributions (S.B.D.), whereas for Contaminated Distributions (C.D.), the two methods perform nearly equivalently. Furthermore, both bandwidth selectors demonstrate relatively strong efficiency across the tested scenarios.

In the second part of this illustration, we examine the robustness of the proposed estimator by comparing it to the LSER-regression defined by

$$\widetilde{RER}={\displaystyle \frac{{\sum }_{i=1}^n{S}_i^{-1}M\left(m^{-1}N(\mathcal{T},T_i)\right)|}{{\sum }_{i=1}^{n}M\left(m^{-1}N(\mathcal{T},T_i)\right)S_i^{-2}}}.$$

It is important to note that robustness can be assessed through various metrics, including the influence function, gross-error sensitivity, breakdown point, local shift sensitivity or bias/variance under contamination, among others (see [27] for a comprehensive overview). While some of these metrics are challenging to apply in the context of functional statistics, they all share the common goal of evaluating the stability of a statistical method and quantifying its sensitivity to deviations from model assumptions. In this experiment, we assess robustness by analyzing the variance under the contaminated data-generating process. Specifically, we construct a contamination model by taking it as a white noise, where

$${\epsilon }_i=t{\epsilon }_{1_i}+(1-t){\epsilon }_{2_i},t\in (0,1),i=1,2,$$

where ${\epsilon }_{1_i}$ is generated from S.N.D and ${\epsilon }_{2_i}$ is drown from C.D. Under this contamination, the empirical variance is defined by

$$V\left(t\right)=Var{\left(\widehat{RER}\left(T_i\right)\right)}_{i=1}^n.$$

In Figure 2 and Figure 3, we plot the function $V\left(t\right)$ for both estimators in both cases (smooth and rough case).

Without surprising Figure 2 and Figure 3 confirm the robustness of $\widehat{RER}$ compared to $\widetilde{RER}$. It shows that the function $V\left(t\right)$ of $\widehat{RER}$ has a small variation compared to the variation in $\widetilde{RER}$. This conclusion is confirmed by computing the range of $V\left(t\right)$, defined as

$$RG=maxV\left(t\right)-minV\left(t\right),$$

which serves as a measure of the function’s variability. This benchmark is computed over 50 discretized points in the interval (0,1). Specifically, in the smooth case, we obtained $0.6$ for the LARE-regression compared to $1.9$ for the LSRE-regression. The same statement for the rough case was observed, and RG =0.3 for the LARE-regressions against $RG=1.45$ for the LSRE-regression case.

5.2. Real Data Application

This section presents a predictive analysis using the LARE-model. Specifically, we conduct a comparative study of several regression predictors. More precisely, we focus in this empirical analysis on prediction of the sulfone content in onion diets using Near-Infrared (NIR) spectral curves. Sulfones is organic compounds highly beneficial for human health. Their biomedical significance is well-documented in the literature, with studies indicating efficacy in treating conditions such as mucous-membrane inflammation and gastrointestinal disorders. Furthermore, research has highlighted their potential to suppress tumor initiation (see [28]). Our objective is to use the NIR spectrometry to quantify the relationship between onion consumption in a rodent model and subsequent sulfone levels in urine. The functional NIR data utilized in this analysis are available in the website http://www.models.life.ku.dk (accessed on 20 October 2025). The sample preparation involved a two-step chemical process: initial spectrometry processing, followed by spectral acquisition in the wavelength region of 9.6 to 0.3 ppm. The resulting dataset comprises 32 observations, which were recolored with sufficient replication. The spectral curves for these observations are plotted in Figure 4. Further details on this dataset are available in [28].

For the prediction task, we compare the performance of three regression estimators:

$$\begin{array}{cc}\hfill {\displaystyle \widehat{RER}\left(\mathcal{T}\right)=arg\underset{s}{min}\frac{{\sum }_{i=1}^{n}M\left(m^{-1}N(\mathcal{T},T_i)\right)\left|\frac{S_i-s}{S_i}\right|}{{\sum }_{i=1}^{n}M\left(m^{-1}N(\mathcal{T},T_i)\right)}}& (\mathrm{LARE}-\mathrm{Regression}),\hfill \\ \hfill \widetilde{RER}\left(\mathcal{T}\right)={\displaystyle \frac{{\sum }_{i=1}^n{S}_i^{-1}M\left(m^{-1}N(\mathcal{T},T_i)\right)|}{{\sum }_{i=1}^{n}M\left(m^{-1}N(\mathcal{T},T_i)\right)S_i^{-2}}}& (\mathrm{LSRE}-\mathrm{Regression}),\hfill \\ \hfill \overline{RE}\left(\mathcal{T}\right)={\displaystyle \frac{{\sum }_{i=1}^n{S}_{i}M\left(m^{-1}N(\mathcal{T},T_i)\right)|}{{\sum }_{i=1}^{n}M\left(m^{-1}N(\mathcal{T},T_i)\right)}}& (\mathrm{NW}-\mathrm{Regression}).\hfill \end{array}$$

To ensure a fair comparison, common computational strategies were employed for the three estimators. This includes the use of the RMT-metric, an identical kernel function M and the application of the rule specified in Equation (4) for selecting the bandwidth parameter for the three methods.

To evaluate the robustness of these estimators, we introduced contamination by multiplying a subset of observations by a factor of 10. The impact of this contamination was assessed by controlling the variation in the Mean Absolute Error (MAE). Specifically, for a given contamination proportion $t\in [0,1]$, we multiplied $\left[tn\right]$ observations (where $[.]$ denotes the integer part) and computed the variability of the prediction error function using

$$\mathrm{AME}\left(t\right)={\displaystyle \frac{1}{35}}\sum _{k=1}^{35}(S_k-\dot{Pred}\left(T_k\right)),$$

where $\dot{Pred}$ denotes one of the three compared estimators. The resulting prediction errors are displayed in Figure 5.

The results indicate that while all three predictors have satisfactory performance under uncontaminated conditions($t=0$), the LARE approaches demonstrates a significant advantage. This advantage is particularly evident in the stability of the LARE estimator ($\widehat{RER}$), which exhibits lower variability in its MAE under contamination. The superior stability of the LARE method is quantitatively confirmed by examining the range of its error function. The range of the LARE-regression was observed to be between $0.02$ and $0.15$, compared to a range of $[0.02,0.6]$ for the LSRE regression and $[0.05,1.4]$ for the standard nonparametric regression. This empirical finding aligns with the theoretical expectation that estimators based on absolute relative error (LARE) are more robust than those based on squared error (LSRE), due to the reduced influence of large residuals.

6. Conclusions and Prospects

This work introduces a robust framework for functional regression estimation by employing the least absolute relative error loss function. We construct a novel estimator and establish its asymptotic properties, specifically, proving its complete consistency. The functional structure is modeled using the small ball probability function that relates the topological structure of the probability measure of the functional explanatory variable. Furthermore, the employment of the LARE rule significantly increases the robustness of the model. The derived asymptotic results are obtained under general assumptions, ensuring their applicability to a broad class of functional data and continuous-time processes. Empirical analyses demonstrate the straightforward implementation of our estimator and confirm its superior robustness and predictive accuracy compared to alternative methods.

This contribution opens several promising questions for future research. A natural extension would be to develop a local linear version of this robust modal regression estimator. It is well-known that local linear modeling can improve the asymptotic convergence rate by reducing boundary bias. Establishing the asymptotic normality of the proposed robust estimator is another critical open question, as it is a prerequisite for constructing confidence intervals and conducting hypothesis tests. Other important directions include exploring linear or partially linear formulations of this regression. A particularly challenging task would be to extend these results to more complex data structures, such as censored data or functional time series. Clearly, the treatment of these complex structures requires new mathematical tools and developments that extend beyond the scope of the present work.