On LSE in regression model for long-range dependent random fields on spheres

TL;DR

This paper investigates the asymptotic behavior of least squares estimators in regression models involving long-range dependent random fields on spheres, deriving their limit distributions and convergence rates.

Contribution

It provides the first derivation of limit distributions and convergence rates for LSEs in long-range dependent random fields on spheres under general conditions.

Findings

Limit distributions can be non-Gaussian.

Convergence rates are explicitly derived.

Simulation studies support theoretical results.

Abstract

We study the asymptotic behaviour of least squares estimators in regression models for long-range dependent random fields observed on spheres. The least squares estimator can be given as a weighted functional of long-range dependent random fields. It is known that in this scenario the limits can be non-Gaussian. We derive the limit distribution and the corresponding rate of convergence for the estimators. The results were obtained under rather general assumptions on the random fields. Simulation studies were conducted to support theoretical findings.

Equations

$$\text{B}(r):=\mathrm{Cov}\left(\eta (x),\eta (y)\right)={\int }_0^{\infty }{Y}_d(rz)\mathrm{d}\Phi (z),x,y\in {\mathbb{R}}^d,$$

$$Y_d(z):=2^{(d-2)/2}\Gamma \left(\frac{d}{2}\right)J_{(d-2)/2}(z)z^{(2-d)/2},z\ge 0,$$

$$\Phi (z)=2{\pi }^{d/2}{\Gamma }^{-1}\left(d/2\right){\int }_0^z{u}^{d-1}f(u)\mathrm{d}u,z\ge 0,u^{d-1}f(u)\in L_1({\mathbb{R}}_{+}).$$

$$\eta (x)={\int }_{{\mathbb{R}}^d}{e}^{i<\lambda ,x>}\sqrt{f\left(\Vert \lambda \Vert \right)}W(\mathrm{d}\lambda ),$$

$$L_2(\mathbb{R},\varphi (w)dw)=\left\{G:{\int }_{\mathbb{R}}{G}^2(w)\varphi (w)\mathrm{d}w<\infty \right\},\varphi (w):=\frac{1}{\sqrt{2\pi }}{e}^{-\frac{w^2}{2}}.$$

$$G(w)=\sum _{j=0}^{\infty }\frac{C_j{H}_j(w)}{j!},C_j:={\int }_{\mathbb{R}}G(w)H_j(w)\varphi (w)\mathrm{d}w.$$

$$\sum _{j=0}^{\infty }\frac{C_j^2}{j!}={\int }_{\mathbb{R}}{G}^2(w)\varphi (w)\mathrm{d}w.$$

$$\mathbf{E}{\xi }_j{\xi }_k=\{\begin{array}{ll}1,& \text{mbox}ifk=j,\\ r_j,& \text{mbox}ifk=j+p\text{mbox}and1\le j\le p,\\ 0,& \text{mbox}otherwise,\end{array}$$

$$\mathbf{E}\prod _{j=1}^p{H}_{k_j}({\xi }_j)H_{m_j}({\xi }_{j+p})=\prod _{j=1}^p{\delta }_{k_j}^{m_j}{k}_j!r_j^{k_j},$$

$$\underset{\partial \Delta (r)}{\int }\underset{\partial \Delta (r)}{\int }G(\Vert x-y\Vert )\sigma (\mathrm{d}x)\sigma (\mathrm{d}y)={\left|\partial \Delta \right|}^2{r}^{2d-2}\mathrm{E}G(\Vert U-V\Vert )=$$

$$={\left|\partial \Delta \right|}^2{r}^{2d-2}{\int }_0^{diam\left\{\partial \Delta (r)\right\}}G(\rho ){\psi }_{\Delta (r)}(\rho )d\rho .$$

$${\psi }_{\Delta (r)}(\rho )=\frac{1}{\sqrt{\pi }}\Gamma \left(\frac{d}{2}\right){\Gamma }^{-1}\left(\frac{d-1}{2}\right)r^{1-d}{\rho }^{d-2}{\left(1-\frac{{\rho }^2}{4u^2}\right)}^{\frac{d-3}{2}},0<\rho <2r.$$

$$\underset{\lambda \to \infty }{\lim }\frac{L(\lambda t)}{L(\lambda )}=1.$$

$$\underset{\lambda \to \infty }{\lim }\frac{g(\lambda t)}{g(\lambda )}=t^{\tau }.$$

$$\forall \lambda >1:\frac{L(\lambda x)}{L(x)}-1\sim k(\lambda )g(x),x\to \infty ,$$

$$h_{\tau }(\lambda )=\{\begin{array}{l}\ln (\lambda ),if\tau =0,\\ \frac{{\lambda }^{\tau }-1}{\tau },if\tau \ne 0.\end{array}$$

$$\frac{L(\lambda x)}{L(x)}-1=\frac{\ln (\lambda )+\ln (x)}{\ln (x)}-1=\ln (\lambda )\cdot \frac{1}{\ln (x)}.$$

$$\rho \left(Y_1,Y_2\right)=\underset{z\in \mathbb{R}}{\sup }\left|P\left(Y_1\le z\right)-P\left(Y_2\le z\right)\right|.$$

$$\rho \left(X+Y,Z\right)\le \rho (X,Z)+\rho \left(Z+\epsilon ,Z\right)+P\left(\left|Y\right|\ge \epsilon \right).$$

$$\xi (x)=a_{r}h(x)+\psi (x),x\in S(r),r>0,$$

$${\hat{a}}_r=\frac{\underset{S(r)}{\int }\xi (x)h(x)\sigma (\mathrm{d}x)}{\underset{S(r)}{\int }{h}^2(x)\sigma (\mathrm{d}x)}=\frac{\underset{S(r)}{\int }\xi (x)h_{rad}(\Vert x\Vert )h_{sp}\left(\frac{x}{\Vert x\Vert }\right)\sigma (\mathrm{d}x)}{\underset{S(r)}{\int }{h}_{rad}^2(\Vert x\Vert )h_{sp}^2\left(\frac{x}{\Vert x\Vert }\right)\sigma (\mathrm{d}x)}$$

$$=\frac{\underset{S(r)}{\int }\xi (x)h_{sp}\left(\frac{x}{\Vert x\Vert }\right)\sigma (\mathrm{d}x)}{h_{rad}(r)\underset{S(r)}{\int }{h}_{sp}^2\left(\frac{x}{\Vert x\Vert }\right)\sigma (\mathrm{d}x)}.$$

$$B(0)=1,B(x)=\mathbf{E}\eta (0)\eta (x)={\Vert x\Vert }^{-\alpha }{L}_0(\Vert x\Vert ),$$

$$f(\Vert \lambda \Vert )=c_2(d,\alpha ){\Vert \lambda \Vert }^{\alpha -d}L\left(\frac{1}{\Vert \lambda \Vert }\right),$$

$$\left|1-\frac{L(tr)}{L(r)}\right|\le Cg(r)h_{\tau }(t),t\ge 1,$$

$$\left|1-\frac{L^{k/2}(tr)}{L^{k/2}(r)}\right|\le Cg(r)h_{\tau }(t)t^{\delta },t\ge 1.$$

$${\hat{a}}_r=a_r+\frac{\underset{S(r)}{\int }G(\eta (x))h_{sp}\left(\frac{x}{\Vert x\Vert }\right)\sigma (\mathrm{d}x)}{h_{rad}(r)\underset{S(r)}{\int }{h}_{sp}^2\left(\frac{x}{\Vert x\Vert }\right)\sigma (\mathrm{d}x)}.$$

$$X_{r,G}:=c_h(r)c_r(d,\alpha )({\hat{a}}_r-a_r),$$

$$c_h(r):=h_{rad}^{-1}(r)\underset{S(r)}{\int }{h}^2(x)\sigma (\mathrm{d}x)\text{ and }{c}_r(d,\alpha )=\frac{\kappa !c_2^{-\kappa /2}(d,\alpha )}{C_{\kappa }{r}^{d-1-\frac{\kappa \alpha }{2}}{L}^{\frac{\kappa }{2}}(r)}.$$

$$X_{r,G}=c_r(d,\alpha )\underset{S(r)}{\int }G(\eta (x))h_{sp}\left(\frac{x}{\Vert x\Vert }\right)\sigma (\mathrm{d}x)$$

Full text

{amscode}

62H11, 62J05, 60G60.

1 Introduction

Studying the asymptotic behavior of the regression parameter estimators is an important topic in Spatial Statistics, see [1], [2], [3]. A potential area of application of the results in this paper is cosmology, and more specifically, models involving the cosmic microwave background (CMB). The recent data obtained by the Wilkinson Microwave Anisotropy Probe and Planck missions, combined with data that will be obtained during the future missions, will allow to predict the behaviour of CMB at increasing distances that are far beyond the reach of any modern telescopes. Another area of application is studying the porosity of mineral resources, see, for example, [4] about data and analysis of the porosity of hydrocarbon reservoirs. Some other potential applications of the model discussed in this article might involve 3D rupture models, soil fertility and paleomagnetism.

One of the first works that developed statistical theory of random fields was the book by [5]. Direct probability techniques for models involving random fields were proposed and various estimators of regression parameters were investigated. The asymptotic behaviour of some of these estimators was studied in [6]. Both Gaussian and non-Gaussian limits were obtained for particular models. No results about rates of convergence were given.

In this paper we consider the regression model $\xi(x)=ah(x)+\psi(x)$, where $h(\cdot)$ is a known deterministic function, $a$ is an unknown parameter and the error $\psi(\cdot)$ is a random field. For results concerning such regression models see [7] and references therein. Since there are numerous situations in practice when random processes and fields are non-Gaussian, we consider a particular case when $\psi(\cdot)$ is a non-linear function of a Gaussian random field. Such random fields are very common in non-Gaussian modeling since they can be analyzed using Wiener chaos expansions and in many cases offer a good data approximation, see [8, 9]. In this article, we consider the underlying Gaussian random fields to be long-range dependent.

Long-range dependence is a well-established phenomenon that can be observed in various fields such as hydrology, agriculture, image analysis, earth sciences, cosmology, just to name a few. For this reason, models involving long-range dependent random fields have been an object of statistical interest for years, see [10, 11, 12, 13].

Various earth science and cosmology applications require studying random fields defined on surfaces, for example, see [14, 15]. In this paper we consider random fields that are defined on expanding spheres. One can find detailed information about such spherical random fields in [16, 5].

The regression model studied in this paper was first considered in [5]. For this model the best linear unbiased estimator (BLUE) of the unknown parameter $a$ was derived. In [17] the least squares estimator (LSE) of $a$ was obtained. Its mean-square efficiency was compared to BLUE and it was shown that, apart from some degenerate cases, LSE is less efficient than BLUE. However, for most functions $h(\cdot)$ BLUE is much more difficult to compute than LSE. For this reason using LSE is preferable to BLUE in practice. No results concerning limit distributions of the regression parameter estimators were obtained in [17, 5] or other literature known to the authors. This article investigates the limit distribution of LSE and its rate of convergence.

For the particular case $h(\cdot)\equiv 1$, some approaches to study non-linear functionals of random fields on surfaces and rates of convergence were proposed in [18]. This article shows how the methodology developed in [18] can be extended to the general case of weighted non-linear functionals of random fields, i.e. the case of arbitrary functions $h(\cdot)$. To reduce repetitions to a minimum, in this article we present only those parts of proofs that are new or require modifications of results in [18].

The main goal of this article is to demonstrate potential applications of the developed methodology to statistical problems. We focus on the LSE for a regression model, however there are numerous other statistical problems that can be studied using weighted non-linear functionals of random fields. For example, by choosing the appropriate integrand functions, various characteristics of random fields, such as moments, can be represented by the considered functionals. In [18] it was shown that Minkowski functionals can be obtained by choosing indicator functions as the integrands. Another important statistical application is assigning different weights to different observation regions. In spatial statistics it is often used to control the edge effects. Also, assigning weights to observations or observation regions is an important case of data tapering that is often used to reduce effects of less reliable observations.

Section 5 presents detailed simulation studies that back up the theory and suggest new research questions. To the best of our knowledge, simulation studies have not been done in the available literature on convergence rates in non-central limit theorems for long-range dependent random fields. The numerical study suggests that the exponential rate of convergence obtained in all known theoretical studies might be improved. For the cases of a sphere and a cube, various observation windows and weight functions were investigated. The detailed methodology for practical simulations and studying weighted functionals of random fields observed on surfaces is provided. The corresponding R codes for one-core and parallel computing are freely available on-line and can be used by researchers in this field.

To obtain the main results some fine properties of the Fourier transforms on surfaces are employed. Note that the Fourier transforms were more frequently studied on solid bodies than on surfaces. The assumptions for obtaining the required rate of decay are weaker in the case of solid bodies. Therefore, all results derived in this article can be analogously obtained for solid bodies.

As mentioned earlier, to study the asymptotic behaviour of LSE we investigate weighted non-linear functionals of Gaussian random fields. Non-linear functionals of random processes were considered in [6, 19]. In [6] only the case of short-range dependent random processes and Gaussian asymptotics for weighted functionals was studied. In [19] non-Gaussian asymptotics were obtained for long-range dependent random processes. The results in [19] can be obtained from the methodology of this article as random processes can be considered as random fields indexed by a one-dimensional Euclidean space. For random fields and the simplest particular case $h(\cdot)\equiv 1$ it was shown that the functionals can produce non-Gaussian limits, see [20, 21, 22]. The more detailed overview of this particular case was given in [23]. In the general case, weighted linear functionals of the Gaussian random fields were studied in [24]. The assumptions on weighted functions used in this article are weaker than those used in [24]. Moreover, we do not apply any conditions on Fourier transforms of weight functions. The known Gaussian asymptotics are a particular case when integrand functions are linear. None of the above studies provided the rate of convergence to the obtained limits. To the best of our knowledge, this paper is the first work that provides the rate of convergence for weighted functionals and LSE of long-range dependent random fields on hypersurfaces.

The article is organized as follows. In Section 2 we recall some basic definitions of the spectral theory of random fields. Section 3 presents the model and states the results. Proofs of the results are provided in Section 4. In Section 5 some simulation studies are presented to confirm theoretical findings. Conclusions and directions for future research are stated in Section 6.

2 Definitions

In this section we provide the main definitions that are used in this article.

In what follows $|\cdot|$ and $\|\cdot\|$ denote the Lebesgue measure and the Euclidean distance in $\mathbb{R}^{d}$, $d\geq 2$, respectively. Let $B(y,s)$ be a $d$-dimensional ball with centre $y$ and radius $s$, and let $S(r)$ be a $(d-1)$-dimensional sphere in $\mathbb{R}^{d}$ with the centre at the origin and radius $r$. We use the symbols $C$ and $\delta$ to denote constants which are not important for our exposition. Moreover, the same symbol may be used for different constants appearing in the same proof.

Let $(\Omega,\mathcal{F},P)$ be a complete probability space and let $T\subseteq\mathbb{R}^{d}$ be a set.

Definition 1*.*

[12] A random field is a function $\eta(w,x):\Omega\times T\rightarrow\mathbb{R}$ such that $\eta(w,x)$ is a random variable for any $x\in T$. It will also be denoted as $\eta(x),\ x\in T.$

Definition 2*.*

[12] A random field satisfying $\mathbf{E}\eta^{2}(x)<\infty$ is called homogeneous in the wide sense if its mathematical expectation $m(x)=\mathbf{E}\eta(x)$ and covariance function ${\rm{B}}(x,y)={\mathbf{E}}\left[\eta(x)-m(x)\right]\left[\eta(y)-m(y)\right]$ are invariant with respect to the group of shifts $\tau$ in $\mathbb{R}^{d}$, that is, $m(x)=m(x+\tau)$, ${\rm{B}}(x,y)={\rm{B}}(x+\tau,y+\tau)$ for any $x,y\in T$.

It means that $\mathbf{E}\eta(x)=c=const$, and the covariance ${\rm{B}}(x,y)$ depends only on the difference $x-y$.

Definition 3*.*

[12] A random field satisfying $\mathbf{E}\eta^{2}(x)<\infty$ is called isotropic in the wide sense if its mathematical expectation $m(x)$ and covariance function ${\rm{B}}(x,y)$ are invariant with respect to rotations.

In this article homogeneity and isotropy in the wide sense will be simply called homogeneity and isotropy.

Let us consider a measurable mean-square continuous zero-mean homogeneous isotropic real-valued random field $\eta(x),\ x\in\mathbb{R}^{d},$ defined on a probability space $(\Omega,\mathcal{F},P),$ with the covariance function

[TABLE]

where $r:=\left\|x-y\right\|,$ $\Phi(\cdot)$ is the isotropic spectral measure, the function $Y_{d}(\cdot)$ is defined by

[TABLE]

and $J_{(d-2)/2}(\cdot)$ is the Bessel function of the first kind of order $(d-2)/2.$

Definition 4*.*

The random field $\eta(x),$ $x\in\mathbb{R}^{d},$ defined above is said to possess an absolutely continuous spectrum if there exists a function $f(\cdot)$ such that

[TABLE]

The function $f(\cdot)$ is called the isotropic spectral density function of the field $\eta(x).$ The field $\eta(x)$ with an absolutely continuous spectrum has the following representation

[TABLE]

where $W(\cdot)$ is the complex Gaussian white noise random measure on $\mathbb{R}^{d}.$

Let $H_{k}(u)$, $k\geq 0$, $u\in\mathbb{R}$, be the Hermite polynomials, see [25]. These polynomials form a complete orthogonal system in the Hilbert space

[TABLE]

An arbitrary function $G\in{L}_{2}(\mathbb{R},\phi(w)\,dw)$ admits the mean-square convergent expansion

[TABLE]

By Parseval’s identity

[TABLE]

It will be shown that studying asymptotics of non-linear functionals defined by a function $G\in{L}_{2}(\mathbb{R},\phi(w)\,dw)$ can be done by investigating leading terms in their expansions. The following definitions and remarks provide basic notations and tools to formulate and prove these asymptotic results.

Definition 5*.*

[21] Let $G\in{L}_{2}(\mathbb{R},\phi(w)\,dw)$ and assume there exists an integer $\kappa\in\mathbb{N}$ such that $C_{j}=0$, for all $0\leq j\leq\kappa-1,$ but $C_{\kappa}\neq 0.$ Then $\kappa$ is called the Hermite rank of $G(\cdot)$ and is denoted by $H\mbox{rank}\,G.$

The following remark gives a well known property of Hermite polynomials of Gaussian vectors, see [12].

Remark 1*.*

If $(\xi_{1},\ldots,\xi_{2p})$ is a $2p$-dimensional zero-mean Gaussian vector with

[TABLE]

then

[TABLE]

where $\delta_{k}^{m}=\begin{cases}1,&\mbox{if }k=m,\\ 0,&\mbox{if }k\neq m,\end{cases}$  is the Dirac delta function.

Let $\Delta$ be a bounded set in $\mathbb{R}^{d},\,d\geq 2,$ with a boundary $\partial\Delta$. Let $\Delta(r)$, $r>0,$ be the homothetic image of the set $\Delta$ with the centre of homothety at the origin and the coefficient $r>0$, that is $|\Delta(r)|=r^{d}|\Delta|$. Let $\sigma(\cdot)$ be the $d-1$-dimensional Lebesgue measure on $\partial\Delta(r)$. Let $U$ and $V$ be two independent and uniformly distributed random vectors on the hypersurface $\partial\Delta(r)$. We denote by $\psi_{\Delta(r)}(\rho),$ $\rho\geq 0,$ the probability density function of the distance $\left\|U-V\right\|$ between $U$ and $V.$ Note that $\psi_{\Delta(r)}(\rho)=0$ if $\rho>diam\left\{\partial\Delta(r)\right\}.$ Using the above notations, we obtain the representation

[TABLE]

[TABLE]

For various hypersurfaces $\partial\Delta$ explicit expressions for $\psi_{\Delta(r)}(\cdot)$ are presented in [12]. In the case of spheres it takes the following form.

Remark 2*.*

If $\partial\Delta(r)=S(r)$, then

[TABLE]

Definition 6*.*

[26] A measurable function $L:(0,\infty)\rightarrow(0,\infty)$ is said to be slowly varying at infinity if for all $t>0$

[TABLE]

If $L(\cdot)$ varies slowly, then $r^{a}L(r)\rightarrow\infty,$ $r^{-a}L(r)\rightarrow 0$ for an arbitrary $a>0$ when $r\rightarrow\infty,$ see Proposition 1.3.6 in [26].

Definition 7*.*

[26] A measurable function $g:(0,\infty)\rightarrow(0,\infty)$ is said to be regularly varying at infinity, denoted $g(\cdot)\in R_{\tau}$, if there exists $\tau$ such that, for all $t>0,$ it holds that

[TABLE]

Definition 8*.*

[26] Let $g:(0,\infty)\rightarrow(0,\infty)$ be a measurable function and $g(x)\rightarrow 0$ as $x\rightarrow\infty$. A slowly varying function $L(\cdot)$ is said to be slowly varying with remainder of type 2, or that it belongs to the class SR2, if

[TABLE]

for some function $k(\cdot)$.

If there exists $\lambda$ such that $k(\lambda)\neq 0$ and $k(\lambda\mu)\neq k(\mu)$ for all $\mu$, then $g(\cdot)\in R_{\tau}$ for some $\tau\leq 0$ and $k(\lambda)=ch_{\tau}(\lambda)$, where

[TABLE]

Remark 3*.*

An example of a function that satisfies Definition 8 for $\tau=0$ is $L(x)=\ln(x).$ Indeed,

[TABLE]

Definition 9*.*

Let $Y_{1}$ and $Y_{2}$ be arbitrary random variables. The uniform (Kolmogorov) metric for the distributions of $Y_{1}$ and $Y_{2}$ is defined by

[TABLE]

The next result follows from Lemma 1.8 in [27].

Lemma 1**.**

If $X,Y$ and $Z$ are arbitrary random variables, then for any $\varepsilon>0$

[TABLE]

3 Model and results

Let us consider the random field

[TABLE]

where $h(\cdot)$ is a known deterministic function, $a_{r}$ is an unknown parameter and $\psi(\cdot)$ is a homogeneous isotropic mean-square random field with mean 0.

Suppose that $h(x):=h_{rad}(\|x\|)h_{sp}\left(\frac{x}{\|x\|}\right),$ where $h_{rad}(\cdot)$ is a radial function such that $h_{rad}(r)\neq 0$ for all $r>0$ and $h_{sp}(\cdot)\not\equiv 0$ is a smooth bounded function defined on the unit sphere $S(1)$.

In this case, the least squares estimator (LSE) of the coefficient $a_{r}$ has the explicit form, see [5],

[TABLE]

[TABLE]

Let $\psi(x)=G(\eta(x)),$ where $G\in{L}_{2}(\mathbb{R},\phi(w)\,dw)$ and $\eta(x)$ is a random field that satisfies the following assumptions.

Assumption 1*.*

Let $\eta(x),$ $x\in\mathbb{R}^{d}$, be a homogeneous isotropic Gaussian random field with $\mathbf{E}\eta(x)=0$ and a covariance function $B(x)$ such that

[TABLE]

where $L_{0}(\left\|\cdot\right\|)$ is a function slowly varying at infinity.

This assumption is a classical way to introduce hyperbolically decaying dependencies between observations, see [20, 11, 22, 12] and references therein. Random fields satisfying this assumption for $\alpha>d$ are weakly-dependent. If $\alpha\in(0,d)$ then the long-range dependence case is considered.

Assumption 2*.*

The random field $\eta(x),$ $x\in\mathbb{R}^{d},$ has the spectral density

[TABLE]

where $c_{2}(d,\alpha):=\frac{\Gamma\left(\frac{d-\alpha}{2}\right)}{2^{\alpha}\pi^{d/2}\Gamma\left(\frac{\alpha}{2}\right)},$ and $L(\left\|\cdot\right\|)$ is a locally bounded function slowly varying at infinity which satisfies for sufficiently large $r$ the condition

[TABLE]

where $g(\cdot)\in R_{\tau},\tau\leq 0$, such that $g(x)\to 0,\ x\to\infty$, and $h_{\tau}(t)$ is defined by (3).

Long-range dependence is usually introduced by requiring the hyperbolic decay of covariance functions or power-type singularity of the corresponding spectral densities. For many real data these two definitions are operationally equivalent, see also Tauberian-Abelian theorems in [28]. However, there are cases when Assumption 2 does not follow from Assumption 1, see [28]. To make all following results rigorous, we require that the both assumptions hold.

Examples of popular classes of random fields that satisfy Assumption 1 and 2 simultaneously are Bessel, Cauchy, and Linnik random fields.

Remark 4*.*

By Tauberian and Abelian theorems, see [28], for $L_{0}(\cdot)$ and $L(\cdot)$ given in Assumptions 1 and 2 it holds $L_{0}(r)\sim L(r),$ $r\to+\infty.$

In [23] the following property of slowly varying functions satisfying conditions (4) was proven.

Remark 5*.*

If $L(\cdot)$ satisfies (4), then for any $k\in\mathbb{N}$, $\delta>0$, and sufficiently large $r$

[TABLE]

Since $\psi(x)=G(\eta(x)),$ the LSE has the form

[TABLE]

Let $H\mathrm{rank}\,G=\kappa$. Our main object of interest is the random variable

[TABLE]

where

[TABLE]

It is straightforward to see that

[TABLE]

[TABLE]

[TABLE]

Theorem 1**.**

Suppose that $H\mathrm{rank}\,G=\kappa\in\mathbb{N}$ and $\eta(x),$ $x\in\mathbb{R}^{d},$ satisfies Assumption 1 for $\alpha\in(0,(d-1)/\kappa)$. If at least one of the following random variables

[TABLE]

has a limit distribution, then the limit distributions of the other random variables also exist and they coincide when $r\rightarrow\infty.$

By Theorem 1 to study limit distributions one can use $X_{r,\kappa}$ instead of $X_{r,G}$.

In Theorem 2 below, we show that the limit distribution of $X_{r,\kappa}$ is a Hermite-type random variable, which can be represented by a Wiener-Itô integral. Let

[TABLE]

be the Fourier transform of the function $h_{sp}(\cdot)$ over the $(d-1)$-dimensional sphere with radius 1. Using the decay properties of the Fourier transform one can prove the following.

Lemma 2**.**

If $\tau_{1},...,\tau_{\kappa},$ $\kappa\geq 1,$ are such positive constants, that $\sum_{i=1}^{\kappa}\tau_{i}<d-1,$ then

[TABLE]

Let

[TABLE]

where ${\int\limits_{\mathbb{R}^{d\kappa}}}^{\prime}$ denotes the multiple Wiener-Itô integral.

Theorem 2**.**

Let $\eta(x),$ $x\in\mathbb{R}^{d},$ be a homogeneous isotropic Gaussian random field with $\mathbf{E}\eta(x)=0.$ If Assumptions 1 and 2 hold, $\alpha\in(0,(d-1)/\kappa),$ and $H\mathrm{rank}\,G=\kappa\in\mathbb{N},$ then for $r\to\infty$ the random variables $X_{r,\kappa}$ converge weakly to $X_{\kappa}$.

To obtain the rate of convergence of $X_{r,G}$ to $X_{\kappa}$ we will use some fine properties of Hermite-type distributions. We will denote the Wiener-Itô integrals of rank $\kappa$ by $I_{\kappa}(f)$, where $f(\cdot)$ is an integrand. For more details about Wiener-Itô integrals and admissible functions $f(\cdot)$ one can refer to [29, 30]. The following result was obtained in [23] for specific form of the integrand. Since the proof does not rely on the form of the integrand, this theorem can be easily generalized as follows.

Theorem 3**.**

[23]* For any $\kappa\in\mathbb{N}$ and an arbitrary positive $\varepsilon$ it holds*

[TABLE]

where $b=1$ if $\kappa<3$ and $b=1/\kappa$ if $\kappa\geq 3$.

Also, we will use the following result.

Lemma 3**.**

[18]* Let $f_{1}(\cdot)$ and $f_{2}(\cdot)$ be symmetric functions in $L_{2}({\mathbb{R}^{d}}),\,d\geq 1$. Then,*

[TABLE]

and

[TABLE]

Now we are ready to formulate the main result.

Theorem 4**.**

Let $H\mathrm{rank}\,G=\kappa\in\mathbb{N}$ and Assumptions 1 and 2 hold for $\alpha\in(0,\frac{d-1}{\kappa})$.

If $\tau\in\left(-\frac{d-\kappa\alpha}{2},0\right)$ then for any $\varkappa<\frac{b}{2+b}\min\left(\frac{\alpha(d-1-\kappa\alpha)}{d-1-(\kappa-1)\alpha},\varkappa_{1}\right)$

[TABLE]

where $\varkappa_{1}:=\min\left(-2\tau,\frac{1}{\frac{1}{d-2\alpha}+\dots+\frac{1}{d-\kappa\alpha}+\frac{1}{d-1-\kappa\alpha}}\right)$ and $b$ is the parameter from Theorem 3.

If $\tau=0$ then

[TABLE]

4 Proofs of the result from section 3

In this section the results stated in Section 2 are proved. In [18] Theorems 1, 2, 4 and Lemma 2 were proved in the particular case $h_{sp}(\cdot)\equiv 1$. Since the proofs for the general setting are very similar to the ones in [18], only parts of the proofs that differ are provided.

Proof of Theorem 1.

By Remark 1 it holds $\mathrm{Var}X_{r,G}=\mathrm{Var}X_{r,\kappa}+\mathrm{Var}V_{r}$.

By (1) and (2) we get

[TABLE]

[TABLE]

[TABLE]

[TABLE]

It follows from $z^{-\alpha}\left|L_{0}\left(z\right)\right|\in[0,1],$ $z\geq 0,$ that

[TABLE]

[TABLE]

Considering $|L_{0}(\cdot)|$ as a new slowly varying function, one can estimate the integral above using bounds (9) and (10) in [18] to obtain

[TABLE]

where $\delta$ is an arbitrary number in $(0,\min(\alpha,d-1-\kappa\alpha))$.

Similar to $\mathrm{Var}\,V_{r}$ we obtain

[TABLE]

Let us consider a constant $c_{h}>\max\limits_{x}|h_{sp}(x)|$ and let $c_{I}:=\int\limits_{{S}(1)}(h_{sp}(x)+c_{h})\mathrm{d}\sigma(x)$. The function $\displaystyle\frac{h_{sp}(x)+c_{h}}{c_{I}}$ can be considered as a probability density function of some random variable $X$ defined on $S(1)$. Therefore,

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where $X_{1},\,X_{2}$ are random variables distributed on $S(1)$ with the probability density $\displaystyle\frac{h_{sp}(x)+c_{h}}{c_{I}}$ and $Y_{1},\,Y_{2}$ are uniformly distributed on $S(1)$ random variables.

Let $Z_{1}:=\left\|X_{1}-X_{2}\right\|$, $Z_{2}:=\left\|X_{1}-Y_{1}\right\|$, $Z_{3}:=\left\|Y_{1}-Y_{2}\right\|$ and denote their corresponding probability densities as $\bar{\psi}_{S(1)}(\cdot)$, $\tilde{\psi}_{S(1)}(\cdot)$ and $\psi_{S(1)}(\cdot)$. Thus

[TABLE]

[TABLE]

If $\alpha\in(0,(d-1)/\kappa)$ then by asymptotic properties of integrals of slowly varying functions (see Theorem 2.7 in [31]) we get

[TABLE]

where $\bar{c}_{1}(\kappa,\alpha):=c_{I}^{2}\int\limits_{0}^{2}\frac{\bar{\psi}_{S(1)}(z)}{z^{\alpha\kappa}}\mathrm{d}z$, $\tilde{c}_{1}(\kappa,\alpha):=-2c_{I}|{S}(1)|c_{h}\int\limits_{0}^{2}\frac{\tilde{\psi}_{S(1)}(z)}{z^{\alpha\kappa}}\mathrm{d}z$, and ${c}_{1}(\kappa,\alpha):={-\displaystyle|{S}(1)|^{2}c_{h}^{2}\int\limits_{0}^{2}\frac{{\psi}_{S(1)}(z)}{z^{\alpha\kappa}}\mathrm{d}z}$.

By (2) it holds

[TABLE]

[TABLE]

[TABLE]

By (7) and (8), one can choose $\beta_{1}=1/2$ and make $\delta$ arbitrary close to 0 to obtain

[TABLE]

Thus

[TABLE]

and

[TABLE]

which completes the proof. ∎

To obtain Lemma 2 a result analogous to the average decay rate of the Fourier transforms of surface measures in [32] is needed. Therefore we prove the following lemma.

Lemma 4**.**

For sufficiently large $r$ it holds

[TABLE]

Proof.

Using Theorem 1 in [33, p. 322] we obtain $|K(x)|\leq C\|x\|^{\frac{1-d}{2}}.$

It follows

[TABLE]

Now, we can prove Lemma 2 using bound (9) analogously to Lemma 2 in [18].

Proof of Lemma 2.

For $\kappa=1$ it holds $d-\tau_{1}>1$. By the properties of the Fourier transform $|K(\lambda)|\leq C|S(1)|$ for all $\lambda\in\mathbb{R}^{d}$. Therefore, using integration formula for polar coordinates we get

[TABLE]

[TABLE]

By Lemma 4 we obtain

[TABLE]

[TABLE]

For $\kappa>1$ one can obtain (6) by the recursive estimation routine and the change of variables $u=\lambda_{\kappa-1}+\lambda_{\kappa}$ and $\tilde{\lambda}_{\kappa-1}={\lambda_{\kappa-1}}/{\left\|u\right\|}$, see [18] for details. ∎

The remaining theorems consider the random variables $X_{r,\kappa}$ and $X_{\kappa}$. As it was mentioned before, in [18] these theorems were proved for the particular case $h_{sp}(\cdot)\equiv 1$. Up to multiplication by a constant, $X_{\kappa}$ in this paper is identical to its counterpart $X_{\kappa}(\Delta)$ in [18] and does not depend on the weight function $h_{sp}(\cdot)$. If $X_{r,\kappa}$ is represented as the multiple Wiener-Itô stochastic integral the function $h_{sp}(\cdot)$ appears only in the function $K(\cdot)$. Note that Lemma 2 and bound (9) for $K(\cdot)$ are the same as the corresponding results for $\mathcal{K}(\cdot)$ in [18]. Since Lemma 2 and bound (9) are the only results concerning $K(\cdot)$ that are required, to prove theorems 2 and 4 one can follow the corresponding proofs in [18]. Therefore, only the key details in the proofs for Theorems 2 and 4 are provided.

Proof of Theorem 2.

Using Itô formula (2.3.1) in [34] we obtain

[TABLE]

[TABLE]

As $h_{sp}(\cdot)$ is bounded and $\prod\limits_{j=1}^{\kappa}\sqrt{f(\|\lambda_{j}\|)}\in L_{2}(\mathbb{R}^{d\kappa})$ then the stochastic Fubini theorem, see Theorem 5.13.1 in [25], can be used to interchange the integrals. It results in

[TABLE]

where

[TABLE]

Using Lemma 2 instead of its counterpart in [18], the rest of the theorem can be proven analogously to Theorem 2 in [18]. ∎

Proof of Theorem 4.

Since $X_{r,G}=X_{r,\kappa}+V_{r}$, applying Chebyshev’s inequality and Lemma 1 to $X=X_{r,\kappa}$, $Y=V_{r}$ and $Z=X_{\kappa},$ we get

[TABLE]

It follows from Theorem 3 that

[TABLE]

By (7) it holds

[TABLE]

Since, by Remark 4, $L_{0}(\cdot)\sim L(\cdot)$, we can replace $L_{0}(\cdot)$ by $L(\cdot)$ in the above estimate. Thus, choosing $\beta_{1}=\frac{\alpha}{d-1-(\kappa-1)\alpha}$ to minimize the upper bound we get

[TABLE]

Choosing $\varepsilon:=r^{-\frac{\alpha(d-1-\kappa\alpha)}{(2+b)(d-1-(\kappa-1)\alpha)}}$ to minimize the second term in estimate (11) we obtain

[TABLE]

By (10)

[TABLE]

Therefore, ${\rho}\left(X_{r,\kappa},X_{\kappa}\right)$ is the Kolmogorov distance between two multiple Wiener-Itô integrals of the rank $\kappa$. Using Lemma 3 we get

[TABLE]

Using Lemma 2 and the bound (9) instead of Lemma 2 and (2) in [18] respectively, the rest of the proof is analogous to Theorem 5 in [18]. ∎

5 Simulation studies

In this section we present some simulation studies to confirm the obtained results. We also investigate random fields on cube surfaces and demonstrate that the rate of convergence for cube surfaces is similar to the rate of convergence for spheres. These simulations support the hypothesis that the results of this paper should hold for more general surfaces.

The theoretical results in Sections 3 and 4 were derived for integral functionals. However, computer simulations are possible only in discrete spaces. Thus, we simulate the functionals for sums instead of integrals. Recent results in [35] show that additive functionals can be used to provide reliable approximations of the integral functionals.

All simulations were done by using parallel computing on the NCI’s high-performance computer Raijin. The codes used for simulations and examples in this article are available in the folder ”Research materials” from https://sites.google.com/site/olenkoandriy/

We consider the Cauchy random field defined on the spheres $S(r),\,r>0$. The covariance function of this field is $\text{\rm{B}}(x)=(1+\|x\|^{2})^{-\frac{\alpha}{2}},$ $x\in S(r)$, where $\alpha$ is a long-range dependence parameter from Assumption 1. By [36, p.293] the spectral density of this field is

[TABLE]

where $K_{\frac{d-\alpha}{2}}(\cdot)$ is the modified Bessel function of the third kind.

By (10.30.2) in [37] $K_{\nu}(z)\sim\frac{1}{2}\Gamma(\nu)(\frac{1}{2}z)^{-\nu},\,z\rightarrow 0$. Thus, one can see that

[TABLE]

is a slowly varying function that satisfies Assumption 2. Hence, Theorem 4 holds for the Cauchy random fields.

Random fields were simulated by using the R package ’RandomFields’, see [38]. Figure 1 provides two examples of the simulated Cauchy random field that were obtained.

To verify the results of Sections 3 and 4 we simulated the random variables

[TABLE]

where $P(r):=\{p\in S(r)\}$ is a discrete set of uniformly distributed points in $S(r)$, $|S(\cdot)|$ is the surface area of the sphere, $|P(\cdot)|$ is the number of points in the set $P(\cdot)$, and $r_{min}$, $r_{max}$ denote the minimum and maximum values of the considered radius.

For the limit random variable $X_{\kappa}$ we considered the random variable $X_{R,G}$, where $R>>r_{max}$.

Let $d=3$, $\alpha=2/3$, $G(y)=H_{2}(y)=y^{2}-1$, $r_{min}=200$, and $r_{max}=3000$. Then the variables in (12) take the form

[TABLE]

Let $h_{sp}^{s}(\theta,\phi):=1.2+0.2\sin(5\theta)\sin(5\phi)$ and $h_{sp}^{c}(\theta,\phi):=2+\cos(3\theta)$, where $\theta$ and $\phi$ are spherical coordinates. These functions embed the sphere of radius 1 and smoothly vary in different directions. We present two examples where we used $h_{sp}^{s}(\cdot)$ as a weighted function in the case of spheres and $h_{sp}^{c}(\cdot)$ as a weighted function in the case of cubes. See Figure 2 for the visualization of these functions.

The sets $P(r),\,r>0$, were created by using the algorithm in [39] and setting the number of points $|P(r)|\sim r^{2}$. The random field $\eta(\cdot)$ was simulated at the points $P(r),r\in[200,3000],$ and at the points $P(R)$. Each random variable (13) was simulated 2000 times. Using the simulated values, the Kolmogorov distance between $X_{r,2}$ and $X_{R,2}$ was computed by the default Kolmogorov-Smirnov function ’ks.test’ in R. We present the results in two figures, obtained by repeating the simulation steps above 100 times. The first figure shows the box plots of the distances between $X_{r,2}$ and $X_{R,2}$. The second figure shows the corresponding box plots of the distance’s logarithms. The first figure displays actual convergence of the Kolmogorov distance for increasing $r$. The second figure displays the changes in the power index of the rate of convergence.

For $R=4000$ and $h_{sp}(\cdot)\equiv 1$ Figure 3 presents the obtained results for the unweighted case on spheres.

In the weighted case $h_{sp}(\cdot)=h_{sp}^{s}(\cdot)$ the obtained results are presented in Figure 4.

As one can see, Figures 3 and 4 are quite similar which indicates that un-weighted and weighted functionals have the same rate of convergence. This observation agrees with Theorem 4 according to which the weighted function has no effect on the rate of convergence. If the upper bound for the convergence rate in Theorem 4 was sharp, one would expect to see the means of Kolmogorov distance’s logarithms converge to a fixed value. Instead, we can see that the means form a declining slope, which might suggest that the real rate of convergence is a power function of the radius $r$. By fitting the linear regression model we obtained the equation

[TABLE]

in the unweighted case and a very similar equation in the weighted case.

Similar simulations were also conducted for cube surfaces. We considered 3-dimensional cube surfaces $SQ(r)$ with the centre at the origin and the side lengths $2r,\,r>0$. Thus, in (13) we replaced $S(r)$ by $SQ(r)$ and $P(r)$ by $PQ(r):=\{p\in SQ(r)\}$. Each $PQ(r)$ was set by creating the equidistant grid for each side of the cube and, therefore, $\left|PQ(r)\right|\sim r^{2}$. To approximate integrals (5) by sums (12) with the equal precision in both sphere and cube cases we created approximately $6/\pi$ times more points on cubes than on spheres of the same radius. Other steps remained unchanged.

In the unweighted case $h_{sp}(\cdot)\equiv 1$ the obtained results for surfaces of cube are presented in Figure 5.

In the weighted case $h_{sp}(\cdot)=h_{sp}^{c}(\cdot)$ we obtained the results presented in Figure 6.

As one can see, in the case of cubes the rate of convergence is also similar in both weighted and unweighted cases. Furthermore, these figures also suggest a power rate of convergence. By fitting the linear regression lines we obtained the equation $\log(\rho\left(X_{r,2},X_{R,2}\right))\approx-1.83-0.000571r$ in the unweighted case and the equation $\log(\rho\left(X_{r,2},X_{R,2}\right))\approx-1.81-0.000584r$ in the weighted case. The slopes of these equations are very similar to the slope in (14), indicating the same rate of convergence in the cases of both surfaces.

6 Directions for future research

This paper discussed the asymptotic behaviour of least squares estimators in regression models for long-range dependent random fields observed on spheres. The results were obtained under rather general assumptions on the random fields.

We expect that the results derived in Sections 3 and 4 also hold in the case of random fields observed on more general surfaces. This conjecture is supported by the results in Section 5, where asymptotic behaviour of LSE was studied for the cases of sphere and cube surfaces. The obtained results were quite similar in both cases.

The upper bound for the rate of convergence presented in the paper is the first result of such type for LSE. It would be interesting to obtain sharp upper bounds and extremal examples. The results of Section 5 suggest that the actual speed of convergence of LSE might be a power function of $r$.

The linear regression model was considered in this paper. It would be interesting to obtain analogous results for non-linear regression models.

Acknowledgements

This research was partially supported under the Australian Research Council’s Discovery Project DP160101366.

This research includes extensive simulation studies using the computational cluster Raijin of the National Computational Infrastructure (NCI), which is supported by the Australian Government and La Trobe University.