Linear Wavelet Estimation for Regression with Additive-Multiplicative Noise for Ergodic Data

Basics on wavelets and BESOV Balls

This section provides a concise overview of wavelets and Besov spaces, which are central to the construction of our estimator. We employ compactly supported wavelets from the Daubechies family. For a comprehensive treatment, we refer to standard texts such as, Daubechies^[14] Mallat^[15]. For any resolution level j ≥ 0, the index set is ${\Lambda }_{1\le i\le n}=\left(0,\mathrm{...},2^{j-1}\right).$ For any $k\in {\Lambda }_j$, in this set, the scaling and wavelet basis functions are generated via translation and dilation of the father and mother wavelets, respectively:

${\varphi }_{j,k}\left(x\right)=2^{j/2}\varphi \left(2^{J}x-k\right),{\psi }_{j,k}\left(x\right)=2^{j/2}\psi \left(2^{J}x-k\right)$

Following the methodology, Cohen et al.^[16] there exists an integer $\tau$ such that, for any integer the collection of functions

$S=\left({\varphi }_{j_0,k},k\in {\Lambda }_{j_0};{\psi }_{j,k};j\in \mathbb{N}-\left(0,\cdots ,j_0-1\right)\right),k\in {\Lambda }_{j_0}$

forms an orthonormal basis of $L^2\left(\right[0,1\left]\right)$. Therefore, for any integer $j_0,\ge \tau$ and any function $f\in L^2\left(\right[0,1\left]\right)$, admits the wavelet expansion:

$f\left(x\right)=\sum _{k\in {\Lambda }_{j_0}}{\alpha }_{j_0,k}{\varphi }_{j_0,k}\left(x\right)+\sum _{j=j_0}^{\infty }\sum _{k\in {\Lambda }_{j_0}}{\beta }_{j_0,k}{\psi }_{j_0,k}\left(x\right),x\in \left(0,1\right)$

where the coefficients are given by the inner products:

${\alpha }_{j_0,k}={\int }_0^{1}f\left(x\right){\varphi }_{j_0,k}\left(x\right)dx,{\beta }_{j_0,k}={\int }_0^{1}f\left(x\right){\psi }_{j_0,k}\left(x\right)dx$

A property that is instrumental in our proofs is the integral of the scaling function over its domain:

${\int }_0^1{\varphi }_{j,k}\left(x\right)dx=2^{-j/2}.$

Let $P_i$ denote the orthogonal projection operator from $L^2\left(\right[0,1\left]\right)$ onto the space with the orthonormal basis $\left({\varphi }_{j,k}\left(\cdot \right)=2^{-j/2}\varphi \left(2^j\cdot -k\right),k\in {\Lambda }_j\right)$. Then, for any $f\in L^2\left(\right[0,1\left]\right)$, we have

$P_{i}f\left(x\right)=\sum _{k\in {\Lambda }_{j_0}}{\alpha }_{j_0,k}{\varphi }_{j_0,k}\left(x\right),x\in \left(0,1\right)$

Besov spaces $B_{p,q}^s\left(\right[0,1\left]\right)$ are a flexible class of function spaces that can characterize a wide range of smoothness properties, including spatially inhomogeneous behavior. For further details, see [Härdle et al.,^[17] Meyer,^[18] Triebel.^[19] Definitions of those spaces are given below. Suppose that $\varphi$ is m regular (i.e. $\varphi \in C^m$ and $\left(D^{\alpha }\varphi \left(x\right)\right)\le c{\left(1+{\left(x\right)}^2\right)}^{-l}$ for each $l\in \mathbb{Z}$, with $\alpha =0,1,\mathrm{...},m$). Let $f\in L^2\left(\right[0,1\left]\right)$, $p,q\in [1,\infty ]$ and $0<s<m$. Then the following assertions are equivalent:

1. $f\in B_{p,q}^s\left(\left(0,1\right)\right)$

2. $\left(2^{js}{P}_{i+1}f-P_i{f}_p\right)\in l_q;$

3. $\left(2^{j\left(s-\frac{1}{p}+\frac{1}{2}\right)}{\beta }_{j,.}{}_p\right)\in l_q.$

The Besov norm of h can be defined by

$f_{B_{p,q}^s}:={\alpha }_{j,\cdot }{}_p+{\left(2^{j\left(s-\frac{1}{p}+\frac{1}{2}\right)}{\beta }_{j,.}{}_p\right)}_{j\ge \tau }{}_q,\text{where}{\beta }_{j,.}{}_p^p=\sum _{k\in {\Lambda }_{j_0}}{\left({\beta }_{j,k}\right)}^p$

Assumptions and Estimators

To facilitate the presentation of our main results, we introduce the following notation. For each integer i, let $F_i$ be the $\sigma -$algebra generated by the collection $\left(X_j:0\le j\le i\right)$. When i ranges from 1 to n, we define the conditional density of $X_i$given $F_{i-1}$ as

$f_{X_i}^{F_{i-1}}\left(\cdot \right)=f^{F_{i-1}}\left(\cdot \right)$

Technical assumptions on the model (1.1) are formulated below.

A.1 We suppose that$f:\left(0,1\right)\to ℝ$ is bounded from above.

A.2 We suppose that $U_1\sim U\left(\left(-\theta ,\theta \right)\right)$ with $\theta >0$ a fixed real number.

A.3 We suppose that $V_1$ has a moment of order 4.

A.4 We suppose that $X_i$ and $V_i$ are independent for any $i\in \left(1,\mathrm{...},n\right)$ and $V_i$ is stationary and ergodic.

A.5 For every x ∈ R, the sequence

$\frac{1}{n}\sum _{i=1}^n{f}^{F_{i-1}}\left(X\right)=f\left(x\right)asn\to \infty$

both almost surely (a.s.) and in the L ² sense.

A.6 The conditional mean of $Y_i^2$ given the σ −field $F_i$ depends only on $X_i$, i.e., for any $i\ge 1$,

$E\left(Y_i^2|F_i\right)=E\left(Y_i^2|X_i\right)=\frac{{\theta }^2}{3}r\left(X_i\right)+E\left(V_1^2\right)$

Comments on the assumptions

A.1 (Boundedness of $f$): This is a standard technical condition that ensures the function $r=f^2$ is also bounded, which is necessary to control the moments of the observations $Y_i$ and guarantee the finiteness of the wavelet coefficients.

A.2 (Uniform Multiplicative Noise): Specifying the uniform distribution for $U_i$ is crucial as it allows us to compute the constant factor ${\theta }^2/3$ that appears in the key identity of A.6. The symmetry about zero is central to canceling the cross-term in $E\left(Y_i^2|X_i\right)$.

A.3 (Additive Noise Moments): Requiring a fourth-order moment for $V_1$ is necessary to control the variance of our coefficient estimator ${\hat{\alpha }}_{j,k}$, which involves the squared observations $Y_i^2$.

A.4 (Independence): The independence between $X_i$ and $V_i$ is crucial to ensure that the additive noise does not confound the relationship between the covariate and the function r. This allows the conditional expectation $E\left(Y_i^2|X_i\right)$ to cleanly separate into a component depending on $r\left(X_i\right)$ and the constant variance of $V_i$.

A.5 (Ergodicity of the Covariate Process): This is a significant generalization of the typical i.i.d. covariate assumption. The ergodic theorem ensures that time averages converge to ensemble expectations, which is the cornerstone for establishing the consistency of our sample-based estimator ${\hat{\alpha }}_{j,k}$[GETIMGFILENAME=D:\ZAZA\SS\JQUS\XML\JQUS_12_2025_1006317\JQUS_12_2025_R21_InDesign\Links\JQUS12_78.wmf]. This relaxation widens the applicability of our model to many time series and spatially dependent data scenarios where mixing conditions may not hold or are difficult to verify Didi and Bouzebda.^[5,6]

A.6 (Conditional Mean Structure): This assumption is fundamental to our estimator construction. It formalizes the fact that the relevant information for predicting $Y_i^2$ is fully captured by the current covariate $X_i$, rather than the entire filtration history. This conditional independence structure is what allows us to derive the pivotal relationship $E\left(Y_i^2|X_i\right)=\frac{{\theta }^2}{3}r\left(X_i\right)+E\left(V_1^2\right)$, which directly isolates $r\left(X_i\right)$ and forms the basis for our estimator in (3.1) and (3.2).

The extension from an i.i.d. to a strictly stationary and ergodic framework for the covariate process is theoretically justified by the ergodic theorem, which ensures that time averages converge to ensemble expectations. This forms the cornerstone for establishing the consistency of our sample-based estimator ${\hat{\alpha }}_{j,k}$. Methodologically, this approach is supported by modern asymptotic theory for dependent data, particularly the convergence results for empirical processes based on ergodic data as developed in Didi and Bouzebda, and Didi et al.^[5,6] This relaxation significantly widens the applicability of our model to include many time series and spatially dependent data scenarios.

The reader may also be referred to Didi S, Bouzebda S.^[6,21,22]

Our linear wavelet estimator for the function r is defined by the projection:

The empirical wavelet coefficients ${\hat{\alpha }}_{j,k}$ are calculated as:

The form of the estimator ${\hat{\alpha }}_{j,k}$ is an original contribution developed by Chesneau et al.^[23] to address the challenges of multiplicative noise in nonparametric regression. While straightforward to compute, the efficacy of this estimator is largely governed by the tuning parameter j₀. This linear wavelet approach is well-established in nonparametric regression [Härdle et al., 1998],^[17] with contemporary extensions found in Chaubey et al.^[24]

The form of the estimator is derived from the mathematical relationship $E\left(Y_i^2|X_i\right)=\frac{{\theta }^2}{3}r\left(X_i\right)+E\left(V_1^2\right)$, which isolates $r\left(X_i\right)$ under the assumptions of the model (e.g., independence between $X_i,U_i,V_i$ (for the same $(i$), uniformity of $U_i$, and known noise moments). The term $2^{-j\text{/}2}$ arises from the integral of the scaling function ${\varphi }_{j,k}\left(x\right)$ over $\left(0,1\right)$.

Simulation study

In order to illustrate the empirical performance of the proposed estimator under ergodic covariates, a numerical study was conducted. To reflect a realistic application context, we implemented an automatic data-driven method for selecting the truncation parameter j₀, which does not require prior knowledge of the regularity of the function to be estimated. Simulations were carried out using R, with wavelet computations supported by the ‘wavethresh’ package Nesson 2023.^[26] This procedure involves randomly partitioning the dataset into two equally sized subsets, constructing an estimator on one subset, and validating its performance on the other. The resolution level j0​ that minimizes the average validation error across both subsets is selected. This well-established routine provides a fully data-driven and computationally efficient method for choosing the tuning parameter without prior knowledge of the function’s smoothness.

In our numerical study, we generated data from model (1.1) with a sample size of n=4096. The covariates $X_i$ were simulated from a stationary AR(1) process, $X_t=\phi X_{t-1}+{\epsilon }_t,$ with $\phi =0.5$ and ${\epsilon }_t\sim N\left(0,0.1\right)$, normalized to the interval [0,1]. This process is known to be ergodic for $\left(\phi \right)<1$, thereby satisfying our core Assumption (A.5).

In accordance with Assumption A.2, the multiplicative noise $U_i$ was drawn from a uniform distribution $U\left(\right[-\theta ,\theta \left]\right).$ For this numerical study, we used $\theta =1$. The additive noise $V_i$ was drawn from a normal distribution $N\left(0,0.01\right)$.

For the wavelet construction, we selected Daubechies’ compactly-supported wavelets with 8 vanishing moments, implemented via the wavethresh package in R. We utilized the Daubechies extremal phase wavelets (filter.number=8) and employed periodic boundary handling to manage the interval [0,1].

We consider the HeaviSine and Bumps test functions for f , and we wish to estimate$r=f^2$.

To implement the estimator with a random design, we sorted the ($X_i,Y_i^2$) pairs by $X_i$ and then applied Mallat’s algorithm to the sequentially ordered $Y^2$ values. The truncation parameter j₀ was selected using a 2-Fold Cross Validation (2FCV) method. Figure 2 displays the observations, sample estimators, and the selection criterion curves. The method successfully identifies a near-optimal truncation level, and the resulting estimator accurately captures the shape of r(x). Similar results were obtained with the Bumps function, confirming the robustness of the procedure.

Close

Figure 1:

The two test (squared) functions to be estimated. (a) HeavySine function, (b) Bumps Function.

Export to PPT

Close

Figure 2:

(a) Mean Integrated Squared Error (log scale) comparing proposed and standard (i.i.d case) methods across increasing sample sizes (n = 512, 2048, 4096) (b) Visual comparison of true HeaviSine function against estimates from both methods in a single simulation trial (c) Average truncation level j₀ with standard deviation bars, showing stable parameter choice across sample sizes.

Export to PPT

The empirical performance of our proposed estimator is further summarized in Figures 3 and 4, which provide a comparative analysis against a standard estimator designed for i.i.d. covariates. For both the HeaviSine [Figure 3] and Bumps [Figure 4] functions, Sub-Figure(a) displays the Mean Integrated Squared Error (MISE) on a log scale across different sample sizes, demonstrating the competitive convergence of our method. Sub-Figure(b) offers a visual comparison from a single simulation trial, illustrating that our estimator accurately captures the underlying function shape. Furthermore, Sub-Figure(c) shows the average data-driven truncation level $j_0$​ selected by the 2FCV method, confirming its stable and adaptive selection across simulations.

Close

Figure 3:

Bumps Function (a) Mean Integrated Squared Error comparison between proposed and standard methods for Bumps function across sample sizes n = 512, 2048, 4096 (b) Visual comparison of true Bumps function against estimates from both methods, demonstrating superior peak preservation by the proposed approach (c) Average truncation level j₀ with standard deviation bars, showing adaptive resolution selection for the spatially complex Bumps function.

Export to PPT

Close

Figure 4:

( a ) Noisy observations (X,Y 2). (b): Sample of the model collection. (c) Graph of the MSE (blue) against j₀ and (re-scaled) 2FCV criterion. (d): Typical estimations from one simulation with n = 4096. Blue lines indicate the true functions; red lines correspond to the estimators r̂_0,n. FCV: Fold cross validation

Export to PPT

Proofs

To prove Theorem 3.1, we use the following lemmas. We establish a new upper bound for the partial sums of unbounded martingale differences. This result is crucial for analyzing the asymptotic properties of a density estimator derived from strictly stationary and ergodic samples. Throughout, C denotes a positive constant that may change from one instance to the next. The corresponding inequality is presented in the following lemmas.

The following well-known inequality is crucial for controlling the growth of martingale sequences in $L^p$ spaces. We use it to bound the variance of our empirical wavelet coefficients.

Lemma 6.1. (Burkholder-Rosenthal inequality) Following Notation 1 in Burkholder.^[27]

Let ${\left(X_i\right)}_{i\ge 1}$ be a stationary martingale adapted to the filtration ${\left(F_i\right)}_{i\ge 1}$, define ${\left(d_i\right)}_{i\ge 1}$ is the sequence of martingale differences adapted to ${\left(F_i\right)}_{i\ge 1}$ and $S_n={\sum }_{i=1}^n{d}_i$, then for any positive integer n

Where as usual the norm ${\left(\cdot \right)}_p={\left(E\left({\left(\cdot \right)}^p\right)\right)}^{1/p}$.

This exponential tail bound for martingale differences is essential for establishing almost-sure convergence rates. It will be used to handle the deviation of the coefficient estimators in the uniform convergence analysis.

Lemma 6.2. (Exponential Inequality for Martingale Differences) Let ${\left(Z_n\right)}_{n\ge 1}$ be a sequence of real-valued martingale differences with respect to the filtration ${\left(F_n=\sigma \left(Z_1,\cdots ,Z_n\right)\right)}_{n\ge 1}$. Define the partial sum:

$S_n={\sum }_{i=1}^n{Z}_i$

Suppose there exist constants $\text{C}>0$ and a sequence ${\left({\text{d}}_{\text{n}}\right)}_{\text{n}\ge 1}$ of nonnegative real numbers such that for all $\text{p}\ge 2$, $\text{n}\ge 1$, the following moment condition holds almost surely:

$E\left(Z_i^p|F_{n-1}\right)\le C^{p-1}p!d_n^2,\text{almost}\text{surely}$

Then, for any $\epsilon >0$, we have the tail bound:

$\mathbb{P}\left(\left(S_n\right)>\epsilon \right)\le 2\exp \left(-\frac{{\epsilon }^2}{2\left(D_n+C\epsilon \right)}\right)$

Where $D_n={\sum }_{i=1}^n{d}_i^2$.

Proof of Lemma 6.2. The proof follows as a particular case of Theorem 8.2.2.^[28]

This lemma establishes that our empirical wavelet coefficient is an unbiased estimator for the true coefficient, which is a direct consequence of the model assumptions and the linearity of the estimator.

Lemma 6.3. Chesneau et al.^[22]

Let $j\ge \tau ,k\in {\Lambda }_j,{\hat{\alpha }}_{j,k}$ be (3.2). Then, under A.1-A.4, we have

$E\left( {\hat{\alpha }}_{j,k}\right)={\alpha }_{j,k}.$

Define

This result strengthens Lemma 6.3 by showing that the conditional expectation of our estimator, given the past, also converges to the true coefficient. This demonstrates the robustness of the estimator under the ergodic dependent structure.

Lemma 6.4. Let $j\ge \tau ,k\in {\Lambda }_j,{\tilde{\alpha }}_{j,k,n}$ be (6.2). Then, under A.1-A.6, we have

${\tilde{\alpha }}_{j,k,n}={\alpha }_{j,k},asn\to \infty$

Proof of Lemma 6.4. Consider the following decomposition

Owing to Lemma 6.3 we have

considering the structure of the $\sigma -$filtration $F_l$ we have

$F_{l-1}\subseteq F_l$

therefore,

$\begin{array}{c}{E}_{j,k,1}={\tilde{\alpha }}_{j,k,n}-E\left({\hat{\alpha }}_{j,k}\right)\\ =\left(\frac{3}{{\theta }^2}\right)\frac{1}{n}\sum _{l=1}^n\left(E\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)|F_{l-1}\right)-E\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)\right)\right)\\ =\left(\frac{3}{{\theta }^2}\right)\frac{1}{n}\sum _{l=1}^n\left(\begin{array}{c}E\left(E\left(Y_l^2|F_l\right){\varphi }_{j,k}\left(X_l\right)|F_{l-1}\right)\\ -E\left(E\left(Y_l^2|F_l\right){\varphi }_{j,k}\left(X_l\right)\right)\end{array}\right)\end{array}$

observe that under assumption (A.6), we have

And

gathering statement (6.5) and (6.6), we have

$\begin{array}{c}{E}_{j,k,1,1}=\frac{1}{n}\sum _{l=1}^n\left(E\left(r\left(X_i\right){\varphi }_{j,k}\left(X_l\right)|F_{l-1}\right)-E\left(r\left(X_i\right){\varphi }_{j,k}\left(X_l\right)\right)\right)\\ +E\left(V_1^2\right)\left(\frac{3}{{\theta }^2}\right)\frac{1}{n}\sum _{l=1}^n\left(E\left({\varphi }_{j,k}\left(X_l\right)|F_{l-1}\right)-E\left({\varphi }_{j,k}\left(X_l\right)\right)\right)\\ =E_{j,k,1,1}+E_{j,k,1,2}\end{array}$

Where, under assumption (A.5), we have

and

using statements (6.7) and (6.8), we deduce

By combining statements (6.3), (6.4) and (6.9), we conclude the proof of Lemma 1.

This lemma provides a bound on the mean-squared error of the empirical wavelet coefficients. It is the key step in the bias-variance decomposition that leads to the final MISE convergence rate in Theorem 3.1. The proof relies on our martingale framework to handle the dependence.

Lemma 6.5. Let $j\ge \tau$ such that $2^j\le n$, $k\in {\text{Λ}}_j, {\hat{\alpha }}_{j,k}$ be (3.2). Then, under A.1-A.4, we have

$E\left({\left({\hat{\alpha }}_{j,k}-{\alpha }_{j,k}\right)}^2\right)\lesssim \frac{1}{n}.$

Proof of Lemma 6.5. Consider first the decomposition

From decomposition (6.10), it follows that

$\begin{array}{c}E\left({\left({\hat{\alpha }}_{j,k}-{\alpha }_{j,k}\right)}^2\right)=E\left({\left({\hat{\alpha }}_{j,k}-{\tilde{\alpha }}_{j,k,n}\right)}^2\right)\\ +2E\left(\left({\hat{\alpha }}_{j,k}-{\alpha }_{j,k}\right)\left({\tilde{\alpha }}_{j,k,n}{}_{j,k}-{\alpha }_{j,k}\right)\right)\\ +E\left({\left({\tilde{\alpha }}_{j,k,n}{}_{j,k}-{\alpha }_{j,k}\right)}^2\right)\\ =A_{j,k,1}+A_{j,k,2}+A_{j,k,3},\end{array}$

Owing to Lemma 6.4 we have

On the other hand, from statements (3.2) and (6.2), we observe that

$\begin{array}{c}{\hat{\alpha }}_{j,k}-{\tilde{\alpha }}_{j,k,n}=\left(\frac{3}{{\theta }^2}\right)\left(\frac{1}{n}\right)\sum _{l=1}^n\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)-E\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)|F_{l-1}\right)\right)\\ = \left(\frac{3}{{\theta }^2}\right)\frac{1}{n}\sum _{l=1}^n{Z}_{j,k,l},\end{array}$

Therefore,

Lastly, notice that ${\left(Z_{j,k,l}\right)}_{0\le l\le n}$ is a sequence of martingale differences with respect to the sequence of $\sigma$-fields ${\left(F_l\right)}_{0\le l\le n}$, Using the inequality provided by Lemma 6.1 we immediately deduce

$A_{j,k,1}=\left(\frac{9}{{\theta }^4}\right)\left(\frac{1}{n^2}\right)E\left({\left(\sum _{l=1}^n{Z}_{j,k,l}\right)}^2\right)$

Furthermore, one can establish:

To analyze these terms, we use a standard decomposition and note that $F_0$ is the trivial $\sigma$-field. Hence, one obtains

following a similar argument as in (6.6), and using the statement $r=f^2$ we obtain

Using the independence assumptions on the random variables, A.1-A.5 with $E\left[U_1\right]=0$, observe that

$\begin{array}{c}E\left(U_1{V}_1^{3}f\left(X_1\right){\varphi }_{j,k}{}^2\left(X_1\right)\right)\\ =E\left(U_1\right)E\left(V_1^3\right)E\left(f\left(X_1\right){\varphi }_{j,k}{}^2\left(X_1\right)\right)=0.\end{array}$

Moreover, $E\left(V_1\right)=0$, observe that

$\begin{array}{c}E\left(U_1^3{V}_1{f}^3\left(X_1\right){\varphi }_{j,k}{}^2\left(X_1\right)\right)\\ =E\left(U_1^3\right)E\left(V_1\right)E\left(f^3\left(X_1\right){\varphi }_{j,k}{}^2\left(X_1\right)\right)=0.\end{array}$

and by assumption (A.3) we have

$\begin{array}{c}E\left(V_1^4{\varphi }_{j,k}{}^2\left(X_1\right)\right)=E\left(V_1^4\right)E\left({\varphi }_{j,k}{}^2\left(X_1\right)\right)\\ =E\left(V_1^4\right)\underset{0}{\overset{1}{\int }}{\varphi }_{j,k}{}^2\left(x\right)dx=E\left(V_1^4\right)\lesssim 1.\end{array}$

Therefore, using similar mathematical arguments with $E\left(U_1^2\right)=\frac{{\theta }^2}{3}$, and assumption (A.1), we have

$\begin{array}{c}E\left(U_1^2{V}_1^2{f}^2\left(X_1\right){\varphi }_{j,k}{}^2\left(X_1\right)\right)\\ =\left(\frac{{\theta }^2}{3}\right)E\left(V_1^2\right)E\left(f^2\left(X_1\right){\varphi }_{j,k}{}^2\left(X_1\right)\right)\le C_f^2\left(\frac{{\theta }^2}{3}\right)E\left(V_1^2\right).\end{array}$

Lastly, by the independence assumption between $U_1$ and $X_1$ and (A.1) and (A.2), we have

$\begin{array}{c}E\left(U_1^4{f}^4\left(X_1\right){\varphi }_{j,k}{}^2\left(X_1\right)\right)\\ =E\left(U_1^4\right)E\left(f^4\left(X_1\right){\varphi }_{j,k}{}^2\left(X_1\right)\right)\lesssim E\left(U_1^4\right)=\frac{{\theta }^4}{5}.\end{array}$

Thus, all the terms of (6.15) are bounded from above

Next, we analyze the second component of the decomposition in equation (6.13). Specifically, we consider

$\begin{array}{c}{\text{Φ}}_{\left(2\right)}={\left(E\left(\sum _{l=1}^{n}E\left(Z_{j,k,l}^2|F_{l-1}\right)\right)\right)}^{1/2}\\ ={\left(\sum _{l=1}^{n}E\left(E\left(Z_{j,k,l}^2|F_{l-1}\right)\right)\right)}^{1/2}=\left(\sum _{l=1}^{n}E\left(Z_{j,k,l}^2\right)\right).\end{array}$

Applying a standard identity, we find

$\begin{array}{c}E\left(Z_{j,k,l}^2\right)=E\left({\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)-E\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)|F_{l-1}\right)\right)}^2\right)\\ \le 2E\left({\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)\right)}^2\right)+2E\left({\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)\right)}^2|F_{l-1}\right)\\ \le 4E\left({\left(Y_l^2{\varphi }_{j,k}\left(X_l\right)\right)}^2\right).\end{array}$

Using equations (6.14) and (6.15), it follows that

Combining the results from equations (6.16) and (6.17), we derive

${\left(E\left({\left(\sum _{l=1}^n{Z}_{j,k,l}\right)}^2\right)\right)}^{1/2}=O\left(n^{1/2}\right).$

Consequently,

$\begin{array}{c}{A}_{j,k,1}=E\left({\left({\hat{\alpha }}_{j,k}-{\tilde{\alpha }}_{j,k,n}\right)}^2\right)\\ =\left(\frac{9}{{\theta }^4}\right)\left(\frac{1}{n^2}\right)E\left({\left(\sum _{l=1}^n{Z}_{j,k,l}\right)}^2\right)=O\left(n^{-1}\right).\end{array}$

This ends the proof of Lemma 2.

Proof of Theorem 3.1. The proof follows the standard bias-variance decomposition paradigm for projection estimators. The novel challenge, addressed by our martingale framework, lies in establishing the variance bound under the general ergodic setting.

We proceed via a conventional bias-variance trade-off analysis [Härdle et al., 2012].^[18] The pivotal step is the application of Lemma 6.5, which provides the necessary variance bound for the coefficient estimators; this lemma’s proof relies critically on our martingale framework. The optimal convergence rate is then achieved by selecting j_ to equilibrate this variance with the projection bias. To this end, we begin with the projector-based decomposition:

From the orthonormality of the wavelet basis, we obtain:

$\begin{array}{c}E\left({\hat{r}}_{j_0,n}-P_{j_0}{r}_2^2\right)=E\left(\sum _{k\in {\Lambda }_{j_0}}\left({\hat{\alpha }}_{j_0,k}-{\alpha }_{j_0,k}\right){\varphi }_{j_0,k}{}_2^2\right)\\ =\sum _{k\in {\Lambda }_{j_0}}E\left({\left({\hat{\alpha }}_{j_0,k}-{\alpha }_{j_0,k}\right)}^2\right)\end{array}$

It follows from Lemma 6.5 that, ${\Lambda }_{j_0}\sim 2^{j_0}$ and $2^{j_0}\sim n^{\frac{1}{2s^{\prime }+1}}$,

Case 1: $p\ge 2,s^{\prime }=s$. By Hölder inequality and $r\in B_{p,q}^s\left(\left(0,1\right)\right)$ we have

$P_{j_0}r-r_2^2\lesssim P_{j_0}r-r_p^2\lesssim 2^{-j_{0}s}\lesssim n^{-\frac{2s^{\prime }}{2s^{\prime }+1}};$

Case 2: $1\le p\left(2,s\right)1/p$. Using the embedding  $B_{p,q}^s\left(\left(0,1\right)\right)\subseteq B_{2,\infty }^{s^{\prime }}\left(\left(0,1\right)\right)$ we have

In both cases,

$P_{j_0}r-r_2^2\lesssim n^{-\frac{2s^{\prime }}{2s^{\prime }+1}};$

By statements (6.18), (6.19) and (6.20), we obtain

$E\left({\int }_0^1{\left({\hat{r}}_{j_0,n}\left(x\right)-r\left(x\right)\right)}^{2}dx\right)\lesssim n^{-\frac{2s^{\prime }}{2s^{\prime }+1}};$

This proves Theorem 3.1.