The exact contents of citations.db main_text.text for this paper — one flattened LaTeX string, title through conclusion, appendix excluded, unmodified except for removing email addresses. This is what our citation measures are computed over.
84,589 characters
Multivariate kernel regression in vector and product metric spaces
\author{Marcia Schafgans\thanks{\textit{Corresponding author}: Economics
Department, London School of Economics, Houghton Street, London WC2A 2AE,
UK. E-mail address: [email removed]} \qquad \qquad Victoria Zinde-Walsh
\thanks{
Economics Department, McGill University, 855 Sherbrooke St. W., Montreal,
Quebec H3A 2T7, Canada. E-mail address: [email removed].} \\
London School of Economics \quad McGill University and CIREQ }
\title{Multivariate kernel regression in vector and product metric spaces}
\date{January 2026}
\maketitle
\begin{abstract}
This paper derives limit properties of nonparametric kernel regression
estimators without requiring existence of density for regressors in $\mathbb{R}^{q}.$
In functional regression limit properties are established for multivariate
functional regression. The rate and asymptotic normality for the
Nadaraya-Watson (NW) estimator is established for distributions of
regressors in $\mathbb{R}^{q}$ that allow for mass points, factor structure,
multicollinearity and nonlinear dependence, as well as fractal distribution;
when bounded density exists we provide statistical guarantees for the
standard rate and the asymptotic normality without requiring smoothness. We
demonstrate faster convergence associated with dimension reducing types of
singularity, such as a fractal distribution or a factor structure in the
regressors. The paper extends asymptotic normality of kernel functional
regression to multivariate regression over a product of any number of metric
spaces. Finite sample evidence confirms rate improvement due to singularity
in regression over $\mathbb{R}^{q}.$ For functional regression the
simulations underline the importance of accounting for multiple functional
regressors. We demonstrate the applicability and advantages of the NW
estimator in our empirical study, which reexamines the job training program
evaluation based on the LaLonde data.\newline
\end{abstract}
\section{Introduction}
This paper extends nonparametric kernel regression to more general regressor
settings than those considered in the literature. The general regression
model
\begin{equation}
Y=m(X)+u,\quad E(u|X)=0, \label{model}
\end{equation}
is free from the difficulty of choosing a parametric specification. We focus
on the Nadaraya-Watson (NW) estimator, introduced by Nadaraya (1965) and
Watson (1964), which recognizes that a continuous regression function can be
estimated pointwise by a weighted average that attaches higher weights to
close-by observations.
Here we emphasize the fact that the regressor $X$ can be a vector in $
\mathbb{R}^{q}$, or alternatively $X$ may belong to a function space, a more
general metric space, or comprise of components from several such metric
spaces. In fact, the components of $X$ do not necessarily have to belong
to spaces of vectors or functions but could be intervals, graphs, or
networks, as long as a metric (or even a semi-metric) can be defined for
each space. $Y$ represents a scalar dependent variable, $u$ denotes an
unobserved error, and the conditional mean function $m$ satisfies some
smoothness assumptions.
The NW estimator has been used extensively with $X\in \mathbb{R}^{q}$ (see
e.g. the textbook Li and Racine, 2007, for discussion and examples) and has
recently been introduced to functional regression by Ferraty and Vieu
(2004). Well known limit distributional results for the NW estimator were
derived in $\mathbb{R}^{q}$ under restrictions requiring existence and
smoothness of the density. For functional regression (where there is no
density) the limit distributional results were derived in a univariate
context only.
In this paper we establish asymptotic normality of the NW estimator for
regression in the presence of a general regressor $X$ that could have a
multivariate {singular distribution} in $\mathbb{R}^{q}$ or is comprised of
any number of functional and vector regressors. A singular distribution does
not admit a density function that integrates to it.
In settings where data has both discrete and continuous components Racine
and Li (2007) obtained asymptotic normality of the NW estimator without
having to deal with the singularity by treating the discrete and
(absolutely) continuous components separately. However, sometimes the
distinction between discrete and continuous variables is not
straightforward; continuous variables could be discretized with different
levels of discretization. When data with both discrete and continuous
components is viewed as a vector in a Euclidean space, $X\in \mathbb{R}^{q},$
the distribution of $X$ is singular. In our simulations we demonstrate that
there may be no gain from avoiding the singularity by considering the
discrete regressors separately.
The presence of latent factors in the continuous regressors, common in
macroeconomic and finance models (e.g., Bai and Ng, 2006, for portfolio,
stock returns and macroeconomic data) could also imply a singular
distribution. For example, if the regressor is a \thinspace $q\times 1$
vector $X\sim N\left( 0,\Sigma \right) $ with $\Sigma $ a singular matrix of
rank $r<q$, the distribution is singular. Similarly, non-linear common
factors, such as in Hotelling's (1929) spatial model of horizontal
differentiation which assumes that each consumer has an `ideal' variety
identified by his location on the unit circle (see also Desmet et al.,
2010), imply a smaller effective dimension for the regressor space,
resulting in a singular distribution over $\mathbb{R}^{q}.$ This also is
true when there exists a functional relation between the regressors (e.g.
with exact collinearity that can arise in production functions, Ackerberg et
al., 2015).
Singularity also originates from a fractal structure in the data; examples
in economics include the daily prices in the cotton market (Mandelbrot,
1997), financial markets, and networks (see Takayasu et al., 2009). Fractals
are common to many geographic features, including coastlines, river networks
and landforms, and have been used in urban growth studies (e.g., Shen, 2002)
and spatial econometrics in general. Furthermore, singularities also result
when continuously distributed variables exhibit mass points (e.g.
Arulampalam et al., 2017, for neonatal mortality and Olson, 1998, for weekly
hours worked).
We demonstrate the benefit of extending the NW estimator to regression with
singular data by applying it to the data from a randomized experiment in a
job training program evaluation study by LaLonde (1986). Following the work
by Rosenbaum and Rubin (1983), Dehejia and Wahba (1999, 2002) applied
propensity score methods to the LaLonde data for estimation of causal
treatment effects in an attempt to generalize the experimental results to
nonexperimental data. The propensity score matching was used instead of a
multivariate nonparametric model with matching on individual characteristics
which was deemed impractical because of the high dimensionality of the
regressors. The benefits of kernel regression were analyzed in Heckman et
al. (1997, 1998) (without allowing for singularity). The LaLonde data and
methodologies were discussed by Angrist and Pischke (2009) and in a recent
review by Imbens and Xu (2024). As shown here the discreteness of most of
the regressors implies reduced dimension of the support of the joint
distribution; for continuous variables existence of continuous density is
not imposed and mass at zero in income is accounted for. Our asymptotic
results provide the validity of the NW estimator for this singular
distribution. The kernel estimators we employ give new insights into the
heterogeneous effects of the program, based on a variety of individual
characteristics and compare quite well with random forest estimates of the
conditional average treatment effect on the treated (CATT) (Wager and Athey,
2018).
Our results also extend to functional regression (see, e.g. Ramsey and
Silverman, 2005) where estimation and inference techniques have been
developed by Ferraty and Vieu (2004) and pointwise asymptotic normality was
established in regression for a Banach or metric space by Ferraty et al. (2007), Ferraty and Vieu (2006), and Geenens (2015) in the i.i.d. case. Masry
(2005) derived the limit distribution for a strongly mixing process.
Recently Kurisu et al. (2025) made a case for extending the univariate
set-up of functional regression by considering jointly a random vector and a
function to obtain an estimate for the propensity score used in evaluating
the average treatment effect. We establish asymptotic normality in
multivariate functional regression with regressors in any number of
heterogeneous metric spaces. This provides a basis for simultaneously
evaluating the impact of the different predictors rather than comparing
their performance in distinct models, as in Caldeira et al. (2020) and
Ferraty and Nagy (2022).\footnote{
E.g., Caldeira et al. (2020) compares the model forecasting aggregate stock
market excess return on a function representing the history of returns with
regression models based on traditional predictors. Ferraty and Nagy (2022)
compare the performance of separate models for predicting adult height with
functional regressors (one being growth velocity profiles from ages 1-10 and
the other for 5-8).} Our simulations show that using multivariate rather
than univariate functional regression can improve the fit of the kernel
estimator.
We derive asymptotic normality results for a random regressor $X$ supported
on some domain in a vector space, $\mathbb{R}^{q},$ or metric, semi-metric
space, $\Xi ^{\left[ 1\right] },$ or a product of such spaces $\Xi ^{\lbrack
q]}\equiv \Xi _{1}^{[1]}\times \cdots \times \Xi _{q}^{[1]}$. The metrics on
$\Xi _{l}^{[1]}$, $\left\Vert .\right\Vert _{l}$, may differ for each of the
$q$ components of function spaces, thus as in Kurisu et al. (2025) one may
be the $\mathbb{R}^{1}$ space and the other one a function space. A key
ingredient in our technical derivations is small cube probability, which
characterizes local properties of $X$ in the general multivariate case in
place of the density. We introduce this concept here.
In the univariate metric space, $\Xi =\Xi ^{\lbrack 1]}$, the probability
measure is characterized by the small ball probability (e.g., see Ferraty
and Vieu, 2006): for the ball $B\left( x,h\right) =\left\{ X:\left\Vert
x-X\right\Vert \leq h\right\} $ {centered at $x$} in $\Xi ^{\left[ 1\right] }
$ {the probability measure} is denoted $P_{X}\left( B\left( x,h\right)
\right) $. Characterizing the measure locally via a ball is insufficient
when we wish to examine heterogeneous regressors in $\mathbb{R}^{q}$ {or, in
general, in product metric spaces $\Xi ^{\left[ q\right] }$.}
Kankanala and Zinde-Walsh (2024) introduced small cube probability for a
cuboid. A cuboid $C\left( x,h\right) $ centered around $x=\left(
x^{1},...,x^{q}\right) \in \mathbb{R}^{q}$ for a vector $h=\left(
h^{1},...h^{q}\right) ^{\prime }$ with positive finite components is defined
as the set $C\left( x,h\right) =$ $\left\{ X\in \mathbb{R}^{q}:\left\vert
X^{l}-x^{l}\right\vert \leq h^{l},l=1,\cdots ,q\right\} .$ With the
distribution function of $X$ given by $F_{X}$ the corresponding probability
measure is $P_{X}\left( C\left( x,h\right) \right) =\int_{C\left( x,h\right)
}dF_{X}.$ The small cube probability permits us to extend the regression on
univariate metric spaces to $\Xi ^{\lbrack q]}$ where the probability
measure $P_{X}$ is defined. For the cuboid
\begin{equation}
C\left( x,h\right) =\left\{ X:\left\Vert X^{l}-x^{l}\right\Vert _{l}\leq
h^{l},l=1,\cdots ,q\right\} =\left\{ X:X^{l}\in B^{l}\left(
x^{l},h^{l}\right) ,l=1,\cdots ,q\right\} . \label{cuboid}
\end{equation}
the corresponding small cube probability is also denoted $P_{X}\left(
C\left( x,h\right) \right) .$\medskip
One of our contributions is the derivation of auxiliary technical results
that express moments for the multivariate kernels and related functions in
terms of the small cube probabilities without appealing to differentiability
on which previous multivariate derivations relied. The moments and moment
bounds are derived for general multivariate local functions under arbitrary
distributions over $\mathbb{R}^{q}$ or probability measures over $\Xi ^{
\left[ q\right] }$. Bounds on a moment functional expressed via power of
small cube probability pinpoint the rate of growth of the functional. These
results generalize the derivations for the univariate kernel used in
functional regression to the multivariate setting. The full details of these
auxiliary technical results are presented in the supplemental material
(Appendix A). The moment expressions could find use in other contexts, for
instance for local linear and local polynomial estimation in $\mathbb{R}^{q}$
or in products of suitable metric spaces, $\Xi ^{\left[ q\right] },$ kernel
estimation of distribution functions and conditional distributions in $
\mathbb{R}^{q}$ as well as to kernel regression of objects in metric spaces
on objects in products of spaces.
Implementation of the NW estimator relies on a tuning bandwidth parameter.
We show that in $\mathbb{R}^{q}$ a popular cross-validation method of
choosing a bandwidth with properties that were worked out for the absolutely
continuous (a.c.) case, has similar properties in some empirically relevant
classes of singular distributions with dimension-reducing singularity. We
also examine adaptive bandwidth selection for regressor distributions that
are represented by a mixture of a continuous distribution with some mass
points.
We provide simulation evidence on some important features of the behavior of
the NW estimator under possible singularity of the distribution of
regressors, $F_X$, in $\mathbb{R}^{q}$, in particular on the pointwise rate
of convergence and specific impact of mass points. We examine the behavior
of the NW estimator for models with dependence on both a functional object
in $\Xi ^{\left[ 1\right] }$ and a random variable.
The structure of the paper is as follows. Section 2 provides the set-up
suitable for the multivariate vector and functional regression highlighting
the probability measure for the regressor. Section 3 gives the asymptotic
normality results under the most general distributional assumptions. Section
4 discusses implementation, in particular, bandwidth selection. Section 5
provides a sketch of the simulation results and Section 6 is devoted to the
empirical study. The supplementary material collects various auxiliary
results and the proofs as well as the details of the Monte Carlo simulations
and the empirical study.
\section{The set-up and assumptions}
This section provides the formula for the Nadaraya-Watson (NW) kernel
estimator over $\Xi ^{\left[ q\right] }$, introduces some useful notation
and gives formal assumptions. The distributional assumptions are very
general in that they do not restrict the distribution over $\mathbb{R}^{q}$
to have absolutely continuous components, and apply to the probability
measure over the multivariate metric space $\Xi ^{\left[ q\right] }$ for an
arbitrary $q.$
\begin{assumption}
\label{A.measure on prod} [Probability Measure] Given the metric measure
spaces $\Xi _{l}^{\left[ 1\right] },$ $l=1,...,q$ with corresponding
sigma-algebras and probability measures $P_{X^{l}}$ assume that the
sigma-algebra for $\Xi ^{\left[ q\right] }=\prod_{l=1}^{q}\Xi _{l}^{\left[ 1
\right] }$ is generated by the products of sets from sigma algebras for $\Xi
_{l}^{\left[ 1\right] }$ and a probability measure $P_{X}$ is defined on
this sigma algebra$;$ the mapping of $X=\left( X^{1},...,X^{q}\right) $ into
each of the components $X^{l}\in \Xi _{l}^{\left[ 1\right] }$ is measurable $
\left( P_{X^{l}}\right) $ with respect to the joint measure.
\end{assumption}
In the product space $\Xi ^{\left[ q\right] }$ we define a vector $w$ as $
\left( w^{1},...,w^{q}\right) ^{T}$ where each component is in the
corresponding space, thus for $\Xi =\mathbb{R}^{q},$ $w$ is a $q$
-dimensional vector of reals, in $\Xi =\Xi ^{\lbrack q]}$ each $w^{l}\in \Xi
_{l}^{[1]},\ l=1,...,q.$ The bandwidth vector is $h=\left(
h^{1},...,h^{q}\right) \in \mathbb{R}^{q}$ with $0<\underline{h}=\min
\left\{ h^{1},...,h^{q}\right\} >0$ and $\bar{h}=\max \left\{
h^{1},...,h^{q}\right\} $. We use the same notation $\left\Vert \cdot
\right\Vert $ for the absolute value of a scalar in $\mathbb{R}^{1}$, the
Euclidean norm for a vector in $\mathbb{R}^{q}$ or norm for a function in $
\Xi =\Xi ^{\lbrack 1]},$ with $\Xi ^{\lbrack 1]}$ a Banach space, or metric
(semi-metric) in a metric space $\Xi ^{\lbrack 1]}$; where the meaning is
not clear from the context we shall specify.
\subsection{The Nadaraya-Watson (NW) estimator}
The NW estimator, $\widehat{m}(x)$, for a sample $\left\{ \left(
Y_{i},X_{i}\right) \right\} _{i=1}^{n}$ generated by (\ref{model}) is
defined below. Generically the argument of the kernel function is
\begin{align}
& W_{X}\left(x\right ) \notag \\
=& \left\{
\begin{array}{lll}
h^{-1}(x-X) & =\left( \left( h^{1}\right) ^{-1}(x^{1}-X^{1}),\cdots ,\left(
h^{q}\right) ^{-1}(x^{q}-X^{q})\right) & \relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{on }\Xi =\mathbb{R}^{q} \\
h^{-1}\left\Vert x-X\right\Vert & =\left( \left( h^{1}\right)
^{-1}\left\Vert x^{1}-X^{1}\right\Vert _{1},\cdots ,\left( h^{q}\right)
^{-1}\left\Vert x^{q}-X^{q}\right\Vert _{q}\right) & \relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{on }\Xi =\Xi ^{
\left[ q\right] }.
\end{array}
\right. \label{W(x)}
\end{align}
The NW estimator is given by
\begin{eqnarray}
\widehat{m}\left( x\right) &=&B_{n}^{-1}\left( x\right) A_{n}\left( x\right)
,\relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{ with} \label{NW} \\
B_{n}\left( x\right) &=&\frac{1}{n}\sum_{i=1}^{n}K\left( W_{i}(x)\right) ;
\relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{ }A_{n}\left( x\right) =\frac{1}{n}\sum_{i=1}^{n}K\left(
W_{i}(x)\right) Y_{i}.\;\;\;\;\;\;\; \label{B, A}
\end{eqnarray}
where $K(W_{i}(x))=K(W_{X_{i}}\left(x\right) )$ is a multivariate
(non-negative) kernel function and $h$ usually depends on $n$; $x$ such that
at least for some $i$ we have that $K\left( W_{i}(x)\right) >0$. The kernel
function $K$ and bandwidth vector $h$ determine the properties for the NW
estimator. In the metric space $\Xi =\Xi ^{\lbrack 1]}$ the kernel function $
K$ is defined for a univariate non-negative argument; in the case $\Xi =\Xi
^{\lbrack q]}$ with $q>1$ different bandwidths could appear for the
different components $W_{X}^{l}\left( x\right) ,l=1,..,q.$ With a symmetric
kernel on $\mathbb{R}^{q}$ we can just write $W_{X}\left( x\right)
=h^{-1}\left\Vert x-X\right\Vert $ for any $\Xi .$
\subsection{The kernel}
We restrict the multivariate kernel functions on $\mathbb{R}^{q}$ to have
bounded support and be suitably differentiable in the interior.
Let $I_{{\Greekmath 0118} }$ denote any subset of the set $\left\{ 1,...,q\right\} $ $\ $
of consecutive non-negative integers; there are $2^{q}$ such subsets
including the empty set $\varnothing ;$ denote by $q\left( {\Greekmath 0118} \right) $ the
cardinality of the set $I_{{\Greekmath 0118} }=\left\{ j_{1},...,j_{q({{\Greekmath 0118} })}\right\} $
with $j_{1}<...<j_{q({{\Greekmath 0118} })}.$ We use $\prod_{j\in I_{{\Greekmath 0118} }}^{{}}\left(
\partial _{j}\right) $ to denote an operator that, when applied to a
differentiable function $g\left( z\right) =g\left( z^{1},...,z^{q}\right) $
at $z$, maps it to its partial derivative for $j_{1}<...<j_{q({{\Greekmath 0118} })}$,
that is
\begin{equation*}
\left( \prod_{j\in I_{{\Greekmath 0118} }}^{{}}\left( \partial _{j}\right) \right) g\left(
z\right) =\frac{\partial ^{q({{\Greekmath 0118} })}}{\partial _{j_{1}}...\partial _{j_{q({
{\Greekmath 0118} })}}}g\left( z\right) .
\end{equation*}
We call a function $g(z)$ \textquotedblleft sufficiently
differentiable\textquotedblright\ if for any set $I_{{\Greekmath 0118} }$ the derivative $
\left( \prod_{j\in I_{{\Greekmath 0118} }}^{{}}\left( \partial _{j}\right) \right) g\left(
z\right) $ exists and is continuous at any point on the interior of its
support.\medskip
The following assumption is made on the kernel function.
\begin{assumption}
\label{A.kernel} [Kernel]
\begin{itemize}
\item[(a)] The kernel function $K\left( w\right) =K\left(
w^{1},...,w^{q}\right) $ is a sufficiently differentiable density function.
\item[(b)] $K\left( w\right) $ is non-negative; $K\left( w\right) $ is
non-increasing for $w:w^{j}\geq 0, \ j=1,...,q$.
\item[(c)] $K\left( w\right) $ is either symmetric (with respect to zero)
with support on $\left[ -1,1\right] ^{q}$ or $K\left( w\right) $ is
supported on $\left[ 0,1\right] ^{q}$.
\item[(d)] $K\left( w\right) $ satisfies $K\left({\Greekmath 0113}\right ) >0$ where $
{\Greekmath 0113}=(1,...,1)^{\prime }$.
\end{itemize}
\end{assumption}
Assumptions \ref{A.kernel}(a--c) are satisfied by the commonly employed
product kernels of Epanechnikov or quartic kernels. Assumption \ref{A.kernel}
(d) is not usual for kernel regression on $\mathbb{R}^{q};$ in the context
of univariate functional regression\ it is satisfied by a Type I kernel
defined in Ferraty and Vieu (2006) as \thinspace $K:C_{1}I_{\left[ 0,1\right]
}\leq K\leq C_{2}I_{\left[ 0,1\right] }$ with some $0<C_{1}\leq C_{2}<\infty
.$ Condition (d) in conjunction with (a-c) provides the same type of
univariate kernel. Extended to a multivariate setting it can be said that a
kernel that satisfies Assumption \ref{A.kernel} (a-d) is a type I kernel.
The uniform kernel is an example. The functional regression literature
demonstrates that with kernels of type I asymptotic normality can be
established in more general univariate settings. As commonly used in $
\mathbb{R}^{q}$ kernels are not of type I, the asymptotic normality results
are given separately to apply under Assumption \ref{A.kernel}(a,b,c) and
under the full Assumption \ref{A.kernel}.
\subsection{Additional Assumptions}
Consider the process $\left\{ \left( X_{i},Y_{i}\right) \right\} _{i\in
\mathbb{N}}.$ An i.i.d sequence would provide the simplest characterization,
but strong mixing makes it possible to extend the results to time series
data. Denote by $\mathcal{F}_{a}^{b}$ the sigma algebra generated by $
\left\{ \left( X_{i},Y_{i}\right) \right\} _{i=a}^{b}.$ Define
\begin{equation*}
{\Greekmath 010B} \left( l\right) =\underset{t}{\sup }\underset{A\in \mathcal{F}
_{-\infty }^{t};B\in \mathcal{F}_{t+l}^{\infty }}{\sup }\left\vert P\left(
AB\right) -P\left( A\right) P\left( B\right) \right\vert .
\end{equation*}
Recall that the process is strong mixing if ${\Greekmath 010B} \left( l\right)
\rightarrow 0$ as $l\rightarrow \infty .$
\begin{assumption}
\label{A.mom} [Data Generating Process and Moments]
\begin{itemize}
\item[(a)] The sequence $\left\lbrace(Y_{i},X_{i})\right\rbrace$ for $
i=1,\cdots,n$ with $Y_{i}\in \mathbb{R};X_{i}\in \Xi ^{\left[ q\right] }$ is
stationary and strong mixing with ${\Greekmath 010B} \left( l\right) $ that satisfies
for some ${\Greekmath 0110} >0$
\begin{equation*}
{\Greekmath 010B} \left( l\right) <C l^{-{\Greekmath 0114} };\ \relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{ }{\Greekmath 0114} >\frac{2\left(
2+{\Greekmath 0110} \right) }{{\Greekmath 0110}}.
\end{equation*}
\item[(b)] $E(u|X=x)=0;$ ${\Greekmath 0116} _{2}\left( x\right) =E\left( u^{2}|X=x\right) $
satisfies $0<L_{{\Greekmath 0116} _{2}}<{\Greekmath 0116} _{2}\left( x\right) <M_{{\Greekmath 0116} _{2}}<\infty ,$ $
{\Greekmath 0116} _{2}\left( x\right) $ is continuous in the neighborhood of $x.$
\item[(c)] $E\left\vert Y_{i}\right\vert ^{2+{\Greekmath 0110} }<\infty $ and $
E(\left\vert u \right\vert^{2+{\Greekmath 0110}}|X=x)<\infty.$
\item[(d)] For $x\in \Xi ^{\left[ q\right] }$ and $i\neq j$ the bivariate
function
\begin{equation*}
{\Greekmath 0116} \left( x_{1},x_{2}\right) =E\left( \left\vert u_{i}u_{j}\right\vert
|X_{i}=x_{1},X_{j}=x_{2}\right)
\end{equation*}
is continuous in a neighborhood of the point $\left( x,x\right) \in \Xi ^{
\left[ q\right] }\times \Xi ^{\left[ q\right] }.$
\item[(e)] The conditional expectation $E\left( \left\vert
Y_{i}Y_{j}\right\vert |X_{i},X_{j}\right) \leq C<\infty $ for all $i,j.$
\end{itemize}
\end{assumption}
\noindent The assumption requires a polynomial bound on the rate of decline
of the mixing coefficient with a link to the moment of $Y;$ it is similar to
those in Masry (2005) and Hong and Linton (2020).
\begin{assumption}
\label{A.mx} [Conditional mean] The function $m\left( x\right) $ on the
space $\Xi ^{\left[ q\right] }$ is such that
\begin{equation*}
\left\vert m\left( x\right) -m\left( z\right) \right\vert \leq M_{\Delta m}
\underset{l}{\max }\left\Vert x^{l}-z^{l}\right\Vert _{l}^{{\Greekmath 010E} };\relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{ }
{\Greekmath 010E} >0.
\end{equation*}
\end{assumption}
Assumption \ref{A.mx} requires Holder continuity of $m\left( x\right) ;$ it
would follow from differentiability or Lipschitz continuity in $\mathbb{R}
^{q}$ with ${\Greekmath 010E} =1$. In the above assumptions, and below, $L$ and $M$
denote lower and upper bounds of functions where the subscript typically
denotes the function whose bounds are provided. The bounds could depend on
the point $x.$
\subsection{The probability measures}
For the probability measure $P_{X}$ on a generic space $\Xi ,$ that could
coincide with $\mathbb{R}^{q},$ $\Xi ^{\left[ 1\right] },$ or $\Xi ^{\left[ q
\right] },$ any point $x\in \Xi $ is a point of support if for $\underline{h}
>0$ {the measure $P_{X}\left( C\left( x,\underline{h}\right) \right) >0.$}
\subsubsection{Measures on $\mathbb{R}^{q}$}
By the Lebesgue decomposition, the distribution $F_{X}$ on $\mathbb{R}^{q}$
can be represented as a mixture of an absolutely continuous distribution, $
F^{a.c.},$ a singular distribution (the distribution function is continuous
but there is no function that integrates to it), $F^{s},$ and a discrete
distribution, $F^{d}:$
\begin{equation*}
F_{X}\left( x\right) ={\Greekmath 010B} _{1}F^{a.c.}\left( x\right) +{\Greekmath 010B}
_{2}F^{s}\left( x\right) +{\Greekmath 010B} _{3}F^{d}\left( x\right) ;{\Greekmath 010B} _{l}\geq
0,\ {l=1,2,3};\ \sum\nolimits_{l=1}^{3}{\Greekmath 010B} _{l}=1.
\end{equation*}
In a multivariate setting as soon as at least one variable is continuously
distributed, mass points do not arise and the joint distribution is a
continuous function, but with some discrete components or mass points in
some of the continuous components the distribution can no longer be
absolutely continuous and is singular. In many applications at least one of
the variables is assumed continuous and in a semiparametric regression often
an index model is assumed (single index in Ichimura, 1993; multiple index in
Donkers and Schafgans, 2008) to avoid singularity as well as to reduce
dimensionality of the model. In a general multivariate distribution the
presence of singularity achieves reduction of dimension (see, e.g. examples
2-4 in Kankanala and Zinde-Walsh, 2024) that will have a similar beneficial
effect on the convergence of the kernel estimator.
\subsubsection{Measures on metric spaces and products}
The discussion in this section applies to the space $\mathbb{R}^{q}$ as a
special case. Particular classes of probability measures considered in
univariate functional regression (e.g. Ferraty et al., 2007, Ferraty and Vieu, 2006) have a small ball probability centered at a point $x$ of support
either with a polynomial (fractal) rate of decline $P_{X}(B(x,h))\sim
C\left( x\right) h^{{\Greekmath 011C} }>0$ (${\Greekmath 011C} >0),$ or with an exponential type rate
of decline $P_{X}(B(x,h))\sim C\left( x\right) \exp (-h^{-{\Greekmath 011C} _{1}}\log
h^{-{\Greekmath 011C} _{2}})$ (${\Greekmath 011C} _{1}>0,{\Greekmath 011C} _{2}>0)$ as $h\rightarrow 0.$ This
characterization can be applied to the multivariate setting by replacing $
B(x,h)$ with the cuboid and the univariate bandwidth in the rate with $\bar{h
}.$ The exponential rate of decay of the small ball probability requires a
type I kernel and leads to slow convergence for the estimators (curse of
dimensionality). There are ways to mitigate the curse of dimensionality
arising from such exponential decay. It is common to apply finite
dimensional approximation of these functionals as suggested in Gasser et al.
(1998). Indeed, the case where functional data can be accurately
approximated in a finite dimensional space is not rare (corresponds to
observation of smooth curves with common shape) as noted by Ferraty and Nagy
(2022).\vspace{0.1in}
Kernels of type I play an important role in establishing pointwise
asymptotic normality in the absence of any restrictions on the decline of
the small cube measure. For kernels that may not be of type I sufficient
conditions on the shrinkage of the probability measure as $h\rightarrow 0$
were proposed in Assumption $H_{3}$ in Ferraty et al. (2007), Ferraty and
Vieu (2006) and were referred to in various subsequent papers on functional
regression, e.g. Hong and Linton (2020). The assumption below generalizes these
conditions to apply to $C\left( x,h\right) $ on $\Xi ^{\left[ q\right] };$
the assumption is both necessary for the conditions to hold (see the
supplementary material, Appendix B) and at the same time sufficient for the convergence
results.
\begin{assumption}
\label{A.dist} [Small ball probability measure] {Given any point $x\in \Xi ^{
\left[ q\right] }$ in the support of the probability measure $P_{X}$} for
all $h$ with $\underline{h}>0$ and for some $0<{\Greekmath 0122} <1,$ there is a
constant $1<C_{{\Greekmath 0122} }\,<\infty $ such that
\begin{equation}
\relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{ }\frac{P_{X}(C(x,h))}{P_{X}(C(x,{\Greekmath 0122} h))}<C_{{\Greekmath 0122}
}<\infty . \label{general band}
\end{equation}
\end{assumption}
\begin{definition}
$\mathcal{D}$ is the class of probability measures that satisfies (\ref
{general band}).\footnote{
Condition (\ref{general band}) is equivalent to the doubling property (e.g.
Vol'berg, Konyagin, 1988) that states that (\ref{general band}) applies with
${\Greekmath 0122} =1/2.$ Indeed for any ${\Greekmath 0122} $ there are positive
integers ${\Greekmath 0114} _{1},{\Greekmath 0114} _{2}:$ ${\Greekmath 0122} \geq 2^{-{\Greekmath 0114} _{1}}$ and $
2^{-1}\geq {\Greekmath 0122} ^{{\Greekmath 0114} _{2}}.$ If the measure is doubling for
constant $C_{1/2}$, then (\ref{general band}) holds for $C_{{\Greekmath 0122}
}=C_{1/2}^{{\Greekmath 0114} _{1}};$ if ( \ref{general band}) holds, then the constant
for doubling is $C_{1/2}=C_{{\Greekmath 0122} }^{{\Greekmath 0114} _{2}}.$ We introduce the
form (\ref{general band}) in case there is a preference for some $
{\Greekmath 0122} .$\medskip}
\end{definition}
A polynomial decay condition places a measure into class $\mathcal{D}.$
Indeed if the small cube probability satisfies
\begin{equation}
0<L_{P}\left( x\right) (2\underline{h})^{s\left( x\right) q}\leq
P_{X}(C(x,h))\leq M_{P}\left( x\right) \left( 2\bar{h}\right) ^{s\left(
x\right) q}<\infty \label{LM F bounds}
\end{equation}
{where for some $c$, $H(x),$ $1\leq c<\infty ,$ $0<H\left(x\right) <\infty $
and $\bar{h}=c\underline{h}<H\left( x\right) ,$ $0\leq s\left( x\right) \leq
1$ and $M_{F}\left( x\right) /L_{F}\left( x\right)<B<\infty $ at all points
of support $x$, (\ref{general band}) holds with $C_{{\Greekmath 0122} }=B\left(
c/{\Greekmath 0122} \right) ^{q}.$}
Condition (\ref{LM F bounds}) applies quite widely and holds for many
distributions of regressors used in econometric models. In $\mathbb{R}^{q}$
it is satisfied by any absolutely continuous distribution with\ a positive
bounded density function $f_X\left( x\right) $ where $M_{P}\left( x\right)
\geq \underset{\tilde{x}\in C\left( x,H\right) }{\sup }f_X\left( \tilde{x}
\right) ;$ $L_{F}\left( x\right) =\underset{\tilde{x}\in C\left(
x,H/c\right) }{\inf }f_X\left( \tilde{x}\right) $ and $s\left( x\right) =1.$
If $x$ is an isolated mass point then (\ref{LM F bounds}) applies with $s=0.$
If $X$ has a linear structure with $r$ common factors, the probability
measure is singular and satisfies (\ref{LM F bounds}) with $s\left( x\right)
=s=\frac{q}{r}.$ For a fractal distribution that is singular with constant $
s,$ $0<s<1$, the bounds also apply.
Condition (\ref{LM F bounds}) is satisfied by the general class of Ahlfors
(1966) regular (A-r) distributions common in statistics, where for this
class $s\left( x\right) =s,$ and $L_{P}\left( x\right) =L$ and $M_{P}\left(
x\right) =M$ are constants, as well as by a finite mixture of such
distributions (as proved in the supplementary material, Appendix B). Thus an absolutely
continuous distribution ($s=1)$ or, more generally, a measure given by a
continuous possibly singular distribution function that satisfies (\ref{LM F
bounds}) contaminated with some mass points ($s=0$) is in $\mathcal{D}$
\textbf{$;$ }this applies to the empirical example examined here, ensuring
the pointwise asymptotic normality of the NW estimator with standard kernels.
\subsubsection{Joint measure}
Consider the product space $\Xi ^{\left[ 2q\right] }=\Xi ^{\left[ q\right]
}\times \Xi ^{\left[ q\right] };$ the measure on this product space has
marginals $P_{X}$ on each $\Xi ^{\left[ q\right] }$ (see, e.g., Pollard,
2001). {The joint measure $P_{s,t}\left( C\left( x,h\right) \times C\left(
x,h\right) \right) ,$ defined as $\Pr \left( X_{t}\in C\left( x,h\right)
,X_{s}\in C\left( x,h\right) \right) ,$ is a product of the measures of the
cuboid in the case of independency.} With dependence an additional
assumption is made on how the joint measure relates to the small cuboid
measure. {We provide the same assumption as in e.g. Masry (2005) and Hong
and Linton (2020) for the small cube probability.}
\begin{assumption}
\label{A.Jointmeasure} [Joint Measure] The joint measure $P_{s,t}\left(
C\left( x,h\right) \times C\left( x,h\right) \right) $ is such that for some
$0<M_{FF}<\infty $
\begin{equation}
\underset{t\neq s}{\sup }P_{s,t}\left( C\left( x,h\right) \times C\left(
x,h\right) \right) \leq M_{FF} \left( \overset{\textcolor{white}.}{P}
_{X}\left( C\left( x,h\right) \right) \right) ^{2}. \label{general product}
\end{equation}
\end{assumption}
\section{Asymptotic normality of the NW estimator}
Consider the NW estimator as given by (\ref{NW}), (\ref{B, A}). As the
sample size increases the bandwidths {are assumed} go to zero. For $\Xi =
\mathbb{R}^{q}$ the denominator, $B_{n}(x)$ is proportional to the usual
kernel density estimator, given by $h^{-q}B_{n}(x),$ at point $x.$ {When the
density, $f_{X}\left( x\right) $, exists and is continuous}, the estimator $
h^{-q}B_{n}(x)$ consistently estimates $f_{X}\left( x\right) ,$ but if the
density does not exist, $h^{-q}B_{n}\left( x\right) $ diverges to infinity.
Consistency of the NW estimator $\widehat{m}\left( x\right) $ over a
univariate metric space was established in Gy\"{o}rfi et al. (2002), the
limit distribution in Ferraty et al. (2007), Masry (2005) and Geenens (2015).
The key to the asymptotic normality result is the derivation of the moments
for multivariate functions of the form $g(X)K^{m}\left( h^{-1}\left\Vert
x-X\right\Vert \right) $ for general probability measures and establishing
lower and upper bounds (derivations in the supplementary material, Appendix B). The
bounds provide expressions in terms of the small cube probability:
\begin{equation}
L_{EgK^{m}}\left( x\right) {P_{X}(C(x,h))}\leq \left\vert E\left[
g(X)K^{m}\left( h^{-1}\left\Vert x-X\right\Vert \right) \right] \right\vert
\leq M_{EgK^{m}}\left( x\right) {P_{X}(C(x,h))} \label{EK^mg}
\end{equation}
with constants $L_{EgK^{m}}\left( x\right) $ and $M_{EgK^{m}}\left( x\right)
$ at $x.$ Most important, (\ref{EK^mg}) provides a lower bound on $
EB_{n}\left( x\right) =EK\left( K^{m}\left( h^{-1}\left\Vert x-X\right\Vert
\right) \right) $, given by $L_{EK}P_{X}\left( C\left( x,h\right) \right) $,
with appropriate conditions for $L_{EK}$ to be strictly positive to ensure
that the denominator of the NW estimator is such that it exists and the
limit does not blow up. Type I kernel automatically entails that $L_{EK}>0$,
but for kernels such as Epanechnikov the bound requires Assumption {\ref
{A.dist}. }With $g\left( X\right) $ that is continuous at $x$
\begin{equation*}
E\left[ g(X)K^{m}\left( h^{-1}\left\Vert x-X\right\Vert \right) \right]
=g\left( x\right) E\left[ K^{m}\left( h^{-1}\left\Vert x-X\right\Vert
\right) \right] \left( 1+o\left( 1\right) \right) .
\end{equation*}
These moment expressions for distributions over $\mathbb{R}^{q}$ hold under
the standard assumptions of existence and continuity of (bounded) density $
f_X,$ and the function $g$, {where}
\begin{equation}
E\left[ g(X)K^{m}\left( h^{-1}\left\Vert x-X\right\Vert \right) \right]
=\prod\limits_{i=1}^{q}\left( -h^{i}\right) g\left( x\right) f_X\left(
x\right) \int K^{m}\left( v\right) dv\left( 1+o\left( 1\right) \right) ,
\label{EgfK}
\end{equation}
with more details about the $o\left( 1\right) $ term under smoothness of $f_X
$ (see, e.g. derivations in Li and Racine, 2007). Once the moments and the
bounds are derived, the proofs of asymptotic normality proceed along similar
lines to those in Masry (2005).
The point-wise limit normality is provided in the next theorem under two
alternative types of conditions: (i) with type I kernel without imposing
further constraints on $F_X$, and (ii) not imposing the type I kernel but
with the distributional Assumption \ref{A.dist}. Denote the bias of the
estimator given $x$, $E\left( \widehat{m}\left( x\right) \right) -m\left(
x\right) ,$ by $bias\left( \widehat{m}\left( x\right) \right) .$ The
difference $\widehat{m}\left( x\right) -m\left( x\right) $ is delivered by $
\frac{A_{n}^{c}(x)}{B_{n}(x)}$ with the \textquotedblleft
centered\textquotedblright\ $A_{n}^{c}(x)=A_{n}(x)-m\left( x\right)
B_{n}(x). $
\setcounter{theorem}{0}
\begin{theorem}
\label{T.1} {Under either of the following sets of assumptions (i)
Assumptions \ref{A.measure on prod}-\ref{A.mx} and \ref{A.Jointmeasure} or
(ii) Assumptions \ref{A.measure on prod}, \ref{A.kernel}(a-c), \ref{A.mom}-
\ref{A.Jointmeasure} } for $h\rightarrow 0$ as $n\rightarrow \infty $ such
that $nP_{X}\left( C\left( x,h\right) \right) \rightarrow \infty $
\begin{itemize}
\item[(a)]
\begin{equation*}
\frac{\sqrt{n}E\left[ K\left( h^{-1}\left\Vert x-X\right\Vert \right)\right]
}{\sqrt{{\Greekmath 0116} _{2}\left( x\right) E\left[ K^{2}\left( h^{-1}\left\Vert
x-X\right\Vert\right) \right] }}\left( \widehat{m}\left( x\right) -m\left(
x\right)-bias(\widehat{m}(x)\right) \rightarrow _{d}Z\sim N\left( 0,1\right)
;
\end{equation*}
\item[(b)] the rates are\
\begin{eqnarray*}
&& bias(\widehat{m}(x) =O(\bar{h}^{{\Greekmath 010E} })+O\left( nP_{X}\left(
C(x,h)\right) \right) ^{-1}; \\
&& \frac{\sqrt{n}E \left[ K\left( h^{-1}\left\Vert x-X\right\Vert \right)
\right] }{ \sqrt{{\Greekmath 0116}_{2}\left( x\right) E\left[ K^{2}\left( h^{-1}\left\Vert
x-X\right\Vert\right)\right] }} \simeq O\left( (nP_{X}\left( C\left(
x,h\right)\right) ^{1/2}\right) .
\end{eqnarray*}
\item[(c)] for $h$ such that $\bar{h}^{2{\Greekmath 010E} }\left( nP_{X}\left(
C(x,h)\right) \right) \rightarrow 0$
\begin{equation*}
\frac{\sqrt{n}E\left[K\left( h^{-1}\left\Vert x-X\right\Vert \right)\right]
}{\sqrt{{\Greekmath 0116} _{2}\left( x\right) E\left[ K^{2}\left( h^{-1}\left\Vert
x-X\right\Vert \right)\right] }}\left( \widehat{m}\left( x\right) -m\left(
x\right) \right) \rightarrow _{d}Z\sim N\left( 0,1\right) .
\end{equation*}
\end{itemize}
\end{theorem}
\noindent \textbf{Remarks.}
\begin{enumerate}
\item A sequence of bandwidths at $x$ that satisfy the conditions of the
theorem always exists. Indeed, whatever the rate of monotonic decline in $
P_{X}(C(x,h))$ as $h\rightarrow 0$ for $n\rightarrow \infty $ a sequence of $
h$ that depends on $n$ such that $nP_{X}\left( C\left( x,h\right) \right)
\rightarrow \infty $ always exists. T{he rate for the bias of $\widehat{m}
\left( x\right) $ in $\Xi ^{\left[ q\right] }$ is established in the theorem
as $O\left( \bar{h}^{{\Greekmath 010E} }\right) +O\left( \left( nP_{X}\left(
C(x,h)\right) \right) ^{-1}\right) .$ For the bias (squared) to disappear in
the limit }$\bar{h}^{2{\Greekmath 010E} }P_{X}\left( C(x,h)\right) n$ needs to go to
zero. If $P_{X}\left( C\left( x,h\right) \right) \rightarrow 0$ a bandwidth
sequence that simultaneously satisfies $nP_{X}\left( C\left( x,h\right)
\right) \rightarrow \infty $ and $\bar{h}^{2{\Greekmath 010E} }P_{X}\left(
C(x,h)\right) n\rightarrow 0$ can always be found; when $x$ is a mass point $
P_{X}\left( C\left( x,h\right) \right) $ will be bounded from below, but
selecting $h=o\left( n^{-1/2{\Greekmath 010E} }\right) $ for such a point makes the
bias term go to zero.
\item The assumptions of Theorem \ref{T.1} and the moment computations in
the supplementary material (Appendix B) imply that $E\left[ K\left( h^{-1}\left\Vert
x-X\right\Vert \right)\right]$ has the same rate as $P_{X}(C(x,h))$ while $
varA_{n}^{c}(x)$ declines at the rate ${P_{X}\left( C\left( x,h\right)
\right) }/{n}.$ The rate for the asymptotic variance for $\widehat{m}\left(
x\right) $ equals $\left( nP_{X}\left( C\left( x,h\right) \right) \right)
^{-1}$ (this goes to zero).
\item The limit result shows that when density exists for a distribution on $
\mathbb{R}^{q},$ the standard convergence rate $n^{1/2}h^{q/2}$ applies
since then $P_{X}\left( C\left( x,h\right) \right) =O\left( h^{q}\right) .$
This rate holds even when the density is discontinuous. Without the usual
smoothness assumptions made in the literature, statistical guarantees for
the rate\textbf{\ }and for asymptotic normality are thus shown to hold.
\item If there is singularity at the point $x$ that satisfies (\ref{LM F
bounds}) with $s<1,$ then the rate is $n^{1/2}h^{sq/2}$, which is faster
than in the absolutely continuous case $\left(
n^{1/2}h^{sq/2}>n^{1/2}h^{q/2}\right) $, mitigating somewhat the
\textquotedblleft curse of dimensionality\textquotedblright . When $x$ is an
isolated mass point then at that point the parametric rate $n^{1/2}$ holds.
\item Under continuous differentiability the rate of the bias can be reduced
by employing a local linear estimator (see, e.g. the standard derivations in
Li and Racine, 2007, and for univariate functional regression in Ferraty and
Nagy, 2022). Establishing the distributional properties of the local linear
estimator with arbitrary probability distributions in $\mathbb{R}^{q}$ and
multivariate probability measures in a metric space can proceed similarly,
but requires stronger assumptions.
\end{enumerate}
The convergence rate in (c) of Theorem \ref{T.1} is $O\left( (nP_{X}\left(
C\left( x,h\right) \right) ^{-1/2}\right).$\footnote{
This convergence rate obtains under $h\rightarrow 0.$ In the presence of an
irrelevant regressor, say $x^{(2)}$, such that $m\left( x\right) =m\left(
x^{\left( 1\right) }\right) $ for all $x=\left( x^{\left( 1\right)
},x^{\left( 2\right) }\right) ,$ this requirement can be restricted to the
function $m\left( x^{\left( 1\right) }\right) $ with the irrelevant $
x^{\left( 2\right) }$ eliminated. For the estimator this elimination can be
achieved by setting the bandwidth on components of $x^{\left( 2\right) }$ to
be larger than the range of those variables, possibly infinite.} Existence
of a limit variance ${\Greekmath 011B} _{\widehat{m}\left( x\right) }^{2}$ requires
that $\left( nP_{X}\left( C\left( x,h\right) \right) \right)\frac{\left[
EK\left( h^{-1}\left\Vert x-X\right\Vert \right) \right]^2}{{{\Greekmath 0116} _{2}\left(
x\right) E\left[ K^{2}\left( h^{-1}\left\Vert x-X\right\Vert\right)\right] }}
$ converges. Without additional assumptions it is possible that the ratio
does not converge; see example in the supplementary material (Appendix B) that provides a
case when convergence does not hold; this happens when the small cube
probability declines very rapidly and the kernel is not uniform. Suitable
additional assumptions on the distribution, such as $H_{3}$ in Ferraty et
al. (2007) and Condition 3(i) in Masry (2005) and similar ones in subsequent
papers provide restrictions on the probability measure on $\Xi ^{\left[ 1
\right] }$ that are sufficient for the convergence. Generally, one needs to
ensure that the limits given below on the expectation of the kernel function
and its square hold.\footnote{
This implies that the extra condition is also required for the Corollary 1
of Hong and Linton (2020).}
\begin{assumption}
\label{A.lim} As $n\rightarrow \infty ,$ $h\rightarrow 0$
\begin{equation*}
\left( P_{X}(C(x,h))\right) ^{-1}E\left[K^{s}\left( h^{-1}\left\Vert
x-X\right\Vert \right) \right] \rightarrow \bar{B}_{s}\left( x\right) ;s=1,2.
\end{equation*}
\end{assumption}
This assumption holds quite widely. From the moment expressions it can
easily be shown that it holds for the uniform kernel without any additional
distributional assumptions. In the case of continuous density it holds by
virtue of (\ref{EgfK}) with
\begin{equation}
\bar{B}_{1}\left( x\right) =f_X\left( x\right) \int K\left( v\right)
dv,\quad \bar{B}_{2}\left( x\right) =f_X\left( x\right) \int K^{2}\left(
v\right) dv. \label{lim for a a c}
\end{equation}
{Suppose that singularity arises, because of combining discrete and
continuous variables in $\mathbb{R}^{q}$ or functional dependence between
the regressors, that restrict the support of the distribution to be in some
subspace of dimension $r<q$, $V\left(r\right) \subset $ $\mathbb{R}^{q}.$}
If the distribution on $V\left( r\right) $ is absolutely continuous with a
continuous density, then derivations provide similar limits to (\ref{lim for
a a c}) with integration over $V\left( r\right) $ and density restricted to $
V\left( r\right) .$
Define now
\begin{equation*}
{\Greekmath 010B} \left( n,h\right) ={nP_{X}\left( C\left( x,h\right) \right) };\quad
{\Greekmath 011B} _{\widehat{m}\left( x\right) }^{2}={\Greekmath 0116} _{2}(x)\bar{B}_{2}(x)/(\bar{B}
_{1}(x))^{2}\ .
\end{equation*}
\begin{theorem}
\label{T.2} Under the conditions of Theorem \ref{T.1} and Assumption \ref
{A.lim} with ${\Greekmath 010B} \left( n,h\right) \rightarrow \infty $ and for $h$ such
that ${\Greekmath 010B} \left( n,h\right) \bar{h}^{2{\Greekmath 010E}}\rightarrow 0$
\begin{equation*}
\sqrt{{\Greekmath 010B} \left( n,h\right) }\left( \widehat{m}\left( x\right)
-m\left(x\right) \right) \rightarrow _{d}N\left( 0,{\Greekmath 011B} _{\widehat{m}
\left( x\right) }^{2}\right) .
\end{equation*}
\end{theorem}
This limit extends the results that were obtained in the literature on
kernel estimation in $\mathbb{R}^{q}$ under smoothness assumptions on the
distribution $F_{X}.$ {For functional regression our assumptions are
comparable to those of Ferraty et al. (2007), Masry (2005), and subsequent
papers while they make the extension to multivariate functional
regression possible. }
\section{\protect\bigskip \textbf{Implementation and bandwidth selection}}
Estimation of $m\left( x\right) $ requires a selection of the kernel, $K,$
and bandwidth, $h.$ As may be clear from the results here and the
literature, type I kernel (such as the uniform) is preferred but other
kernels can also deliver asymptotic rates provided the small cube
probability does not decline exponentially fast. Aside from the estimator of
the conditional mean, estimators of variance and mean squared error are
needed to evaluate the performance of the estimator. While in the literature
on kernel regression on $\mathbb{R}^{q},$ the leading term of the limit
variance is expressed via the density function, often in the actual
implementation the corresponding estimators do not make use of plug-in
expressions, instead estimating the variance directly from the data and
possibly with bootstrap (see Hall and Horowitz, 2013).
Cross-validation procedures in popular statistical packages (such as R)
provide a single bandwidth (vector) that was shown to be consistent for the
\textquotedblleft optimal\textquotedblright\ bandwidth: minimizer of
weighted integrated mean squared error, WIMSE, (e.g. Li and Racine, 2007). The
proofs of consistency relied on absolute continuity of the regressors. The
consistency results extend to some classes of singular distributions.
WIMSE is defined for an absolutely continuous distribution with density
function $f_X\left( x\right) $ as
\begin{equation*}
\int E\left( \widehat{m}\left( x\right) -m\left( x\right) \right)
^{2}M\left( x\right) f_X\left( x\right) dx
\end{equation*}
with some weighting function $M\left( x\right) $ chosen to mitigate boundary
effects. The expression can be written with $dF_X$ replacing $f_X\left(
x\right) dx$ (valid in the case of singularity):
\begin{equation}
\int E\left( \widehat{m}\left( x\right) -m\left( x\right) \right)
^{2}M\left( x\right) dF_X=\int \left[ var\left( \widehat{m}\left(
x\right)\right) +bias^{2}\left( \widehat{m}\left( x\right) \right) \right]
M\left( x\right) dF_X. \label{WIMSE}
\end{equation}
This function depends on the bandwidth vector $h$ used in the estimator (see
the review of bandwidth selection methods, including cross-validation and
plug-in in K\"{o}hler et al., 2014). The \textquotedblleft
optimal\textquotedblright\ bandwidth vector $h^{0}$ is a minimizer of the
WIMSE criterion function based on a trade-off between the variance and bias
of the NW estimator.
In the cross-validation procedure the finite sample analogue of WIMSE
replaces the expectation by
\begin{equation*}
CV=n^{-1}\sum_{i=1}^{n}\left( Y_{i}-\widehat{m}_{-i}\left( X_{i}\right)
\right) ^{2}M\left( X_{i}\right)
\end{equation*}
employing the leave-one-out kernel estimator, $\widehat{m}_{-i},$ and
provides the bandwidth vector $h_{cv}$ by minimizing the CV criterion.
Hall et al. (2007) gave a general result about consistency of the
cross-validated bandwidth for regression over $\mathbb{R}^{q}$ with discrete
and continuous regressors, with some of the regressors possibly being
irrelevant. Their general result in Theorem 2.1 was obtained under a set of
assumptions that required independent identically distributed observations,
restrictions on the support of the probability measure, two continuous
derivatives for density, the regression function, and the conditional
variance of the error; in addition, for the $d$ continuous relevant
regressors $h^{o}=n^{-\frac{1}{4+rd}}a^{o}$ holds with the vector $a^{o}$
having unique, positive and finite components. This result was extended to
weakly dependent data by Li et al. (2009) under assumptions that replaced
the i.i.d. assumption by requiring strict stationarity and ${\Greekmath 010C} -$mixing
in the process for $\left\{ x,y\right\} $ and martingale difference error,
with suitable restrictions on the mixing parameters.
The result on the cross-validated bandwidth applies more widely. For instance, consider a singular distribution of $X\in \mathbb{R}^{q}\,$\ where there is a functional dependence among the continuous variables in the presence of possibly some discrete covariates such that the support of the distribution is restricted to a subspace $V\left( r\right) \subset \mathbb{R}^{q}$ of dimension $r<q$ represented by a union of affine subspaces. If, restricted to $V\left(
r\right) ,$ the distribution function is such that the conditions of Theorem
2.1 of Hall et al. (2007) or Theorem 1 of Li et al. (2009) are satisfied
(Assumption CV) then the conclusions of those theorems are valid and the
consistency of the bandwidth and automatic dimension reduction by smoothing
out irrelevant regressors hold for this singular distribution. More details
are provided in the supplementary material (Appendix B).
Importantly, no knowledge of $V\left( r\right) $ or $r$ is required. This
implies that for functionally dependent continuous regressors the knowledge
of the number of factors is not required for the consistency of the
cross-validated bandwidth or the automatic dimension reduction. We
conjecture that in many other cases with possible singularity the
cross-validation procedure will facilitate dimension reduction by smoothing
out irrelevant variables.
Bandwidth selection could benefit from adaptation to different types of
singularity. The treatment of adaptive bandwidth selection in the literature
(Fan and Gijbels, 1996, Sain, 1994, Demir et al., 2010) typically focuses on
adjusting the smoothing parameter to accommodate the varying data density,
but not dealing with singularity or mass points. Adaptive bandwidths can
provide a better fit of the criterion function by increasing the number of
observations used to estimate the function at a point of sparsity.\footnote{
Given some initial bandwidth $\tilde{h}$ and density estimate at this
bandwidth, $\hat{f}_{X},$ an adaptive bandwidth is defined for each point as
$h\left( X_{i}\right) =\tilde{h}\left( \frac{\hat{f}_{X}\left( X_{i}\right)
}{G}\right) ^{-{\Greekmath 010B} }$ where $G=\left( \prod \hat{f}_{X}\left(
X_{j}\right) \right) ^{1/n}$ is the geometric mean of the densities and $
{\Greekmath 010B} $ is typically selected to be $1/2.$ {One could construct $\hat{f}
_{X}(x)$ $\ $\ with a uniform kernel in which case it is identical to an
estimate of $P(C(x,\tilde{h}))$ by the proportion of observations in the $
\tilde{h}$ cuboid around $x.$}} Such bandwidths can similarly be constructed
for cases of singular distributions. But these adaptation procedures still
need to be investigated in the case of general mixtures of singular
distributions. However, singularity adaptation simplifies considerably for
the empirically important case of a mixture of an absolutely continuous
distribution with mass points, where the two levels of singularity can be
separated. The approach is detailed in the supplementary material (Appendix B).
\section{Simulations}
This section provides the highlights of various simulations that show
features of the finite sample performance of the NW estimator under
singularity. Additional details and features are in the supplementary
material (Appendix C).
\subsection{Univariate (Point mass example)}
In this example we consider the regression distribution with mass points.
Alongside we examine the trinormal mixture considered in Kotlyarova et al.
(2016), an a.c. distribution which represents features (high density
derivatives) that makes it comparable to a singular distribution.
The distribution with mass points, following Jun and Song (2019), is given
by
\begin{equation*}
F_{X}(x)=pF^{d}(x)+(1-p)\Phi (x)\quad \relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{with }p=0.2,
\end{equation*}
where $F^{d}$ is the discrete uniform distribution function with $
D=\{-1,0,1\}$ the set of mass points; $\Phi $ is the standard Gaussian
distribution function.
We simulated 500 random samples $\{(Y_{i},X_{i})\}_{i=1}^{n}$ using the
model
\begin{equation*}
Y_i=\sin (2.5X_i)+{\Greekmath 011B} {\Greekmath 0122}_i ,
\end{equation*}
for different sample sizes. The error $\{{\Greekmath 0122} _{i}\}_{i=1}^{n}$ is
drawn independently of the regressor and has a standard Gaussian
distribution; ${\Greekmath 011B} $ is selected to yield a given signal to noise ratio, $
snr$, here selected to equal one. We use the Epanechnikov kernel $K(u)=\frac{
3}{4}(1-{u^{2}})1(u^{2}\leq 1)$ and obtain the leave-one-out cross-validated
bandwidth.
We analyze the pointwise RMSE at a coarse grid of points across samples of
size $n$ equal to 50, 100, 200, 400, 800, 1600, 3200 based on 500
replications from the above DGP. To obtain empirical rates of convergence we
regress $\log (RMSE)$ on $\log (n)$ and a constant. The coefficient on $\log
(n)$ is the \textquotedblleft realized\textquotedblright\ rate of
convergence; for example if $RMSE\propto n^{-2/5}$ (univariate kernel
regression with smooth density and second order kernel) then $\log
(RMSE)={\Greekmath 010B} _{0}+{\Greekmath 010B} _{1}\log (n)$ and ${\Greekmath 010B} _{1}$ should be close to
-0.4.\footnote{
The authors thank Jeff Racine for suggesting this insightful exercise. See
also Hall and Racine (2015).}
In Table 1, illustrative results are provided for the regressor distribution
with mass points and the trinormal distribution on a set of support points.
\begin{table}[t]
\caption{Empirical rate of convergence (i.e., $-\protect{\Greekmath 010B}_1$ for $O(n^{-
\protect{\Greekmath 010B}_1})$) in the mass point and high derivative setting.}
\begin{tabular}{p{.5cm}p{1.25cm}p{1.25cm}p{1.25cm}p{1.25cm}}
\multicolumn{5}{c}{$F_X(x)=0.2F^d(x)+0.8\Phi(x)$} \\
\cmidrule(lr){1-5} & X & NW & NW$_a$ & \\
& 0.00 & -0.455 & -0.515 & \\
& 0.10 & -0.186 & -0.443 & \\
& 0.20 & -0.411 & -0.438 & \\
& 0.30 & -0.465 & -0.431 & \\
& 0.40 & -0.461 & -0.425 & \\
& & & &
\end{tabular}
\begin{tabular}{lp{1.25cm}p{1.25cm}}
\multicolumn{3}{c}{$F_X(x)=$ trinormal $(x)$} \\
\cmidrule(lr){1-3} & X & NW \\
& 0.00 & -0.449 \\
& 0.50 & -0.381 \\
& 0.75 & -0.416 \\
& 1.00 & -0.413 \\
\ & & \\
& &
\end{tabular}
\newline
\begin{minipage}{1.0\textwidth}{Note: The column labeled NW$_a$ contains the results implementing the adaptive bandwidth selection procedure in the presence of masspoints.}
\end{minipage}
\end{table}
For the distribution with mass points, the NW estimator with cross-validated
bandwidth performs remarkably well at points sufficiently far from our mass
points (faster than the expected rate of -0.4). The empirical rate at mass
points is close to -0.5 when the bandwidth is set equal to zero. The
empirical convergence rate is slow for points close to the mass points
(within the small ball probability measure under cross validated bandwidth)
due to the boundary weight associated with mass in the neighborhood. {
Bandwidth adaptive to masspoints improves the rate.} The convergence rates
for the trinormal distribution, are reflective of usual smooth nonparametric
regression although are somewhat faster at points with high derivatives.
\subsection{Bivariate (with effective dimension 1)}
We consider a model where $m(X)=\log (X_{1})+\log (X_{2})$ with regressors $
X_{1}$ and $X_{2}$ satisfying $X_{1}+X_{2}=d(k),$ with fixed $d(k)$
corresponding to $k=1,2,3$.\footnote{
An example could be where $X_{1}$ and $X_{2}$ represent earnings of the
husband and wife and, for tax purposes, their combined income is set at some
$d(k)$.} This is equivalent to a model with one continuous and one discrete
regressor $m(X)=\log (X_{1})+\log (D-X_{1}),$ with $D=d(k)$.
We simulated 500 random samples $\{(Y_{i},X_{1i},X_{2i}\}_{i=1}^{n}$ using
the model
\begin{equation*}
Y_i=\log (X_{1i})+\log (X_{2i})+{\Greekmath 011B} {\Greekmath 0122}_i
\end{equation*}
for different sample sizes with the additive error chosen as in the previous
simulation. The probability of an observation belonging to a sub-population
with $k=1,2,3$ is set equal to $0.5$, $0.3$, and $0.2$ respectively and $
d(1)=4,d(2)=6,d(3)=7$; $X_{1}$ is drawn from the uniform distribution: $
U[1,3]$.
We implement the NW estimator first using $X_{1}$ and $X_{2}$ as regressors
(NW.c) and second using $X_{1}$ and $D$ as regressors (NW.d) and obtain the leave-one-out cross-validated bandwidths. For the
discrete regressor $D$ we use special discrete kernel weights proposed by
Wang and van Ryzin (1981) in accordance with Racine and Li (2004).\
In Table 2 we provide illustrative results comparing the empirical rate of
convergence of the NW.c and NW.d at a grid of points.
\begin{table}[t]
\caption{Empirical Rates of the NW.c and NW.d estimators.}
\begin{tabular}{llllp{1.25cm}p{1.25cm}p{1.25cm}}
\ & & & & & & \\
& $X_1$ & $X_2$ & $d(k)$ & \multicolumn{1}{c}{NW.c} & \multicolumn{2}{c}{NW.d
} \\
& & & & $(X_1,X_2)$ & \multicolumn{2}{c}{$(X_1,d(k))$} \\
\cmidrule(lr){6-7} & & & & & ordered & unordered \\
& 1.5 & 2.5 & 4 & -0.451 & -0.422 & -0.419 \\
& 2.0 & 2.0 & 4 & -0.436 & -0.429 & -0.426 \\
& 2.5 & 1.5 & 4 & -0.429 & -0.406 & -0.404 \\
& 1.5 & 4.5 & 6 & -0.457 & -0.410 & -0.422 \\
& 1.5 & 5.5 & 7 & -0.445 & -0.392 & -0.410 \\
\ & & & & & & \\
& & & & & &
\end{tabular}
\newline
\begin{minipage}{1.0\textwidth}{Note: The column labeled ``ordered'' contains the NW.d estimator where the discrete kernel is used for the discrete regressor; the column labeled ``unordered'' uses the Epanechnikov kernel.}
\end{minipage}
\end{table}
The reduced dimensionality is reflected in the estimates of the pointwise
rate of convergence which are around $-0.40$ rather than the slower rate of $
-0.33$ the presence of two continuous regressors would suggest ($q=2$). The
estimate of the empirical rate for NW.c is slightly faster than NW.d,
moreover, indicating that there is no gain from separate treatment of
discrete regressors. With the reduced dimension structure here therefore one
gets the rate corresponding to the Hausdorf dimension of the regressor space
automatically without the need to recognize that it is possible to transform
the regressors to one discrete, and one continuous variable.
\subsection{Bivariate (in the presence of a functional regressor)}
Here we examine a functional regressor in a multivariate setting. Consider a
bivariate conditional mean function $m(X)=m(X_{1},X_{2})$, where $X_{1}$ is
a functional regressor and $X_{2}\in \mathbb{R}$ \ may be correlated with
some $m_{1}(X_{1})$. Let
\begin{equation*}
Y_i=m_{1}(X_{1i})+X_{2i}+{\Greekmath 011B} {\Greekmath 0122}_i.
\end{equation*}
Following Ferraty et al. (2007), the functional regressor is defined as
\begin{equation*}
X_{1i}(t)=\sin (w_it)+(a_i+2{\Greekmath 0119} )t+b_i,\quad t\in (-1,1)
\end{equation*}
with $a_i$ and $b_i$ drawn from $U(-1,1)$, $w_i$ drawn from $U(-{\Greekmath 0119} ,{\Greekmath 0119} )$
and
\begin{equation*}
m_{1}(X_{1i})=\int_{-1}^{1}|X_{1i}^{\prime }(t)|(1-\cos ({\Greekmath 0119} t))dt.
\end{equation*}
For $X_{2}$ we consider two possibilities: (a) a N(0,1) random variable
independent of $X_{1}$; (b) $X_{2}=m_{1}(Z)$ where $Z(t)$ is a functional
regressor similar to $X_{1}(t)$ with $(a_i,b_i,w_i)$ replaced by $
(a_i^{\prime },b_i^{\prime },w_i^{\prime })$ where the correlation between $
(a_i^{\prime },b_i^{\prime },w_i^{\prime })$ and $(a_i,b_i,w_i)$ is given by
${\Greekmath 011A} $ (and set equal to either $0$ or $0.8$).
For the functional regressor $X_{1}$ we use the same metric as in Ferraty et
al. (2007), that is $\left\Vert x_{1}-X_{1}\right\Vert _{1}=\sqrt{
\int_{-1}^{1}\left( x_{1}^{\prime }(t)-X_{1}^{\prime }(t)\right) ^{2}dt}.$
We use a product kernel with kernel $K(u)=1-u^{2}$ defined on $[0,1]$ for
the functional regressor and the Epanechnikov kernel defined on $[-1,1]$ for
$X_{2}$.
Table 3 shows RMSE of the NW estimator at the cross-validated bandwidths as
well as RMSE where either the functional or scalar regressor is dropped. The
loss from misspecifying the functional regression as univariate can be
substantial.
\begin{table}[t]
\caption{RMSE of the NW estimator in the presence of functional regressor $
X_1$ at cross validated bandwidth, $n=250$.}
\begin{tabular}{lccc}
\ & & & \\
& \multicolumn{1}{c}{$X_2 = N(0,1)$} & \multicolumn{2}{c}{$X_2 = m_1(Z)$} \\
\cmidrule(lr){3-4} & & ${\Greekmath 011A} =0.0$ & ${\Greekmath 011A}=0.8$ \\
\textbf{In-sample} & & & \\
\ & & & \\
RMSE & 0.746 & 0.915 & 1.058 \\
\ & & & \\
\quad \textbf{Misspecification:} & & & \\
\quad RMSE$_1$ & 1.140 & 1.918 & 2.099 \\
\quad RMSE$_2$ & 1.833 & 1.854 & 1.746 \\
\ & & & \\
\textbf{Out-of-sample} & & & \\
\ & & & \\
RMSE & 0.915(4) & 1.026(19) & 1.210(18) \\
\ & & & \\
& & &
\end{tabular}
\begin{minipage}{1.0\textwidth}{Note: RMSE$_1$ stands for the RMSE where the $X_2$ regressor is excluded and RMSE$_2$ stands for the RMSE when ignoring the functional regressor. The number in brackets indicates the number of simulations (out of 500) where at the cross-validation bandwidth no neighbor to the out-of-sample observation exists.}
\end{minipage}
\end{table}
\section{Empirical study}
The causal inference literature has made extensive use of the LaLonde (1986)
data on the National Supported Work Demonstration (NSW) program following
the release of that data by Dehejia and Wahba (1999, 2002). Their finding,
that propensity score-based methods provide a way to generalize the
experimental results on the impact of training to nonexperimental data, was
influential and led to significant methodological advances and practical
changes as discussed in the review by Imbens and Xu (2024). Here, we
consider the experimental sample to analyze potential heterogeneous
treatment effects using multivariate kernel estimation. Kernel-based
matching on individual characteristics, advocated in Heckman et al. (1997,
1998), was not considered due to the claimed high dimensionality of the
regressors. We show that it is both feasible and insightful for this data
due to the dimension reduction implied by the presence of several discrete,
discretized and categorical regressors. The mass point of the regressor on
pre-treatment earnings at zero further contributes to the regressor
singularity and we do not require continuity or indeed existence of density
over positive values, thus kinks or mass at positive values are not
excluded. The kernel estimator is applicable to such singular distributions.
We focus here on the full LaLonde NSW male sample which contains 297 treated
individuals and 425 controls where the pre-intervention variables are
well-matched.\footnote{
In the sub-sample with 1974 earnings data in Dehejia and Wahba (1999) the
distribution of 1975 earnings exhibits a significantly different mass at
zero between the treated (68\%) and untreated (60\%).}
With $Y$ denoting the post-treatment outcome, $T$ the treatment and $X$ the
individual pre-treatment characteristic(s), we use the nonparametric
regression model
\begin{equation*}
m(x,j)=E(Y|X=x,T=j)\relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{ for }j=0,1
\end{equation*}
to evaluate the heterogeneous effects of the treatment as
\begin{eqnarray*}
{\Greekmath 011C} (x)&=& m(x,1)-m(x,0).
\end{eqnarray*}
The heterogeneous effect of treatment on the treated, also known as the
conditional average treatment effect, CATT, is given by
\begin{eqnarray*}
{\Greekmath 011C}_T (x)&=& E(m(x,1)-m(x,0) |T=1)
\end{eqnarray*}
Focusing on the latter, we use the NW estimates to evaluate
\begin{equation*}
\hat{{\Greekmath 011C}}_T(x_i)= \widehat{m}(x_i,1) - \widehat{m}(x_i,0) \quad
i=1,\cdots,n_T
\end{equation*}
for all treated individuals $n_T$ (i.e., we use both the actual and the
counterfactual treatment for our estimates).
First, we consider a bivariate kernel regression model where we only use the
pre-treatment earnings (re75) as regressor $X$. Following that, we estimate
the multivariate model with the full set of variables $X$, where in addition
to the pre-treatment earnings we include years of education, high school
``no degree'' status, race, age, marital status, and pre-treatment
unemployment status, u75. It is not unreasonable to attempt nonparametric
estimation for this problem where the only truly continuous regressor is
earnings (and possibly age and education) as singularity provides dimension
reduction.
As with our simulations, we use the np package in R for the nonparametric
estimation where we consider the Epanechnikov (e), Uniform (u) and discrete
(d) kernel.\footnote{
We use the discrete kernel proposed by Aitchison and Aitkin, 1976, where $
K((d-d_{i})/h)=1-h$ if $d=d_{i}$, else $h$ where $h\in \lbrack 0,1/2]$.}
Bandwidth selection is based on cross validation and we consider the
adaptive bandwidth selection approach that accounts for the masspoint. As
was shown in our simulations the rate improvement associated with
singularities does not require special attention to discrete variables to
benefit from it.
Estimation results are reported in detail in the supplementary material
(Appendix D). Below the main findings are summarized.
For the bivariate regression model, the cross validated bandwidths confirm
that we should not smooth across treated and untreated observations and that
local heterogeneous treatment effects as related to pre-treatment earnings
are present. Figure 1, displays estimates of the conditional expectation
using the Epanechnikov kernel by treatment status and pre-treatment earnings
together with the bootstrapped confidence bounds. It suggests that treatment
for individuals at low levels of pre-treatment earnings, in particular, is
beneficial.
\begin{figure}[H]
\caption{Nonparametric fit of the conditional expectation by pre-treatment
earnings and treatment status (cross validated bandwidth, Epanechnikov
kernel)}\includegraphics[scale=.7]{Emp_fig1.jpeg}\newline
\begin{minipage}{1\textwidth}{Note: All graphs related to the empirical application are rescaled with all numbers denoted in '000\$s. }
\end{minipage}
\end{figure}
The adaptive bandwidth results in a slightly better in-sample correlation
between the post-treatment outcome, $re78$, and its fit (increasing from
0.2098 to 0.2116 (for OLS the correlation is 0.1697)); bandwidths obtained
using non-masspoint-observations only are quite similar to those obtained
when including the masspoints in this case. The NW estimates with the
adaptive bandwidth provide values of CATT that on average equal $\$920$
(76), \$920 (76), and \$906 (80) (standard error in brackets) for the (e,e), (d,e), (d,u) kernels on $(T,X)$,
respectively.\footnote{
As discussed in the supplemental material (Appendix D), we denote the kernel with two
arguments: the first argument denotes the kernel applied to all binary
regressors (treat, u75, nodegree, black, hispanic, and married) and the
second argument denotes the kernel applied to the other regressors (re75,
educ, and age).} For comparison, the local linear kernel based estimates on
average equal \$822 (49) with the (e,e) kernel, while the average of the
CATT estimates based on random forest (RF) equal \$848 (52). The CATT
results of the kernel regression based approach are more variable than those
obtained using the random forest approach. For observations at mass points,
CATT estimates using adaptive bandwidth are closer to those obtained using
the random forest based approach.
For the multivariate model the cross-validated bandwidths provide important
insights. Firstly, even though pre-treatment earnings is still relevant, the
bandwidth is much larger than in the baseline model for all kernels,
suggesting a reduction of the heterogeneous impact with individual's
pre-treatment earnings. The bandwidths selected for nodegree, hispanic and
married are large, signaling that these variables are not relevant (these
regressors are automatically smoothed out from the regression function). At
the same time, the bandwidths for education and age imply a heterogeneous
impact associated with those characteristics, although the size of the
bandwidth for age is fairly large.
The inclusion of additional controls yields an improvement in the in-sample
correlation between the post-treatment outcome and its fit. For the (e,e)
kernel we see an increase in correlation from 0.210 in the bivariate model
to 0.338 (for comparison, for OLS the correlation equals 0.209 when age
squared is included as well); the results for the (d,e) and (d,u) kernel are
comparable.
To highlight the heterogeneity of the treatment effect of education and its
interplay with race, we display in Figure 2 estimates of the conditional
expectation by treatment status, years of education, and race for an
individual with median age and pre-treatment earnings together with the
bootstrapped confidence bounds.
\begin{figure}[H]
\caption{Nonparametric fit of the conditional expectation by years of
education, race, and treatment status with median pre-treatment earnings and
age) (cross validated bandwidth, Epanechnikov kernel)}
\medskip \relax\protect\ifmmode\expandafter\text@\else\expandafter\mbox\fi{\hspace{.61in} Black =0 (N=144) \hspace{.9in} Black=1 (N=579)}
\newline
\includegraphics[scale=.9]{Emp_fig2.jpeg}\newline
\begin{minipage}{1\textwidth}{Note: The median pre-treatment earnings equals $\$936$ and the median age is $23$. The estimates are rescaled and are denoted in '000\$. }
\end{minipage}
\end{figure}
The graph reflects a heterogeneity of the impact of treatment whereby the
more educated individuals identified as black appear to benefit more from
treatment than their nonblack counterparts. Gains of treatment arise where
the confidence band around the estimated nonparametric fit $\widehat{m}(x,1)$
lies above that of $\widehat{m}(x,0)$; for non-black individuals this is at
the middle range of education, for black individuals this starts around 10
years of education and is rising over that range. These results are further
supported when evaluating the average CATT for black individuals across
different levels of education (see supplemental material, Appendix D).
Box-plots of the CATT estimates for the multivariate model using the NW
regression estimate and the RF estimates are presented in Figure 3. The
limit distributional results of Wager and Athey (2018) do not apply here as
many components of $X$ are not continuously distributed.
\begin{figure}[]
\caption{Box-plots of the CATT estimates (NW and RF)}
\medskip \includegraphics[scale=.7]{Emp_fig3.jpeg}\newline
\begin{minipage}{1\textwidth}{Note: The estimates are rescaled and are denoted in '000\$.
}
\end{minipage}
\end{figure}
The kernel based regression CATT results remain more variable than those
provided by the random forest approach, but their interquartile range is
comparable. The NW kernel based estimates of the CATT on average exceed the
RF based estimates: $\$1,045$ (107), \$1,018 (108), and \$1,019 (104) for
the (e,e), (d,e) and (d,u) kernel on $(T,X)$ against \$794 (54) based on the
random forest. \medskip
The NW based results are stable across kernel, give interpretable insights
and with cross-validation make it possible to detect irrelevant regressors.
\newpage
\noindent \textbf{Acknowledgements:} The authors thank the participants at
the Econometric Study Group conference in Bristol, the Canadian Econometric
Study Group conference, and Saraswata Chaudhuri for their comments. We thank
Jeffrey Racine for his discussion and valuable suggestions at the CESG 2023
and Sid Kankanala for insightful comments on earlier versions of the paper.
We thank the Associate Editor and three anonymous referees for their careful
reading of the paper and very helpful comments and suggestions. \vspace{.2in}
\noindent \textbf{Funding:} Victoria Zinde-Walsh gratefully acknowledges
financial support from the Natural Sciences and Engineering Research Council
of Canada (NSERC) grant 253139.
\begin{thebibliography}{99}
\bibitem{} Ackerberg, D.A., K. Caves, and G. Frazer (2015) ``Identification
properties of recent production function estimators,'' \textit{Econometrica}
, \textbf{83}, 2411--2451.
\bibitem{} Ahlfors, Lars (1966) \textit{Lectures on quasiconformal mappings}
,\ Princeton University Press.
\bibitem{} Aitchison, J. and C.G.G. Aitkin (1976) ``Multivariate binary
discrimination by the kernel method,'' \textit{Biometrika}, \textbf{63},
413--420.
\bibitem{} Angrist, J.D. and J.-S. Pischke (2009) \textit{Mostly harmless
econometrics: An empiricist's companion}, Princeton University Press.
\bibitem{} Arulampalam, W., V. Corradi, and D. Gutknecht (2017) ``Modeling
heaped duration data: An application to neonatal mortality,'' \textit{
Journal of Econometrics}, \textbf{200}, 363--377.
\bibitem{} Bai, J., S. Ng (2006) ``Confidence intervals for diffusion index forecarsts and inference with factor-augmented regressors,'' \textit{Econometrica}, \textbf{74}, 1133--1150.
\bibitem{} Caldeira J.F., R. Gupta, H.S. Torrent (2020) ``Forecasting U.S.
aggregate stock market excess return: Do functional data analysis add
economic value?,'' \textit{Mathematics}, \textbf{8}, 2042.
https://doi.org/10.3390/math8112042 .
\bibitem{} Dehejia R.H. and S.Wahba (1999) ``Causal effect in
nonexperimental studies: Reevaluating the evaluation of training programs,''
\textit{Journal of the American Statistical Association}, \textbf{94},
1053--1062.
\bibitem{} Dehejia R.H. and S.Wahba (2002) ``Propensity score-matching
methods for nonexperimental causal studies,'' \textit{The Review of
Economics and Statistics}, \textbf{84}, 151--161.
\bibitem{} Demir, S. and \"{O}. Toktamis (2010) ``On the adaptive
Nadaraya-Watson kernel regression estimators,'' \textit{Hacettepe Journal of
Mathematics and Statistics}, \textbf{39}, 429--437.
\bibitem{} Desmet, K. and S.L. Parente (2010) ``Bigger is better: Market
size, demand elasticity, and innovation,'' \textit{International Economic
Review}, \textbf{51}, 319--333.
\bibitem{} Donkers, A.C. and M.M.A. Schafgans (2008) ``Estimation and
specification of semiparametric index models,'' \textit{Econometric Theory},
\textbf{24}, 1584--1606.
\bibitem{} Fan, J. and I. Gijbels (1996) \textit{Local polynomial modelling
and its applications}, Chapman and Hall.
\bibitem{} Ferraty F., A. Mas, and P. Vieu (2007) ``Nonparametric regression
on functional data: Inference and practical aspects,'' \textit{Australian
and New Zealand Journal of Statistics}, \textbf{49}, 267--286.
\bibitem{} Ferraty, F. and S. Nagy (2022) ``Scalar-on-function local linear
regression and beyond,''\ \textit{Biometrika}, \textbf{109}, 439--455.
\bibitem{} Ferraty, F. and P. Vieu (2004) ``Nonparametric models for
functional data, with application in regression, time series prediction and
curve discrimination,'' \textit{Journal of Nonparametric Statistics},
\textbf{16}, 111--125.
\bibitem{} Ferraty F. and P. Vieu (2006) \textit{Nonparametric functional
data analysis: Theory and Practice}, Springer, New York.
\bibitem{} Gasser, T., P. Hall, and B. Presnell (1998) ``Nonparametric
estimation of the mode of a distribution of random curves,'' \textit{Journal
of Royal Statistical Society, Series B}, \textbf{60}, 681--691.
\bibitem{} Geenens, G. (2015) ``Moments, errors, asymptotic normality and
large deviation principle in nonparametric functional regression,'' \textit{
Statistics and Probability Letters}, \textbf{107}, 369--377.
\bibitem{} Gy$\ddot{o}$rfi, L., M. Kohler, A. Krzyzak, and H. Walk (2002)
\textit{A distribution-free theory of nonparametric regression}, Springer,
New York.
\bibitem{} Hall, P. and J. Horowitz (2013) ``A simple bootstrap method for
constructing nonparametric confidence bands for functions,'' \textit{Annals
of Statistics}, \textbf{41}, 1892--1921.
\bibitem{} Hall, P., Q. Li, and J.S. Racine (2007) ``Nonparametric
estimation of regression functions in the presence of irrelevant
regressors,'' \textit{The Review of Economics and Statistics}, \textbf{89},
784--789.
\bibitem{} Hall, P. and J.S. Racine (2015) ``Infinite order cross-validated
local polynomial regression,'' \textit{Journal of Econometrics}, \textbf{185}
, 510--525.
\bibitem{} Heckman, J.J., H. Ichimura, and P.E. Todd (1997) ``Matching as an
econometric evaluation estimator: Evidence from evaluating a job training
programme,'' \textit{Review of Economic Studies}, \textbf{64}, 605--654.
\bibitem{} Heckman, J.J., H. Ichimura, and P.E. Todd (1998) ``Matching as an
econometric evaluation estimator,'' \textit{Review of Economic Studies},
\textbf{65}, 261--294.
\bibitem{} Hong,S. and O. Linton (2020) ``Nonparametric estimation of
infinite order regression and its application to the risk-return tradeoff,''
\textit{Journal of Econometrics}, \textbf{219}, 389--424.
\bibitem{} Hotelling, H. (1929) ``Stability in competition,'' \textit{The
Economic Journal}, \textbf{39}, 41--57.
\bibitem{} Ichimura, H. (1993) ``Semiparametric least squares (SLS) and
weighted SLS estimation of single index models,'' \textit{Journal of
Econometrics}, \textbf{58}, 71--120.
\bibitem{} Imbens, G. and Y. Xu (2024) ``LaLonde (1986) after nearly four
decades: Lessons learned,'' arXiv 2406.00827 (econ.EM).
\bibitem{} Jun, B and H. Song (2019) ``Tests for detecting probability mass
points,'' \textit{Korean Economic Review}, \textbf{35}, 205--248.
\bibitem{} Kankanala, S. and V. Zinde-Walsh (2024) ``Kernel-weighted
specification testing under general distributions,'' \textit{Bernoulli},
\textbf{30}, 1921--1944.
\bibitem{} K\"{o}hler, M., A. Schindler, and S. Sperlich (2014) ``A review
and comparison of bandwidth selection methods for kernel regression,''
\textit{International Statistical Review / Revue Internationale de
Statistique}, \textbf{82}, 243--274.
\bibitem{} Kotlyarova, Y, M. Schafgans, and V. Zinde-Walsh (2016)
``Smoothness: Bias and efficiency of non-parametric kernel estimators,'' in
\textit{Advances in Econometrics: Essays in Honor of Aman Ullah}, Vol. 36,
G. Gonzales-Rivera, R.C. Hill and T.-H. Lee, eds. 561--589.
\bibitem{} Kurisu, D., Otsu, T., and M. Xu (2025) ``Nonparametric Causal
Inference with Functional Covariates,'' \textit{Journal of Business and
Economic Statistics}, 1--14. https://doi.org/10.1080/07350015.2025.2501563
\bibitem{} LaLonde, R. (1986) ``Evaluation the Econometric Evaluations of
Training Programs with Experimental Data,'' \textit{American Economic Review}
, \textbf{76}, 604--620.
\bibitem{} Li, C., D. Ouyang, and J.S. Racine (2009) ``Nonparametric
regression with weakly dependent data: the discrete and continuous regressor
case,'' \textit{Journal of Nonparametric Statistics}, \textbf{21}, 697--711.
\bibitem{} Li, Q. and J.S. Racine (2007) \textit{Nonparametric econometrics:
Theory and practice}, Princeton University Press.
\bibitem{} Mandelbrot, B. (1997) \textquotedblleft Fractals and scaling in
finance, discontinuity, concentration, risk,\textquotedblright\ \textit{
Selecta}, Volume E, Springer.
\bibitem{} Masry, E. (2005) ``Nonparametric regression estimation for
dependent functional data: asymptotic normality,'' \textit{Stochastic
Processes and their Applications}, \textbf{115}, 155--177.
\bibitem{} Nadaraya, E. (1965) ``On non-parametric estimates of density
functions and regression curves,'' \textit{Theory of Probability and its
Applications}, \textbf{10}, 186--190.
\bibitem{} Olson, C.A. (1998) ``A comparison of parametric and
semiparametric estimates of the effect of spousal health insurance coverage
on weekly hours worked by wives,'' \textit{Journal of Applied Econometrics},
\textbf{13}, 543--565.
\bibitem{} Pollard, D. (2001) ``A User's Guide to Measure Theoretic
Probability,'' Cambridge series in Statistical and Probabilistic Mathematics.
\bibitem{} Racine, J. and Q. Li (2004) ``Nonparametric estimation of
regression function with both categorical and continuous data,'' \textit{
Journal of Econometrics}, \textbf{119}, 99--130.
\bibitem{} Ramsay, J. O. and B.W. Silverman (2005) ``Functional Data
Analysis,'' Springer, New York.
\bibitem{} Rosenbaum, P.R. and D.B. Rubin (1983) ``The central role of the
propensity score in observational studies for causal effects,'' \textit{
Biometrika}, \textbf{70}, 41--55.
\bibitem{} Sain, S.R. (1994) ``Adaptive kernel density estimation,'' \textit{
Computational Statistics and Data Analysis}, \textbf{39}, 165--186.
\bibitem{} Shen, G. (2002) ``Fractal dimension and fractal growth of
urbanized areas,'' \textit{International Journal of Geographical Information
Science}, \textbf{16}, 419---437.
\bibitem{} Takayasu, M. and H. Takaysu (2009) ``Fractals and Economics'', in
Encyclopedia of Complex Systems in Finance and Econometrics, Meyers, R.A.,
Eds., 444-463, Springer.
\bibitem{} Tibshirani, J. and S.Athey (2024) ``Package `grf' '',
\url{https://cran.r-project.org/web/packages/grf/grf.pdf}.
\bibitem{} Vol'berg, A.L. and S.V. Konyagin (1988) ``On measures with the
doubling condition,'' \textit{Mathematics of the USSR-Izvestiya}, \textbf{30}
, 629--638.
\bibitem{} Wager, S. and S. Athey (2018) ``Estimation and inference of
heterogeneous treatment effects using random forests,'' \textit{Journal of
the American Statistical Association}, \textbf{113}, 1228--1242.
\bibitem{} Wang, M.-C. and J. van Ryzin (1981) ``A class of smooth
estimators for discrete distributions,'' \textit{Biometrika}, \textbf{68},
301--309.
\bibitem{} Watson, G. S. (1964) ``Smooth regression analysis,'' \textit{
Sankhya: The Indian Journal of Statistics}, Series (1961-2002), \textbf{26},
359--372.\vspace{.4in}
\end{thebibliography}
\pagebreak