1 Introduction.
Consider the problem of estimating a periodic twodimensional function, , based on noisy convolutions that are described by the equations
(1) 
Here, , is the convolution kernel and it is assumed to be known, is an anisotropic twodimensional fractional Brownian sheet (fBs), , and , , are the parameters of the longmemory in the direction of and , respectively. A twodimensional fractional Brownian sheet is defined by the formula
(2) 
where is a twodimensional standard Brownian sheet, is some explicit constant and
(3) 
with . In addition, the twoparameter fBs is characterized by a covariance function of the form
(4) 
for some , , , where are the Hurst parameters in the direction of and , respectively (see Kamont (1996)).
The discrete, and the more realistic, version of model (1) is given by
(5) 
where , ,
is a positive variance constant, and
are zeromean secondorder stationary Gaussian random variables satisfying the autocovariance structure
(6) 
Assumption A.1. The error structure in model (5) is a zeromean secondorder stationary process satisfying
(7) 
where
(8) 
and is a twoparameter fractional Brownian sheet.
Assumption A.1 is valid under autocovariance structure (6) and it appeared in Adu and Richardson (2018). It guarantees that, with the calibration , models (1) and (5) are asymptotically equivalent. Therefore, for sufficiently large sample size , model (1) can be used to approximate model (5).
Deconvolution model has been the subject of a great deal of papers since late 1980s, but the most significant contribution was that of Donoho (1995) who was the first to devise a wavelet solution to the problem. Other attempts include, Abramovich and Silverman (1998), Walter and Shen (1999), Johnstone et al. (2004), Donoho and Raimondo (2004), among others. In the case of functional deconvolution model with , Pensky and sapatinas (2009, 2010, 2011) pioneered into the formulation and further development of the problem.
Functional deconvolution problem of type (1) with
corresponds to the white noise case and it was investigated under the
risk in Benhaddou et al. (2013), and under risk, , in Benhaddou (2017), where they constructed an adaptive hardthresholding wavelet estimator, and showed that it is asymptotically quasioptimal within a logarithmic factor of over a wide range of Besov balls of mixed regularity. This model is motivated by experiments in which one needs to recover a twodimensional function using observations of its convolutions along profiles . This situation occurs, for example, in seismic inversions (see Robinson (1999)). In these articles, it is assumed that the error terms are white noise processes or i.i.d noise. However, empirical evidence has shown that, even at large lags, the correlation structure in the errors takes the powerlike form (6). This phenomenon is referred to as longmemory (LM) or longrange dependence (LRD). The presence of long memory in oceanic seismic data was pointed out in Wood et al. (2014) for instance, and it may be due for example to sea floor temperature variations.Longmemory has been investigated quite considerably in the standard deconvolution model. One can list; Wang (1996, 1997), Wishart (2013), Benhaddou et al. (2014) and Kulik et al. (2015). In a few other relevant contexts, LM was also investigated in density deconvolution in Comte et al. (2008), and in the Laplace deconvolution in Benhaddou (2018) where the unknown response function is nonperiodic and defined on the entire positive real halfline.
The objective of the paper is to extend the work of Benhaddou et al. (2013) to the case when the noise is an anisotropic multiparameter fractional Brownian sheet. Following a standard procedure, we derive minimax lower bounds for the risk, with , when belongs to an anisotropic Besov ball and the blurring function is regular smooth. Taking advantage of the Riesz polypotential operator, and inspired by the work of Wang (1996, 1997), we apply the waveletvaguelette expansion (WVD) to decorrelate the anisotropic fBs in two directions. In addition, we take advantage of the flexibility of the Meyer wavelet basis in the Fourier domain to construct a wavelet hardthresholding estimator that is adaptive and asymptotically quasioptimal within a logarithmic factor of over a large array of Besov balls. Furthermore, the proposed estimator attains convergence rates that depend on a delicate balance between the parameters of the Besov ball, the degree of illposedness of the convolution, and the parameters and associated with the fBs. Similar to the white noise case studied in Benhaddou et al. (2013), the proposed approach is easily extended to recovering an dimensional function, , and the corresponding convergence rates turn out to be dimensionfree. Finally, a simulation study is carried out to further confirm the results of our theoretical findings.
The rest of the paper is organized as follows. Section 2 introduces some notation as well as the estimation algorithm. Section 3 describes the derivation of the lower bounds for the risk of estimators of as well as the upper bounds and establishes the asymptotic optimality of the estimator. Section 4 presents a limited simulations study to complement the theoretical results from Section 3. Section 5 extends the results in Sections 2 and 3 to the general dimensional case. Finally, Section 6 contains the proofs of the theoretical results.
2 Estimation Algorithm.
In what follows, denote the complex conjugate of by .
Let , , and
be the twodimensional Fourier coefficients of functions , , , and , respectively.
Consider a bandlimited wavelet basis (e.g., Meyertype) on , and form the product wavelet basis
(9) 
where , and
(10) 
Let be the lowest resolution level for the Meyer basis and denote the scaling function for the wavelet by , . Since the functions (9) form an orthonormal basis of the space, the function can be expanded over these bases with coefficients into a wavelet series
(11) 
where
. Denote the twodimensional Fourier transform of
by . It is wellknown (see, e.g, Johnstone et al (2004), section 3.1) that under the Fourier domain and for any , , one has(12) 
Apply the twodimensional Fourier transform to equation (1) to obtain
(13) 
Now, applying Plancherel formula in both directions, the wavelet coefficients in (11) can be rewritten as
(14) 
Therefore, an unbiased estimator for
is given by(15) 
Finally, consider the hard thresholding estimator for
(16) 
where
(17) 
and the quantities , and are to be specified.
Assumption A.2. In the Fourier domain, the convolution kernel , for some positive constants , and and , independent of and , is such that
(18) 
Let us now evaluate the variance of (15).
Lemma 1
According to Lemma 1, choose the thresholds of the form
(21) 
where is some positive constant independent of . Furthermore, the finest resolution levels and are such that
(22) 
3 Convergence rates and asymptotic optimality.
Denote
(23)  
(24)  
(25) 
where and .
Assumption A.3. The function belongs to an anisotropic twodimensional Besov space. In particular, its wavelet coefficients satisfy
(26) 
To construct minimax lower bounds for the risk, we define the risk over the set as
(27) 
where is the norm of a function and the infimum is taken over all possible estimators of .
Proof of Theorem 1. In order to prove the theorem, first we establish the equivalence between a process that involves the fractional Wiener sheet (eq. (1)) and another that involves the standard Wiener sheet by applying an appropriate integral operator, taking advantage of relation (2). Then, we consider two cases, the case
when is dense in both dimensions (the densedense case) and the case when is dense in but sparse in . Lemma of Bunea et al. (2007) is then applied to find such lower bounds using conditions (18) and (26), combined with some useful properties of Meyer wavelet basis. To complete the proof, we choose the highest of the lower bounds.
The derivation of upper bounds of the risk relies on the following two lemmas.
Lemma 2
Under condition (26), one has
(29) 
Lemma 3
Theorem 2
The proof of Theorem 2. The proof is very similar to that of Theorem 2 in Benhaddou. (2017).
Remark 1
Remark 2
Remark 3
The convergence rates depend on a delicate balance between the parameters of the Besov ball, smoothness of the convolution kernel and the parameters of the anisotropic fractional Brownian sheet, and . These rates deteriorate as the values of and get closer and closer to zero.
Remark 4
To decorrelate the twodimensional fractional Gaussian noise a waveletvaguelette expansion based on Meyer wavelet bases has been applied to the fractional Brownian sheet. The validity of a waveletbased series representation for the fractional Brownian motion has been established in Wang (1997), Meyer, Sellan and Taqqu (1999) and Ayache and Taqqu (2002). In particular, Ayache and Taqqu (2002) established the optimality of waveletbased series representation of the fractional Brownian sheets.
Remark 5
For , the rates of convergence match exactly those in Benhaddou (2017), and, with , those in Benhaddou et al. (2013) in their bivariate white noise case.
Remark 6
If we hold the second dimension fixed, our rates are comparable to those in Wang (1997) and Wishart (2013) in their standard deconvolution with the oneparameter fractional Gaussian noise case.
4 Simulation Study
In order to investigate the effect of the longmemory on the performance of our estimator, we carried out a limited simulation study. We implemented an estimation algorithm for the model (5) using a modification of WaveD method of Raimondo and Stewart (2007). We evaluated mean integrated square error (MISE) of the functional deconvolution estimator.

We generated the data using equation (5) with convolution kernel , and . In particular, we chose to be a LIDAR or Doppler signal over the grid with , and we chose over the grid with , and . The LIDAR and Doppler signals were generated from the waved package. All test functions were scaled to have a unit norm. Bear in mind that, though is a product of two univariate functions, the method does not “recognize” this and, therefore, cannot take advantage of this information. Also, notice that with our choice of , the degree of illposedness (DIP) of the convolution is DIP=0.5.

We simulated the LM error in the direction of both and using the fracdiff package of R available from CRAN. In particular, the fracdiff.sim command was used to simulate two onedimensional fractionally differenced ARIMA (fARIMA) sequences which behave similar to a fractional Gaussian noise, and then multiply them together to generate a twodimensional error structure. For the LM parameters , we used different combinations of the levels , , and . To create a dependence structure that is anisotropic, we only used combinations such that . Note that the fractional differencing parameter is obtained from by the relation .

The choice of in (5) was determined by the blurred signaltonoise ratio (SNR), where We considered three choices, SNR=10dB (high noise), 20dB (medium noise) and 30dB (low noise).

To compute the estimator, we used the function FWaveD and IWaveD from the waved package. First, we applied the Meyer wavelet transform to the data in direction by FWaveD, with no thresholding by setting thr
to be a zero vector. Then, we applied the second thresholded Meyer wavelet transform in
direction by FWaveD, however, the level dependent threshold were modified according to (21). For the tuning parameter , we tried different values ranging from to , but the performance of the estimator reached its best at the latter value. Therefore, our default choice of was . The finest resolution level was estimated using the default method in the waved package, together with the choice .
Table 1 reports the averages of the errors over 2000 simulation runs using a sample size of . Figure (1) illustrates the function and its estimator and shows a relatively good precision that deteriorates as the anisotropic pair of LM parameters decreases from 0.8 to 0.2. Notice here that according to Table 1, the influence of LM is a little more pronounced in the direction of the dimension than it is in the direction of . Table 1 complements Figure (1) and confirms our theoretical results in the sense that as the level of LM increases ( gets smaller), the mean integrated squared error increases (the performance deteriorates).
DIP=0.5  (,)  

signal  (0.6,0.8)  (0.4,0.8)  (0.2,0.8)  (0.4,0.6)  (0.2,0.6)  (0.2,0.4)  
Lidar (10dB)  0.20  0.0760 (3)  0.0851 (3)  0.0956 (3)  0.1038 (3)  0.1293 (3)  0.2251 (3) 
Lidar (20dB)  0.06  0.0423 (4)  0.0453 (4)  0.0486 (4)  0.0536 (4)  0.0612 (4)  0.0791 (4) 
Lidar (30dB)  0.02  0.0120 (5)  0.0132 (5)  0.0149 (5)  0.0180 (5)  0.0210 (5)  0.0331 (5) 
Doppler (10dB)  0.15  0.2601 (3)  0.2640 (3)  0.2749 (3)  0.2836 (3)  0.3090 (3)  0.3922 (3) 
Doppler (20dB)  0.05  0.1247 (4)  0.1264 (4)  0.1292 (4)  0.1331 (4)  0.1376 (4)  0.1578 (4) 
Doppler (30dB)  0.015  0.0513 (5)  0.0510 (5)  0.0522 (5)  0.0527 (5)  0.0541 (5)  0.0591 (5) 
5 The dimensional case: estimation, convergence rates and asymptotic optimality.
Consider the dimensional case of model (1)
(35) 
Here, , , is the convolution kernel, is an anisotropic dimensional fBs, , and , , is the parameter of the longmemory in the direction of the variable. An parameter fBs is defined as a zeromean Gaussian process whose covariance function is of the form
(36) 
for some , , , where is the Hurst parameter in the direction of the dimension.
Consider a bandlimited wavelet basis (e.g., Meyertype) on , and form the product wavelet basis
(37) 
where , and
(38) 
Then, the function can be expanded over these bases with coefficients into a wavelet series
(39) 
where . Apply the dimensional Fourier transform to equation (35) to obtain
(40) 
where . Then, applying Plancherel formula in all directions, the wavelet coefficients in (39) can be rewritten as
(41) 
Therefore, an unbiased estimator for is given by
(42) 
Finally, consider the hard thresholding estimator for
(43) 
where
(44) 
with , and the quantities , , , and are to be determined.
Assumption A.4. In the Fourier domain, the convolution kernel , for some positive constants , and and , independent of , is such that
(45) 
Lemma 4
Consequently, choose the thresholds of the from
(48) 
where is some positive constant that is independent of ,
Comments
There are no comments yet.