# Local polynomial regression for spatial data on $\mathbb{R}^{d}$

Daisuke Kurisu  and  Yasumasa Matsuda  Graduate School of International Social
Sciences, Yokohama National University
79-4, Tokiwadai, Hodogaya-ku, Yokohama 240-8501, Japan.  kurisu-daisuke-
<EMAIL_ADDRESS>Graduate School of Economics and Management, Tohoku University
Sendai 980-8576, Japan<EMAIL_ADDRESS>
(Date: First version: July 18, 2022. This version: )

###### Abstract.

This paper develops a general asymptotic theory of local polynomial (LP)
regression for spatial data observed at irregularly spaced locations in a
sampling region $R_{n}\subset\mathbb{R}^{d}$. We adopt a stochastic sampling
design that can generate irregularly spaced sampling sites in a flexible
manner including both pure increasing and mixed increasing domain frameworks.
We first introduce a nonparametric regression model for spatial data defined
on $\mathbb{R}^{d}$ and then establish the asymptotic normality of LP
estimators with general order $p\geq 1$. We also propose methods for
constructing confidence intervals and establish uniform convergence rates of
LP estimators. Our dependence structure conditions on the underlying processes
cover a wide class of random fields such as Lévy-driven continuous
autoregressive moving average random fields. As an application of our main
results, we discuss a two-sample testing problem for mean functions and their
partial derivatives.

###### Key words and phrases:

irregularly spaced spatial data, Lévy-driven moving average random field,
local polynomial regression, two-sample test
MSC2020 subject classifications: 62M30, 62G08, 62G20

D. Kurisu is partially supported by JSPS KAKENHI Grant Number 20K13468. The
authors would like to thank Takuya Ishihara, Taisuke Otsu, Peter Robinson,
Masayuki Sawada, and Yoshihiro Yajima for their helpful comments and
suggestions.

##  1\. Introduction

The goal of this paper is to develop a general asymptotic theory for local
polynomial (LP) estimators of any order $p\geq 1$ for spatial data under
irregular sampling on $\mathbb{R}^{d}$. We propose a nonparametric regression
model for spatial data $\\{Y(\bm{x}_{n,i})\\}_{i=1}^{n}$ observed at
irregularly spaced sampling sites $\\{\bm{x}_{n,i}\\}_{i=1}^{n}$ over a
sampling region $R_{n}\subset\mathbb{R}^{d}$ ($d\geq 1$). Precisely, each
$Y(\bm{x}_{n,i})$ is explained by the sum of a deterministic spatial trend
function (i.e. mean function), a random field on $\mathbb{R}^{d}$ that
represents spatial dependence, and a location specific measurement error (see
Section 2.1 for details). In many scientific fields, such as ecology, geology,
meteorology, and seismology, spatial samples are often collected over
irregularly spaced points from continuous random fields because of physical
constraints. To cope with irregularly spaced sampling sites, we adopt the
stochastic sampling scheme of Lahiri, 2003a , which allows the sampling sites
to have a non-uniform density in the sampling region and allows the number of
sampling sites $n$ to grow at a different rate from the volume of the sampling
region $A_{n}$. We design this scheme to accommodates both the pure increasing
domain case ($\lim_{n\to\infty}A_{n}/n=\kappa\in(0,\infty)$) and the mixed
increasing domain case ($\lim_{n\to\infty}A_{n}/n=0$). We note that this
scheme covers possible asymptotic regimes that would validate asymptotic
inference for spatial data. Although the infill asymptotics is excluded from
our regime, our sampling design is general enough as it is known that the
infill asymptotics does not work for that of even sample mean (cf.Lahiri,
(1996)). Refer to Lahiri, 2003b , Lahiri and Zhu, (2006), Matsuda and Yajima,
(2009), Bandyopadhyay et al., (2015), Kurisu et al., (2021), and Kurisu,
(2022) for discussions on the stochastic spatial sampling design. Further, our
model can be seen as a spatial extension of locally stationary time series
introduced in Dahlhaus, (1997).

The contributions of this paper are as follows. First, we (i) establish the
asymptotic normality of LP estimators of the mean function of the proposed
model, (ii) construct consistent estimators of their asymptotic variances, and
(iii) derive uniform convergence rates of LP estimators over a compact set.
The results (i) and (ii) enable us to evaluate the bias and
variance/covariance matrix (of the asymptotic distribution) of LP estimators
and, as a result, to construct confidence intervals of LP estimators, which
would work for a hypothesis testing on the mean function. We discuss a two-
sample test for the partial derivatives as well as the mean function as an
application of our results. Additionally, in the literature of causal
inference, local polynomial fitting is know as important tools to analyze
average treatment effect of interventions, an example of which is the
regression discontinuity desings (RDDs) (cf. Hahn et al., (2001) and Calonico
et al., (2014)). Existing methods for RDDs often assume i.i.d. even for
spatial data (cf. Keele and Titiunik, (2015) and Ehrlich and Seidel, (2018)).
We claim our results pave the way for a new framework of RDDs for spatially
dependent data. To establish the result (iii), we first consider general
kernel estimators and derive their uniform convergence rates. The uniform
convergence rates of LP estimators can be given as special cases of the
results. Since the general estimators include many kernel-based estimators
such as, kernel density, local constant (LC), local linear (LL), and LP
estimators for random fields on $\mathbb{R}^{d}$ with irregularly spaced
sampling sites, the results are of independent theoretical interest. We note
that the general results are also useful for evaluating both the bias and
variance terms of LP estimators. Particularly, the results on uniform
convergence rates enable us to predict the values of the mean function
uniformly on a spatial region that does not contain sampling sites.

Second, we provide examples of random fields that satisfy the mixing
assumptions under which the asymptotic normality of LP estimators will be
established. Specifically, we show that a broad class of Lévy-driven moving
average (MA) random fields, which include continuous autoregressive moving
average (CARMA) random fields (cf. Brockwell and Matsuda, (2017)), satisfies
our assumptions. The CARMA random fields are known as a rich class of models
for spatial data that can represent non-Gaussian random fields by introducing
non-Gaussian Lévy random measures (cf. Brockwell and Matsuda, (2017), Matsuda
and Yajima, (2018), and Kurisu, (2022)). However, mixing properties of Lévy-
driven MA random fields have not been investigated since it is often difficult
to check mixing conditions in the ways considered by Lahiri and Zhu, (2006)
and Bandyopadhyay et al., (2015) for general (possibly non-Gaussian) random
fields on $\mathbb{R}^{d}$, which will be discussed later from the viewpoint
of our theoretical analysis. We show that a wide class of Lévy-driven MA
random fields can be approximated by $m_{n}$-dependent random fields with
$m_{n}\to\infty$ as $n\to\infty$. We claim that the approximation will work
for the flexible modeling of nonparametric, nonstationary and possibly non-
Gaussian spatial data on $\mathbb{R}^{d}$ by addressing an open question on
dependence structure of statistical models built on Lévy-driven MA random
fields.

### Connections to the literature

There is fairly extensive literature on LC, LL, and LP estimators for
dependent data. For stationary and regularly spaced time series (this case
corresponds to stationary random fields with regular sampling on
$\mathbb{Z}$), we refer to Hansen, (2008) and Zhao and Wu, (2008) for LC
estimators and Masry, 1996a ; Masry, 1996b , and Masry and Fan, (1997) for LP
estimators. For nonstationary and regularly spaced time series, we refer to
Kristensen, (2009) and Vogt, (2012) for LC estimators, and Zhou and Wu, (2009)
and Zhang and Wu, (2015) for LL estimators of quantile curves and conditional
mean functions, respectively. For stationary spatial data with regular
sampling on $\mathbb{Z}^{d}$, refer to El Machkouri and Stoica, (2010) for LC
estimators and Hallin et al., (2004) for LL estimators. For stationary spatial
data with irregular sampling on $\mathbb{Z}^{d}$, we refer to El Machkouri et
al., (2017) for LL estimators. For nonstationary spatial data with (possibly)
irregular sampling on $\mathbb{Z}^{d}$, we refer to Robinson, (2011) for LC
estimators and Jenish, (2012) for LL estimators. For spatial data with
irregular sampling on $\mathbb{R}^{d}$, we refer to Kurisu, (2019) and Kurisu,
(2022) who investigate LC estimators for the stationary and nonstationary
case, respectively. Existing results on LP estimators are available only for
stationary random fields under regular sampling on $\mathbb{Z}$, i.e.,
regularly spaced stationary time series, while no studies on LL and LP
estimators have been known under irregular sampling on $\mathbb{R}^{d}$ with
$d\geq 1$.

To the best of our knowledge, our work is the first attempt to establish an
asymptotic theory on local polynomial fitting for a mean function of spatial
data on $\mathbb{R}^{d}$ by (i) establishing the asymptotic normality and
uniform convergence rates of LP estimators, (ii) providing a way to construct
confidence intervals of LP estimators, and (iii) showing the applicability of
our theoretical results to a wide class of Lévy-driven MA random fields. From
a theoretical point of view, this paper has advantages over the existing
studies of Lahiri, 2003a  and Lahiri and Zhu, (2006) in the fields of
irregularly spaced data analysis. Specifically, (i) we extend the coupling
technique used in Yu, (1994) for time series to that for irregularly spatial
data to establish uniform convergence rates of LP estimators. The difficulties
in the extension come from no natural ordering for spatial data and the number
of observations in each block constructed is random, and hence our approach to
blocking construction for establishing uniform rates is quite different from
those in Lahiri, 2003a  and Lahiri and Zhu, (2006) whose proofs essentially
rely on approximating the characteristic function of the weighted sample mean
by that of independent blocks. (ii) We have confirmed concrete examples of
random fields that satisfy our assumptions in detail. Verification of our
regularity conditions to Lévy-driven MA fields is indeed non-trivial and
relies on several probabilistic techniques from Lévy process theory and theory
of infinitely divisible random measures (cf. Bertoin, (1996), Sato, (1999),
and Rajput and Rosinski, (1989)).

The rest of the paper is organized as follows. In Section 2, we introduce our
nonparametric regression model for spatial data with irregularly spaced
sampling sites. In Section 3, we define local polynomial estimators as
solutions of a multivariate weighted least squares problem. In Section 4, we
establish the asymptotic normality of LP estimators, construct estimators of
their asymptotic variances, and discuss a two-sample test for the mean
functions and their partial derivatives. In Section 5, we provide the uniform
convergence rates of LP estimators. In Section 6, we provide examples of the
random fields that satisfies our assumptions. All proofs are included in
Appendix.

###  1.1. Notation

For any vector $\bm{x}=(x_{1},\dots,x_{q})^{\prime}\in\mathbb{R}^{q}$, let
$|\bm{x}|=\sum_{j=1}^{q}|x_{j}|$ and
$\|\bm{x}\|=\sqrt{\sum_{j=1}^{q}x_{j}^{2}}$ denote the $\ell^{1}$-norm and
$\ell^{2}$-norms of $\bm{x}$, respectively. For any set
$A\subset\mathbb{R}^{d}$ and any vector
$\bm{a}=(a_{1},\dots,a_{d})^{\prime}\in(0,\infty)^{d}$, let $|A|$ denote the
Lebesgue measure of $A$, let $[\\![A]\\!]$ denote the number of elements in
$A$, and let $\bm{a}A=\\{(a_{1}x_{1},\dots
a_{d}x_{d}):\bm{x}=(x_{1},\dots,x_{d})\in A\\}$. For any positive sequences
$a_{n},b_{n}$, we write $a_{n}\lesssim b_{n}$ if there is a constant $C>0$
independent of $n$ such that $a_{n}\leq Cb_{n}$ for all $n$, $a_{n}\sim b_{n}$
if $a_{n}\lesssim b_{n}$ and $b_{n}\lesssim a_{n}$. For a sequence of random
variables $\\{\bm{X}_{i}\\}_{i\geq 1}$, let $\sigma(\\{\bm{X}_{i}\\}_{i\geq
1})$ denote the $\sigma$-field generated by $\\{\bm{X}_{i}\\}_{i\geq 1}$. Let
$E_{\bm{X}}$ denote the expectation with respect to a sequence of random
variables $\\{\bm{X}_{i}\\}_{i\geq 1}$ and let $P_{\cdot\mid\bm{X}}$ and
$E_{\cdot\mid\bm{X}}$ denote the conditional probability and expectation given
$\sigma(\\{\bm{X}_{i}\\}_{i\geq 1})$, respectively. For any real-valued random
variable $X$ and $\tau\in(0,1)$, let $q_{1-\tau}=\inf\\{x\in\mathbb{R}:P(X\leq
x)\geq 1-\tau\\}$ be the $(1-\tau)$-quantile of $X$. For $a\in\mathbb{R}$ and
$b>0$, we use the shorthand notation $[a\pm b]=[a-b,a+b]$.

##  2\. Settings

In this section, we discuss the mathematical settings of our model (Section
2.1), sampling design (Section 2.2), and spatial dependence structure (Section
2.3).

###  2.1. Model

Usually it is impossible to estimate consistently a model for nonstationary
processes, since the domain of functions to be estimated gets larger. Dahlhaus
avoids the difficulty by desigining a function over a fixed interval in the
following way. Dahlhaus, (1997) introduced a locally stationary process with a
time-varying mean function for the modeling of nonstationary time series:

$Y_{T}(t)=m\left({t\over T}\right)+\xi_{T}(t),\ t=1,\dots,T,$

where $m:[0,1]\to\mathbb{R}$ is a (time-varying) mean function and
$\\{\xi_{T}(t)\\}$ is a sequence of zero-mean locally stationary time series
with a time-varying transfer function (see Definition 2.1 in Dahlhaus, (1997)
for details). The model setting of $m(t/T)$ instead of $m(t)$ makes the mean
function have the fixed domain of $[0,1]$, which provides the asymptotic
scheme on which consistent estimation is available. We extend his framework to
spatial data with irregular sampling on $\mathbb{R}^{d}$.

In particular, consider the following nonparametric regression model:

$\displaystyle Y(\bm{x}_{n,i})$ $\displaystyle=m\left({\bm{x}_{n,i}\over
A_{n}}\right)+\eta\left({\bm{x}_{n,i}\over
A_{n}}\right)e(\bm{x}_{n,i})+\sigma_{\varepsilon}\left({\bm{x}_{n,i}\over
A_{n}}\right)\varepsilon_{i},$ (2.1) $\displaystyle:=m\left({\bm{x}_{n,i}\over
A_{n}}\right)+e_{n,i}+\varepsilon_{n,i},\
\bm{x}_{n,i}=(x_{ni,1},\dots,x_{ni,d})^{\prime}\in R_{n},\ i=1,\dots,n,$

where $R_{n}=\prod_{j=1}^{d}[-A_{n,j}/2,A_{n,j}/2]^{d}$,
$A_{n}=\prod_{j=1}^{d}A_{n,j}$, ${\bm{x}_{n,i}\over
A_{n}}=\left({x_{ni,1}\over A_{n,1}},\dots,{x_{ni,d}\over
A_{n,d}}\right)^{\prime}$ with $A_{n,j}\to\infty$ as $n\to\infty$,
$m:[-1/2,1/2]^{d}\to\mathbb{R}$ is the mean function,
$\bm{e}=\\{e(\bm{x}):\bm{x}\in\mathbb{R}^{d}\\}$ is a stationary random field
defined on $\mathbb{R}^{d}$ with $E[e(\bm{x})]=0$ and $E[e^{2}(\bm{x})]=1$ for
any $\bm{x}\in\mathbb{R}^{d}$, $\eta:[-1/2,1/2]^{d}\to(0,\infty)$ is the
variance function of spatially dependent random variables $\\{e_{n,i}\\}$,
$\\{\varepsilon_{i}\\}$ is a sequence of i.i.d. random variables such that
$E[\varepsilon_{i}]=0$ and $E[\varepsilon_{i}^{2}]=1$, and
$\sigma_{\varepsilon}:[-1/2,1/2]^{d}\to(0,\infty)$ is the variance function of
random variables $\\{\varepsilon_{n,i}\\}$. The mean function $m$ represents
deterministic spatial trend, the random field $\bm{e}$ represents spatial
correlation, and the random variables $\\{\varepsilon_{n,i}\\}$ can represent
location specific measurement error.

We assume the following conditions on the mean function $m$, the variance
function $\eta$, and $\\{\varepsilon_{n,j}\\}$:

######  Assumption 2.1.

Let $U_{\bm{z}}$ be a neighborhood of
$\bm{z}=(z_{1},\dots,z_{d})\in(-1/2,1/2)^{d}$.

  * (i)

The mean function $m$ is $(p+1)$-times continuously partial differentiable on
$U_{\bm{z}}$ and define $\partial_{j_{1}\dots j_{L}}m(\bm{z}):=\partial
m(\bm{z})/\partial z_{j_{1}}\dots\partial z_{j_{L}}$, $1\leq
j_{1},\dots,j_{L}\leq d$, $0\leq L\leq p+1$. When $L=0$, we set
$\partial_{j_{1}\dots j_{L}}m(\bm{z})=\partial_{j_{0}}m(\bm{z})=m(\bm{z})$.

  * (ii)

The function $\eta$ is continuous over $U_{\bm{z}}$ and $\eta(\bm{z})>0$.

  * (iii)

The random variables $\\{\varepsilon_{i}\\}_{i=1}^{n}$ are i.i.d. with
$E[\varepsilon_{1}]=0$, $E[\varepsilon_{1}^{2}]=1$,
$E[|\varepsilon_{1}|^{q_{1}}]<\infty$ for some integer $q_{1}>4$, and the
function $\sigma_{\varepsilon}(\cdot)$ is continuous over $U_{\bm{z}}$ with
$\sigma_{\varepsilon}(\bm{z})>0$.

###  2.2. Sampling design

To account for irregularly spaced data, we consider a stochastic sampling
design. First, we define the sampling region $R_{n}$. For $j=1,\dots,d$, let
$\\{A_{n,j}\\}_{n\geq 1}$ be a sequence of positive numbers such that
$A_{n,j}\to\infty$ as $n\to\infty$. We consider the following set as the
sampling region.

$\displaystyle R_{n}=\prod_{j=1}^{d}[-A_{n,j}/2,A_{n,j}/2].$ (2.2)

Next, we introduce our (stochastic) sampling designs. Let
$g(\bm{z})=g(z_{1},\dots,z_{d})$ be a probability density function on
$R_{0}=[-1/2,1/2]^{d}$, and let $\\{\bm{X}_{n,i}\\}_{i\geq 1}$ be a sequence
of i.i.d. random vectors with probability density
$A_{n}^{-d}g(\bm{x}/A_{n})=A_{n}^{-d}g(x_{1}/A_{n,1},\dots,x_{d}/A_{n,d})$
where $A_{n}=\prod_{j=1}^{d}A_{n,j}$. We assume that the sampling sites
$\bm{x}_{n,1},\ldots,\bm{x}_{n,n}$ are obtained from the realizations of
random vectors $\bm{X}_{n,1},\ldots,\bm{X}_{n,n}$. To simplify the notation,
we will write $\bm{x}_{n,i}$ and $\bm{X}_{n,i}$ as
$\bm{x}_{i}=(x_{i,1},\dots,x_{i,d})^{\prime}$ and
$\bm{X}_{i}=(X_{i,1},\dots,X_{i,d})^{\prime}$, respectively.

We summarize conditions on the stochastic sampling design as follows:

######  Assumption 2.2.

Recall that $U_{\bm{z}}$ is a neighborhood of $\bm{z}\in(-1/2,1/2)^{d}$. Let
$g$ be a probability density function with support $R_{0}=[-1/2,1/2]^{d}$.

  * (i)

$A_{n}/n\to\kappa\in[0,\infty)$ as $n\to\infty$,

  * (ii)

$\\{\bm{X}_{i}=(X_{i,1},\dots,X_{i,d})^{\prime}\\}_{i=1}^{n}$ is a sequence of
i.i.d. random vectors with density $A_{n}^{-d}g(\cdot/A_{n})$ and $g$ is
continuous over $U_{\bm{z}}$ and $g(\bm{z})>0$.

  * (iii)

$\\{\bm{X}_{i}\\}_{i=1}^{n}$,
$\bm{e}=\\{e(\bm{x}):\bm{x}\in\mathbb{R}^{d}\\}$, and
$\\{\varepsilon_{i}\\}_{i=1}^{n}$ are mutually independent.

Condition (i) implies that our sampling design allows both the pure increasing
domain case ($\lim_{n\to\infty}A_{n}/n=\kappa\in(0,\infty)$) and the mixed
increasing domain case ($\lim_{n\to\infty}A_{n}/n=0$). Condition (ii) implies
that the sampling density can be nonuniformly distributed over the sampling
region $R_{n}=\prod_{j=1}^{d}[-A_{n,j}/2,A_{n,j}/2]$. The definition (2.2) is
only for convenience, since it is possible to consider sampling regions of
various shapes including non-standard shapes (e.g., ellipsoids, polyhedrons,
and non-convex sets) by adjusting the support of the density $g$.

###  2.3. Dependence structure

We assume that random field $\bm{e}$ satisfies a mixing condition. First, we
define the $\alpha$\- and $\beta$-mixing coefficients for the random field
$\bm{e}$. Let $\mathcal{F}_{\bm{e}}(T)=\sigma(\\{e(\bm{x}):\bm{x}\in T\\})$ be
the $\sigma$-field generated by the variables $\\{e(\bm{x}):\bm{x}\in T\\}$,
$T\subset\mathbb{R}^{d}$. For any two subsets $T_{1}$ and $T_{2}$ of
$\mathbb{R}^{d}$, let

$\displaystyle\bar{\alpha}(T_{1},T_{2})$ $\displaystyle=\sup\\{|P(A\cap
B)-P(A)P(B)|:A\in\mathcal{F}_{\bm{e}}(T_{1}),B\in\mathcal{F}_{\bm{e}}(T_{2})\\},$
$\displaystyle\bar{\beta}(T_{1},T_{2})$ $\displaystyle=\sup{1\over
2}\sum_{j=1}^{J}\sum_{k=1}^{K}|P(A_{j}\cap B_{k})-P(A_{j})P(B_{k})|,$

where the supremum for $\bar{\beta}(T_{1},T_{2})$ is taken over all pairs of
(finite) partitions $\\{A_{1},\ldots,A_{J}\\}$ and $\\{B_{1},\ldots,B_{K}\\}$
of $\mathbb{R}^{d}$ such that $A_{j}\in\mathcal{F}_{\bm{e}}(T_{1})$ and
$B_{k}\in\mathcal{F}_{\bm{e}}(T_{2})$. The $\alpha$\- and $\beta$-mixing
coefficients of the random field $\bm{e}$ are defined as

$\displaystyle\alpha(a;b)$
$\displaystyle=\sup\\{\bar{\alpha}(T_{1},T_{2}):d(T_{1},T_{2})\geq
a,T_{1},T_{2}\in\mathcal{R}(b)\\},$ $\displaystyle\beta(a;b)$
$\displaystyle=\sup\\{\bar{\beta}(T_{1},T_{2}):d(T_{1},T_{2})\geq
a,T_{1},T_{2}\in\mathcal{R}(b)\\}.$

where $a,b>0$, $d(T_{1},T_{2})=\inf\\{|\bm{x}-\bm{y}|:\bm{x}\in
T_{1},\bm{y}\in T_{2}\\}$, and $\mathcal{R}(b)$ is the collection of all the
finite disjoint unions of cubes in $\mathbb{R}^{d}$ with a total volume not
exceeding $b$. Moreover, we assume that there exist a non-increasing functions
$\alpha_{1}$ and $\beta_{1}$ with $\alpha_{1}(a),\beta_{1}(a)\to 0$ as
$a\to\infty$ and a non-decreasing functions $\varpi_{1}$ and $\varpi_{2}$
(that may be unbounded) such that

$\displaystyle\alpha(a;b)\leq\alpha_{1}(a)\varpi_{1}(b),\
\beta(a;b)\leq\beta_{1}(a)\varpi_{2}(b).$

######  Remark 2.1.

The definitions of the $\alpha$\- and $\beta$-mixing coefficients are based on
the argument in Bradley, (1989). It is crucial to restrict the size of the
index sets $T_{1}$ and $T_{2}$ in the definition of $\alpha$\- (or $\beta$-)
mixing coefficients since no restrictions on $T_{1}$ and $T_{2}$ make the
$\alpha$\- and $\beta$-mixing be equivalent to $m$ dependent for a fixed
$m>0$, which would not work for our asymptotic inference. Let us define the
$\beta$-mixing coefficient of a random field $\bm{e}$ similarly to the time
series as follows: For any subsets $T_{1}$ and $T_{2}$ of $\mathbb{R}^{d}$,
the $\beta$-mixing coefficient between $\mathcal{F}_{\bm{e}}(T_{1})$ and
$\mathcal{F}_{\bm{e}}(T_{2})$ is defined by
$\widetilde{\beta}(T_{1},T_{2})=\sup{\textstyle\sum_{j=1}^{J}\sum_{k=1}^{K}}|P(A_{j}\cap
B_{k})-P(A_{j})P(B_{k})|/2$, where the supremum is taken over all partitions
$\\{A_{j}\\}_{j=1}^{J}\subset\mathcal{F}_{\bm{e}}(T_{1})$ and
$\\{B_{k}\\}_{k=1}^{K}\subset\mathcal{F}_{\bm{e}}(T_{2})$ of $\mathbb{R}^{d}$.
Let $\mathcal{O}_{1}$ and $\mathcal{O}_{2}$ be half-planes with boundaries
$L_{1}$ and $L_{2}$, respectively. For each $a>0$, define
$\beta(a)=\sup\\{\widetilde{\beta}(\mathcal{O}_{1},\mathcal{O}_{2}):d(\mathcal{O}_{1},\mathcal{O}_{2})\geq
a\\}$. According to Theorem 1 in Bradley, (1989), if
$\\{e(\bm{x}):\bm{x}\in\mathbb{R}^{2}\\}$ is strictly stationary, then
$\beta(a)=0$ or $1$ for $a>0$. This implies that if a random field $\bm{e}$ is
$\beta$-mixing ($\lim_{a\to\infty}\beta(a)=0$), then it is automatically $m$
dependent, that is, $\beta(a)=0$ for some $a>m$, where $m$ is a positive
constant. To allow a certain flexibility, we restrict the size of $T_{1}$ and
$T_{2}$ in the definitions of $\alpha(a;b)$ and $\beta(a;b)$. We refer to
Bradley, (1993) and Doukhan, (1994) for more details on mixing coefficients
for random fields.

For the asymptotic normality of LP estimators, we assume the following
conditions for the random field $\bm{e}$:

######  Assumption 2.3.

For $j=1,\dots,d$, let $\\{A_{n1,j}\\}_{n\geq 1}$ and $\\{A_{n2,j}\\}_{n\geq
1}$ be sequences of positive numbers such that
$\min\left\\{A_{n2,j},{A_{n1,j}\over A_{n2,j}}\right\\}\to\infty$ as
$n\to\infty$.

  * (i)

The random field $\bm{e}$ is stationary and
$E[|\bm{e}(\bm{0})|^{q_{2}}]<\infty$ for some integer $q_{2}>4$.

  * (ii)

Define $\sigma_{\bm{e}}(\bm{x})=E[e(\bm{0})e(\bm{x})]$. Assume that
$\sigma_{\bm{e}}(\bm{0})=1$ and
$\int_{\mathbb{R}^{d}}|\sigma_{\bm{e}}(\bm{v})|d\bm{v}<\infty$.

  * (iii)

The random field $\bm{e}$ is $\alpha$-mixing with mixing coefficients
$\alpha(a;b)$ such that as $n\to\infty$,

$\displaystyle
A_{n}^{(1)}\left(\alpha_{1}^{1-2/q}(\underline{A}_{n2})+\sum_{k=\underline{A}_{n1}}^{\infty}k^{d-1}\alpha_{1}^{1-2/q}(k)\right)\varpi_{1}^{1-2/q}(A_{n}^{(1)})\to
0,$

where $q=\min\\{q_{1},q_{2}\\}$, $A_{n}^{(1)}=\prod_{j=1}^{d}A_{n1,j}$,
$\underline{A}_{n1}=\min_{1\leq j\leq d}A_{n1,j}$, and
$\underline{A}_{n2}=\min_{1\leq j\leq d}A_{n2,j}$.

The sequences $\\{A_{n1,j}\\}$ and $\\{A_{n2,j}\\}$ will be used in the large-
block-small-block argument, which is commonly used in proving CLTs for sums of
mixing random variables. Specifically, $A_{n1,j}$ corresponds to the side
length of large blocks, while $A_{n2,j}$ corresponds to the side length of
small blocks. In Section 6, we provide examples of random fields that satisfy
Assumptions 2.3 and 4.1 below. In particular, a wide class of Lévy-driven
moving average (MA) random fields that includes continuous autoregressive and
moving average (CARMA) random fields (cf. Brockwell and Matsuda, (2017))
satisfies our assumptions.

##  3\. Local polynomial regression of order $p$

In this section, we introduce local polynomial (LP) estimators of order $p\geq
1$ for the estimation of the mean function $m$ of the model (2.1) and their
partial derivatives.

Define

$\displaystyle D$ $\displaystyle=[\\![\\{(j_{1},\dots,j_{L}):1\leq
j_{1}\leq\dots\leq j_{L}\leq d,0\leq L\leq p\\}]\\!],$ $\displaystyle\bar{D}$
$\displaystyle=[\\![\\{(j_{1},\dots,j_{p+1}):1\leq j_{1}\leq\dots\leq
j_{p+1}\leq d\\}]\\!],$

$(s_{j_{1}\dots j_{L}1},\dots,s_{j_{1}\dots j_{L}d})\in\mathbb{Z}_{\geq
0}^{d}$ such that $s_{j_{1}\dots
j_{L}k}=[\\![\\{j_{\ell}:j_{\ell}=k,1\leq\ell\leq L\\}]\\!]$, and define

$\bm{s}_{j_{1}\dots j_{L}}!=\prod_{k=1}^{d}s_{j_{1}\dots j_{L}k}!.$

When $L=0$, we set $(j_{1},\dots,j_{L})=j_{0}=0$ and $\bm{s}_{j_{1}\dots
j_{L}}!=1$. Note that $\sum_{k=1}^{d}s_{j_{1}\dots j_{L}k}=L$. Further, for
$p\geq 1$ and $\bm{z}\in[-1/2,1/2]^{d}$, define

$\displaystyle\bm{M}(\bm{z})$
$\displaystyle:=\left(m(\bm{z}),\partial_{1}m(\bm{z}),\dots,\partial_{d}m(\bm{z}),{\partial_{11}m(\bm{z})\over
2!},{\partial_{12}m(\bm{z})\over 1!1!},\dots,{\partial_{dd}m(\bm{z})\over
2!},\right.$ $\displaystyle\left.\quad\quad\dots,{\partial_{1\dots
1}m(\bm{z})\over p!},{\partial_{1\dots
2}m(\bm{z})\over(p-1)!1!}\dots,{\partial_{d\dots d}m(\bm{z})\over
p!}\right)^{\prime}$ $\displaystyle=\left({1\over\bm{s}_{j_{1}\dots
j_{L}}!}\partial_{j_{1},\dots j_{L}}m(\bm{z})\right)^{\prime}_{1\leq
j_{1}\leq\dots\leq j_{L}\leq d,0\leq L\leq p}\in\mathbb{R}^{D}.$

We define the local polynomial regression estimator of order $p$ for
$M(\bm{z})$ as a solution of the following problem:

$\displaystyle\widehat{\bm{\beta}}(\bm{z})$ $\displaystyle:=\mathop{\rm
arg~{}min}\limits_{\bm{\beta}\in\mathbb{R}^{D}}\sum_{i=1}^{n}\left(Y(\bm{X}_{i})-\sum_{L=0}^{p}\sum_{1\leq
j_{1}\leq\dots\leq j_{L}\leq d}\beta_{j_{1}\dots
j_{L}}\prod_{\ell=1}^{L}\left({X_{i,j_{\ell}}-A_{n,j_{\ell}}z_{j_{\ell}}\over
A_{n,j_{\ell}}}\right)\right)^{2}K_{Ah}\left(\bm{X}_{i}-A_{n}\bm{z}\right)$
(3.1)
$\displaystyle=(\widehat{\beta}_{0}(\bm{z}),\widehat{\beta}_{1}(\bm{z}),\dots,\widehat{\beta}_{d}(\bm{z}),\widehat{\beta}_{11}(\bm{z}),\dots,\widehat{\beta}_{dd}(\bm{z}),\dots,\widehat{\beta}_{1\dots
1}(\bm{z}),\dots,\widehat{\beta}_{d\dots d}(\bm{z}))^{\prime}$
$\displaystyle=(\widehat{\beta}_{j_{1}\dots j_{L}}(\bm{z}))^{\prime}_{1\leq
j_{1}\leq\dots\leq j_{L}\leq d,0\leq L\leq p},$

where $\bm{\beta}=(\beta_{j_{1}\dots j_{L}})^{\prime}_{1\leq
j_{1}\leq\dots\leq j_{L}\leq d,0\leq L\leq p}$,
$K:\mathbb{R}^{d}\to\mathbb{R}$ is a kernel function, and each $h_{j}$ is a
sequence of positive constants (bandwidths) such that $h_{j}\to 0$ as
$n\to\infty$, and where

$K_{Ah}(\bm{X}_{i}-A_{n}\bm{z})=K\left({X_{i,1}-A_{n,1}z_{1}\over
A_{n,1}h_{1}},\dots,{X_{i,d}-A_{n,d}z_{d}\over A_{n,d}h_{d}}\right)$

and $\sum_{1\leq j_{1}\leq\dots\leq j_{L}\leq d}\beta_{j_{1}\dots
j_{L}}\prod_{\ell=1}^{L}(X_{i,j_{\ell}}-A_{n,j_{\ell}}z_{j_{\ell}})/A_{n,j_{\ell}}=\beta_{0}$
when $L=0$.

To compute LP estimators, we introduce some notations:
$\bm{Y}:=(Y(\bm{X}_{1}),\dots,Y(\bm{X}_{n}))^{\prime}$,

$\displaystyle\bm{X}$
$\displaystyle:=(\widetilde{\bm{X}}_{1},\dots,\widetilde{\bm{X}}_{n})=\left(\begin{array}[]{ccc}1&\dots&1\\\
{\left(\bm{X}_{1}-A_{n}\bm{z}\right)_{1}\over
A_{n}}&\dots&{\left(\bm{X}_{n}-A_{n}\bm{z}\right)_{1}\over A_{n}}\\\
\vdots&\dots&\vdots\\\ {\left(\bm{X}_{1}-A_{n}\bm{z}\right)_{p}\over
A_{n}}&\dots&{\left(\bm{X}_{n}-A_{n}\bm{z}\right)_{p}\over
A_{n}}\end{array}\right)=\left(\begin{array}[]{ccc}1&\dots&1\\\
\check{(\bm{X}_{1}-A_{n}\bm{z})}&\dots&\check{(\bm{X}_{n}-A_{n}\bm{Z})}\end{array}\right),$
$\displaystyle\bm{W}$ $\displaystyle:={\rm
diag}\left(K_{Ah}\left(\bm{X}_{1}-A_{n}\bm{z}\right),\dots,K_{Ah}\left(\bm{X}_{n}-A_{n}\bm{z}\right)\right),$

where

$\displaystyle{\left(\bm{X}_{i}-A_{n}\bm{z}\right)_{L}\over A_{n}}$
$\displaystyle=\left(\prod_{\ell=1}^{L}\left({X_{i,j_{\ell}}-A_{n,j_{\ell}}z_{j_{\ell}}\over
A_{n,j_{\ell}}}\right)\right)^{\prime}_{1\leq j_{1}\leq\dots\leq j_{L}\leq
d}.$

The minimization problem (3.1) can be rewritten as

$\displaystyle\widehat{\bm{\beta}}(\bm{z})$ $\displaystyle=\mathop{\rm
arg~{}min}\limits_{\bm{\beta}\in\mathbb{R}^{D}}(\bm{Y}-\bm{X}^{\prime}\bm{\beta})^{\prime}\bm{W}(\bm{Y}-\bm{X}^{\prime}\bm{\beta})=:\mathop{\rm
arg~{}min}\limits_{\bm{\beta}\in\mathbb{R}^{D}}Q_{n}(\bm{\beta}).$

Then the first order condition of the problem (3.1) is given by

$\displaystyle{\partial\over\partial\bm{\beta}}Q_{n}(\bm{\beta})$
$\displaystyle=-2\bm{X}\bm{W}\bm{Y}+2\bm{X}\bm{W}\bm{X}^{\prime}\bm{\beta}=0.$

Hence the solution of the problem (3.1) is given by

$\displaystyle\widehat{\bm{\beta}}(\bm{z})$
$\displaystyle=(\bm{X}\bm{W}\bm{X}^{\prime})^{-1}\bm{X}\bm{W}\bm{Y}$
$\displaystyle=\left[\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}-A_{n}\bm{z}\right)\widetilde{\bm{X}}_{i}\widetilde{\bm{X}}^{\prime}_{i}\right]^{-1}\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}-A_{n}\bm{z}\right)\widetilde{\bm{X}}_{i}Y(\bm{X}_{i}).$

We assume the following conditions on the kernel function $K$:

######  Assumption 3.1.

Let $K:\mathbb{R}^{d}\to\mathbb{R}$ be a kernel function such that

  * (i)

$\int K(\bm{z})d\bm{z}=1$.

  * (ii)

The kernel function $K$ is bounded and supported on
$S_{K}\subset[-1/2,1/2]^{d}$ with $U_{\bm{z}}\subset S_{K}$.

  * (iii)

Define $\kappa_{0}^{(r)}:=\int K^{r}(\bm{z})d\bm{z}$, $\kappa_{j_{1}\dots
j_{M}}^{(r)}:=\int\prod_{\ell=1}^{M}z_{j_{\ell}}K^{r}(\bm{z})d\bm{z}$, and

$\check{\bm{z}}:=(1,(\bm{z})^{\prime}_{1},\dots,(\bm{z})^{\prime}_{p})^{\prime},\
(\bm{z})_{L}=\left(\prod_{\ell=1}^{L}z_{j_{\ell}}\right)^{\prime}_{1\leq
j_{1}\leq\dots\leq j_{L}\leq d},\ 1\leq L\leq p.$

The matrix $S=\int\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)(1\ \check{\bm{z}}^{\prime})K(\bm{z})d\bm{z}$
is non-singular.

##  4\. Main results

In this section, we discuss asymptotic properties of LP estimators defined in
Section 3. In particular, we establish the asymptotic normality of LP
estimator (Section 4.1) and construct an estimator of the asymptotic variance
of LP estimators (Section 4.2). As an application of our main results, we
discuss a two-sample test for the mean functions and their partial derivatives
(Section 4.3).

###  4.1. Asymptotic normality of local polynomial estimators

We assume the following conditions for the sample size $n$, bandwidths
$h_{j}$, constants $A_{n,j}$, $A_{n1,j}$, and $A_{n2,j}$, and mixing
coefficients $\alpha(a;b)$:

######  Assumption 4.1.

Recall $q=\min\\{q_{1},q_{2}\\}$, $A_{n}^{(1)}=\prod_{j=1}^{d}A_{n1,j}$, and
$\underline{A}_{n1}=\min_{1\leq j\leq d}A_{n1,j}$. Define
$\overline{A}_{n1}=\max_{1\leq j\leq d}A_{n1,j}$,
$\overline{A}_{n2}=\max_{1\leq j\leq d}A_{n2,j}$, and
$\overline{A_{n}h}=\max_{1\leq j\leq d}A_{n,j}h_{j}$. As $n\to\infty$,

  * (i)

$h_{j}\to 0$, ${A_{n,j}h_{j}\over A_{n1,j}}\to\infty$ for $1\leq j\leq d$.

  * (ii)

$nh_{1}\dots h_{d}\to\infty$.

  * (iii)

$A_{n}h_{1}\dots h_{d}\times h_{j_{1}}^{2}\dots h_{j_{p}}^{2}\to\infty$ for
$1\leq j_{1}\leq\dots\leq j_{p}\leq d$.

  * (iv)

$A_{n}h_{1}\dots h_{d}\times h_{j_{1}}^{2}\dots
h_{j_{p}}^{2}h_{j_{p+1}}^{2}\to c_{j_{1}\dots j_{p+1}}\in[0,\infty)$ for
$1\leq j_{1}\leq\dots\leq j_{p+1}\leq d$.

  * (v)

$\displaystyle\left({A_{n}h_{1}\dots h_{d}\over
A_{n}^{(1)}}\right)\alpha_{1}(\underline{A}_{n2})\varpi_{1}(A_{n}h_{1}\dots
h_{d})\to 0,$ (4.1) $\displaystyle\left({A_{n}^{(1)}\over A_{n}h_{1}\dots
h_{d}}\right)\sum_{k=1}^{\overline{A}_{n1}}k^{2d-1}\alpha_{1}^{1-4/q}(k)\to
0,$ (4.2)
$\displaystyle\left\\{\left({\overline{A}_{n1}\over\underline{A}_{n1}}\right)^{d}\left({\overline{A}_{n2}\over\overline{A}_{n1}}\right)+\left({A_{n}^{(1)}\over\underline{A}_{n1}^{d}}\right)\left({\left(\overline{A_{n}h}\right)^{d}\over
A_{n}h_{1}\dots
h_{d}}\right)\left({\overline{A}_{n1}\over\overline{A_{n}h}}\right)\right\\}\sum_{k=1}^{\overline{A}_{n1}}k^{d-1}\alpha_{1}^{1-2/q}(k)\to
0.$ (4.3)

We need Condition (ii) to compute the asymptotic variances of LP estimators.
Conditions (iii) and (iv) are concerned with the rates of convergence of
variance and bias terms of LP estimators, respectively. Condition (v) is
concerned with the large-block-small-block argument to show the asymptotic
normality of LP estimators. Indeed, we use the condition (4.1) to approximate
a weighted sum of spatially dependent data of the form

$\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{i}\end{array}\right)(e_{n,i}+\varepsilon_{n,i})$

by a sum of independent large blocks where

$H:={\rm
diag}(1,h_{1},\dots,h_{d},h_{1}^{2},h_{1}h_{2},\dots,h_{d}^{2},\dots,h_{1}^{p},h_{1}^{p-1}h_{2},\dots,h_{d}^{p})\in\mathbb{R}^{D\times
D}.$

The condition (4.2) is used to apply Lyapunov’s central limit theorem to the
sum of independent blocks. The condition (4.3) is used to show the asymptotic
negligibility of a sum of small blocks. See the proof of Theorem 4.1 for
detailed definitions of large and small blocks.

Throughout Sections 4.1, 4.2, and 4.3, we set $\bm{z}=\bm{0}$ without loss of
generality. Extending the results in this section to the case
$\bm{z}\in(-1/2,1/2)^{d}$ is straightforward.

######  Theorem 4.1 (Asymptotic normality of local polynomial estimators).

Suppose Assumptions 2.1, 2.2, 2.3, 3.1, and 4.1 hold. Then, as $n\to\infty$,
the following result holds:

$\displaystyle\sqrt{A_{n}h_{1}\dots
h_{d}}\left(H\left(\widehat{\bm{\beta}}(\bm{0})-\bm{M}(\bm{0})\right)-S^{-1}B^{(d,p)}\bm{M}_{n}^{(d,p)}(\bm{0})\right)$
$\displaystyle\quad\stackrel{{\scriptstyle
d}}{{\to}}N\left(\left(\begin{array}[]{c}0\\\ \vdots\\\
0\end{array}\right),\left\\{{\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))\over
g(\bm{0})}+\eta^{2}(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right\\}S^{-1}\mathcal{K}S^{-1}\right),$

where

$\displaystyle B^{(d,p)}$ $\displaystyle=\int\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)(\bm{z})^{\prime}_{p+1}K(\bm{z})d\bm{z}\in\mathbb{R}^{D\times\bar{D}},\
\mathcal{K}=\int\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)(1\
\check{\bm{z}}^{\prime})K^{2}(\bm{z})d\bm{z}\in\mathbb{R}^{D\times D},$
$\displaystyle\bm{M}_{n}^{(d,p)}(\bm{z})$
$\displaystyle=\left({\partial_{j_{1}\dots
j_{p+1}}m(\bm{z})\over\bm{s}_{j_{1}\dots
j_{p+1}}!}\prod_{\ell=1}^{p+1}h_{j_{\ell}}\right)^{\prime}_{1\leq
j_{1}\leq\dots\leq j_{p+1}\leq d}$ $\displaystyle=\left({\partial_{1\dots
1}m(\bm{z})\over(p+1)!}h_{1}^{p+1},{\partial_{1\dots 2}m(\bm{z})\over
p!}h_{1}^{p}h_{2},\dots,{\partial_{d\dots
d}m(\bm{z})\over(p+1)!}h_{d}^{p}\right)^{\prime}\in\mathbb{R}^{\bar{D}}.$

Theorem 4.1 differs from the asymptotic normality of LP estimators under
i.i.d. observations in several points. First, the convergence rates of LP
estimators depends not on the sample size $n$ explicitly but on the volume of
the sampling region $A_{n}$. Second, the asymptotic variance is represented as
a sum of two components
$\\{\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))\\}S^{-1}\mathcal{K}S^{-1}/g(\bm{0})$
and
$\eta^{2}(\bm{0})\left(\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right)S^{-1}\mathcal{K}S^{-1}$.
When the sampling design satisfies the mixed increasing domain asymptotics,
that is, $\kappa=0$, then the asymptotic variance depends only on the second
term, which represents the effect of the spatial dependence, and does not
includes $\sigma_{\varepsilon}^{2}(\bm{0})$, the effect of the measurement
error $\\{\varepsilon_{n,j}\\}$. This is completely different from i.i.d.
case. We also note that the form of the asymptotic variance in Theorem 4.1 is
different from that of Theorem 4 in Masry, 1996b  who investigates asymptotic
properties of LP estimators for equidistant time series. Indeed, in his
result, the variance term that corresponds to the second term of the
asymptotic variance in our result does not appear. When the sampling design
satisfies the pure increasing domain asymptotics, that is,
$\kappa\in(0,\infty)$, then the asymptotic variance depends on both first and
second terms. In this case, the asymptotic variance includes the effect of the
sampling design $1/g(\bm{0})$, which implies that the more likely the sampling
sites are distributed around $\bm{0}$, the more accurate the estimation of
$M(\bm{0})$. Moreover, if $\eta(\cdot)\equiv 0$, then the asymptotic variance
coincides with that of i.i.d. case.

######  Remark 4.1 (General form of the mean squared error of
$\partial_{j_{1}\dots j_{L}}\widehat{m}(\bm{0})$).

Define

$\displaystyle\bm{b}_{n}^{(d,p)}(\bm{x})$
$\displaystyle:=B^{(d,p)}\bm{M}_{n}^{(d,p)}(\bm{x})$
$\displaystyle=\left(b_{n,0}(\bm{x}),b_{n,1}(\bm{x}),\dots,b_{n,d}(\bm{x}),\right.$
$\displaystyle\left.\quad\quad\quad
b_{n,11}(\bm{x}),b_{n,12}(\bm{x}),\dots,b_{n,dd}(\bm{x}),\dots,b_{n,1\dots,1}(\bm{x}),b_{n,1\dots
2}(\bm{x}),\dots,b_{n,d\dots d}(\bm{x})\right)^{\prime}$

and let $e_{j_{1}\dots j_{L}}=(0,\dots,0,1,0,\dots,0)^{\prime}$ be a
$D$-dimensional vector such that $e_{j_{1}\dots
j_{L}}^{\prime}\bm{b}_{n}^{(d,p)}(\bm{x})=b_{j_{1}\dots j_{L}}(\bm{x})$.
Theorem 4.1 yields that

$\displaystyle b_{n,j_{1},\dots,j_{L}}(\bm{0})$ $\displaystyle=\sum_{1\leq
j_{1,1}\leq\dots\leq j_{1,p+1}\leq d}{\partial_{j_{1,1}\dots
j_{1,p+1}}m(\bm{0})\over\bm{s}_{j_{1,1}\dots
j_{1,p+1}}!}\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\kappa_{j_{1}\dots
j_{L}j_{1,1}\dots j_{1,p+1}}^{(1)},$

for $1\leq j_{1}\leq\dots\leq j_{L}\leq d$, $0\leq L\leq p$, and the mean
squared error (MSE) of LP estimator $\partial_{j_{1}\dots
j_{L}}\widehat{m}(\bm{0})$ is given as follows:

$\displaystyle\text{MSE}(\partial_{j_{1}\dots j_{L}}\widehat{m}(\bm{0}))$
$\displaystyle=E\left[\left(\partial_{j_{1}\dots
j_{L}}m(\bm{0})-\partial_{j_{1}\dots
j_{L}}\widehat{m}(\bm{0})\right)^{2}\right]$
$\displaystyle=\left\\{\bm{s}_{j_{1}\dots j_{L}}!{(S^{-1}e_{j_{1}\dots
j_{L}})^{\prime}B^{(d,p)}\bm{M}_{n}^{(d,p)}(\bm{0})\over\prod_{\ell=1}^{L}h_{j_{\ell}}}\right\\}^{2}$
$\displaystyle\quad+\left({\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))\over
g(\bm{0})}+\eta^{2}(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right)\left(\bm{s}_{j_{1}\dots
j_{L}}!\right)^{2}{e^{\prime}_{j_{1}\dots
j_{L}}S^{-1}\mathcal{K}S^{-1}e_{j_{1}\dots j_{L}}\over A_{n}h_{1}\dots
h_{d}\times\left(\prod_{\ell=1}^{L}h_{j_{\ell}}\right)^{2}}.$ (4.4)

###  4.2. Estimation of asymptotic variances

An estimator of the asymptotic variance of the statistics
$\widehat{\bm{\beta}}(\bm{0})$ can be constructed by using leave-one (or
two)-out estimators. For $\bm{z}\in(-1/2,1/2)^{d}$, let
$\widehat{m}_{-I}(\bm{z})$ be the LP estimator (of order $p$) of $m(\bm{z})$
computed without $\\{(Y(\bm{X}_{i}),\bm{X}_{i})\\}_{i\in I}$,
$I\subset\\{1,\dots,n\\}$.

Define

$\displaystyle\widehat{g}(\bm{0})$ $\displaystyle={1\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}(\bm{X}_{i}),$
$\displaystyle\widehat{V}_{n,1}(\bm{0})$ $\displaystyle={1\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}(\bm{X}_{i})\left(Y(\bm{X}_{i})-\widehat{m}_{-\\{i\\}}(\bm{X}_{i}/A_{n})\right)^{2},$
$\displaystyle\widehat{V}_{n,2}(\bm{0})$ $\displaystyle={A_{n}\over
nh_{1}\dots h_{d}}\sum_{i=1}^{n-1}K_{Ah}(\bm{X}_{i})K_{Ah}(\bm{X}_{i+1})$
$\displaystyle\quad\times\left(Y(\bm{X}_{i})-\widehat{m}_{-\\{i,i+1\\}}(\bm{X}_{i}/A_{n})\right)\left(Y(\bm{X}_{i+1})-\widehat{m}_{-\\{i,i+1\\}}(\bm{X}_{i+1}/A_{n})\right),$

Note that $\widehat{m}_{-\\{i\\}}(\bm{z})$ and
$\widehat{m}_{-\\{i,i+1\\}}(\bm{z})$ are leave-$i$-out and leave-$(i,i+1)$-out
version of $\widehat{m}(\bm{z})$, respectively and then
$\widehat{m}_{-\\{i\\}}(\bm{z})$ and $\bm{X}_{i}$ (or
$\widehat{m}_{-\\{i,i+1\\}}(\bm{z})$ and $\\{\bm{X}_{i},\bm{X}_{i+1}\\}$) are
independent under Assumption 2.2.

######  Proposition 4.1.

Under the assumptions of Theorem 4.1, as $n\to\infty$,

$\displaystyle\widehat{V}_{n}(\bm{0}):={(A_{n}/n)\widehat{V}_{n,1}(\bm{0})\over\widehat{g}^{2}(\bm{0})}+{(\kappa_{0}^{(2)})^{-1}\widehat{V}_{n,2}(\bm{0})\over\widehat{g}^{2}(\bm{0})}\stackrel{{\scriptstyle
p}}{{\to}}{\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))\over
g(\bm{0})}+\eta^{2}(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}.$

Theorem 4.1 and Proposition 4.1 enable us to construct confidence intervals of
$\partial_{j_{1}\dots j_{L}}m(\bm{0})$. Consider a confidence interval of the
form

$\displaystyle C_{n,j_{1}\dots j_{L}}(1-\tau)=\left[\partial_{j_{1}\dots
j_{L}}\widehat{m}(\bm{0})\pm\sqrt{{\widehat{V}_{n}(\bm{0})\left(\bm{s}_{j_{1}\dots
j_{L}}!\right)^{2}\left(e^{\prime}_{j_{1}\dots
j_{L}}S^{-1}\mathcal{K}S^{-1}e_{j_{1}\dots j_{L}}\right)\over A_{n}h_{1}\dots
h_{d}\left(\prod_{\ell=1}^{L}h_{j_{\ell}}\right)^{2}}}q_{1-\tau/2}\right],$

where $q_{1-\tau}$ is the $(1-\tau)$-quantile of the standard normal random
variable. Then we can show the asymptotic validity of the confidence interval
as follows:

######  Corollary 4.1.

Let $\tau\in(0,1)$. Under the assumptions of Theorem 4.1 with

$A_{n}h_{1}\dots h_{d}\left((S^{-1}e_{j_{1}\dots
j_{L}})^{\prime}B^{(d,p)}M_{n}^{(d,p)}(\bm{0})\right)^{2}\to 0$

as $n\to\infty$. Then, $\lim_{n\to\infty}P(\partial_{j_{1}\dots
j_{L}}m(\bm{0})\in C_{n,j_{1}\dots j_{L}}(1-\tau))=1-\tau$.

###  4.3. Two-sample test for spatially dependent data

In this section, we discuss a two-sample test for the partial derivatives of
the mean function as an application of our main results.

Consider the following nonparametric regression model:

$\displaystyle Y_{1}(\bm{x}_{1,\ell_{1}})$
$\displaystyle=m_{1}\left({\bm{x}_{1,\ell_{1}}\over
A_{n}}\right)+\eta_{1}\left({\bm{x}_{1,\ell_{1}}\over
A_{n}}\right)e_{1}(\bm{x}_{1,\ell_{1}})+\sigma_{\varepsilon,1}\left({\bm{x}_{1,\ell_{1}}\over
A_{n}}\right)\varepsilon_{1,\ell_{1}},\ \ell_{1}=1,\dots,n_{1}$ $\displaystyle
Y_{2}(\bm{x}_{2,\ell_{2}})$
$\displaystyle=m_{2}\left({\bm{x}_{2,\ell_{2}}\over
A_{n}}\right)+\eta_{2}\left({\bm{x}_{2,\ell_{2}}\over
A_{n}}\right)e_{2}(\bm{x}_{2,\ell_{2}})+\sigma_{\varepsilon,2}\left({\bm{x}_{2,\ell_{2}}\over
A_{n}}\right)\varepsilon_{2,\ell_{2}},\ \ell_{2}=1,\dots,n_{2},$

where $\bm{x}_{1,\ell_{1}},\bm{x}_{2,\ell_{2}}\in R_{n}$,
$\bm{e}=\\{e(\bm{x})=(e_{1}(\bm{x}),e_{2}(\bm{x}))^{\prime}:\bm{x}\in\mathbb{R}^{d}\\}$
is a bivariate stationary random field such that $E[e_{k}(\bm{0})]=0$,
$E[e_{k}^{2}(\bm{0})]=1$, and $\\{\varepsilon_{k,\ell_{k}}\\}$ is a sequence
of i.i.d. random variables such that $E[\varepsilon_{k,\ell_{k}}]=0$, $k=1,2$.

Assume that $\\{\bm{x}_{k,\ell_{k}}\\}$ are realizations of a sequence of
random variables $\\{\bm{X}_{k,\ell_{k}}\\}$ with density
$A_{n}^{-1}g_{k}(\cdot/A_{n})$ where $g_{k}(\cdot)$ is a probability density
function with support $[-1/2,1/2]^{d}$, $k=1,2$. This allows the sampling
sites $\\{\bm{x}_{1,\ell_{1}}\\}$ and $\\{\bm{x}_{2,\ell_{2}}\\}$ to be
different.

######  Assumption 4.2.

The bivariate random field $\bm{e}$ satisfies the following conditions:

  * (i)

$E[|e_{k}(\bm{0})|^{q_{2}}]<\infty$, $k=1,2$ for some integer $q_{2}>4$.

  * (ii)

Define $\Sigma_{\bm{e}}(\bm{x})=(\sigma_{\bm{e},jk}(\bm{x}))_{1\leq j,k\leq
2}$ where $\sigma_{\bm{e},jk}(\bm{x})=E[e_{j}(\bm{0})e_{k}(\bm{x})]$,
$j,k=1,2$. Assume that $\sigma_{\bm{e},kk}(\bm{0})=1$, $k=1,2$ and
$\int_{\mathbb{R}^{d}}|\sigma_{\bm{e},jk}(\bm{v})|d\bm{v}<\infty$, $j,k=1,2$.

  * (iii)

The random field $\bm{e}$ is $\alpha$-mixing with mixing coefficients
$\alpha(a;b)\leq\alpha_{1}(a)\varpi_{1}(b)$ such that as $n\to\infty$,

$\displaystyle
A_{n}^{(1)}\left(\alpha_{1}^{1-2/q}(\underline{A}_{n2})+\sum_{k=\underline{A}_{n1}}^{\infty}k^{d-1}\alpha_{1}^{1-2/q}(k)\right)\varpi_{1}^{1-2/q}(A_{n}^{(1)})\to
0,$

where $q=\min\\{q_{1},q_{2}\\}$, $A_{n}^{(1)}=\prod_{j=1}^{d}A_{n1,j}$,
$\underline{A}_{n1}=\min_{1\leq j\leq d}A_{n1,j}$, and
$\underline{A}_{n2}=\min_{1\leq j\leq d}A_{n2,j}$. Here,
$\\{A_{n1,j}\\}_{n\geq 1}$ and $\\{A_{n2,j}\\}_{n\geq 1}$ are sequences of
constants with $\min\left\\{A_{n2,j},{A_{n1,j}\over
A_{n2,j}}\right\\}\to\infty$ as $n\to\infty$, and $q_{1}$ is the integer that
appear in Assumption 2.1.

  * (iv)

$\\{\bm{X}_{1,\ell_{1}}\\}_{\ell_{1}=1}^{n_{1}}$,
$\\{\bm{X}_{2,\ell_{2}}\\}_{\ell_{2}=1}^{n_{2}}$, $\bm{e}$,
$\\{\varepsilon_{1,\ell_{1}}\\}_{\ell_{1}=1}^{n_{1}}$, and
$\\{\varepsilon_{2,\ell_{2}}\\}_{\ell_{2}=1}^{n_{2}}$ are mutually
independent.

In Section 6, we give examples of bivariate random fields that satisfy
Assumptions 4.1 and 4.2. We note that a wide class of bivariate Lévy-driven MA
random fields satisfies our assumptions.

We are interested in testing the null hypothesis

$\displaystyle\mathbb{H}_{0,j_{1}\dots j_{L}}:\partial_{j_{1}\dots
j_{L}}m_{1}(\bm{0})-\partial_{j_{1}\dots j_{L}}m_{2}(\bm{0})=0$ (4.5)

against the alternative $\mathbb{H}_{1,j_{1}\dots j_{L}}:\partial_{j_{1}\dots
j_{L}}m_{1}(\bm{0})-\partial_{j_{1}\dots j_{L}}m_{2}(\bm{0})\neq 0$.

Define $\bm{M}_{k}(\bm{0})$ as $\bm{M}(\bm{0})$ with $m=m_{k}$ and
$\overline{\bm{\beta}}_{k}(\bm{0})$ as LP estimators of order $p$ for
$\bm{M}_{k}(\bm{0})$ computed by using
$\\{(Y_{k}(\bm{x}_{k,\ell_{k}}),\bm{x}_{k,\ell_{k}})\\}$, bandwidths
$h_{1},\dots,h_{d}$, and a common kernel function $K$, $k=1,2$, respectively.
The next theorem is a building block of the two-sample test (4.5).

######  Proposition 4.2.

Suppose Assumptions 2.1, 2.2 (i), 3.1, 4.1, and 4.2 hold with $m=m_{k}$,
$\eta=\eta_{k}$, $\sigma_{\varepsilon}=\sigma_{\varepsilon,k}$,
$\\{\varepsilon_{j}\\}=\\{\varepsilon_{k,\ell_{k}}\\}$, $g=g_{k}$, $k=1,2$.
Moreover, assume that $n=n_{1}$, $n_{1}/n_{2}\to\theta\in(0,\infty)$ as
$n_{1}\to\infty$ and
$(\eta_{1}(\bm{0}),-\eta_{2}(\bm{0}))\left(\int\Sigma_{\bm{e}}(\bm{v})d\bm{v}\right)(\eta_{1}(\bm{0}),-\eta_{2}(\bm{0}))^{\prime}\geq
0$. Then, as $n\to\infty$,

$\displaystyle\sqrt{A_{n}h_{1}\dots
h_{d}}\left\\{H\left((\overline{\bm{\beta}}_{1}(\bm{0})-\overline{\bm{\beta}}_{2}(\bm{0}))-(\bm{M}_{1}(\bm{0})-\bm{M}_{2}(\bm{0}))\right)-(\overline{B}_{n1}(\bm{0})-\overline{B}_{n2}(\bm{0}))\right\\}$
$\displaystyle\quad\stackrel{{\scriptstyle
d}}{{\to}}N\left(\left(\begin{array}[]{c}0\\\ \vdots\\\
0\end{array}\right),\left(\overline{V}_{1}(\bm{0})+\overline{V}_{2}(\bm{0})-2\overline{V}_{3}(\bm{0})\right)S^{-1}\mathcal{K}S^{-1}\right),$

where

$\displaystyle\overline{B}_{n1}(\bm{0})$
$\displaystyle=S^{-1}B^{(d,p)}\bm{M}_{n1}^{(d,p)}(\bm{0}),\
\overline{B}_{n2}(\bm{0})=S^{-1}B^{(d,p)}\bm{M}_{n2}^{(d,p)}(\bm{0}),$
$\displaystyle\overline{V}_{1}(\bm{0})$
$\displaystyle=\left({\kappa(\eta_{1}^{2}(\bm{0})+\sigma_{\varepsilon,1}^{2}(\bm{0}))\over
g_{1}(\bm{0})}+\eta_{1}^{2}(\bm{0})\int\sigma_{\bm{e},11}(\bm{v})d\bm{v}\right),$
$\displaystyle\overline{V}_{2}(\bm{0})$
$\displaystyle=\left({\theta\kappa(\eta_{2}^{2}(\bm{0})+\sigma_{\varepsilon,2}^{2}(\bm{0}))\over
g_{2}(\bm{0})}+\eta_{2}^{2}(\bm{0})\int\sigma_{\bm{e},22}(\bm{v})d\bm{v}\right),$
$\displaystyle\overline{V}_{3}(\bm{0})$
$\displaystyle=\eta_{1}(\bm{0})\eta_{2}(\bm{0})\int\sigma_{\bm{e},12}(\bm{v})d\bm{v},$

where $\bm{M}_{nk}^{(d,p)}(\bm{0})$ are defined as
$\bm{M}_{n}^{(d,p)}(\bm{0})$ with $m=m_{k}$.

An estimator of the asymptotic variance of the statistics
$\overline{\bm{\beta}}_{1}(\bm{0})-\overline{\bm{\beta}}_{2}(\bm{0})$ can be
constructed as follows. For $\bm{z}\in(-1/2,1/2)^{d}$, let
$\widehat{m}_{k,-I_{k}}(\bm{z})$ be the LP estimator (of order $p$) of
$m_{k}(\bm{z})$ computed without
$\\{(Y_{k}(\bm{X}_{k,\ell_{k}}),\bm{X}_{k,\ell_{k}})\\}_{\ell_{k}\in I_{k}}$,
$I_{k}\subset\\{1,\dots,n_{k}\\}$, $k=1,2$.

Define

$\displaystyle\overline{g}_{n_{k}}(\bm{0})$ $\displaystyle={1\over
n_{k}h_{1}\dots h_{d}}\sum_{\ell_{k}=1}^{n_{k}}K_{Ah}(\bm{X}_{k,\ell_{k}}),$
$\displaystyle\overline{V}_{n,1k}(\bm{0})$ $\displaystyle={1\over
n_{k}h_{1}\dots
h_{d}}\sum_{\ell_{k}=1}^{n_{k}}K_{Ah}(\bm{X}_{k,\ell_{k}})\left(Y_{k}(\bm{X}_{k,\ell_{k}})-\widehat{m}_{k,-\\{\ell_{k}\\}}(\bm{X}_{k,\ell_{k}}/A_{n})\right)^{2},\
k=1,2,$ $\displaystyle\overline{V}_{n,2k}(\bm{0})$ $\displaystyle={A_{n}\over
n_{k}h_{1}\dots
h_{d}}\sum_{\ell_{k}=1}^{n_{k}-1}K_{Ah}(\bm{X}_{k,\ell_{k}})K_{Ah}(\bm{X}_{k,\ell_{k}+1})$
$\displaystyle\quad\times\left(Y_{k}(\bm{X}_{k,\ell_{k}})-\widehat{m}_{k,-\\{\ell_{k},\ell_{k}+1\\}}(\bm{X}_{k,\ell_{k}}/A_{n})\right)$
$\displaystyle\quad\times\left(Y_{k}(\bm{X}_{k,\ell_{k}+1})-\widehat{m}_{k,-\\{\ell_{k},\ell_{k}+1\\}}(\bm{X}_{k,\ell_{k}+1}/A_{n})\right),\
k=1,2,$ $\displaystyle\overline{V}_{n,3}(\bm{0})$ $\displaystyle={A_{n}\over
n_{1}n_{2}h_{1}\dots
h_{d}}\sum_{\ell_{1}=1}^{n_{1}}\sum_{\ell_{2}=1}^{n_{2}}K_{Ah}(\bm{X}_{1,\ell_{1}})K_{Ah}(\bm{X}_{2,\ell_{2}})$
$\displaystyle\quad\times\left(Y_{1}(\bm{X}_{1,\ell_{1}})-\widehat{m}_{1,-\\{\ell_{1}\\}}(\bm{X}_{1,\ell_{1}}/A_{n})\right)\left(Y_{2}(\bm{X}_{2,\ell_{2}})-\widehat{m}_{2,-\\{\ell_{2}\\}}(\bm{X}_{2,\ell_{2}}/A_{n})\right).$

######  Proposition 4.3.

Under the assumptions of Theorem 4.2, as $n\to\infty$,

$\displaystyle\check{V}_{n}(\bm{0})$
$\displaystyle:=\left\\{(A_{n}/n_{1})\overline{V}_{n,11}(\bm{0})+(\widehat{V}_{n,21}(\bm{0})/\kappa_{0}^{(2)})\over\overline{g}_{n_{1}}^{2}(\bm{0})\right\\}+\left\\{(A_{n}/n_{2})\overline{V}_{n,12}(\bm{0})+(\widehat{V}_{n,22}(\bm{0})/\kappa_{0}^{(2)})\over\overline{g}_{n_{2}}^{2}(\bm{0})\right\\}$
$\displaystyle\quad-2{(\overline{V}_{n,3}(\bm{0})/\kappa_{0}^{(2)})\over\overline{g}_{n_{1}}(\bm{0})\overline{g}_{n_{2}}(\bm{0})}\stackrel{{\scriptstyle
p}}{{\to}}\overline{V}_{1}(\bm{0})+\overline{V}_{2}(\bm{0})-2\overline{V}_{3}(\bm{0}).$

Define the test statistics

$\displaystyle T_{n,j_{1}\dots j_{L}}:={\sqrt{A_{n}h_{1}\dots
h_{d}\left(\prod_{\ell=1}^{L}h_{j_{\ell}}\right)^{2}}\left(\partial_{j_{1}\dots
j_{L}}\widehat{m}_{1}(\bm{0})-\partial_{j_{1}\dots
j_{L}}\widehat{m}_{2}(\bm{0})\right)\over\sqrt{\overline{V}_{n}(\bm{0})\left(\bm{s}_{j_{1}\dots
j_{L}}!\right)^{2}\left(e^{\prime}_{j_{1}\dots
j_{L}}S^{-1}\mathcal{K}S^{-1}e_{j_{1}\dots j_{L}}\right)}}.$

The asymptotic properties of the test statistics under both null and
alternative hypotheses are given as follows:

######  Corollary 4.2.

Let $\tau\in(0,1/2)$. Under the assumptions of Theorem 4.2 with

$A_{n}h_{1}\dots h_{d}\left((S^{-1}e_{j_{1}\dots
j_{L}})^{\prime}B^{(d,p)}M_{n1}^{(d,p)}(\bm{0})\right)^{2}\to 0,\ n\to\infty,$

then $\lim_{n\to\infty}P(|T_{n,j_{1}\dots j_{L}}|\geq q_{1-\tau/2})=\tau$
under $\mathbb{H}_{0,j_{1}\dots j_{L}}$ and
$\lim_{n\to\infty}P(|T_{n,j_{1}\dots j_{L}}|\geq q_{1-\tau/2})=1$ under
$\mathbb{H}_{1,j_{1}\dots j_{L}}$, where $q_{1-\tau}$ is the
$(1-\tau)$-quantile of the standard normal random variable.

##  5\. Uniform convergence of local polynomial estimators

In this section, we derive the uniform convergence rates of LP estimators for
the mean function of the model (2.1) and their partial derivatives. We note
that these results can be derived as special cases of the results on the
uniform convergence rates of more general kernel estimators provided in
Appendix. We assume the following conditions on the mean function $m$, the
variance function $\eta$, and $\\{\varepsilon_{n,i}\\}$:

######  Assumption 5.1.

  * Recall $R_{0}=[-1/2,1/2]^{d}$.

  * (i)

The mean function $m$ is $(p+1)$-times continuously partial differentiable on
$R_{0}$ and define $\partial_{j_{1}\dots j_{L}}m(\bm{z}):=\partial
m(\bm{z})/\partial z_{j_{1}}\dots\partial z_{j_{L}}$, $1\leq
j_{1},\dots,j_{L}\leq d$, $0\leq L\leq p+1$. When $L=0$, we set
$\partial_{j_{1}\dots j_{L}}m(\bm{z})=\partial_{j_{0}}m(\bm{z})=m(\bm{z})$.

  * (ii)

The function $\eta$ is continuous over $R_{0}$ and $\inf_{\bm{z}\in
R_{0}}\eta(\bm{z})>0$.

  * (iii)

The sequence of random variables $\\{\varepsilon_{i}\\}_{i=1}^{n}$ are i.i.d.
with $E[\varepsilon_{1}]=0$, $E[\varepsilon_{1}^{2}]=1$,
$E[|\varepsilon_{1}|^{q_{1}}]<\infty$ for some integer $q_{1}>4$ and the
function $\sigma_{\varepsilon}(\cdot)$ is continuous over $R_{0}$ and
$\inf_{\bm{z}\in R_{0}}\sigma_{\varepsilon}(\bm{z})>0$.

For the sampling sites $\\{\bm{X}_{i}\\}_{i=1}^{n}$, we assume the following
conditions:

######  Assumption 5.2.

Let $g$ be a probability density function with support $R_{0}=[-1/2,1/2]^{d}$.

  * (i)

$A_{n}/n\to\kappa\in[0,\infty)$ as $n\to\infty$,

  * (ii)

$\\{\bm{X}_{i}=(X_{i,1},\dots,X_{i,d})^{\prime}\\}_{i=1}^{n}$ is a sequence of
i.i.d. random vectors with density $A_{n}^{-d}g(\cdot/A_{n})$ and $g$ is
continuous and positive on $R_{0}$.

  * (iii)

$\\{\bm{X}_{i}\\}_{i=1}^{n}$,
$\bm{e}=\\{e(\bm{x}):\bm{x}\in\mathbb{R}^{d}\\}$, and
$\\{\varepsilon_{i}\\}_{i=1}^{n}$ are mutually independent.

We also assume the following conditions on the bandwidth $h_{j}$ and the
random field $\bm{e}=\\{e(\bm{x}):\bm{x}\in\mathbb{R}^{d}\\}$:

######  Assumption 5.3.

For $j=1,\dots,d$, let $\\{A_{n1,j}\\}_{n\geq 1}$, $\\{A_{n2,j}\\}_{n\geq 1}$
be sequence of positive numbers.

  * (i)

The random field $\bm{e}$ is stationary and $E[|e(\bm{0})|^{q_{2}}]<\infty$
for some integer $q_{2}>4$.

  * (ii)

Define $\sigma_{\bm{e}}(\bm{x})=E[e(\bm{0})e(\bm{x})]$. Assume that
$\int_{\mathbb{R}^{d}}|\sigma_{\bm{e}}(\bm{v})|d\bm{v}<\infty$.

  * (iii)

$\min\left\\{A_{n2,j},{A_{n1,j}\over A_{n2,j}},{A_{n,j}h_{j}\over
A_{n1,j}}\right\\}\to\infty$ as $n\to\infty$.

  * (iv)

The random field $\bm{e}$ is $\beta$-mixing with mixing coefficients
$\beta(a;b)\leq\beta_{1}(a)\varpi_{2}(b)$ such that as $n\to\infty$, $h_{j}\to
0$, $1\leq j\leq d$,

$\displaystyle{A_{n}^{(1)}\over(\overline{A}_{n1})^{d}}\sim 1,\
{A_{n}^{{1\over 2}}(h_{1}\dots h_{d})^{{1\over 2}}\over
n^{1/q_{2}}(\overline{A}_{n1})^{d}(\log n)^{{1\over 2}+\iota}}\gtrsim 1\
\text{for some $\iota\in(0,\infty)$},$ (5.1)
$\displaystyle\sqrt{{n^{2}A_{n}h_{1}\dots h_{d}\over(A_{n}^{(1)})^{2}\log
n}}\beta_{1}(\underline{A}_{n2})\varpi_{2}(A_{n}h_{1}\dots h_{d})\to 0,$ (5.2)

where

$\displaystyle A_{n}^{(1)}$ $\displaystyle=\prod_{j=1}^{d}A_{n1,j},\
\overline{A}_{n1}=\max_{1\leq j\leq d}A_{n1,j},\
\underline{A}_{n1}=\min_{1\leq j\leq d}A_{n1,j},$
$\displaystyle\overline{A}_{n2}$ $\displaystyle=\max_{1\leq j\leq d}A_{n2,j},\
\underline{A}_{n2}=\min_{1\leq j\leq d}A_{n2,j}.$

Condition (5.2) is concerned with large-block-small-block argument for
$\beta$-mixing sequences. In order to derive uniform convergence rates of LP
estimators, more careful arguments on the effects of non-equidistant sampling
sites are necessary than those for proving asymptotic normality and this also
requires additional works in comparison with the equidistant time series or
spatial data. We first approximate LP estimators excluding bias terms, which
can be written as a sum of spatially dependent data, by a sum of independent
blocks by extending the blocking technique in Yu, (1994)(Corollary 2.7) that
does not require regularly spaced sampling sites. Then we derive the uniform
convergence rates of LP estimators by applying maximum inequalities for
independent and possibly not identically distributed random variables to the
independent blocks. In Section 6, we will show that a wide class of Lévy-
driven MA random fields satisfies our $\beta$-mixing conditions.

######  Remark 5.1 (Discussion on $\beta$-mixing conditions).

Lahiri, 2003b  established central limit theorems for weighted sample means of
bounded spatial data under $\alpha$-mixing conditions. Lahiri’s proof relies
essentially on approximating the characteristic function of the weighted
sample mean by that of independent blocks using the Volkonskii-Rozanov
inequality (cf. Proposition 2.6 in Fan and Yao, (2003)) and then showing that
the characteristic function corresponding to the independent blocks converges
to the characteristic function of its Gaussian limit. However, characteristic
functions are difficult to capture the uniform behavior of LP estimators over
a compact set so we rely on a different argument from that of Lahiri, 2003b .
Indeed, we use a blocking argument tailored to $\beta$-mixing sequences (cf.
Corollary 2.7 in Yu, (1994)) and this enables us to compare the uniform
convergence rates of LP estimators with that of a sum of independent blocks
that approximates LP estimators. Another approach for handling spatial
dependence is $m$-dependent approximation under a physical dependence
structure (cf. El Machkouri et al., (2013)), but this approach is designed for
regularly spaced spatial data on $\mathbb{Z}^{d}$ and does not work in our
framework. We also note that it is not known that the results corresponding to
Corollary 2.7 in Yu, (1994) hold for $\alpha$-mixing sequences; see Remark
(ii) right after the proof of Lemma 4.1 in Yu, (1994).

We assume the following conditions on the kernel function $K$:

######  Assumption 5.4.

Let $K:\mathbb{R}^{d}\to\mathbb{R}$ be a kernel function such that

  * (i)

$\int K(\bm{z})d\bm{z}=1$.

  * (ii)

The kernel function $K$ is bounded and supported on
$[-C_{K},C_{K}]^{d}\subset[-1/2,1/2]^{d}$ for some $C_{K}>0$. Moreover, $K$ is
Lipschitz continuous on $\mathbb{R}^{d}$, i.e.,
$|K(\bm{v}_{1})-K(\bm{v}_{2})|\leq L_{K}|\bm{v}_{1}-\bm{v}_{2}|$ for some
$L_{K}\in(0,\infty)$ and all $\bm{v}_{1},\bm{v}_{2}\in\mathbb{R}^{d}$.

  * (iii)

Define $\kappa_{0}^{(r)}:=\int K^{r}(\bm{z})d\bm{z}$,
$\kappa_{j_{1},\dots,j_{M}}^{(r)}:=\int\prod_{\ell=1}^{M}z_{j_{\ell}}K^{r}(\bm{z})d\bm{z}$,
and

$\check{\bm{z}}:=(1,(\bm{z})^{\prime}_{1},\dots,(\bm{z})^{\prime}_{p})^{\prime},\
(\bm{z})_{L}=\left(\prod_{\ell=1}^{L}z_{j_{\ell}}\right)^{\prime}_{1\leq
j_{1}\leq\dots\leq j_{L}\leq d},\ 1\leq L\leq p.$

The matrix $S=\int\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)(1\ \check{\bm{z}}^{\prime})K(\bm{z})d\bm{z}$
is non-singular.

The next result provides uniform convergence rates of LP estimators
$\partial_{j_{1}\dots j_{L}}\widehat{m}(\bm{z})$.

######  Theorem 5.1.

Define $\mathrm{T}_{n}=\prod_{j=1}^{d}[-1/2+C_{K}h_{j},1/2-C_{K}h_{j}]$.
Suppose that Assumptions 5.1, 5.2, 5.3, and 5.4 hold with $q_{1}\geq q_{2}$.
Then for $1\leq j_{1}\leq\dots\leq j_{L}\leq d$, $0\leq L\leq p$, as
$n\to\infty$, we have

$\displaystyle\sup_{\bm{z}\in\mathrm{T}_{n}}\left|\partial_{j_{1}\dots
j_{L}}\widehat{m}(\bm{z})-\partial_{j_{1}\dots
j_{L}}\widehat{m}(\bm{z})\right|$ $\displaystyle\quad=O_{p}\left({\sum_{1\leq
j_{1}\leq\dots\leq j_{p+1}\leq
d}\prod_{\ell=1}^{p+1}h_{j_{\ell}}\over\prod_{\ell=1}^{L}h_{j_{\ell}}}+\sqrt{{\log
n\over A_{n}h_{1}\dots
h_{d}\left(\prod_{\ell=1}^{L}h_{j_{\ell}}\right)^{2}}}\right).$

##  6\. Examples

In this section, we discuss examples of random fields to which our theoretical
results can be applied. To this end, we consider Lévy-driven moving average
(MA) random fields and discuss their dependence structure. Lévy-driven MA
random fields include many Gaussian and non-Gaussian random fields and
constitute a flexible class of models for spatial data. We refer to Bertoin,
(1996) and Sato, (1999) for standard references on Lévy processes, and Rajput
and Rosinski, (1989) and Kurisu, (2022) for details on the theory of
infinitely divisible measures and fields. In particular, we show that a broad
class of Lévy-driven MA random fields, which includes continuous
autoregressive and moving average (CARMA) random fields as special cases (cf.
Brockwell and Matsuda, (2017)), satisfies our assumptions.

For the two-sample test discussed in Section 4, we considered nonparametric
regression models for spatial data
$\\{Y_{1}(\bm{x}_{1,\ell_{1}}),Y_{2}(\bm{x}_{2,\ell_{2}})\\}$ with bivariate
random field
$\bm{e}=\\{e(\bm{x})=(e_{1}(\bm{x}),e_{2}(\bm{x}))^{\prime}:\bm{x}\in\mathbb{R}^{d}\\}$.
Hence, we give examples of bivariate random fields that satisfy Assumptions
4.1 and 4.2. The examples of univariate random fields that satisfy Assumptions
2.3 and 4.1 (for Theorem 4.1, Proposition 4.1, and Corollary 4.1), and
Assumption 5.3 (for Theorem 5.1) can be given as special class of bivariate
cases. Indeed, for univariate cases, it is sufficient to consider the first
component of the examples of bivariate random fields.

Let
$\bm{L}=\\{\bm{L}(A)=(L_{1}(A),L_{2}(A))^{\prime}:A\in\mathcal{B}(\mathbb{R}^{d})\\}$
be an $\mathbb{R}^{2}$-valued random measure on the Borel subsets
$\mathcal{B}(\mathbb{R}^{d})$ that satisfies the following conditions:

  * 1.

For each sequence $\\{A_{m}\\}_{m\geq 1}$ of disjoint sets in
$\mathbb{R}^{d}$,

    * (a)

$\bm{L}(\cup_{m\geq 1}A_{m})=\sum_{m\geq 1}\bm{L}(A_{m})$ a.s. whenever
$\cup_{m\geq 1}A_{m}\in\mathcal{B}(\mathbb{R}^{d})$,

    * (b)

$\\{\bm{L}(A_{m})\\}_{m\geq 1}$ is a sequence of independent random variables.

  * 2.

For every Borel subset $A$ of $\mathbb{R}^{d}$ with finite Lebesgue measure
$|A|$, $\bm{L}(A)$ has an infinitely divisible distribution, that is,

$\displaystyle E[\exp(\mathrm{i}\bm{\theta}^{\prime}\bm{L}(A))]$
$\displaystyle=\exp(|A|\psi(\bm{\theta})),\ \bm{\theta}\in\mathbb{R}^{2},$
(6.1)

where $\mathrm{i}=\sqrt{-1}$ and $\psi$ is the logarithm of the characteristic
function of an $\mathbb{R}^{2}$-valued infinitely divisible distribution,
which is given by

$\displaystyle\psi(\bm{\theta})$
$\displaystyle=\mathrm{i}\bm{\theta}^{\prime}\bm{\gamma}_{0}-{1\over
2}\bm{\theta}^{\prime}\Sigma_{0}\bm{\theta}+\int_{\mathbb{R}^{2}}\left\\{e^{\mathrm{i}\bm{\theta}^{\prime}\bm{x}}-1-\mathrm{i}\bm{\theta}^{\prime}\bm{x}1_{\\{\|\bm{x}\|\leq
1\\}}\right\\}\nu_{0}(d\bm{x}),$

where $\bm{\gamma}_{0}=(\gamma_{0,1},\gamma_{0,2})^{\prime}\in\mathbb{R}^{2}$,
$\Sigma_{0}=(\sigma_{0,jk})_{1\leq j,k\leq 2}$ is a $2\times 2$ positive semi-
definite matrix, and $\nu_{0}$ is a Lévy measure with
$\int_{\mathbb{R}^{2}}\min\\{1,\|\bm{x}\|^{2}\\}\nu_{0}(d\bm{x})<\infty$. If
$\nu_{0}(d\bm{x})$ has a Lebesgue density, i.e.,
$\nu_{0}(d\bm{x})=\nu_{0}(\bm{x})d\bm{x}$, we call $\nu_{0}(\bm{x})$ as the
Lévy density. The triplet $(\bm{\gamma}_{0},\Sigma_{0},\nu_{0})$ is called the
Lévy characteristic of $\bm{L}$ and uniquely determines the distribution of
$\bm{L}$.

By equation (6.1), the first and second moments of the random measure $L$ are
determined by

$E[L_{j}(A)]=\mu_{j}^{(\bm{L})}|A|,\
\operatorname{Cov}(L_{j}(A),L_{k}(A))=\sigma_{j,k}^{(\bm{L})}|A|,$

where
$\mu_{j}^{(\bm{L})}=-\mathrm{i}{\partial\psi(\bm{0})\over\partial\theta_{j}}$
and
$\sigma_{j,k}^{(\bm{L})}=-{\partial^{2}\psi(\bm{0})\over\partial\theta_{j}\partial\theta_{k}}$.

The following are a couple of examples of Lévy random measures.

  * •

If $\psi(\bm{\theta})=-\bm{\theta}^{\prime}\Sigma_{0}^{2}\bm{\theta}/2$ with a
$2\times 2$ positive semi-definite matrix $\Sigma_{0}$, then $\bm{L}$ is a
Gaussian random measure.

  * •

If
$\psi(\bm{\theta})=\lambda\int_{\mathbb{R}^{2}}(\exp(\mathrm{i}\bm{\theta}^{\prime}\bm{x})-1)F(d\bm{x})$,
where $\lambda>0$ and $F$ is a probability distribution function with no jump
at the origin, then $\bm{L}$ is a compound Poisson random measure with
intensity $\lambda$ and jump size distribution $F$. More specifically,

$\displaystyle\bm{L}(A)$
$\displaystyle=\sum_{i=1}^{\infty}\bm{J}_{i}1_{\\{\bm{s}_{i}\\}}(A),\
A\in\mathcal{B}(\mathbb{R}^{d}),$

where $\bm{s}_{i}$ denotes the location of the $i$th unit point mass of a
Poisson random measure on $\mathbb{R}^{d}$ with intensity $\lambda>0$ and
$\\{\bm{J}_{i}\\}$ is a sequence of i.i.d. random vectors in $\mathbb{R}^{2}$
with distribution function $F$ independent of $\\{\bm{s}_{i}\\}$.

Let $\bm{\phi}=(\phi_{j,k})_{1\leq j,k\leq 2}$ be a measurable function on
$\mathbb{R}^{d}$ with $\phi_{j,k}\in L^{1}(\mathbb{R}^{d})\cap
L^{\infty}(\mathbb{R}^{d})$. A bivariate Lévy-driven MA random field with
kernel $\bm{\phi}$ driven by a Lévy random measure $\bm{L}$ is defined by

$\displaystyle\bm{e}(\bm{x})$
$\displaystyle=\int_{\mathbb{R}^{d}}\bm{\phi}(\bm{x}-\bm{u})\bm{L}(d\bm{u}),\
\bm{x}\in\mathbb{R}^{d}.$ (6.2)

Define $\bm{\mu}_{\bm{L}}=(\mu_{1}^{(\bm{L})},\mu_{2}^{(\bm{L})})^{\prime}$
and $\Sigma_{\bm{L}}=(\sigma_{j,k}^{(\bm{L})})_{1\leq j,k\leq 2}$. The first
and second moments of $\bm{e}(\bm{x})$ satisfy

$\displaystyle E[\bm{e}(\bm{0})]$
$\displaystyle=\bm{\mu}_{(L)}\int_{\mathbb{R}^{d}}\bm{\phi}(\bm{u})d\bm{u},\
\operatorname{Cov}(\bm{e}(\bm{0}),\bm{e}(\bm{x}))=\int_{\mathbb{R}^{d}}\bm{\phi}(\bm{x}-\bm{u})\Sigma_{\bm{L}}\bm{\phi}(\bm{u})d\bm{u}.$

We refer to Brockwell and Matsuda, (2017) for more details on the computation
of moments of Lévy-driven MA processes.

Before discussing theoretical results, we look at some examples of univariate
random fields defined by (6.2). Let
$a_{*}(z)=z^{p_{0}}+a_{1}z^{p_{0}-1}+\cdots+a_{p_{0}}=\prod_{i=1}^{p_{0}}(z-\lambda_{i})$
be a polynomial of degree $p_{0}$ with real coefficients and distinct negative
zeros $\lambda_{1},\ldots,\lambda_{p_{0}}$, and let
$b_{*}(z)=b_{0}+b_{1}z+\cdots+b_{q_{0}}z^{q_{0}}=\prod_{i=1}^{q_{0}}(z-\xi_{i})$
be a polynomial of degree $q_{0}$ with real coefficients and real zeros
$\xi_{1},\ldots,\xi_{q_{0}}$ such that $b_{q_{0}}=1$ and $0\leq q_{0}<p_{0}$
and $\lambda_{i}^{2}\neq\xi_{j}^{2}$ for all $i$ and $j$. Define
$a(z)=\prod_{i=1}^{p_{0}}(z^{2}-\lambda_{i}^{2})$ and
$b(z)=\prod_{i=1}^{q_{0}}(z^{2}-\xi_{i}^{2})$. Then, the Lévy-driven MA random
field driven by an infinitely divisible random measure $L$ with

$\displaystyle\phi(\bm{x})=\sum_{i=1}^{p_{0}}{b(\lambda_{i})\over
a^{\prime}(\lambda_{i})}e^{\lambda_{i}\|\bm{x}\|},$

where $a^{\prime}$ denotes the derivative of the polynomial $a$, is called a
univariate (isotropic) CARMA($p_{0},q_{0}$) random field. For example, if the
Lévy random measure of a CARMA random field is compound Poisson, then the
resulting random field is called a compound Poisson-driven CARMA random field.
In particular, when

$\phi(\bm{x})=(1-\varsigma)\exp(\lambda_{1}\|\bm{x}\|)+\varsigma\exp(\lambda_{2}\|\bm{x}\|),$

where $\varsigma$ is a parameter that satisfies

$-{\lambda_{2}^{2}-\xi^{2}\lambda_{1}\over\lambda_{1}^{2}-\xi^{2}\lambda_{2}}={\varsigma\over
1-\varsigma},\ \lambda_{1}<\lambda_{2}<0,\ \xi\leq 0,$

then the random field (6.2) is called a CARMA($2,1$) random field. This random
field includes normalized CAR($1$) (when $\varsigma=0$) and CAR($2$) (when
$\varsigma=-\lambda_{1}/(\lambda_{2}-\lambda_{1})$) as special cases. See
Brockwell and Matsuda, (2017) for more details. We note that although we focus
on isotropic case, it is possible to extend the results in this section to
anisotropic Lévy-driven MA random fields.

######  Remark 6.1 (Connections to Matérn covariance functions).

In spatial statistics, Gaussian random fields with the following Matérn
covariance functions play an important role (cf. Matérn,, 1986; Stein,, 1999;
Guttorp and Gneiting,, 2006):

$M(\bm{x};\nu,a,\sigma)=\sigma^{2}\|a\bm{x}\|^{\nu}F_{\nu}(\|a\bm{x}\|),\
\nu>0,a>0,\sigma>0,$

where $F_{\nu}$ denotes the modified Bessel function of the second kind of
order $\nu$ (we call $\nu$ the index of Matérn covariance function). Brockwell
and Matsuda, (2017) showed that in the univariate case, when the kernel
function is $\phi(\bm{x})=\|a\bm{x}\|^{\nu}F_{\nu}(\|a\bm{x}\|)$, which they
call a Matérn kernel with index $\nu$, then the Levy-driven MA random field
has a Matérn covariance function with index $d/2+\nu$. For example, a
normalized CAR($1$) random field has a Matérn covariance function since its
kernel function is given by
$\phi(\bm{x})=\exp(-\|\lambda_{1}\bm{x}\|)=\sqrt{(2/\pi)}\|\lambda_{1}\bm{x}\|^{1/2}F_{1/2}(\|\lambda_{1}\bm{x}\|)$
for some $\lambda_{1}<0$.

In general, if $\bm{\phi}$ depends only on $\|\bm{x}\|$, i.e.,
$\bm{\phi}(\bm{x})=\bm{\phi}(\|\bm{x}\|)$, then $\bm{e}$ is a strictly
stationary isotropic random field and the second moment of $\bm{e}(\bm{x})$
satisfies

$\displaystyle\operatorname{Cov}(\bm{e}(\bm{0}),\bm{e}(\bm{x}))=\int_{\mathbb{R}^{d}}\bm{\phi}(\|\bm{x}-\bm{u}\|)\Sigma_{\bm{L}}\bm{\phi}(\|\bm{u}\|)d\bm{u}.$

Consider the following decomposition:

$\displaystyle\bm{e}(\bm{x})$
$\displaystyle=\int_{\mathbb{R}^{d}}\bm{\phi}(\bm{x}-\bm{u})\psi_{0}\left(\|\bm{x}-\bm{u}\|:m_{n}\right)\bm{L}(d\bm{u})+\int_{\mathbb{R}^{d}}\bm{\phi}(\bm{x}-\bm{u})\left(1-\psi_{0}\left(\|\bm{x}-\bm{u}\|:m_{n}\right)\right)\bm{L}(d\bm{u})$
$\displaystyle=:\bm{e}_{1,m_{n}}(\bm{x})+\bm{e}_{2,m_{n}}(\bm{x}),$

where $m_{n}$ is a sequence of positive constants with $m_{n}\to\infty$ as
$n\to\infty$ and $\psi_{0}(\cdot:c):\mathbb{R}\to[0,1]$ is a truncation
function defined by

$\displaystyle\psi_{0}(x:c)=\begin{cases}1&\text{if $|x|\leq c/4$},\\\
-{4\over c}\left(x-{c\over 2}\right)&\text{if $c/4<|x|\leq c/2$},\\\
0&\text{if $x>c/2$}.\end{cases}$

The random field
$\bm{e}_{1,m_{n}}=\\{\bm{e}_{1,m_{n}}(\bm{x})=(e_{11,m_{n}}(\bm{x}),e_{12,m_{n}}(\bm{x}))^{\prime}:\bm{x}\in\mathbb{R}^{d}\\}$
is $m_{n}$-dependent (with respect to the $\ell^{2}$-norm), i.e.,
$\bm{e}_{1,m_{n}}(\bm{x}_{1})$ and $\bm{e}_{1,m_{n}}(\bm{x}_{2})$ are
independent if $\|\bm{x}_{1}-\bm{x}_{2}\|\geq m_{n}$. Also, if the tail of the
kernel function $\bm{\phi}(\cdot)$ decays sufficiently fast, then the random
field
$\bm{e}_{2,m_{n}}=\\{\bm{e}_{2,m_{n}}(\bm{x})=(e_{21,m_{n}}(\bm{x}),e_{22,m_{n}}(\bm{x}))^{\prime}:\bm{x}\in\mathbb{R}^{d}\\}$
is asymptotically negligible. In such cases, we can approximate $\bm{e}$ by
the $m_{n}$-dependent process $\bm{e}_{1,m_{n}}$ and verify conditions on
mixing coefficients in Assymptions 2.3, 4.1, and 4.2 as shown in the following
proposition.

######  Proposition 6.1.

Consider a Lévy-driven MA random field $\bm{e}$ defined by (6.2). Assume that
$\phi_{j,k}(\bm{x})=r_{0,jk}e^{-r_{1,jk}\|\bm{x}\|}$ where $|r_{0,jk}|>0$ and
$r_{1,jk}>0$, $j,k=1,2$. Additionally, assume that

  * (a)

the random measure $\bm{L}(\cdot)$ is Gaussian with triplet $(0,\Sigma_{0},0)$
or

  * (b)

the random measure $\bm{L}(\cdot)$ is non-Gaussian with triplet
$(\bm{\gamma}_{0},0,\nu_{0})$, $\bm{\mu}_{(L)}=(0,0)^{\prime}$, and the
marginal Lévy density $\nu_{0,j}(x)$ of $L_{j}(\cdot)$ is given by

$\displaystyle\nu_{0,j}(x)$
$\displaystyle={1\over|x|^{1+\beta_{0,j}}}\left(C_{0,j}e^{-c_{0,j}|x|^{\alpha_{0,j}}}+{C_{1,j}\over
1+|x|^{\beta_{1,j}}}\right)1_{\mathbb{R}\backslash\\{0\\}}(x),\ $ (6.3)

where $\alpha_{0,j}>0$, $\beta_{0,j}\geq-1$, $\beta_{1,j}>0$,
$\beta_{0,j}+\beta_{1,j}>6$, $c_{0,j}>0$, $C_{0,j}\geq 0$, $C_{1,j}\geq 0$,
and $C_{0,j}+C_{1,j}>0$, $j=1,2$.

Then $\bm{e}_{2,m_{n}}$ is asymptotically negligible, that is, we can replace
$\bm{e}$ with $\bm{e}_{1,m_{n}}$ in the results in Section 4. Further,
$\bm{e}_{1,m_{n}}$ satisfies Assumptions 2.3, 4.1, and 4.2 with $A_{n,j}\sim
n^{\zeta_{0}/d}$, $A_{n1,j}=A_{n,j}^{\zeta_{1}}$,
$A_{n2,j}=A_{n1,j}^{\zeta_{2}}$, $m_{n}=\underline{A}_{n2}^{1/2}$, and
$h_{j}\sim n^{-\zeta_{3}/d}$ where $\zeta_{0}$, $\zeta_{1}$, $\zeta_{2}$, and
$\zeta_{3}$ are positive constants such that

$\displaystyle\zeta_{0}$ $\displaystyle\in\left(0,{2p+2\over d}\right),\
\zeta_{1}\in\left({\zeta_{0}d\over d+2p+2},{2p+2\over d+2p+2}\right),$
$\displaystyle\zeta_{2}$ $\displaystyle\in\left(0,\min\left\\{{2\over
2+d\max\\{1,\zeta_{0}\\}},1-{\zeta_{0}d\over\zeta_{1}(d+2p+2)},{2p+2\over\zeta_{1}(d+2p+2)}-1\right\\}\right),$
$\displaystyle\zeta_{3}$ $\displaystyle\in\left({d\zeta_{0}\over
2p+d+2},\min\left\\{{d\zeta_{0}\over
2p+d},\zeta_{0}\left(1-\zeta_{1}(1+\zeta_{2})\right),\zeta_{1}\left(1-\left(1+{d\over
2}\zeta_{0}\right)\zeta_{2}\right)\right\\}\right).$

######  Remark 6.2.

When $d=2$ and $p\geq 1$, the conditions on $\\{\zeta_{j}\\}_{j=0}^{3}$ are
typically satisfied when $\zeta_{0}=1$, $\zeta_{1}={3\over 2p+4}$,
$\zeta_{2}\in\left(0,{1\over 6}\right)$. The Lévy density of the form (6.3)
corresponds to a compound Poisson random measure if $\beta_{0,j}\in[-1,0)$, a
Variance Gamma random measure if $\alpha_{0,j}=1$, $\beta_{0,j}=0$,
$C_{1,j}=0$, and a tempered stable random measure if $\beta_{0,j}\in(0,1)$,
$C_{1,j}=0$ (cf. Section 5 in Kato and Kurisu, (2020)). It is straight forward
to extend Proposition 6.1 to the case that $\phi$ is a finite sum of kernel
functions with exponential decay. Therefore, our results in Section 4 can be
applied to a wide class of CARMA($p_{0},q_{0}$) random fields and extending
the results to anisotropic CARMA random fields (cf. Brockwell and Matsuda,
(2017)) is straightforward.

The next result provides examples of Lévy-driven MA random fields that satisfy
assumptions in Theorem 5.1.

######  Proposition 6.2.

Consider a univariate Lévy-driven MA random field $\bm{e}$ defined by (6.2).
Assume that $\phi(\bm{x})=r_{0}e^{-r_{1}\|\bm{x}\|}$ where $|r_{0}|>0$ and
$r_{1}>0$. Additionally, assume Conditions (a) or (b) in Proposition 6.1. Then
$\bm{e}_{2,m_{n}}$ is asymptotically negligible, that is, we can replace
$\bm{e}$ with $\bm{e}_{1,m_{n}}$ in Theorem 5.1. Further, $\bm{e}_{1,m_{n}}$
satisfies Assumption 5.3 with $A_{n,j}\sim n^{\zeta_{0}/d}$,
$A_{n1,j}=A_{n,j}^{\zeta_{1}}$, $A_{n2,j}=A_{n1,j}^{\zeta_{2}}$,
$m_{n}=\underline{A}_{n2}^{1/2}$, and $h_{j}\sim n^{-\zeta_{3}/d}$ where
$\zeta_{0}$, $\zeta_{1}$, $\zeta_{2}$, and $\zeta_{3}$ are positive constants
such that $\zeta_{0}>{2\over q_{2}}$, $\zeta_{1}\in\left(0,{1\over
2}-{1\over\zeta_{0}q_{2}}\right)$, $\zeta_{2}\in(0,1)$, and
$\zeta_{3}\in\left(0,\min\\{1,\zeta_{0}(1-2\zeta_{1})-{2\over
q_{2}}\\}\right)$.

##  7\. Conclusion

In this paper, we have advanced statistical theory of nonparametric regression
for irregularly spaced spatial data. For this, we introduced a nonparametric
regression model defined on a sampling region $R_{n}\subset\mathbb{R}^{d}$ and
derived asymptotic normality and uniform convergence rates of the local
polynomial estimators of order $p\geq 1$ for the mean function of the model
under a stochastic sampling design. As an application of our main results, we
discussed a two-sample test for the mean functions and their partial
derivatives. We also provided examples of random fields that satisfy our
assumptions. In particular, our assumptions hold for a wide class of random
fields that includes Lévy-driven moving average random fields and popular
Gaussian random fields as special cases.

##  Appendix A Proofs for Section 4

###  A.1. Proof of Theorem 4.1

###### Proof.

Define $\bm{h}:=(h_{1},\dots,h_{d})^{\prime}$ and for
$\bm{x},\bm{y}\in\mathbb{R}^{d}$, let
$\bm{x}\circ\bm{y}=(x_{1}y_{1},\dots,x_{d}y_{d})^{\prime}$ be the Hadamard
product. Considering Taylor’s expansion of $m(\bm{z})$ around $\bm{z}$,

$\displaystyle m(\bm{X}_{i}/A_{n})$
$\displaystyle=(1,\check{{\text{\boldmath$X$}}}^{\prime}_{i})M(\bm{z})+{1\over(p+1)!}\sum_{1\leq
j_{1}\leq\dots\leq j_{p+1}\leq d}{(p+1)!\over\bm{s}_{j_{1}\dots
j_{p+1}}!}\partial_{j_{1},\dots,j_{p+1}}m(\dot{{\text{\boldmath$X$}}}_{i}/A_{n})\prod_{\ell=1}^{p+1}{X_{i,j_{\ell}}\over
A_{n,j_{\ell}}},$

where $\dot{\bm{X}}_{i}=\bm{z}+\theta_{i}\bm{X}_{i}$ for some
$\theta_{i}\in[0,1)$. Then we have

$\displaystyle\widehat{\bm{\beta}}(\bm{0})-\bm{M}(\bm{0})$
$\displaystyle=(\bm{X}\bm{W}\bm{X}^{\prime})^{-1}\bm{X}\bm{W}(\bm{Y}-\bm{X}^{\prime}\bm{M}(\bm{0}))$
$\displaystyle=\left[\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{i}\end{array}\right)(1\
\check{\bm{X}}^{\prime}_{i})\right]^{-1}\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)\left(\begin{array}[]{c}1\\\
\check{{\text{\boldmath$X$}}}_{i}\end{array}\right)$
$\displaystyle\quad\times\left(e_{n,i}+\varepsilon_{n,i}+\sum_{1\leq
j_{1}\leq\dots\leq j_{p+1}\leq d}{1\over\bm{s}_{j_{1}\dots
j_{p+1}}!}\partial_{j_{1},\dots,j_{p+1}}m(\dot{\bm{X}}_{i}/A_{n})\prod_{\ell=1}^{p+1}{X_{i,j_{\ell}}\over
A_{n,j_{\ell}}}\right).$

This yields

$\displaystyle\sqrt{A_{n}h_{1}\dots
h_{d}}H(\widehat{\bm{\beta}}(\bm{0})-\bm{M}(\bm{0}))$
$\displaystyle=S_{n}^{-1}(\bm{0})(V_{n}(\bm{0})+B_{n}(\bm{0})),$

where

$\displaystyle S_{n}(\bm{0})$ $\displaystyle={1\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{i}\end{array}\right)(1\ \check{\bm{X}}^{\prime}_{i})H^{-1},$
$\displaystyle V_{n}(\bm{0})$ $\displaystyle={\sqrt{A_{n}h_{1}\dots
h_{d}}\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{i}\end{array}\right)(e_{n,i}+\varepsilon_{n,i})$
$\displaystyle=:(V_{n,j_{1}\dots j_{L}}(\bm{0}))^{\prime}_{1\leq
j_{1}\leq\dots\leq j_{L}\leq d,0\leq L\leq p},$ $\displaystyle B_{n}(\bm{0})$
$\displaystyle={\sqrt{A_{n}h_{1}\dots h_{d}}\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{i}\end{array}\right)$ $\displaystyle\quad\times\sum_{1\leq
j_{1}\leq\dots\leq j_{p+1}\leq d}{1\over\bm{s}_{j_{1}\dots
j_{p+1}}!}\partial_{j_{1},\dots,j_{p+1}}m(\dot{\bm{X}}_{i}/A_{n})\prod_{\ell=1}^{p+1}{X_{i,j_{\ell}}\over
A_{n,j_{\ell}}}$ $\displaystyle=:(B_{n,j_{1}\dots
j_{L}}(\dot{\bm{X}}))^{\prime}_{1\leq j_{1}\leq\dots\leq j_{L}\leq d,0\leq
L\leq p}.$

(Step 1) Now we evaluate $S_{n}(\bm{0})$. By a change of variables and the
dominated convergence theorem, we have

$\displaystyle E[S_{n}(\bm{0})]$ $\displaystyle={A_{n}^{-1}\over h_{1}\dots
h_{d}}\int K_{Ah}(\bm{x})H^{-1}\left(\begin{array}[]{c}1\\\
\check{(\bm{x}/A_{n})}\end{array}\right)(1\
\check{(\bm{x}/A_{n})}^{\prime})H^{-1}g(\bm{x}/A_{n})d\bm{x}$
$\displaystyle={A_{n}^{-1}\over h_{1}\dots h_{d}}A_{n}h_{1}\dots h_{d}\int
K(\bm{w})\left(\begin{array}[]{c}1\\\ \check{\bm{w}}\end{array}\right)(1\
\check{\bm{w}}^{\prime})g(\bm{w}\circ\bm{h})d\bm{w}$
$\displaystyle=\left(g(\bm{0})\int K(\bm{w})\left(\begin{array}[]{c}1\\\
\check{\bm{w}}\end{array}\right)(1\
\check{\bm{w}}^{\prime})d\bm{w}\right)(1+o(1)).$

For $1\leq j_{1,1}\leq\dots\leq j_{1,L_{1}}\leq d$, $1\leq
j_{2,1}\leq\dots\leq j_{2,L_{2}}\leq d,0\leq L_{1},L_{2}\leq p$, we define

$\displaystyle I_{n,j_{1,1}\dots j_{1,L_{1}},j_{2,1}\dots j_{2,L_{2}}}$
$\displaystyle:={1\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)\prod_{\ell_{1}=1}^{L_{1}}\left({X_{i,j_{1,\ell_{1}}}\over
A_{n,j_{1,\ell_{1}}}h_{j_{1,\ell_{1}}}}\right)\prod_{\ell_{2}=1}^{L_{2}}\left({X_{i,j_{2,\ell_{2}}}\over
A_{n,j_{2,\ell_{2}}}h_{j_{2,\ell_{2}}}}\right).$

Then, by a change of variables and the dominated convergence theorem, we have

$\displaystyle\operatorname{Var}(I_{n,j_{1,1}\dots j_{1,L_{1}},j_{2,1}\dots
j_{2,L_{2}}})$ $\displaystyle\quad={1\over n(h_{1}\dots
h_{d})^{2}}\operatorname{Var}\left(K_{Ah}\left(\bm{X}_{1}\right)\prod_{\ell_{1}=1}^{L_{1}}\left({X_{i,j_{1,\ell_{1}}}\over
A_{n,j_{1,\ell_{1}}}h_{j_{1,\ell_{1}}}}\right)\prod_{\ell_{2}=1}^{L_{2}}\left({X_{i,j_{2,\ell_{2}}}\over
A_{n,j_{2,\ell_{2}}}h_{j_{2,\ell_{2}}}}\right)\right)$
$\displaystyle\quad={1\over nh_{1}\dots
h_{d}}\left\\{\int\prod_{\ell_{1}=1}^{L_{1}}z_{j_{1,\ell_{1}}}^{2}\prod_{\ell_{2}=1}^{L_{2}}z_{j_{2,\ell_{2}}}^{2}K^{2}(\bm{z})g(\bm{z}\circ\bm{h})d\bm{z}\right.$
$\displaystyle\left.\quad\quad-h_{1}\dots
h_{d}\left(\int\prod_{\ell_{1}=1}^{L_{1}}z_{j_{1,\ell_{1}}}\prod_{\ell_{2}=1}^{L_{2}}z_{j_{2,\ell_{2}}}K(\bm{z})g(\bm{z}\circ\bm{h})d\bm{z}\right)^{2}\right\\}$
$\displaystyle\quad={1\over nh_{1}\dots
h_{d}}\left(g(\bm{0})\kappa_{j_{1,1}\dots j_{1,L_{1}}j_{2,1}\dots
j_{2,L_{2}}j_{1,1}\dots j_{1,L_{1}}j_{2,1}\dots
j_{2,L_{2}}}^{(2)}+o(1)\right)$ $\displaystyle\quad\quad-{1\over
n}(g(\bm{0})\kappa_{j_{1,1}\dots j_{1,L_{1}}j_{2,1}\dots
j_{2,L_{2}}}^{(1)}+o(1))^{2}$
$\displaystyle\quad={g(\bm{0})\kappa_{j_{1,1}\dots j_{1,L_{1}}j_{2,1}\dots
j_{2,L_{2}}j_{1,1}\dots j_{1,L_{1}}j_{2,1}\dots j_{2,L_{2}}}^{(2)}\over
nh_{1}\dots h_{d}}+o\left({1\over nh_{1}\dots h_{d}}\right).$

Then for any $\rho>0$,

$\displaystyle P\left(|I_{n,j_{1,1}\dots j_{1,L_{1}},j_{2,1}\dots
j_{2,L_{2}}}-g(\bm{0})\kappa_{j_{1,1}\dots j_{1,L_{1}}j_{2,1}\dots
j_{2,L_{2}}}^{(1)}|>\rho\right)$
$\displaystyle\quad\leq\rho^{-1}\left\\{\operatorname{Var}(I_{n,j_{1,1}\dots
j_{1,L_{1}},j_{2,1}\dots j_{2,L_{2}}})+\left(E[I_{n,j_{1,1}\dots
j_{1,L_{1}},j_{2,1}\dots j_{2,L_{2}}}]-g(\bm{0})\kappa_{j_{1,1}\dots
j_{1,L_{1}}j_{2,1}\dots j_{2,L_{2}}}^{(1)}\right)^{2}\right\\}$
$\displaystyle=O\left({1\over nh_{1}\dots h_{d}}\right)+o(1)=o(1).$

This yields $I_{n,j_{1,1}\dots j_{1,L_{1}},j_{2,1}\dots
j_{2,L_{2}}}\stackrel{{\scriptstyle p}}{{\to}}g(\bm{0})\kappa_{j_{1,1}\dots
j_{1,L_{1}}j_{2,1}\dots j_{2,L_{2}}}^{(1)}$. Hence we have

$\displaystyle S_{n}(\bm{0})$ $\displaystyle\stackrel{{\scriptstyle
p}}{{\to}}g(\bm{0})S.$

(Step 2) Now we evaluate $V_{n}(\bm{0})$. For any
$\bm{t}=(t_{0},t_{1},\dots,t_{d},t_{11},\dots,t_{dd},\dots,t_{1\dots
1},\dots,t_{d\dots d})^{\prime}\in\mathbb{R}^{D}$, we define

$\displaystyle\widetilde{V}_{n}(\bm{0})$ $\displaystyle:={nh_{1}\dots
h_{d}\over\sqrt{A_{n}h_{1}\dots
h_{d}}}\bm{t}^{\prime}V_{n}(\bm{0})=\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{i}\end{array}\right)\right](e_{n,i}+\varepsilon_{n,i}).$

In this step, we will show that

$\displaystyle\bm{t}^{\prime}V_{n}(\bm{0})$
$\displaystyle\stackrel{{\scriptstyle
d}}{{\to}}N\left(\bm{0},g(\bm{0})\left\\{\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))+\eta^{2}(\bm{0})g(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right\\}\int
K^{2}(\bm{z})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]^{2}d\bm{z}\right).$ (A.3)

Before we show (A.3), we introduce some notations. For
$\bm{z}_{0}=(z_{0,1},\dots z_{0,d})^{\prime}\in\mathbb{R}^{d}$ and
$\bm{\ell}=(\ell_{1},\dots,\ell_{d})^{\prime}\in\mathbb{Z}^{d}$, let

$\Gamma_{n,\bm{z}_{0}}(\bm{\ell};\bm{0})=\prod_{j=1}^{d}(A_{n,j}z_{0,j}+(\ell_{j}-1/2)A_{n3,j},A_{n,j}z_{0,j}+(\ell_{j}+1/2)A_{n3,j}]$

with $A_{n3,j}=A_{n1,j}+A_{n2,j}$, and define the following hypercubes,

$\Gamma_{n,\bm{z}_{0}}(\bm{\ell};\bm{\Delta})=\prod_{j=1}^{d}I_{j,\bm{z}_{0}}(\Delta_{j}),\
\bm{\Delta}=(\Delta_{1},\dots,\Delta_{d})^{\prime}\in\\{1,2\\}^{d},$

where

$\displaystyle I_{j,\bm{z}_{0}}(\Delta_{j})$
$\displaystyle=\begin{cases}(A_{n,j}z_{0,j}+(\ell_{j}-1/2)A_{n3,j},A_{n,j}z_{0,j}+(\ell_{j}-1/2)A_{n3,j}+A_{n1,j}]&\text{if
$\Delta_{j}=1$},\\\
(A_{n,j}z_{0,j}+(\ell_{j}-1/2)A_{n3,j}+A_{n1,j},A_{n,j}z_{0,j}+(\ell_{j}+1/2)A_{n3,j}]&\text{if
$\Delta_{j}=2$}.\end{cases}$

Let $\bm{\Delta}_{0}=(1,\dots,1)^{\prime}$. The partitions
$\Gamma_{n,\bm{z}_{0}}(\bm{\ell};\bm{\Delta}_{0})$ correspond to “large
blocks” and the partitions $\Gamma_{n,\bm{z}_{0}}(\bm{\ell};\bm{\Delta})$ for
$\bm{\Delta}\neq\bm{\Delta}_{0}$ correspond to “small blocks”. Let
$L_{n1}(\bm{z}_{0})=\\{\bm{\ell}\in\mathbb{Z}^{d}:\Gamma_{n,\bm{z}_{0}}(\bm{\ell};\bm{0})\subset
R_{n}\cap(\bm{h}R_{n}+A_{n}\bm{z}_{0})\\}$ denote the index set of all
hypercubes $\Gamma_{n,\bm{z}_{0}}(\bm{\ell};\bm{0})$ that are contained in
$R_{n}\cap(\bm{h}R_{n}+A_{n}\bm{z}_{0})$, and let
$L_{n2}(\bm{z}_{0})=\\{\bm{\ell}\in\mathbb{Z}^{d}:\Gamma_{n,\bm{z}_{0}}(\bm{\ell};\bm{0})\cap
R_{n}\cap(\bm{h}R_{n}+A_{n}\bm{z}_{0})\neq
0,\Gamma_{n}(\bm{\ell};\bm{0})\cap(R_{n}\cap(\bm{h}R_{n}+A_{n}\bm{z}_{0}))^{c}\neq\emptyset\\}$
be the index set of boundary hypercubes. Define
$\Gamma_{n}(\bm{\ell};\bm{\Delta})=\Gamma_{n,\bm{0}}(\bm{\ell};\bm{\Delta})$,
$L_{n1}=L_{n1}(\bm{0})$, $L_{n2}=L_{n2}(\bm{0})$, and

$\widetilde{V}_{n}(\bm{\ell};\bm{\Delta})=\sum_{i:\bm{X}_{i}\in\Gamma_{n}(\bm{\ell};\bm{\Delta})\cap\bm{h}R_{n}}K_{Ah}\left(\bm{X}_{i}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{i}\end{array}\right)\right](e_{n,i}+\varepsilon_{n,i}).$

Note that by our summation convention, ${V}_{n}(\bm{\ell};\bm{\Delta})=0$ if
the set
$\\{i:\bm{X}_{i}\in\Gamma_{n}(\bm{\ell};\bm{\Delta})\cap\bm{h}R_{n}\\}$ is
empty for some $\bm{\ell}$. Then we have

$\displaystyle\widetilde{V}_{n}(\bm{0})$ $\displaystyle=\sum_{\bm{\ell}\in
L_{n1}}\widetilde{V}_{n}(\bm{\ell};\bm{\Delta}_{0})+\sum_{\bm{\Delta}\neq\bm{\Delta}_{0}}\sum_{\bm{\ell}\in
L_{n1}}\widetilde{V}_{n}(\bm{\ell};\bm{\Delta})+\sum_{\bm{\Delta}\in\\{1,2\\}^{d}}\sum_{\bm{\ell}\in
L_{n2}}\widetilde{V}_{n}(\bm{\ell};\bm{\Delta})$
$\displaystyle=:\widetilde{V}_{n1}+\widetilde{V}_{n2}+\widetilde{V}_{n3}.$

Note that for $\bm{\ell}_{1},\bm{\ell}_{2}\in L_{n1}$,

$\displaystyle
d\left(\Gamma_{n}(\bm{\ell}_{1};\bm{\Delta}_{0}),\Gamma_{n}(\bm{\ell}_{2};\bm{\Delta}_{0})\right)$
$\displaystyle\geq\min\\{|\bm{\ell}_{1}-\bm{\ell}_{2}|,0\\}\underline{A}_{n3}+\underline{A}_{n2},$
(A.4)

where $\underline{A}_{n3}=\min_{1\leq j\leq d}A_{n3,j}$ and
$\underline{A}_{n2}=\min_{1\leq j\leq d}A_{n2,j}$.

Hence, by the Volkonskii-Rozanov inequality (cf. Proposition 2.6 in Fan and
Yao, (2003)), we have

$\displaystyle\left|E[\exp(\mathrm{i}u\widetilde{V}_{n1})]-\prod_{\ell\in
L_{n1}}E[\exp(\mathrm{i}u\widetilde{V}_{n}(\bm{\ell};\bm{\Delta}_{0}))]\right|$
$\displaystyle\lesssim\left({A_{n}h_{1}\dots h_{d}\over
A_{n}^{(1)}}\right)\alpha(\underline{A}_{n2};A_{n}h_{1}\dots h_{d}).$ (A.5)

From Lyapounov’s CLT, it is sufficient to verify the following conditions to
show (A.3): As $n\to\infty$,

$\displaystyle{A_{n}\over n^{2}h_{1}\dots
h_{d}}E[\widetilde{V}_{n}^{2}(\bm{0})]$ $\displaystyle\to
g(\bm{0})\left\\{\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))+\eta^{2}(\bm{0})g(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right\\}$
$\displaystyle\quad\quad\times\int
K^{2}(\bm{z})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]^{2}d\bm{z},$ (A.8)
$\displaystyle\sum_{\bm{\ell}\in
L_{n1}}E[\widetilde{V}_{n}^{2}(\bm{\ell};\bm{\Delta}_{0})]-E[\widetilde{V}_{n}^{2}(\bm{0})]$
$\displaystyle=o\left(n^{2}A_{n}^{-1}h_{1}\dots h_{d}\right),$ (A.9)
$\displaystyle\sum_{\bm{\ell}\in
L_{n1}}E[\widetilde{V}_{n}^{4}(\bm{\ell};\bm{\Delta}_{0})]$
$\displaystyle=o\left(\left(n^{2}A_{n}^{-1}h_{1}\dots
h_{d}\right)^{2}\right),$ (A.10)
$\displaystyle\operatorname{Var}(\widetilde{V}_{n2})$
$\displaystyle=o\left(n^{2}A_{n}^{-1}h_{1}\dots h_{d}\right),$ (A.11)
$\displaystyle\operatorname{Var}(\widetilde{V}_{n3})$
$\displaystyle=o\left(n^{2}A_{n}^{-1}h_{1}\dots h_{d}\right).$ (A.12)

In the following steps, we show (A.8) (Step 2-1), (A.10) (Step 2-2), (A.11)
and (A.12) (Step 2-3), and (A.9) (Step 2-4).

(Step 2-1) Now we show (A.8). Let $\delta_{ij}$ be a function such that
$\delta_{ij}=1$ if $i=j$ and $\delta_{ij}=0$ if $i\neq j$. Observe that

$\displaystyle\sigma_{n}^{2}(\bm{0})$
$\displaystyle:=E_{\cdot\mid\bm{X}}\left(\widetilde{V}_{n}^{2}(\bm{0})\right)$
$\displaystyle=\sum_{i,j=1}^{n}\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}_{i}}\end{array}\right)\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{j}\end{array}\right)K_{Ah}\left(\bm{X}_{i}\right)K_{Ah}\left(\bm{X}_{j}\right)$
$\displaystyle\quad\times\left\\{\eta(\bm{X}_{i}/A_{n})\eta(\bm{X}_{j}/A_{n})\sigma_{\bm{e}}(\bm{X}_{i}-\bm{X}_{j})+\sigma_{\varepsilon}^{2}(\bm{X}_{i}/A_{n})\delta_{ij}\right\\}.$

Thus we have

$\displaystyle E_{\bm{X}}\left[\sigma_{n}^{2}(\bm{0})\right]$
$\displaystyle\quad=nA_{n}^{-1}\int\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{(\bm{x}/A_{n})}\end{array}\right)\right]^{2}K^{2}_{Ah}(\bm{x})\left\\{\eta^{2}(\bm{x}/A_{n})+\sigma_{\varepsilon}^{2}(\bm{x}/A_{n})\right\\}g(\bm{x}/A_{n})d\bm{x}$
$\displaystyle\quad\quad+n(n-1)A_{n}^{-2}\int\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{(\bm{x}_{1}/A_{n})}\end{array}\right)\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{(\bm{x}_{2}/A_{n})}\end{array}\right)K_{Ah}(\bm{x}_{1})K_{Ah}(\bm{x}_{2})$
$\displaystyle\quad\quad\quad\times\eta(\bm{x}_{1}/A_{n})\eta(\bm{x}_{2}/A_{n})\sigma_{\bm{e}}(\bm{x}_{1}-\bm{x}_{2})g(\bm{x}_{1}/A_{n})g(\bm{x}_{2}/A_{n})d\bm{x}_{1}d\bm{x}_{2}$
$\displaystyle\quad=:\sigma^{2}_{n,1}+\sigma^{2}_{n,2}.$

For $\sigma^{2}_{n,1}$, we have

$\displaystyle\sigma^{2}_{n,1}$ $\displaystyle=nh_{1}\dots
h_{d}\int\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]^{2}K^{2}(\bm{z})\left\\{\eta^{2}(\bm{z}\circ\bm{h})+\sigma^{2}_{\varepsilon}(\bm{z}\circ\bm{h})\right\\}g(\bm{z}\circ\bm{h})d\bm{z}$
(A.15) $\displaystyle=nh_{1}\dots
h_{d}(\eta^{2}(\bm{0})+\sigma^{2}_{\varepsilon}(\bm{0}))g(\bm{0})\left(\int
K^{2}(\bm{z})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]^{2}d\bm{z}\right)(1+o(1)).$ (A.18)

For $\sigma^{2}_{n,2}$, we have

$\displaystyle\sigma^{2}_{n,2}$
$\displaystyle=n(n-1)\int_{R_{0}^{2}}\sigma_{\bm{e}}(A_{n}(\bm{y}_{1}-\bm{y}_{2}))\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{y}}_{1}\end{array}\right)\right]\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{y}}_{2}\end{array}\right)\right]$ $\displaystyle\quad\quad\times
K_{h}(\bm{y}_{1})K_{h}(\bm{y}_{2})\eta(\bm{y}_{1})\eta(\bm{y}_{2})g(\bm{y}_{1})g(\bm{y}_{2})d\bm{y}_{1}d\bm{y}_{2}$
$\displaystyle=n(n-1)(h_{1}\dots
h_{d})^{2}\int_{\bm{h}^{-1}R_{0}^{2}}\sigma_{\bm{e}}(A_{n}(\bm{z}_{1}-\bm{z}_{2})\circ\bm{h})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}_{1}\end{array}\right)\right]\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}_{2}\end{array}\right)\right]$ $\displaystyle\quad\quad\times
K(\bm{z}_{1})K(\bm{z}_{2})\eta(\bm{z}_{1}\circ\bm{h})\eta(\bm{z}_{2}\circ\bm{h})g(\bm{z}_{1}\circ\bm{h})g(\bm{z}_{2}\circ\bm{h})d\bm{z}_{1}d\bm{z}_{2}$
$\displaystyle=n(n-1)(h_{1}\dots
h_{d})^{2}\int_{R^{\prime}_{\bm{h},0}}\sigma_{\bm{e}}(A_{n}\bm{w}\circ\bm{h})\left(\int_{R_{\bm{h},0}(\bm{w})}\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{(\bm{z}_{2}+\bm{w})}\end{array}\right)\right]\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}_{2}\end{array}\right)\right]\right.$
$\displaystyle\left.\quad\quad\times
K(\bm{z}_{2}+\bm{w})K(\bm{z}_{2})\eta((\bm{z}_{2}+\bm{w})\circ\bm{h})\eta(\bm{z}_{2}\circ\bm{h})g((\bm{z}_{2}+\bm{w})\circ\bm{h})g(\bm{z}_{2}\circ\bm{h})d\bm{z}_{2}\right)d\bm{w}$
$\displaystyle=n(n-1)h_{1}\dots
h_{d}\int_{\bm{h}R^{\prime}_{\bm{h},0}}\sigma_{\bm{e}}(A_{n}\bm{u})\left(\int_{R_{\bm{h},0}(\bm{u}/\bm{h})}\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{(\bm{z}_{2}+\bm{u}\circ\bm{h}^{-1})}\end{array}\right)\right]\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}_{2}\end{array}\right)\right]\right.$
$\displaystyle\left.\quad\quad\times
K(\bm{z}_{2}+\bm{u}\circ\bm{h}^{-1})K(\bm{z}_{2})\eta(\bm{z}_{2}\circ\bm{h}+\bm{u})\eta(\bm{z}_{2}\circ\bm{h})g((\bm{z}_{2}\circ\bm{h}+\bm{u})g(\bm{z}_{2}\circ\bm{h})d\bm{z}_{2}\right)d\bm{u}$
$\displaystyle=n(n-1)A_{n}^{-1}h_{1}\dots
h_{d}\int_{A_{n}\bm{h}R^{\prime}_{\bm{h},0}}\sigma_{\bm{e}}(\bm{v})\left(\int_{R_{\bm{h},0}((\bm{v}\circ\bm{h}^{-1})/A_{n})}\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\left(\bm{z}_{2}+{\bm{v}\circ\bm{h}^{-1}\over
A_{n}}\right)}\end{array}\right)\right]\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}_{2}\end{array}\right)\right]\right.$
$\displaystyle\left.\quad\quad\times
K\left(\bm{z}_{2}+{\bm{v}\circ\bm{h}^{-1}\over
A_{n}}\right)K(\bm{z}_{2})\eta\left(\bm{z}_{2}\circ\bm{h}+{\bm{v}\over
A_{n}}\right)\eta(\bm{z}_{2}\circ\bm{h})g\left(\bm{z}_{2}\circ\bm{h}+{\bm{v}\over
A_{n}}\right)g(\bm{z}_{2}\circ\bm{h})d\bm{z}_{2}\right)d\bm{v}$

where

$\displaystyle R_{\bm{h},0}^{\prime}$
$\displaystyle=\\{\bm{w}=\bm{z}_{1}-\bm{z}_{2}:\bm{z}_{1},\bm{z}_{2}\in\bm{h}^{-1}R_{0}\\},\
R_{\bm{h},0}(\bm{w})=\\{\bm{z}_{2}:\bm{z}_{2}\in\bm{h}^{-1}R_{0}\cap(\bm{h}^{-1}R_{0}+\bm{w})\\},$
$\displaystyle A_{n}\bm{h}R_{\bm{h},0}^{\prime}$
$\displaystyle=\\{(A_{n,1}x_{1},\dots,A_{n,d}x_{d}):\bm{x}=(x_{1},\dots,x_{d})^{\prime}\in\bm{h}R_{\bm{h},0}^{\prime}\\}.$

We divide the integral $\int_{A_{n}\bm{h}R_{\bm{h},0}^{\prime}}$ into two
parts $\int_{A_{n}\bm{h}R_{\bm{h},0}^{\prime}\cap\\{|\bm{v}|\leq M\\}}$ and
$\int_{A_{n}\bm{h}R_{\bm{h},0}^{\prime}\cap\\{|\bm{v}|>M\\}}$ for some $M>0$
and define these as $\sigma_{n,21}^{2}$ and $\sigma_{n,22}^{2}$, respectively.
Observe that as $n\to\infty$

$\displaystyle|\sigma_{n,22}^{2}|$
$\displaystyle\lesssim\int_{\\{|\bm{v}|>M\\}}|\sigma_{\bm{e}}(\bm{v})|d\bm{v}$

which can be made arbitrary small by choosing a large $M$. Further, observe
that as $n\to\infty$

$\displaystyle 1\\{A_{n}\bm{h}R_{\bm{h},0}^{\prime}\cap\\{|\bm{v}|\leq
M\\}\\}\int_{R_{\bm{h},0}(\bm{v}/(A_{n}\bm{h}))}\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\left(\bm{z}_{2}+{\bm{v}\circ\bm{h}^{-1}\over
A_{n}}\right)}\end{array}\right)\right]\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}_{2}\end{array}\right)\right]$
$\displaystyle\quad\quad\quad\times
K\left(\bm{z}_{2}+{\bm{v}\circ\bm{h}^{-1}\over
A_{n}}\right)K(\bm{z}_{2})\eta\left(\bm{z}_{2}\circ\bm{h}+{\bm{v}\over
A_{n}}\right)\eta(\bm{z}_{2}\circ\bm{h})g\left(\bm{z}_{2}\circ\bm{h}+{\bm{v}\over
A_{n}}\right)g(\bm{z}_{2}\circ\bm{h})d\bm{z}_{2}$
$\displaystyle=1\\{|\bm{v}|\leq M\\}\eta^{2}(\bm{0})g^{2}(\bm{0})\left(\int
K^{2}(\bm{z}_{2})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}_{2}\end{array}\right)\right]^{2}d\bm{z}_{2}\right)(1+o(1)).$

Then as $n\to\infty$, we have

$\displaystyle\sigma_{n,21}^{2}$
$\displaystyle=\eta^{2}(\bm{0})g^{2}(\bm{0})\left(\int_{\\{|\bm{v}|\leq
M\\}}\sigma_{\bm{e}}(\bm{v})d\bm{v}\right)\left(\int
K^{2}(\bm{z}_{2})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}_{2}\end{array}\right)\right]^{2}d\bm{z}_{2}\right)(1+o(1)).$

Therefore, we have

$\displaystyle\sigma^{2}_{n,2}$ $\displaystyle=n^{2}A_{n}^{-1}h_{1}\dots
h_{d}\eta^{2}(\bm{0})g^{2}(\bm{0})\left(\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right)\left(\int
K^{2}(\bm{z})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]^{2}d\bm{z}\right)(1+o(1)).$ (A.21)

By (A.15) and (A.21), we have

$\displaystyle\operatorname{Var}(\bm{t}^{\prime}V_{n}(\bm{0}))$
$\displaystyle=g(\bm{0})\left\\{\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))+\eta^{2}(\bm{0})g(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right\\}\left(\int
K^{2}(\bm{z})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]^{2}d\bm{z}\right)(1+o(1)).$

(Step 2-2) Now we show (A.10). Define
$I_{n}(\bm{\ell})=\\{\bm{i}\in\mathbb{Z}^{d}:\bm{i}+(-1/2,1/2]^{d}\subset\Gamma_{n}(\bm{\ell};\bm{\Delta}_{0})\\}$
for $\bm{\ell}\in L_{n1}$ and

$\widetilde{V}_{n}(\bm{i})=\sum_{i=1}^{n}K_{Ah}\left(\bm{X}_{i}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{i}\end{array}\right)\right](e_{n,i}+\varepsilon_{n,i})1\\{\bm{X}_{i}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}.$

Observe that

$\displaystyle E[\widetilde{V}_{n}^{4}(\bm{\ell};\bm{\Delta}_{0})]$
$\displaystyle=E\left[\left(\sum_{\bm{i}\in
I_{n}(\bm{\ell})}\widetilde{V}_{n}(\bm{i})\right)^{4}\right]$
$\displaystyle=\sum_{\bm{i}\in
I_{n}(\bm{\ell})}E\left[\widetilde{V}_{n}^{4}(\bm{i})\right]+\sum_{\bm{i},\bm{j}\in
I_{n}(\bm{\ell}),\bm{i}\neq\bm{j}}E\left[\widetilde{V}_{n}^{3}(\bm{i})\widetilde{V}_{n}(\bm{j})\right]+\sum_{\bm{i},\bm{j}\in
I_{n}(\bm{\ell}),\bm{i}\neq\bm{j}}E\left[\widetilde{V}_{n}^{2}(\bm{i})\widetilde{V}_{n}^{2}(\bm{j})\right]$
$\displaystyle\quad\quad+\sum_{\bm{i},\bm{j},\bm{k}\in
I_{n}(\bm{\ell}),\bm{i}\neq\bm{j}\neq\bm{k}}E\left[\widetilde{V}_{n}^{2}(\bm{i})\widetilde{V}_{n}(\bm{j})\widetilde{V}_{n}(\bm{k})\right]+\sum_{\bm{i},\bm{j},\bm{k},\bm{p}\in
I_{n}(\bm{\ell}),\bm{i}\neq\bm{j}\neq\bm{k}\neq\bm{p}}E\left[\widetilde{V}_{n}(\bm{i})\widetilde{V}_{n}(\bm{j})\widetilde{V}_{n}(\bm{k})\widetilde{V}_{n}(\bm{p})\right]$
$\displaystyle=:Q_{n1}+Q_{n2}+Q_{n3}+Q_{n4}+Q_{n5}.$

For $Q_{n1}$, we have

$\displaystyle E[\widetilde{V}_{n}^{4}(\bm{i})]$
$\displaystyle=E_{\bm{X}}[E_{\cdot\mid\bm{X}}[\widetilde{V}_{n}^{4}(\bm{i})]]$
$\displaystyle=\sum_{j_{1},j_{2},j_{3},j_{4}=1}^{n}E\left[\prod_{k=1}^{4}K_{Ah}\left(\bm{X}_{j_{k}}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{j_{k}}\end{array}\right)\right]1\\{\bm{X}_{j_{k}}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$ $\displaystyle\left.\quad\quad\quad\quad\times
E_{\cdot\mid\bm{X}}[e_{n,j_{k}}+\varepsilon_{n,j_{k}}]\right]$
$\displaystyle\lesssim\sum_{j_{1},j_{2},j_{3},j_{4}=1}^{n}\\!\\!\\!\\!\\!E\left[\prod_{k=1}^{4}\left|K_{Ah}\left(\bm{X}_{j_{k}}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{j_{k}}\end{array}\right)\right]\right|1\\{\bm{X}_{j_{k}}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\eta(\bm{X}_{j_{k}}/A_{n})\right]$
$\displaystyle\quad+\sum_{j_{1},j_{2},j_{3},j_{4}=1}^{n}\\!\\!\\!\\!\\!\\!E\left[\prod_{k=1}^{4}\left|K_{Ah}\left(\bm{X}_{j_{k}}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{j_{k}}\end{array}\right)\right]\right|1\\{\bm{X}_{j_{k}}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\sigma_{\varepsilon}(\bm{X}_{j_{k}}/A_{n})\right]$
$\displaystyle=:Q_{n11}+Q_{n12}.$

For $Q_{n11}$, we have

$\displaystyle Q_{n11}$ $\displaystyle\lesssim
nE\left[\left|K_{Ah}\left(\bm{X}_{1}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1}\end{array}\right)\right]\right|^{4}1\\{\bm{X}_{1}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\eta^{4}(\bm{X}_{1}/A_{n})\right]$
$\displaystyle\quad+n^{2}E\left[\left|K_{Ah}\left(\bm{X}_{1}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1}\end{array}\right)\right]\right|^{3}1\\{\bm{X}_{1}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\left|K_{Ah}\left(\bm{X}_{2}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{2}\end{array}\right)\right]\right|1\\{\bm{X}_{2}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\eta^{3}(\bm{X}_{1}/A_{n})\eta(\bm{X}_{2}/A_{n})\right]$
$\displaystyle\quad+n^{2}E\left[\left|K_{Ah}\left(\bm{X}_{1}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1}\end{array}\right)\right]\right|^{2}1\\{\bm{X}_{1}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\left|K_{Ah}\left(\bm{X}_{2}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{2}\end{array}\right)\right]\right|^{2}1\\{\bm{X}_{2}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\eta^{2}(\bm{X}_{1}/A_{n})\eta^{2}(\bm{X}_{2}/A_{n})\right]$
$\displaystyle\quad+n^{3}E\left[\left|K_{Ah}\left(\bm{X}_{1}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1}\end{array}\right)\right]\right|^{2}1\\{\bm{X}_{1}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\left|K_{Ah}\left(\bm{X}_{2}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{2}\end{array}\right)\right]\right|1\\{\bm{X}_{2}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\left|K_{Ah}\left(\bm{X}_{3}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{3}\end{array}\right)\right]\right|1\\{\bm{X}_{3}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\eta^{2}(\bm{X}_{1}/A_{n})\eta(\bm{X}_{2}/A_{n})\eta(\bm{X}_{3}/A_{n})\right]$
$\displaystyle\quad+n^{4}E\left[\left|K_{Ah}\left(\bm{X}_{1}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1}\end{array}\right)\right]\right|1\\{\bm{X}_{1}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\left|K_{Ah}\left(\bm{X}_{2}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{2}\end{array}\right)\right]\right|1\\{\bm{X}_{2}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\left|K_{Ah}\left(\bm{X}_{3}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{3}\end{array}\right)\right]\right|1\\{\bm{X}_{3}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\left|K_{Ah}\left(\bm{X}_{4}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{4}\end{array}\right)\right]\right|1\\{\bm{X}_{4}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}\right.$
$\displaystyle\left.\quad\times\eta^{2}(\bm{X}_{1}/A_{n})\eta(\bm{X}_{2}/A_{n})\eta(\bm{X}_{3}/A_{n})\eta(\bm{X}_{4}/A_{n})\right]$
$\displaystyle=:Q_{n111}+Q_{n112}+Q_{n113}+Q_{n114}.$

For $Q_{n111}$, we have

$\displaystyle Q_{n111}$
$\displaystyle=nA_{n}^{-1}\int\left|K_{Ah}(\bm{x})\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{(\bm{x}/A_{n})}\end{array}\right)\right]\right|^{4}1\\{\bm{x}\in[\bm{i}+(-1/2,1/2]^{d}]\cap
R_{n}\\}$
$\displaystyle\quad\times\eta^{4}(\bm{x}/A_{n})g(\bm{x}/A_{n})d\bm{x}$
$\displaystyle=nA_{n}^{-1}A_{n}h_{1}\dots
h_{d}\int\left|K(\bm{z})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]\right|^{4}1\\{\bm{z}\circ\bm{h}\in[\bm{i}+(-1/2,1/2]^{d}]/A_{n}\cap[-1/2,1/2]^{d}\\}$
$\displaystyle\quad\times\eta^{4}(\bm{z}\circ\bm{h})g(\bm{z}\circ\bm{h})d\bm{z}$
$\displaystyle=O\left(nA_{n}^{-1}\right).$

Likewise, $Q_{n112}=O(n^{2}A_{n}^{-2})$, $Q_{n113}=O(n^{3}A_{n}^{-3})$, and
$Q_{n114}=O(n^{4}A_{n}^{-4})$. Then we have $Q_{n11}=O(n^{4}A_{n}^{-4})$. We
can also show that $Q_{n12}=O(n^{4}A_{n}^{-4})$. Therefore, we have

$\displaystyle Q_{n1}$
$\displaystyle\lesssim[\\![I_{n}(\bm{\ell})]\\!]n^{4}A_{n}^{-4}\lesssim
A_{n}^{(1)}(nA_{n}^{-1})^{4}.$ (A.22)

For $Q_{n2}$, by the $\alpha$-mixing property of $\bm{e}$ and Proposition 2.5
in Fan and Yao, (2003), we have

$\displaystyle Q_{n2}$
$\displaystyle\lesssim\sum_{k=1}^{\overline{A}_{n1}}\sum_{\bm{i},\bm{j}\in
I_{n}(\bm{\ell}),|\bm{i}-\bm{j}|=k}\alpha^{1-4/q}(\min\\{k-d,0\\};1)E[|\widetilde{V}_{n}(\bm{i})|^{q}]^{3/q}E[|\widetilde{V}_{n}(\bm{j})|^{q}]^{1/q}$
$\displaystyle\lesssim
A_{n}^{(1)}(nA_{n}^{-1})^{4}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{d-1}\alpha_{1}^{1-4/q}(k)\right).$
(A.23)

where $\overline{A}_{n1}=\max_{1\leq j\leq d}A_{n1,j}$. Likewise,

$\displaystyle Q_{n3}$ $\displaystyle\lesssim
A_{n}^{(1)}(nA_{n}^{-1})^{4}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{d-1}\alpha_{1}^{1-4/q}(k)\right).$
(A.24)

Now we evaluate $Q_{n4}$ and $Q_{n5}$. For distinct indices
$\bm{i},\bm{j},\bm{k},\bm{p}\in I_{n}(\bm{\ell})$, let

$\displaystyle d_{1}(\bm{i},\bm{j},\bm{k})$
$\displaystyle=\max\\{d(\\{\bm{i}\\},\\{\bm{j},\bm{k}\\}),d(\\{\bm{k}\\},\\{\bm{i},\bm{j}\\})\\},$
$\displaystyle d_{2}(\bm{i},\bm{j},\bm{k},\bm{p})$
$\displaystyle=\max\\{d(J,\\{\bm{i},\bm{j},\bm{k},\bm{p}\\}):J\subset\\{\bm{i},\bm{j},\bm{k},\bm{p}\\},[\\![J]\\!]=1\\},$
$\displaystyle d_{3}(\bm{i},\bm{j},\bm{k},\bm{p})$
$\displaystyle=\max\\{d(J,\\{\bm{i},\bm{j},\bm{k},\bm{p}\\}):J\subset\\{\bm{i},\bm{j},\bm{k},\bm{p}\\},[\\![J]\\!]=2\\}.$

Here, $d_{1}$ denotes the maximal gap in the set of integer-indices
$\\{\bm{i},\bm{j},\bm{k}\\}$ from either $\bm{j}$ or $\bm{k}$ which
corresponds to
$E\left[\widetilde{V}_{n}^{2}(\bm{i})\widetilde{V}_{n}(\bm{j})\widetilde{V}_{n}(\bm{k})\right]$.
Similarly, $d_{2}$ and $d_{3}$ are the maximal gap in the index set
$\\{\bm{i},\bm{j},\bm{k},\bm{p}\\}$ from any of its single index-subsets or
two-index subsets, respectively. Applying the argument in the proof of Lemma
4.1 of Lahiri, (1999), for any given values $1\leq
d_{01},d_{02},d_{03}<[\\![I_{n}(\bm{\ell})]\\!]$, we have

$\displaystyle[\\![\\{(\bm{i},\bm{j},\bm{k})\in
I_{n}^{3}(\bm{\ell}):\bm{i}\neq\bm{j}\neq\bm{k}\ \text{and}\
d_{1}(\bm{i},\bm{j},\bm{k})=d_{01}\\}]\\!]\lesssim
d_{01}^{2d-1}[\\![I_{n}(\bm{\ell})]\\!],$ (A.25)
$\displaystyle[\\![\\{(\bm{i},\bm{j},\bm{k},\bm{p})\in
I_{n}^{4}(\bm{\ell}):\bm{i}\neq\bm{j}\neq\bm{k}\neq\bm{p},\
d_{2}(\bm{i},\bm{j},\bm{k},\bm{p})=d_{02},\ \text{and}\
d_{3}(\bm{i},\bm{j},\bm{k},\bm{p})=d_{03}\\}]\\!]$
$\displaystyle\quad\lesssim(d_{02}+d_{03})^{3d-1}[\\![I_{n}(\bm{\ell})]\\!].$
(A.26)

For $Q_{n4}$, by (A.25) and applying the same argument to show (A.23), we have

$\displaystyle Q_{n4}$ $\displaystyle\lesssim
A_{n}^{(1)}\sum_{k=1}^{\overline{A}_{n1}}k^{2d-1}\alpha^{1-4/q}(\min\\{k-d,0\\};2)E[|\widetilde{V}_{n}(\bm{i})|^{q}]^{2/q}E[|\widetilde{V}_{n}(\bm{j})|^{q}]^{1/q}E[|\widetilde{V}_{n}(\bm{k})|^{q}]^{1/q}$
$\displaystyle\lesssim
A_{n}^{(1)}(nA_{n}^{-1})^{4}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{2d-1}\alpha_{1}^{1-4/q}(k)\right).$
(A.27)

Define

$\displaystyle I_{n1}(\bm{\ell})$
$\displaystyle=\\{(\bm{i},\bm{j},\bm{k},\bm{p})\in
I_{n}^{4}(\bm{\ell}):\bm{i}\neq\bm{j}\neq\bm{k}\neq\bm{p},\
d_{2}(\bm{i},\bm{j},\bm{k},\bm{p})\geq d_{3}(\bm{i},\bm{j},\bm{k},\bm{p})\\},$
$\displaystyle I_{n2}(\bm{\ell})$
$\displaystyle=\\{(\bm{i},\bm{j},\bm{k},\bm{p})\in
I_{n}^{4}(\bm{\ell}):\bm{i}\neq\bm{j}\neq\bm{k}\neq\bm{p},\
d_{2}(\bm{i},\bm{j},\bm{k},\bm{p})<d_{3}(\bm{i},\bm{j},\bm{k},\bm{p})\\}.$

For $Q_{n5}$, by (A.26) and applying the same argument to show (A.23), we have

$\displaystyle Q_{n5}$ $\displaystyle=\sum_{(\bm{i},\bm{j},\bm{k},\bm{p})\in
I_{n1}(\bm{\ell})}E\left[\widetilde{V}_{n}(\bm{i})\widetilde{V}_{n}(\bm{j})\widetilde{V}_{n}(\bm{k})\widetilde{V}_{n}(\bm{p})\right]+\sum_{(\bm{i},\bm{j},\bm{k},\bm{p})\in
I_{n2}(\bm{\ell})}E\left[\widetilde{V}_{n}(\bm{i})\widetilde{V}_{n}(\bm{j})\widetilde{V}_{n}(\bm{k})\widetilde{V}_{n}(\bm{p})\right]$
$\displaystyle\lesssim
A_{n}^{(1)}\sum_{k=1}^{\overline{A}_{n1}}k^{3d-1}\alpha^{1-4/q}(\min\\{k-d,0\\};3)$
$\displaystyle\quad\times
E[|\widetilde{V}_{n}(\bm{i})|^{q}]^{1/q}E[|\widetilde{V}_{n}(\bm{j})|^{q}]^{1/q}E[|\widetilde{V}_{n}(\bm{k})|^{q}]^{1/q}E[|\widetilde{V}_{n}(\bm{p})|^{q}]^{1/q}$
$\displaystyle\quad+\left(\sum_{\bm{i},\bm{j}\in
I_{n}(\bm{\ell}),\bm{i}\neq\bm{j}}\left|E[\widetilde{V}_{n}(\bm{i})\widetilde{V}_{n}(\bm{j})]\right|\right)^{2}$
$\displaystyle\quad+A_{n}^{(1)}\sum_{k=1}^{\overline{A}_{n1}}k^{3d-1}\alpha^{1-4/q}(\min\\{k-d,0\\};2)$
$\displaystyle\quad\times
E[|\widetilde{V}_{n}(\bm{i})|^{q}]^{1/q}E[|\widetilde{V}_{n}(\bm{j})|^{q}]^{1/q}E[|\widetilde{V}_{n}(\bm{k})|^{q}]^{1/q}E[|\widetilde{V}_{n}(\bm{p})|^{q}]^{1/q}$
$\displaystyle\lesssim(A_{n}^{(1)})^{2}(nA_{n}^{-1})^{4}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{2d-1}\alpha_{1}^{1-4/q}(k)\right).$
(A.28)

Combining (A.22)-(A.28), we have

$\displaystyle\sum_{\bm{\ell}\in
L_{n1}}E[\widetilde{V}_{n}^{4}(\bm{\ell};\bm{\Delta}_{0})]$
$\displaystyle=\sum_{\bm{\ell}\in L_{n1}}E\left[\left(\sum_{\bm{i}\in
I_{n}(\bm{\ell})}\widetilde{V}_{n}(\bm{i})\right)^{4}\right]$
$\displaystyle\lesssim[\\![L_{n1}]\\!](A_{n}^{(1)})^{2}(nA_{n}^{-1})^{4}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{2d-1}\alpha_{1}^{1-4/q}(k)\right)$
$\displaystyle\lesssim\left({A_{n}h_{1}\dots h_{d}\over
A_{n}^{(1)}}\right)(A_{n}^{(1)})^{2}(nA_{n}^{-1})^{4}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{2d-1}\alpha_{1}^{1-4/q}(k)\right)$
$\displaystyle=o\left((n^{2}A_{n}^{-1}h_{1}\dots h_{d})^{2}\right).$

(Step 2-3) Now we show (A.11) and (A.12). Define

$\displaystyle J_{n}$
$\displaystyle=\\{\bm{i}\in\mathbb{Z}^{d}:(\bm{i}+(-1/2,1/2]^{d})\cap\bm{h}R_{n}\neq\emptyset\\},$
$\displaystyle J_{n1}$ $\displaystyle=\cup_{\bm{\ell}\in
L_{n1}}I_{n}(\bm{\ell}),$ $\displaystyle J_{n2}$ $\displaystyle=\\{\bm{i}\in
J_{n}:\bm{i}+(-1/2,1/2]^{d}\subset\Gamma_{n}(\bm{\ell};\bm{\Delta})\ \text{for
some}\ \bm{\ell}\in L_{n1},\bm{\Delta}\neq\bm{\Delta}_{0}\\},$ $\displaystyle
J_{n3}$ $\displaystyle=J_{n}\backslash(J_{n1}\cup J_{n2}).$

Note that
$[\\![J_{n2}]\\!]\lesssim(\overline{A}_{n1})^{d-1}\overline{A}_{n2}\left({A_{n}h_{1}\dots
h_{d}\over A_{n}^{(1)}}\right)$ and $[\\![J_{n3}]\\!]\lesssim
A_{n}^{(1)}\left({\overline{A_{n}h}\over\underline{A}_{n1}}\right)^{d-1}$.
Then, applying the same argument to show (A.23), we have

$\displaystyle\operatorname{Var}(\widetilde{V}_{n2})$
$\displaystyle\lesssim[\\![J_{n2}]\\!](nA_{n}^{-1})^{2}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{d-1}\alpha_{1}^{1-2/q}(k)\right)$
$\displaystyle\lesssim\left({\overline{A}_{n1}\over\underline{A}_{n1}}\right)^{d}\left({\overline{A}_{n2}\over\overline{A}_{n1}}\right)A_{n}h_{1}\dots
h_{d}(nA_{n}^{-1})^{2}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{d-1}\alpha_{1}^{1-2/q}(k)\right)$
$\displaystyle=o\left(n^{2}A_{n}^{-1}h_{1}\dots h_{d}\right).$
$\displaystyle\operatorname{Var}(\widetilde{V}_{n3})$
$\displaystyle\lesssim[\\![J_{n3}]\\!](nA_{n}^{-1})^{2}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{d-1}\alpha_{1}^{1-2/q}(k)\right)$
$\displaystyle\lesssim\left({A_{n}^{(1)}\over\underline{A}_{n1}^{d}}\right)\left({\left(\overline{A_{n}h}\right)^{d}\over
A_{n}h_{1}\dots
h_{d}}\right)\left(\underline{A}_{n1}\over\overline{A_{n}h}\right)A_{n}h_{1}\dots
h_{d}(nA_{n}^{-1})^{2}\left(1+\sum_{k=1}^{\overline{A}_{n1}}k^{d-1}\alpha_{1}^{1-2/q}(k)\right)$
$\displaystyle=o\left(n^{2}A_{n}^{-1}h_{1}\dots h_{d}\right).$

(Step 2-4) Now we show (A.9). By (A.11) and (A.12), we have for sufficiently
large $n$,

$\displaystyle E[\widetilde{V}_{n1}^{2}]$
$\displaystyle=E[(\widetilde{V}_{n}(\bm{0})-(\widetilde{V}_{n2}+\widetilde{V}_{n3}))^{2}]\leq
2\left(E[(\widetilde{V}_{n}(\bm{0}))^{2}]+E[(\widetilde{V}_{n2}+\widetilde{V}_{n3})^{2}]\right)\leq
4E[\widetilde{V}_{n}^{2}(\bm{0})].$

Thus, by (A.4), (A.11), and (A.12), we have

$\displaystyle\left|\sum_{\bm{\ell}\in
L_{n1}}E[\widetilde{V}_{n}^{2}(\bm{\ell};\bm{\Delta}_{0})]-E[\widetilde{V}_{n}^{2}(\bm{0})]\right|$
$\displaystyle\leq\left|\sum_{\bm{\ell}\in
L_{n1}}E[\widetilde{V}_{n}^{2}(\bm{\ell};\bm{\Delta}_{0})]-E[\widetilde{V}_{n1}^{2}]\right|+2E[(\widetilde{V}_{n2}+\widetilde{V}_{n3})^{2}]^{1/2}E[\widetilde{V}_{n1}^{2}]^{1/2}+E[(\widetilde{V}_{n2}+\widetilde{V}_{n3})^{2}]$
$\displaystyle\lesssim\left(A_{n}^{(1)}nA_{n}^{-1}\right)^{2}\sum_{\bm{\ell}_{1}\neq\bm{\ell}_{2}}\alpha^{1-2/q}(\min\\{|\bm{\ell}_{1}-\bm{\ell}_{2}|-d,0\\}\underline{A}_{n3}+\underline{A}_{n2};A_{n}^{(1)})$
$\displaystyle\quad+o\left(n^{2}A_{n}^{-1}h_{1}\dots h_{d}\right)$
$\displaystyle\lesssim\left(A_{n}^{(1)}nA_{n}^{-1}\right)^{2}\left({A_{n}h_{1}\dots
h_{d}\over A_{n}^{(1)}}\right)$
$\displaystyle\quad\times\left(\alpha^{1-2/q}(\underline{A}_{n2};A_{n}^{(1)})+\sum_{k=1}^{\overline{A}_{n}/\underline{A}_{n1}}k^{d-1}\alpha^{1-2/q}(\min\\{|\bm{\ell}_{1}-\bm{\ell}_{2}|-d,0\\}\underline{A}_{n3}+\underline{A}_{n2};A_{n}^{(1)})\right)$
$\displaystyle\quad+o\left(n^{2}A_{n}^{-1}h_{1}\dots h_{d}\right)$
$\displaystyle=o\left(n^{2}A_{n}^{-1}h_{1}\dots h_{d}\right),$

where $\overline{A}_{n}=\max_{1\leq j\leq d}A_{n,j}$.

(Step 3) Now we evaluate $B_{n}(\bm{0})$. Decompose

$\displaystyle B_{n,j_{1}\dots j_{L}}(\dot{\bm{X}})$
$\displaystyle=\left\\{B_{n,j_{1}\dots j_{L}}(\dot{\bm{X}})-B_{n,j_{1}\dots
j_{L}}(\bm{0})-E\left[B_{n,j_{1}\dots j_{L}}(\dot{\bm{X}})-B_{n,j_{1}\dots
j_{L}}(\bm{0})\right]\right\\}$ $\displaystyle\quad+E\left[B_{n,j_{1}\dots
j_{L}}(\dot{\bm{X}})-B_{n,j_{1}\dots j_{L}}(\bm{0})\right]$
$\displaystyle\quad+\left\\{B_{n,j_{1}\dots
j_{L}}(\bm{0})-E\left[B_{n,j_{1}\dots j_{L}}(\bm{0})\right]\right\\}$
$\displaystyle\quad+E\left[B_{n,j_{1}\dots j_{L}}(\bm{0})\right]$
$\displaystyle=:\sum_{\ell=1}^{4}B_{n,j_{1}\dots j_{L}\ell}.$

Define $N_{\bm{x}}(h):=\prod_{j=1}^{d}[x_{j}-h_{j},x_{j}+h_{j}]$ and
$\bm{x}=(x_{1},\dots,x_{d})\in(-1/2,1/2)^{d}$. For $B_{n,j_{1}\dots j_{L}1}$,
by a change of variables and the dominated convergence theorem, we have

$\displaystyle\operatorname{Var}(B_{n,j_{1}\dots j_{L}1})$
$\displaystyle\leq{A_{n}\over\\{(p+1)!\\}^{2}nh_{1}\dots
h_{d}}E\left[K_{Ah}^{2}\left(\bm{X}_{i}\right)\prod_{\ell=1}^{L}\left({X_{i,j_{\ell}}\over
A_{n,j_{\ell}}h_{j_{\ell}}}\right)^{2}\right.$
$\displaystyle\left.\quad\times\sum_{1\leq j_{1,1}\leq\dots\leq j_{1,p+1}\leq
d,1\leq j_{2,1}\leq\dots\leq j_{2,p+1}\leq d}{1\over\bm{s}_{j_{1,1}\dots
j_{1,p+1}}!}{1\over\bm{s}_{j_{2,1}\dots j_{2,p+1}}!}\right.$
$\displaystyle\left.\quad\times(\partial_{j_{1,1}\dots
j_{1,p+1}}m(\dot{\bm{X}}_{i}/A_{n})-\partial_{j_{1,1}\dots
j_{1,p+1}}m(\bm{0}))(\partial_{j_{2,1}\dots
j_{2,p+1}}m(\dot{\bm{X}}_{i}/A_{n})-\partial_{j_{2,1}\dots
j_{2,p+1}}m(\bm{0}))\right.$
$\displaystyle\left.\quad\times\prod_{\ell_{1}=1}^{p+1}{X_{i,j_{1,\ell_{1}}}\over
A_{n,j_{\ell_{1}}}}\prod_{\ell_{2}=1}^{p+1}{X_{i,j_{2,\ell_{2}}}\over
A_{n,j_{\ell_{2}}}}\right]$
$\displaystyle\leq{A_{n}\over\\{(p+1)!\\}^{2}n}\max_{1\leq j_{1}\leq\dots\leq
j_{p+1}\leq d}\sup_{\bm{y}\in N_{\bm{0}}(h)}|\partial_{j_{1}\dots
j_{p+1}}m(\bm{y})-\partial_{j_{1}\dots j_{p+1}}m(\bm{0})|^{2}$
$\displaystyle\quad\quad\times\sum_{1\leq j_{1,1}\leq\dots\leq j_{1,p+1}\leq
d,1\leq j_{2,1}\leq\dots\leq j_{2,p+1}\leq
d}\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\prod_{\ell_{2}=1}^{p+1}h_{j_{2,\ell_{2}}}$
$\displaystyle\quad\quad\times\int\left(\prod_{\ell=1}^{L}z_{j_{\ell}}^{2}\prod_{\ell_{1}=1}^{p+1}|z_{j_{1,\ell_{1}}}|\prod_{\ell_{2}=1}^{p+1}|z_{j_{2,\ell_{2}}}|\right)K^{2}(\bm{z})g(\bm{z}\circ\bm{h})d\bm{z}$
$\displaystyle=o\left({A_{n}\over n}\sum_{1\leq j_{1,1}\leq\dots\leq
j_{1,p+1}\leq d,1\leq j_{2,1}\leq\dots\leq j_{2,p+1}\leq
d}\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\prod_{\ell_{2}=1}^{p+1}h_{j_{2,\ell_{2}}}\right)$
$\displaystyle=o(1).$ (A.29)

Then we have $B_{n,j_{1}\dots j_{L}1}=o_{p}(1)$.

For $B_{n,j_{1}\dots j_{L}2}$,

$\displaystyle|B_{n,j_{1}\dots j_{L}2}|$
$\displaystyle\leq{1\over(p+1)!}\max_{1\leq j_{1},\dots,j_{p+1}\leq
d}\sup_{\bm{y}\in N_{\bm{0}}(h)}|\partial_{j_{1}\dots
j_{p+1}}m(\bm{y})-\partial_{j_{1}\dots j_{p+1}}m(\bm{0})|$
$\displaystyle\quad\times\sqrt{A_{n}h_{1}\dots h_{d}}\sum_{1\leq
j_{1,1}\leq\dots\leq j_{1,p+1}\leq
d}\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\int\left(\prod_{\ell=1}^{L}|z_{j_{\ell}}|\prod_{\ell_{1}=1}^{p+1}|z_{j_{1,\ell_{1}}}|\right)|K(\bm{z})|g(\bm{z}\circ\bm{h})d\bm{z}$
$\displaystyle=o(1).$ (A.30)

For $B_{n,j_{1}\dots j_{L}3}$,

$\displaystyle\operatorname{Var}(B_{n,j_{1}\dots j_{L}3})$
$\displaystyle\leq{A_{n}h_{1}\dots h_{d}\over\\{(p+1)!\\}^{2}nh_{1}\dots
h_{d}}\sum_{1\leq j_{1,1}\leq\dots\leq j_{1,p+1}\leq d,1\leq
j_{2,1}\leq\dots\leq j_{2,p+1}\leq d}\partial_{j_{1,1}\dots
j_{1,p+1}}m(\bm{0})\partial_{j_{2,1}\dots j_{2,p+1}}m(\bm{0})$
$\displaystyle\quad\quad\times\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\prod_{\ell_{2}=1}^{p+1}h_{j_{2,\ell_{2}}}\int\left(\prod_{\ell=1}^{L}z_{j_{\ell}}^{2}\prod_{\ell_{1}}^{p+1}|z_{j_{1,\ell_{1}}}|\prod_{\ell_{2}=1}^{p+1}|z_{j_{2,\ell_{2}}}|\right)K^{2}(\bm{z})g(\bm{z}\circ\bm{h})d\bm{z}$
$\displaystyle=O\left({A_{n}\over n}\sum_{1\leq j_{1,1}\leq\dots\leq
j_{1,p+1}\leq d,1\leq j_{2,1}\leq\dots\leq j_{2,p+1}\leq
d}\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\prod_{\ell_{2}=1}^{p+1}h_{j_{2,\ell_{2}}}\right).$
(A.31)

Then we have $B_{n,j_{1}\dots j_{L}3}=o_{p}(1)$.

For $B_{n,j_{1}\dots j_{L}4}$,

$\displaystyle B_{n,j_{1}\dots j_{L}4}$ $\displaystyle=\sqrt{A_{n}h_{1}\dots
h_{d}}\sum_{1\leq j_{1,1}\leq\dots\leq j_{1,p+1}\leq d}{\partial_{j_{1,1}\dots
j_{1,p+1}}m(\bm{0})\over\bm{s}_{j_{1,1}\dots j_{1,p+1}}!}$
$\displaystyle\quad\times\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\int\left(\prod_{\ell=1}^{L}z_{j_{\ell}}\prod_{\ell_{1}=1}^{p+1}z_{j_{1,\ell_{1}}}\right)K(\bm{z})g(\bm{z}\circ\bm{h})d\bm{z}$
$\displaystyle=g(\bm{0})\sqrt{A_{n}h_{1}\dots h_{d}}\sum_{1\leq
j_{1,1}\leq\dots\leq j_{1,p+1}\leq d}{\partial_{j_{1,1}\dots
j_{1,p+1}}m(\bm{0})\over\bm{s}_{j_{1,1}\dots
j_{1,p+1}}!}\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\kappa_{j_{1}\dots
j_{L}j_{1,1}\dots j_{1,p+1}}^{(1)}+o(1).$ (A.32)

Combining (A.1)-(A.1),

$\displaystyle B_{n,j_{1}\dots j_{L}}(\dot{\bm{X}})$
$\displaystyle=g(\bm{0})\sqrt{A_{n}h_{1}\dots h_{d}}\sum_{1\leq
j_{1,1}\leq\dots\leq j_{1,p+1}\leq d}{\partial_{j_{1,1}\dots
j_{1,p+1}}m(\bm{0})\over\bm{s}_{j_{1,1}\dots j_{1,p+1}}!}$
$\displaystyle\quad\times\prod_{\ell_{1}=1}^{p+1}h_{j_{1,\ell_{1}}}\kappa_{j_{1}\dots
j_{L}j_{1,1}\dots j_{1,p+1}}^{(1)}+o_{p}(1)$
$\displaystyle=g(\bm{0})\sqrt{A_{n}h_{1}\dots
h_{d}}(B^{(d,p)}\bm{M}_{n}^{(d,p)}(\bm{0}))_{j_{1}\dots j_{L}}+o_{p}(1).$

(Step 4) Combining the results in Steps 2 and 3, we have

$\displaystyle A_{n}(\bm{0})$
$\displaystyle:=V_{n}(\bm{0})+\left(B_{n}(\bm{0})-g(\bm{0})\sqrt{A_{n}h_{1}\dots
h_{d}}B^{(d,p)}M_{n}^{(d,p)}(\bm{0})\right)$
$\displaystyle\stackrel{{\scriptstyle
d}}{{\to}}N\left(\left(\begin{array}[]{c}0\\\ \vdots\\\
0\end{array}\right),g(\bm{0})\left\\{\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))+\eta^{2}(\bm{0})g(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right\\}\mathcal{K}\right).$

This and the result in Step 1 yield the desired result. ∎

###  A.2. Proof of Proposition 4.1

###### Proof.

It is easy to see that $\widehat{g}(\bm{0})\stackrel{{\scriptstyle
p}}{{\to}}g(\bm{0})$ as $n\to\infty$. For $\widehat{V}_{n,1}(\bm{0})$, observe
that

$\displaystyle
E\left[\left.\left(Y(\bm{X}_{i})-\widehat{m}_{-\\{i\\}}(\bm{X}_{i}/A_{n})\right)^{2}\right|\bm{X}_{i}\right]$
$\displaystyle=E\left[\left.\left(m(\bm{X}_{i}/A_{n})-\widehat{m}_{-\\{i\\}}(\bm{X}_{i}/A_{n})\right)^{2}\right|\bm{X}_{i}\right]$
$\displaystyle\quad+2E\left[\left.\left(m(\bm{X}_{i}/A_{n})-\widehat{m}_{-\\{i\\}}(\bm{X}_{i}/A_{n})\right)(e_{n,i}+\varepsilon_{n,i})\right|\bm{X}_{i}\right]$
$\displaystyle\quad+E\left[\left.(e_{n,i}+\varepsilon_{n,i})^{2}\right|\bm{X}_{i}\right]$
$\displaystyle=:\widehat{V}_{n1,i}+\widehat{V}_{n2,i}+E\left[\left.(e_{n,i}+\varepsilon_{n,i})^{2}\right|\bm{X}_{i}\right].$

The representation of the MSE of $\partial_{j_{1}\dots
j_{L}}\widehat{m}(\bm{0})$ (4.1) implies that for $\bm{z}\in(-1/2,1/2)^{d}$,

$\displaystyle\text{MSE}(\widehat{m}(\bm{z}))$
$\displaystyle=\left\\{(S^{-1}e_{0})^{\prime}B^{(d,p)}M_{n}^{(d,p)}(\bm{z})\right\\}^{2}$
$\displaystyle\quad+\left({\kappa(\eta^{2}(\bm{z})+\sigma_{\varepsilon}^{2}(\bm{z}))\over
g(\bm{z})}+\eta^{2}(\bm{z})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}\right){e^{\prime}_{0}S^{-1}\mathcal{K}S^{-1}e_{0}\over
A_{n}h_{1}\dots h_{d}}$ $\displaystyle=O\left(\left(\sum_{1\leq
j_{1}\leq\dots\leq j_{p+1}\leq
d}\prod_{\ell=1}^{p+1}h_{j_{\ell}}\right)^{2}+{1\over A_{n}h_{1}\dots
h_{d}}\right),$

where $e_{0}=(1,0,\dots,0)\in\mathbb{R}^{D}$. Then by Cauchy-Schwarz
inequality, we have

$\displaystyle\max_{1\leq i\leq
n}K_{Ah}(\bm{X}_{i})\left(\widehat{V}_{n1,i}+\widehat{V}_{n2,i}\right)$
$\displaystyle\quad=O\left(\left(\sum_{1\leq j_{1}\leq\dots\leq j_{p+1}\leq
d}\prod_{\ell=1}^{p+1}h_{j_{\ell}}\right)+\sqrt{{1\over A_{n}h_{1}\dots
h_{d}}}\right)=o(1)\ a.s.,\ n\to\infty.$

Applying a similar argument in the proof of Theorem 4.1, this implies that

$\displaystyle\widehat{V}_{n,1}(\bm{0})$ $\displaystyle={1\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}(\bm{X}_{i})(e_{n,i}+\varepsilon_{n,i})^{2}+o_{p}(1)$
$\displaystyle={1\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}(\bm{X}_{i})\left(\eta^{2}(\bm{X}_{i}/A_{n})e^{2}(\bm{X}_{i})+\sigma^{2}_{\varepsilon}(\bm{X}_{i}/A_{n})\varepsilon^{2}_{i}\right)+o_{p}(1)$
$\displaystyle=\left(\eta^{2}(\bm{0})+\sigma^{2}_{\varepsilon}(\bm{0})\right)g(\bm{0})+o_{p}(1),\
n\to\infty.$

Likewise,

$\displaystyle\widehat{V}_{n,2}(\bm{0})$ $\displaystyle={A_{n}\over
nh_{1}\dots
h_{d}}\sum_{i=1}^{n-1}K_{Ah}(\bm{X}_{i})K_{Ah}(\bm{X}_{i+1})(e_{n,i}+\varepsilon_{n,i})(e_{n,i+1}+\varepsilon_{n,i+1})+o_{p}(1)$
$\displaystyle={A_{n}\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n-1}K_{Ah}(\bm{X}_{i})K_{Ah}(\bm{X}_{i+1})\eta(\bm{X}_{i}/A_{n})\eta(\bm{X}_{i+1}/A_{n})e(\bm{X}_{i})e(\bm{X}_{i+1})+o_{p}(1)$
$\displaystyle=\kappa_{0}^{(2)}\eta^{2}(\bm{0})g^{2}(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}+o_{p}(1),\
n\to\infty.$

Therefore, as $n\to\infty$,

$\displaystyle\widehat{V}_{n}(\bm{0}):={(A_{n}/n)\widehat{V}_{n,1}(\bm{0})+(\widehat{V}_{n,2}(\bm{0})/\kappa_{0}^{(2)})\over\widehat{g}^{2}(\bm{0})}\stackrel{{\scriptstyle
p}}{{\to}}{\kappa(\eta^{2}(\bm{0})+\sigma_{\varepsilon}^{2}(\bm{0}))\over
g(\bm{0})}+\eta^{2}(\bm{0})\int\sigma_{\bm{e}}(\bm{v})d\bm{v}.$

∎

###  A.3. Proof of Corollary 4.1

###### Proof.

Corollary 4.1 follows immediately from Theorem 4.1 and Proposition 4.1. ∎

###  A.4. Proof of Proposition 4.2

###### Proof.

For any $\bm{t}=(t_{0},t_{1},\dots,t_{d},t_{11},\dots,t_{dd},\dots,t_{1\dots
1},\dots,t_{d\dots d})^{\prime}\in\mathbb{R}^{D}$, we define

$\displaystyle\overline{W}_{n1}(\bm{0})$
$\displaystyle:=\sum_{\ell_{1}=1}^{n_{1}}K_{Ah}\left(\bm{X}_{1,\ell_{1}}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1,\ell_{1}}\end{array}\right)\right]\underbrace{\left(\eta_{1}\left({\bm{X}_{1,\ell_{1}}\over
A_{n}}\right)e_{1}(\bm{X}_{1,\ell_{1}})+\sigma_{\varepsilon,1}\left({\bm{X}_{1,\ell_{1}}\over
A_{n}}\right)\varepsilon_{1,\ell_{1}}\right)}_{=:\overline{e}_{n1,\ell_{1}}+\overline{\varepsilon}_{n1,\ell_{1}}},$
$\displaystyle\overline{W}_{n2}(\bm{0})$
$\displaystyle:=\sum_{\ell_{2}=1}^{n_{2}}K_{Ah}\left(\bm{X}_{2,\ell_{2}}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{2,\ell_{2}}\end{array}\right)\right]\underbrace{\left(\eta_{2}\left({\bm{X}_{2,\ell_{2}}\over
A_{n}}\right)e_{2}(\bm{X}_{2,\ell_{2}})+\sigma_{\varepsilon,2}\left({\bm{X}_{2,\ell_{2}}\over
A_{n}}\right)\varepsilon_{2,\ell_{2}}\right)}_{=:\overline{e}_{n2,\ell_{2}}+\overline{\varepsilon}_{n2,\ell_{2}}}.$

By inspecting the proof of Theorem 4.1, to show Proposition 4.2, it is
sufficient to verify

$\displaystyle
E\left[\left(\overline{W}_{n1}(\bm{0})-\overline{W}_{n1}(\bm{0})\right)^{2}\right]/(h_{1}\dots
h_{d})$
$\displaystyle=\left(n_{1}\left\\{(\eta_{1}^{2}(\bm{0})+\sigma_{\varepsilon,1}^{2}(\bm{0}))g_{1}(\bm{0})+n_{1}A_{n}^{-1}\eta_{1}^{2}(\bm{0})g_{1}^{2}(\bm{0})\int\sigma_{\bm{e},11}(\bm{v})d\bm{v}\right\\}\right.$
$\displaystyle\left.\quad+n_{2}\left\\{(\eta_{2}^{2}(\bm{0})+\sigma_{\varepsilon,2}^{2}(\bm{0}))g_{2}(\bm{0})+n_{2}A_{n}^{-1}\eta_{2}^{2}(\bm{0})g_{2}^{2}(\bm{0})\int\sigma_{\bm{e},22}(\bm{v})d\bm{v}\right\\}\right.$
$\displaystyle\left.\quad-2n_{1}n_{2}A_{n}^{-1}\eta_{1}(\bm{0})\eta_{2}(\bm{0})g_{1}(\bm{0})g_{2}(\bm{0})\int\sigma_{\bm{e},12}(\bm{v})d\bm{v}\right)$
$\displaystyle\quad\times\left(\int
K^{2}(\bm{z})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]^{2}d\bm{z}\right)(1+o(1)),\
n\to\infty.$

Let $E_{\bm{X}_{12}}$ denote the expectation with respect to
$\\{\bm{X}_{1,\ell_{1}}\\}$ and $\\{\bm{X}_{2,\ell_{2}}\\}$ and let
$E_{\cdot\mid\bm{X}_{12}}$ denote the conditional expectation given
$\sigma(\\{\bm{X}_{1,\ell_{1}}\\}\cup\\{\bm{X}_{2,\ell_{2}}\\})$. Observe that

$\displaystyle
E_{\cdot\mid\bm{X}_{12}}\left[\left(\overline{W}_{n1}(\bm{0})-\overline{W}_{n2}(\bm{0})\right)^{2}\right]$
$\displaystyle=\sum_{\ell_{11},\ell_{12}=1}^{n_{1}}E_{\cdot\mid\bm{X}_{12}}\left[K_{Ah}\left(\bm{X}_{1,\ell_{11}}\right)K_{Ah}\left(\bm{X}_{1,\ell_{12}}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1,\ell_{11}}\end{array}\right)\right]\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1,\ell_{12}}\end{array}\right)\right]\right.$
$\displaystyle\left.\quad\times(\overline{e}_{n1,\ell_{11}}+\overline{\varepsilon}_{n1,\ell_{11}})(\overline{e}_{n1,\ell_{12}}+\overline{\varepsilon}_{n1,\ell_{12}})\right]$
$\displaystyle\quad+\sum_{\ell_{21},\ell_{22}=1}^{n_{2}}E_{\cdot\mid\bm{X}_{12}}\left[K_{Ah}\left(\bm{X}_{2,\ell_{21}}\right)K_{Ah}\left(\bm{X}_{2,\ell_{22}}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{2,\ell_{21}}\end{array}\right)\right]\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{2,\ell_{22}}\end{array}\right)\right]\right.$
$\displaystyle\left.\quad\times(\overline{e}_{n2,\ell_{21}}+\overline{\varepsilon}_{n2,\ell_{21}})(\overline{e}_{n2,\ell_{22}}+\overline{\varepsilon}_{n2,\ell_{22}})\right]$
$\displaystyle\quad-2\sum_{\ell_{1}=1}^{n_{1}}\sum_{\ell_{2}=1}^{n_{2}}E_{\cdot\mid\bm{X}_{12}}\left[K_{Ah}\left(\bm{X}_{1,\ell_{1}}\right)K_{Ah}\left(\bm{X}_{2,\ell_{2}}\right)\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{1,\ell_{1}}\end{array}\right)\right]\left[\bm{t}^{\prime}H^{-1}\left(\begin{array}[]{c}1\\\
\check{\bm{X}}_{2,\ell_{2}}\end{array}\right)\right]\right.$
$\displaystyle\left.\quad\times(\overline{e}_{n1,\ell_{1}}+\overline{\varepsilon}_{n1,\ell_{1}})(\overline{e}_{n2,\ell_{2}}+\overline{\varepsilon}_{n2,\ell_{2}})\right]$
$\displaystyle=:\overline{W}_{n11}+\overline{W}_{n12}-2\overline{W}_{n13}.$

Applying the same argument in Step 2 of the proof of Theorem 4.1, we have

$\displaystyle E_{\bm{X}_{12}}[\overline{W}_{n1\ell}]$
$\displaystyle=n_{\ell}h_{1}\dots
h_{d}g_{\ell}(\bm{0})\left\\{(\eta_{\ell}^{2}(\bm{0})+\sigma_{\varepsilon,\ell}^{2}(\bm{0}))+n_{\ell}A_{n}^{-1}\eta_{\ell}^{2}(\bm{0})g_{\ell}(\bm{0})\int\sigma_{\bm{e},\ell\ell}(\bm{v})d\bm{v}\right\\}$
$\displaystyle\quad\times\left(\int
K^{2}(\bm{z})\left[\bm{t}^{\prime}\left(\begin{array}[]{c}1\\\
\check{\bm{z}}\end{array}\right)\right]^{2}d\bm{z}\right)(1+o(1)),\ \ell=1,2,$
$\displaystyle E_{\bm{X}_{12}}[\overline{W}_{n13}]$
$\displaystyle=n_{1}n_{2}A_{n}^{-1}h_{1}\dots
h_{d}\left(\eta_{1}(\bm{0})\eta_{2}(\bm{0})g_{1}(\bm{0})g_{2}(\bm{0})\int\sigma_{\bm{e},12}(\bm{v})d\bm{v}\right)(1+o(1))$

as $n\to\infty$. Therefore, we obtain the desired result. ∎

###  A.5. Proof of Proposition 4.3

###### Proof.

Applying the same argument in the proof of Proposition 4.1, we have that as
$n\to\infty,$

$\displaystyle\overline{g}_{n_{k}}(\bm{0})$
$\displaystyle=g_{k}(\bm{0})+o_{p}(1),\ k=1,2,$
$\displaystyle\overline{V}_{n,1k}(\bm{0})$
$\displaystyle=(\eta_{k}^{2}(\bm{0})+\sigma_{\varepsilon,k}^{2}(\bm{0}))g_{k}(\bm{0})+o_{p}(1),\
k=1,2,$ $\displaystyle\overline{V}_{n,2k}(\bm{0})$
$\displaystyle=\kappa_{0}^{(2)}\eta_{k}^{2}(\bm{0})g_{k}^{2}(\bm{0})\int\sigma_{\bm{e},kk}(\bm{v})d\bm{v}+o_{p}(1),\
k=1,2,$ $\displaystyle\overline{V}_{n,3}(\bm{0})$
$\displaystyle=\kappa_{0}^{(2)}\eta_{1}(\bm{0})\eta_{2}(\bm{0})g_{1}(\bm{0})g_{2}(\bm{0})\int\sigma_{\bm{e},12}(\bm{v})d\bm{v}+o_{p}(1).$

Therefore, $\check{V}_{n}(\bm{0})\stackrel{{\scriptstyle
p}}{{\to}}\overline{V}_{1}(\bm{0})+\overline{V}_{2}(\bm{0})-2\overline{V}_{3}(\bm{0})$
as $n\to\infty$. ∎

###  A.6. Proof of Corollary 4.2

###### Proof.

Corollary 4.2 follows immediately from Propositions 4.2 and 4.3. ∎

##  Appendix B Proofs for Section 5

Before we prove Theorem 5.1, we consider general kernel estimators and derive
their uniform convergence rates (Section B.1). Since the estimators include
many kernel-based estimators such as, kernel density, LC, LL, and LP
estimators for random fields on $\mathbb{R}^{d}$ with irregularly spaced
sampling sites, the results are of independent theoretical interest. As
applications of the results, we derive uniform convergence rates of LP
estimators (Section B.2).

###  B.1. Uniform convergence rates for general kernel estimators

For $j=1,2,3$, let $f_{j}:\mathbb{R}^{d}\to\mathbb{R}$ be functions such that
$f_{j}$ is continuous on $R_{0,\delta}:=(-1/2-\delta,1/2+\delta)^{d}$ for some
$\delta>0$. Define

$\displaystyle\widehat{\Psi}_{\mathrm{I}}(\bm{z})$ $\displaystyle={1\over
n^{2}A_{n}^{-1}h_{1}\dots h_{d}}\sum_{i=1}^{n}K_{Ah}(\bm{X}_{i}-A_{n}\bm{z})$
$\displaystyle\quad\times
f_{1,Ah}\left(\bm{X}_{i}-A_{n}\bm{z}\right)f_{2,A}\left(\bm{X}_{i}-A_{n}\bm{z}\right)f_{3,A}\left(\bm{X}_{i}\right)Z_{\bm{X}_{i}},$
(B.1) $\displaystyle\widehat{\Psi}_{\mathrm{II}}(\bm{z})$
$\displaystyle={1\over nh_{1}\dots
h_{d}}\sum_{i=1}^{n}K_{Ah}(\bm{X}_{i}-A_{n}\bm{z})$ $\displaystyle\quad\times
f_{1,Ah}\left(\bm{X}_{i}-A_{n}\bm{z}\right)f_{2,A}\left(\bm{X}_{i}-A_{n}\bm{z}\right)f_{3,A}\left(\bm{X}_{i}\right),$
(B.2)

where $f_{j,Aa}(\bm{x})=f_{j}\left({x_{1}\over A_{n,1}a_{1}},\dots,{x_{d}\over
A_{n,d}a_{d}}\right)$ for $\bm{a}=(a_{1},\dots
a_{d})^{\prime}\in(0,\infty)^{d}$ and $\\{Z_{\bm{X}_{i}}\\}_{i=1}^{n}$ is a
sequence of real-valued random variables. Many kernel estimators, such as
kernel density, Nadaraya-Watson, and LP estimators, can be represented by
combining special cases of estimators (B.1) or (B.2). In this study, we use
the uniform convergence rates of these estimators with

$\displaystyle f_{1}$ $\displaystyle\in\left\\{e^{\prime}_{j_{1}\dots
j_{L}}\left(\begin{array}[]{c}1\\\
\check{\bm{x}}\end{array}\right),e^{\prime}_{j_{1,1}\dots
j_{1,L_{1}}}\left(\begin{array}[]{c}1\\\ \check{\bm{x}}\end{array}\right)(1\
\check{\bm{x}}^{\prime})e_{j_{2,1}\dots j_{2,L_{2}}}\right\\},$ $\displaystyle
f_{2}$ $\displaystyle\in\left\\{1,\prod_{\ell=1}^{L}x_{j_{\ell}}\right\\},\
f_{3}\in\left\\{1,\eta,\sigma_{\varepsilon},\\{\partial_{j_{1}\dots
j_{p+1}}m\\}_{1\leq j_{1}\leq\dots\leq j_{p+1}\leq d}\right\\},\
Z_{\bm{X}_{i}}\in\left\\{e(\bm{X}_{i}),\varepsilon_{i}\right\\}.$

We assume the following conditions for the sampling sites
$\\{\bm{X}_{i}\\}_{i=1}^{n}$:

######  Assumption B.1.

Let $g$ be a probability density function with support $R_{0}=[-1/2,1/2]^{d}$.

  * (i)