EconBase
← Back to paper

Semiparametric Learning of Integral Functionals on Submanifolds

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

124,194 characters · 22 sections · 47 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Thin Sets Are Not Equally Thin: Minimax Learning of Submanifold Integrals

\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\def\mc#1{\mathscr{#1}} \global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\global\long\def\abs#1{\left|#1\right|} \global\long\def\norm#1{\left\Vert #1\right\Vert } \global\long\def\rest#1{\left.#1\right|} \global\long\def\bracket#1#2{\left\langle #1\middle\vert#2\right\rangle } \global\long\def\sandvich#1#2#3{\left\langle #1\middle\vert#2\middle\vert#3\right\rangle } \global\long\def\third#1{\frac{#1}{3}} \global\long\global\long\def\sand#1{\left\lceil #1\right\vert } \global\long\def\wich#1{\left\vert #1\right\rfloor } \global\long\def\sandwich#1#2#3{\left\lceil #1\middle\vert#2\middle\vert#3\right\rfloor } \global\long\def\inprod#1{\left\langle #1\right\rangle } \global\long\def\ol#1{\overline{#1}} \global\long\def\ul#1{#1} \global\long\def\td#1{\tilde{#1}} \global\long\def\bs#1{\boldsymbol{#1}} \global\long\global\long\global\long\global\long\global\long{6pt} {6pt}

\sloppy

abstractMany economic parameters are identified by “thin sets” (submanifolds with Lebesgue measure zero) and hence difficult to recover from data in an ambient space. This paper provides a unified theory for estimation and inference of such “thin-set” identified functionals. We show that thin sets are not equally thin: their intrinsic dimensionality $m$ matters in a precise manner. For a nonparametric regression $h_0$ with H\"{o}lder smoothness $s$ and $d$-dimensional covariates in the ambient space, we show that $n^{-\frac{s}{2s+d-m}}$ is the minimax optimal rate of estimating linear and nonlinear (e.g., quadratic, upper contour) integrals of $h_0$ on an $m$-dimensional submanifold ($0\leq m < d$), which is the fastest possible attainable rate among all estimators. The minimax lower bound rate result is generalized to estimating submanifold integrals when $h_0$ is a nonparametric density and a nonparametric instrumental variable function. The asymptotic normality of t statistics is established via sieve Riesz representation, and the corresponding inference is computed using Sobol points. \\ Keywords: submanifold, contour set, level set, Hausdorff measure, differential geometry, minimax rate, sieve Riesz representer, asymptotic normality.

Introduction

Many parameters of interest in economics are identified by information contained in lower-dimensional thin sets, i.e., subsets of the covariate space that have Lebesgue measure zero yet may still carry economic meanings. \citet*{khan2010irregular} coined the term thin-set identification to describe settings in which identifying information is concentrated on such measure-zero sets, and showed that the resulting parameters are irregular in the sense that they cannot be estimated at the parametric $n^{-1/2}$-rate, where $n$ denotes the sample size.

In this paper, we provide a more nuanced and, in some sense, refined view about thin sets. While all thin sets possess Lebesgue measure zero in their ambient space, we show that they may differ substantively in terms of their intrinsic dimensions and geometric structures. These refined differences lead to quantitatively different convergence rate of estimation and different forms of first-order expansions for inference on aggregate parameters (integrals) defined on such thin sets.

Specifically, this paper considers semiparametric estimation and inference for general integral functionals on submanifolds of the following form:

equation[equation omitted — 117 chars of source]

where $\phi:\mathbb{R} \times \mathcal{X} \mapsto \mathbb{R}$ is a known transformation of an unknown function $h_0$ that can be estimated nonparametrically from data, $w:\mathcal{X} \mapsto \mathbb{R}$ is a known (weight) function, and $\mathcal{X}$ is a known closed convex subset (with positive Lebesgue measure) in $\mathbb{R}^d$. Here ${\cal H}^{m}$ denotes the $m$-dimensional Hausdorff measure in $\mathbb{R}^{d}$, which coincides with the Lebesgue measure in $\mathbb{R}^m$ (and hence has Lebesgue measure zero in $\mathbb{R}^d$ whenever $m<d$).\footnote{See Section (ref) for a precise definition of ${\cal H}^{m}$. For $m<d$, the $m$-dimensional Hausdorff measure ${\cal H}^{m}$ can be thought as a generalization of a variety of “uniform” measures for lower-dimensional subsets in $\mathbb{R}^d$ such as point, line, area, surface, and volume measures.} In this paper,

equation[equation omitted — 97 chars of source]

denotes an $m$-dimensional ($0\leq m <d$) submanifold, where $g:\mathcal{X} \mapsto \mathbb{R}^{d-m}$ may be known, or unknown but can be estimated parametrically, semiparametrically or nonparametrically from data. By definition, the $m$-dimensional submanifold ${\cal M}$ has positive $m$-dimensional Hausdorff volume but zero $d$-dimensional Lebesgue volume, and hence is a thin set in the covariate ambient space $\mathcal{X}\subset \mathbb{R}^{d}$.

Submanifold integrals emerge naturally in many economic settings where we take derivatives of a standard Lebesgue integral ($dx$ in $\mathbb{R}^d$), such as an expectation, with respect to a changing domain of integration:

equation[equation omitted — 219 chars of source]

Mathematically, the derivative of an integral with respect to its domain $\Omega_{t}$ translates into an integral over the boundary $\partial\Omega_t$ of its domain, and the boundary is often a lower-dimensional submanifold of the original domain. As we show in Section (ref), many parameters of interest are identified by solutions to first-order condition for optimization over subpopulation, which often take the form (ref).

Differentiation of an integral with respect to its domain also shows up in asymptotic analysis. For example, consider estimation of the following integral functional on upper contour set of the unknown $h_0$

equation[equation omitted — 187 chars of source]

where $w(x)$ is the marginal density of $x$ with support $\mathcal{X}\subset \mathbb{R}^d$. When $h_{0}\left(x\right)=\text{CATE}\left(x\right):=\mathbb{E}\left[\rest{Y_{i}\left(1\right)-Y_{i}\left(0\right)}X_{i}=x\right]$ is the conditional average treatment effect (CATE) for type $x \in \mathcal{X}$, we call $V\left(h_{0}\right)$ a value functional, which is a policy parameter of interest in many applications. If $h_{0}$ is nonparametrically estimated by $\hat{h}$, the impact of the estimation error on the plug-in estimation of $V(h_{0})$ can be analyzed using the pathwise derivative of $V(h)-V(h_{0})$ with respect to $h_0$ in the direction $[h-h_{0}]$, which becomes a submanifold integral of the form: \[ DV(h_0)[h-h_0]:=\rest{\frac{dV\left(h_0 +t[h-h_0]\right)}{dt}}_{t=0}=\int_{\mathcal{M}_0}\frac{h\left(x\right)-h_{0}\left(x\right)}{\left\| \nabla_{x}h_{0}\left(x\right) \right\|}w\left(x\right)d{\cal H}^{d-1}\left(x\right), \] where the level set\footnote{In this paper we use ${\cal M}_0$ to highlight the level set when it is an unknown function $h_0$.} ${\cal M}_0:=\{x\in \mathcal{X}:~ h_{0}\left(x\right)=0\}$ is a $m=(d-1)$-dimensional submanifold in $\mathbb{R}^{d}$.

In this paper, we assume that a data set $\mathcal{D}_n:=\{(Y_i,X_i)\}_{i=1}^n$ is a random sample of size $n$ drawn from an unknown probability distribution of $(Y,X)$, in which the unknown density of $X$ is bounded away from zero and infinity on its support $\mathcal{X}$ in $\mathbb{R}^d$. We also assume that the unknown function $h_0:\mathcal{X} \mapsto \mathbb{R}$ and/or the unknown submanifold mapping $g:\mathcal{X} \mapsto \mathbb{R}^{d-m}$ can be identified and estimated from the data $\mathcal{D}_n$ in the ambient space.\footnote{Our framework is different from the literature that assumes the data is sampled directly from lower dimensional submanifolds.} Our parameters of interest are the $m$-dimensional submanifold integrals $\Gamma(h_0)$ and the upper contour integrals $V(h_0)$. In Section (ref) we show that many economic functionals of interest can be represented as integrals of the forms $\Gamma(h_0)$ and $V(h_0)$.

We first establish the minimax lower bound rates of estimation for the linear integral functional on a submanifold:

equation[equation omitted — 120 chars of source]

a core special case of $\Gamma(h_0)$ in (ref) with $\phi=h_0$ and ${\cal M}$ a known $m$-dimensional submanifold ($m < d$). For the sake of concreteness we assume that the unknown function $h_0:\mathcal{X} \mapsto \mathbb{R}$ belongs to a H\"{o}lder class of finite smoothness $s>0$. When $h_{0}$ is respectively a nonparametric regression $\mathbb{E}[Y|X]$, a nonparametric density of $X$ and a nonparametric instrumental variables (NPIV) regression $\mathbb{E}[Y-h_0 (X)|Z]=0$, we establish the minimax lower bound rates for estimating $L(h_0)$, which are the fastest possible convergence rates among all estimators for $L(h_0)$. These lower bound rates are all slower than the parametric convergence rate of $n^{-1/2}$ whenever $m<d$, confirming that $L(h_0)$ is an irregular functional whenever $m<d$. Specifically, we show that $r^*_n:=n^{-\frac{s}{2s+d-m}}$ is the minimax lower bound rate for estimating $L(h_0)$ when $h_0$ is a nonparametric regression and a density. Interestingly, this rate coincides with the famous lower bound rate of stone1980optimal for the pointwise estimation of a nonparametric regression with $(d-m)$-dimensional covariates. In Subsection (ref) we also establish that $r_{NPIV,n}$ is the minimax lower bound rate for estimating $L(h_0)$ when $h_0$ satisfies a NPIV restriction but does not need to be point-identified. Importantly, the rate $r_{NPIV,n}$ for $L(h_0)$ coincides with the previously established minimax lower bound rate for estimating a $(d-m)$-dimensional point-identified NPIV function (chen2011rate).

Under some mild regularity conditions, we show that the minimax lower bound rates for $L(h_0)$ are also the minimax lower bounds rates for nonlinear integrals $\Gamma(h_0)$ with possibly unknown submanifolds, upper contour integrals $V(h_0)$, and surface integrals of the form

equation[equation omitted — 147 chars of source]

For example, when $h_{0}$ is a nonparametric regression or a density, the rate $r^*_n=n^{-\frac{s}{2s+d-m}}$ is the minimax lower bound rate for estimating nonlinear integral $\Gamma(h_0)$ on $m$-dimensional possibly unknown submanifolds. In particular, when $m=d-1$ the minimax rate $r^*_n$ becomes $n^{-\frac{s}{2s+1}}$, which is the fastest possible convergence rate for estimating $V(h_0)$, $S(h_0)$ and $\Gamma(h_0)$ on $m=(d-1)$-dimensional $\mathcal{M}$ among all estimators. This result recovers the minimax lower bound of horowitz1993et for smooth maximum score estimation that corresponds to an $m=(d-1)$-dimensional parametric submanifold $\{x\in\mathcal{X}:~x'\beta_0=0\}$.\footnote{Our minimax lower bound proof follows the nonparametric literature such as Tsybakov2009, which differs from the proof of horowitz1993et for his smooth maximum score model.}

We then show that the above lower bound rates are attainable, and hence minimax-optimal, by presenting sieve-based estimators for $\theta_0=L(h_0),~\Gamma(h_0)$ and $V(h_0)$ for the case when $h_0$ is a nonparametric regression with H\"{o}lder smoothness $s$ and $d$-dimensional covariates. For the linear integral $L(h_0)$, we consider the plug-in sieve estimator. For the nonlinear integrals $\Gamma(h_0)$ and $V(h_0)$, we consider plug-in, split-sample and leave-one-out sieve estimators. We provide low level sufficient conditions under which they achieve the optimal rate $r^*_n=n^{-\frac{s}{2s+d-m}}$, with the smoothness requirement for the split-sample and leave-one-out debiased estimators weaker than that for the plug-in estimators for $\Gamma(h_0)$ and $V(h_0)$.

Given the irregularity of submanifold integral functionals, they do not admit well-defined Riesz representers; however, the sieve Riesz representers remain well-defined and computable in closed form. Following \citet*{chen2014sieve,chen2014sieveM} and chenpouzo2015sieve, we construct valid confidence intervals via sieve student-$t$ statistics. By exploiting the submanifold structure, we characterize the growth rate of the sieve Riesz representer norm and obtain tighter control of nonlinear remainders. For the upper contour integral $V(h_0)$, its pathwise derivatives are computed using the calculus of moving submanifolds.

Monte Carlo simulations confirm our theoretical results: the sieve estimators produce reasonably small RMSEs that shrink with sample sizes, and the realized confidence interval coverage is close to the nominal 95% level. The submanifold integrals in the simulations are numerically computed using Sobol quasi-random sequences sobol1967distribution for its better numerical performance than uniform random sampling.

In a companion paper CCG2025, we apply the theory developed here to inference on value functionals of the CATE under first-best nonparametric treatment assignment. Using the Job Training Partnership Act (JTPA) data set, we compute confidence intervals for the nonparametric first-best welfare and the treatment share for the JTPA job training program, two parameters estimated in kitagawa2018should without reported confidence intervals.

Related Literature

Our paper contributes to the literature on semiparametric estimation and inference on irregular integral functionals. To our best knowledge, our paper is the first to provide a unified theory on general submanifold integral functionals of an unknown function $h_0$ with H\"{o}lder smoothness $s>0$. We establish the minimax-optimal estimation rate for linear and nonlinear submanifold integral functionals of $h_0$ and for integrals of upper contour set on $h_0$, in which $h_0$ could be a regression, a density, and a NPIV function. In addition, we provide simple asymptotic normality based confidence intervals for these irregular integral functionals using sieve Riesz representation.

Our minimax lower bound estimation rate results can be viewed as quantitative refinements of the famous singular semiparametric information bound results of chamberlain1986asymptotic and \citet*{khan2010irregular} on “thin-set” identified parameters. In Section (ref) we show many semiparametric “thin-set” identified parameters can be viewed as integral functionals on $m$-dimensional submanifolds (for $m<d$). Previously, \citet*{kim1990cube} shows that the first-order condition for the maximum score criterion of \citet*{manski1975maximum} is a submanifold integral with a $m=(d-1)$-dimensional hyperplane, and derives the famous $n^{-1/3}$ convergence rate using empirical process theory.\footnote{horowitz1992smoothed,horowitz1993et obtains the upper and lower bound rates of $n^{-s/(2s+1)}$ for the smoothed maximum score estimator without mentioning $(d-1)$-dimensional submanifold.} \citet*{sasaki2015quantile} also notes that differentiation with respect to the domain of integration produces a $m=(d-1)$-dimensional submanifold integrals in his identification paper on quantile functions in nonseparable structural models.

We view our minimax lower bound results for the large class of submanifold functionals $L(h_0)$, $\Gamma (h_0)$, $V(h_0)$ and $S(h_0)$ important, as they provide a minimax informational criterion to compare many different machine learning estimators. Our paper presents various sieve estimators for $L(h_0)$, $\Gamma(h_0)$ and $V(h_0)$ and establishes that they can attain the minimax lower bound rate of $r^*_n=n^{-\frac{s}{2s+d-m}}$ when $h_0$ is a nonparametric regression. Other estimators for these functionals can also be presented and be verified if they are minimax rate optimal. For example, \citet*{qiao2021nonparametric} proposes a kernel plug-in density estimator of the surface integral $S(h_0)$ of the form (ref) and establishes a convergence rate $n^{-s/(2s+1)}$ when $s\geq d+1$. \citet*{qiao2021nonparametric} concludes his paper by stating that “Another open problem is the minimax rates of estimating the surface integrals on level sets”. Notice that the surface integral $S(h_0)$ is a $m=(d-1)$-dimensional nonlinear submanifold functional, his kernel estimator matches our minimax lower bound rate $r^*_n$ provided $s\geq d+1$. In works that are concurrent to ours, \citet*{cattaneo2025dist,cattaneo2025loc} consider estimation of linear integrals over $m=1$-dimensional known submanifolds that arise in the boundary discontinuity designs using local polynomial regressions and obtain a minimax optimal convergence rate of $n^{-\frac{s}{2s+d-1}}$ in their contexts (they mostly consider $d=2,m=1$). Due to the lack of space, we leave it to future work to compare finite sample performance of different minimax rate-optimal estimators for these submanifold integrals.

Technically, \citet*{qiao2021nonparametric}, \citet*{cattaneo2025dist,cattaneo2025loc} and our paper all utilize some mathematical tools in differential geometry and geometric measure theory to establish the asymptotic properties of different estimators of different irregular submanifold integrals. Differential geometry tools have also been used in asymptotic analysis of regular functionals (i.e., the ones that can be estimated at the parametric convergence rate of $n^{-1/2}$). \citet*{chernozhukov2018sorted} uses integrals on level sets to study sorted partial effects in heterogeneous coefficient models and establishes the convergence rate of $n^{-1/2}$ for their regular functionals. \citet*{feng2024statistical} shows how Hausdorff integrals can be used in the analysis of regular integral functionals.

\paragraph{Organization of the Paper} Section (ref) provides motivating examples for submanifold integrals in econometrics. Section (ref) establishes minimax lower bounds rates for estimating linear and nonlinear submanifold integrals $L(h_0),~\Gamma(h_0)$ when $h_0$ is a nonparametric regression, a nonparametric density and a NPIV function respectively. Section (ref) shows that the minimax lower bound rate is attainable by sieve estimators for estimating $L(h_0),~\Gamma(h_0)$ and $V(h_0)$ when $h_0$ is a nonparametric regression. Section (ref) provides the asymptotic normality of the sieve estimators proposed in Section (ref), along with consistency sieve variance estimators, for inference on both linear and nonlinear submanifold integrals. Section (ref) presents Monte Carlo simulation results. The Appendix contains sections on mathematical tools in differential geometry used in this paper, technical lemmas, as well as all the proofs.

Setup and Examples

Setup

Throughout this paper, we let $\left\{Y_{i},X_{i}\right\}_{i=1}^{n}$ be a random sample of size $n$ drawn from a unknown joint probability distribution $P_{\left(Y,X\right)}$, where $Y_{i}$ is a scalar-valued outcome variable, and $X_{i}$ is a vector of observed covariates with a convex and compact support ${\cal X}\subseteq\mathbb{R}^{d}$. Let $h_{0}:{\cal X} \to \mathbb{R}$ be a nonparametric function\footnote{More generally, $h_0$ may be a vector of nonparametric functions.} that is directly identified from the data and can be estimated using standard nonparametric estimation methods. A leading example of $h_0$ is the conditional expectation (nonparametric regression) function $h_0(x) = \mathbb{E}[\rest{Y_i}X_i=x]$, which will be our main focus in the paper. That said, $h_{0}$ may also take the form of density functions, conditional quantiles and structural regression functions in NPIV models.

We consider submanifolds ${\cal M}:=\left\{ x\in{\cal X}:g\left(x\right)={\bf 0}\right\}$ that take the form of level sets of functions $g$. Depending on the problem setup, $g$ may be known or unknown, parametric or nonparametric, and it may be taken to be different from or the same as $h_0$, which will be illustrated in the examples below and treated in subsequent sections. We maintain the following standard regularity condition in this paper:

assumption[Regular Level Set] (i) $g:{\cal X}\to\mathbb{R}^{d-m}$ is a continuously differentiable function with ${\bf 0}\in\text{int}\left(g\left({\cal X}\right)\right)\subseteq\mathbb{R}^{d-m}$; (ii) $\nabla_{x}g\left(x\right)$ has full rank $d-m$ for every $x\in{\cal M}:=\left\{ x\in{\cal X}:g\left(x\right)={\bf 0}\right\}$.\footnote{Note that ${\bf 0}$ may be replaced with any other constant vector ${\bf c} \in\text{int}\left(g\left({\cal X}\right)\right)\subseteq\mathbb{R}^{d-m}$ without affecting the results in this paper.}

Let $\mathcal{J}_g\left(x\right)$ denote the Jacobian of $g:\mathbb{R}^{d}\to\mathbb{R}^{d-m}$ defined by \[ \mathcal{J}_g\left(x\right):=\sqrt{\sum_{B(x)}\text{det}\left(B(x)\right)^{2}}=\sqrt{\det(\nabla_x g(x)\nabla_x g(x)')}, \] where $B$ indexes all $\left(d-m\right)\times\left(d-m\right)$ minors of $\nabla_x g\left(x\right)$. Under Assumption (ref) and compactness of ${\cal X}$, we have $\mathcal{J}_g\left(x\right)\geq Const.>0$ for all $x\in\mathcal{M}$.

Under Assumption (ref), the set ${\cal M}$ defined in (ref) is an $m$-dimensional submanifold of $\mathbb{R}^{d}$ (e.g. by Theorem 12.1 of loomis2014advanced). We note that ${\cal M}$ has zero Lebesgue measure on $\mathbb{R}^{d}$ for all $0\leq m <d$, although it has positive $m$-dimensional Hausdorff measure. In this paper we use ${\cal H}^{m}\left(x\right)$ to denote the $m$-dimensional Hausdorff measure on $\mathbb{R}^{d}$, which is defined as follows: For a set $A\subseteq \mathbb{R}^d$, define ${\cal H}^m(A) := \lim_{\delta \to 0} {\cal H}^m_\delta$, where, for any $\delta\in(0,\infty)$, $${\cal H}^m_\delta(A) := \inf\left\{\sum_{j=1}^\infty \alpha(m) \left(\frac{\text{diam}(C_j)}{2}\right)^m: A \subseteq \cup_{j=1}^\infty C_j, \text{diam}(C_j) \leq \delta\right\},$$ with $\alpha(m)=\frac{\pi^{m/2}}{\text{Gamma}((m/2)+1)}=\frac{\pi^{m/2}}{\int_{0}^{\infty} e^{-x} x^{m/2}dx}$. The $m$-dimensional Hausdorff measure becomes the standard $m$-dimensional Lebesgue measure in the lower dimensional $\mathbb{R}^m$; see, for example, evans2015measure for more details on the Hausdorff measure.

Motivating Examples

In this subsection, we provide some econometric examples for submanifold integrals, which roughly belong to two categories.

\paragraph{Category 1:} Examples in which researchers are interested in some aggregate parameters of estimated or optimized subpopulations. Then, either the first-order expansion in asymptotic analysis (of estimators), or the first-order condition for optimality, often takes the form of the time derivative of an integral with a changing region of integration $\Omega_{t}$, which produces a submanifold integral term by the generalized Leibniz rule:\footnote{See, for example, Theorem 4.2 in Chapter 9 of delfour2001shapes.} under regularity conditions, (ref) can be expressed as

align[align omitted — 444 chars of source]

where term (I) captures the effect of the change in the integrand $w_{t}\left(x\right)$ with the region of integration $\Omega_{t}$ held fixed, while term (II) captures the effect of the change in the region of integration $\Omega_{t}$ with the integrand $w_{t}\left(x\right)$ held fixed. Importantly, $\partial\Omega_{t}$, the boundary of $\Omega_{t}$, is often a submanifold of dimension $d-1$, and thus term (II) takes the form of an integral over the submanifold $\partial\Omega_t$ with respect to the $(d-1)$-dimensional Hausdorff measure, which can also be viewed as the surface measure on the boundary submanifold $\partial\Omega_t$. For the integrand terms, ${\bf n}_{t}\left(x\right)$ is the outward-pointing unit normal vector, ${\bf v}_{t}\left(x\right)$ is the velocity vector associated with the time movement in the $\partial\Omega_{t}$. For example, consider a simple domain $\Omega_{t}=\left[a_{t},b_{t}\right]\subset\mathbb{R}$, its boundary $\partial\Omega_t$ is a set of two points $\{a_t,b_t\}$, which is the $0$-dimensional submanifold, and the $0$-dimensional Hausdorff measure is simply the point counting measure. When $x$ is one-dimensional, (ref) becomes the standard Leibniz rule for univariate calculus: \[ \rest{\frac{d}{dt}\int_{a_{t}}^{b_{t}}\omega_{t}\left(x\right)dx}_{t=0}=\rest{\int_{a_{t}}^{b_{t}}\frac{\partial}{\partial t}\omega_{t}\left(x\right)dx}_{t=0}+\rest{\left(\omega_{t}\left(b_{t}\right)\frac{d}{dt}b_{t}-\omega_{t}\left(a_{t}\right)\frac{d}{dt}a_{t}\right)}_{t=0}. \] Throughout this paper we focus on situations where term (II) in (ref) is not vanishing. \paragraph{Category 2:} Researchers are interested (for some other reasons than above) in some aggregate parameters of certain boundary or marginal subpopulations with the boundary or margin characterized by a lower-dimensional submanifold.

Roughly speaking, Examples (ref)-(ref) are of Category 1, Examples (ref)-(ref) are of Category 2, while Examples (ref)-(ref) can be of both categories. We emphasize that we do not intend the categorization above to be exact nor exhaustive, but more to provide a high-level summary of the origins of submanifold integrals.

example[Maximum Score Estimation of Binary Choice Models] Consider any model that satisfies the following sign alignment restriction \begin{equation} h_{0}\left(x\right)\gtrless0\ \Leftrightarrow\ x^{'}\beta_{0}\gtrless0 \end{equation} where $h_{0}$ is a nonparametrically identified and estimable function of $x$ and $\beta_{0}$ is a $d$-dimensional parameter normalized to lie on the unit sphere, i.e., $\beta_{0}\in\mathbb{S}^{d-1}:=\left\{ \beta\in\mathbb{R}^{d}:\left\| \beta \right\|=1\right\} $. For example, in the following binary choice model with a conditional median independence restriction as in \citet*{manski1975maximum}, \[ y_{i}=\mathbb 1\left\{ X_{i}^{'}\beta_{0}+\epsilon_{i}\geq0\right\} ,\quad\text{med}\left(\rest{\epsilon_{i}}X_{i}\right)=0, \] the sign alignment restriction (ref) is satisfied with $h_{0}\left(x\right)=\mathbb{E}\left[\rest{Y_{i}}X_{i}=x\right]-\frac{1}{2}$. The population criterion function for maximum score estimator can then be written as \begin{align} W\left(\beta\right) & :=\int h_{0}\left(x\right)\mathbb 1\left\{ x^{'}\beta\geq0\right\} p\left(x\right)dx, \end{align} where $p\left(\cdot\right)$ denotes the density of $X_{i}$. Under appropriate conditions, $\beta_{0}$ can be point identified under a scale normalization, $\beta_{0}=\arg\max_{\beta\in\mathbb{S}^{d-1}}W\left(\beta\right).$ By the generalized Leibniz rule, the first order condition (FOC) for the optimality of $\beta_{0}$ is given by \begin{equation} {\bf 0}=\nabla_{\beta}W\left(\beta_{0}\right):=\int_{\mathcal{M}_0} h_{0}\left(x\right)xp\left(x\right)d{\cal H}^{d-1}\left(x\right), \end{equation} where $\mathcal{M}=\left\{x\in\mathcal{X}:~ x^{'}\beta_{0}=0\right\}$ is a $m=\left(d-1\right)$-dimensional hyperplane (see, e.g., \citet*{kim1990cube}). We note that the hyperplane $\mathcal{M}_0$ has Lebesgue measure $0$ in $\mathbb{R}^{d}$, which is why the identification of $\beta_{0}$ is referred to as a type of “thin-set identification” in \citet*{khan2010irregular}. Consequently, the right-hand side of (ref) cannot be represented by a regular Lebesgue integral, but instead by a Hausdorff integral over a $m=\left(d-1\right)$-dimensional manifold (hyperplane) in $\mathbb{R}^{d}$.
example[Optimal Linear Treatment Assignment] It has been recognized, say in \citet*{kitagawa2018should}, that optimal linear treatment assignment problem shares some similarity with maximum score estimation. Specifically, consider the problem of optimizing over a parametric family of treatment assignment rules that assigns the treatment status $0/1$ according to $\mathbb 1\left\{ x^{'}\beta\geq0\right\} $ for a given observed type $x$, where $\beta$ is a $d$-dimensional choice parameter. Then the welfare function of the assignment rule parameter $\beta$ is given by \begin{equation} W\left(\beta\right):=\int\mathbb 1\left\{ x^{'}\beta\geq0\right\} h_{0}\left(x\right)p\left(x\right)dx, \end{equation} where $h_{0}\left(x\right)=\mathbb{E}\left[\rest{Y_{i}\left(1\right)-Y_{i}\left(0\right)}X_{i}=x\right]$ is the conditional average treatment effect (CATE) for type $x$. Hence (ref) is of exactly the same form as the maximum score criterion function (ref). As a result, the FOC of (ref) for the optimal $\beta_{0}$ (under scale normalization) is again given by the submanifold integral (ref). Often times, researchers conduct (costly) experiments on and estimate CATE from a sample of moderate sample size, but the target population on which the treatment in question might be implemented can be of a much larger scale. In such settings $p(\cdot)$ may be known or estimable using a much larger sample size than that used to estimate $h_{0}$, and thus we may focus on the estimation error for $h_{0}$ in the optimization of $W\left(\beta\right)$.
example[Aggregate Parameter over Estimated Subpopulation] Consider the estimation of the following general welfare or value functional parameter \[ W\left(h_{0},f_{0}\right):=\int\mathbb 1\left\{ h_{0}\left(x\right)\geq0\right\}f_{0}\left(x\right)p\left(x\right)dx, \] where $h_{0}$ and $f_{0}$ may be both unknown but nonparametrically estimable. For example, if we set $h_{0}(x)=\text{CATE}(x)$ and $f_{0}(x)=1$, then $W\left(h_{0},f_{0}\right)$ becomes the value functional $V(h_0)$ defined in (ref). If we set $h_{0}(x)=\text{CATE}(x)$ and $f_{0}(x)=\lambda^{'}x$ for some fixed $\lambda \neq \textbf{0}$, then $W\left(h_{0},f_{0}\right)$ becomes the average characteristics of the subpopulation with nonnegative CATE. If we set $h_{0}\left(x\right)=f_{0}\left(x\right)=\text{CATE}\left(x\right)$, then $W\left(h_{0},f_{0}\right)$ becomes $W\left(\text{CATE},\text{CATE}\right)$---the welfare under the “first-best” treatment assignment. Alternatively, we may take $f_{0}\left(x\right)$ to be any other value/cost function associated with type $x$ (see \citep*{CCG2025}). When $h_{0}$ and $f_{0}$ are real-valued functions and are nonparametrically estimated by $\hat{h}$ and $\hat{f}$, the pathwise derivative of $W\left(h_{0},f_{0}\right)$ with respect to $\left(h_{0},f_0\right)$ is a key object in the characterization of the asymptotic behaviors of the plug-in estimator $W\left(\hat{h},\hat{f}\right)$. Generally, the pathwise derivative of $W$ in the direction of $\left[h-h_{0},f-f_{0}\right]$ is given by \begin{align} & DW\left(h_{0},f_{0}\right)\left[h-h_{0},f-f_{0}\right]\nonumber \\ :=\ & \lim_{t\searrow0}\frac{W\left(h_{0}+t\left(h-h_{0}\right),f_{0}+t\left(f-f_{0}\right)\right)-W\left(h_{0},f_{0}\right)}{t}\nonumber \\ =\ & \int_{\left\{ h_{0}\left(x\right)\geq0\right\} }\left(f\left(x\right)-f_{0}\left(x\right)\right)p\left(x\right)dx+\int_{\left\{ h_{0}\left(x\right)=0\right\} }\left(h\left(x\right)-h_{0}\left(x\right)\right)\frac{f_{0}\left(x\right)p\left(x\right)}{\left\| \nabla_{x}h_{0}\left(x\right) \right\|}d{\cal H}^{d-1}\left(x\right) \end{align} which consists of two terms: the first is the perturbation of the integrand $f_{0}$ in the direction of $[f-f_{0}]$ over the true region of integral $\left\{ h_{0}\left(x\right)\geq0\right\} $, while the second is the perturbation of (the boundary) region of integration induced by the perturbation of $h_{0}$ in the direction of $[h-h_{0}]$. While the first term is standard in the semiparametric estimation literature, the second term takes the nonstandard form of a Hausdorff integral over the $m=\left(d-1\right)$-dimensional submanifold ${\cal M}_0:=\left\{x\in \mathcal{X}:~ h_{0}\left(x\right)=0\right\}$.\footnote{In particular, $\left\| \nabla_{x}h_{0}\left(x\right) \right\|$, which measures the “thinness” of the level set, enters into the derivative formula explicitly here. In the previous examples with hyperplane boundaries, we took $g(x)=x^{'}\beta$ with $\beta\in\mathbb{S}^{d-1}$, and thus $\left\| \nabla_{x}g(x) \right\|=\left\| \beta \right\|=1$, which is why this term becomes implicit in formula (ref).} Hence, as long as $f_0(x)p(x)$ does not vanish on the submanifold ${\cal M}_0:=\left\{x\in \mathcal{X}:~ h_{0}\left(x\right)=0\right\}$, the second term in (ref) will remain as the leading term in the asymptotic behavior of $W\left(\hat{h},\hat{f}\right)-W(h_0,f_0)$. See our companion paper \citep*{CCG2025} for a detailed analysis of this example.
example[Average Treatment Effects under Propensity Score or Density Trimming]\footnote{We thank Tim Armstrong for suggesting this example and for sending us his note on this.} \citet*[CHIM thereafter]{crump2009dealing} proposes as a systematic approach to deal with limited overlap problems in the estimation of ATEs, and shows that the optimal subpopulation that minimizes the asymptotic variance of ATE under homoskedastic errors takes the form of propensity score trimming: \[ \text{ATE}_{p\text{-trimmed}}:=\mathbb{E}\left[\rest{\text{CATE}\left(X_{i}\right)}\alpha\leq p_{0}\left(X_{i}\right)\leq1-\alpha\right]. \] In practice, the propensity score $p_{0}\left(x\right)$ might require nonparametric estimation. CHIM did not provide theoretical results on the asymptotic distribution of their proposed estimators with the first-stage nonparametric estimation error of $p_{0}\left(x\right)$ taken into account. Clearly $\text{ATE}_{p\text{-trimmed}}$ is a two-sided version of Example (ref) with the region of integration defined by the two-sided inequality $\alpha\leq p_{0}\left(X_{i}\right)\leq1-\alpha$ on the propensity score function, and the directional derivative of $\text{ATE}_{p\text{-trimmed}}$ with respect to $p_{0}$ will again features submanifold integrals on the level sets of $\left\{ x:p_{0}\left(x\right)=\alpha\right\} $ and $\left\{ x:p_{0}\left(x\right)=1-\alpha\right\} $.
example[ReLU-Based Generalized Maximum Score Estimator under Multi-Index Single-Crossing Conditions] In a follow-up paper \citep*{CGW2025}, we investigate a ReLU-based generalization of the maximum score estimator for multi-index single-crossing condition (MISC) models gao2019robust. Formally, consider a random sample $\left(Y_{i},X_{i}\right)_{i=1}^{n}$ where $Y_{i}$ is an outcome with support ${\cal Y}\subseteq\mathbb{R}^{d_{y}}$ and $X_{i}:=\left(X_{i1},...,X_{iJ}\right)\in\mathbb{R}^{d\times J}$ is a $d\times J$ random matrix with support ${\cal X}\subset\mathbb{R}^{d\times J}$. Let $h_{0}:{\cal X}\to\mathbb{R}$ be a real-valued functional of the conditional distribution of $Y_{i}$ given $X_{i}$ that is directly identified and nonparametrically estimable from the data, and $\beta_{0}\in \mathbb{S}^{d-1}$ be a finite-dimensional parameter. We say that $\left(h_{0},\beta_{0}\right)$ satisfies the (multi-index single-crossing) MISC condition if, for all $x=\left(x_{1},...,x_{J}\right)\in{\cal X}$, \begin{align} x_{j}^{'}\beta_{0}>0,\ \forall j=1,...,J & \quad\Rightarrow\quad h_{0}\left(x\right)\geq0,\nonumber \\ x_{j}^{'}\beta_{0}<0,\ \forall j=1,...,J & \quad\Rightarrow\quad h_{0}\left(x\right)\leq0. \end{align} The condition is said to be strict if the inequalities on the right-hand side of (ref) are strict, i.e., $h_{0}\left(x\right)>0$ whenever $x_{j}^{'}\beta_{0}>0$ for all $j$, and $h_{0}\left(x\right)<0$ whenever $x_{j}^{'}\beta_{0}<0$ for all $j$. Note that, when $J=1$, (ref) reduces exactly to the sign-alignment restriction in the original maximum score formulation manski1975maximum,manski1985semiparametric. However, the maximum score criterion function in manski1975maximum,manski1985semiparametric does not generalize to multi-index settings with $J\geq2$. CGW2025 proposes a novel criterion function the encodes the sign restrictions in the MISC condition framework using ReLU functions (instead of indicator functions as Manski's maximum score criterion function). Specifically, given $h:{\cal X}\to\mathbb{R}$ and $\beta\in \mathbb{S}^{d-1}$, define \begin{align} g_{+,\beta,h}\left(x\right) & :=\left[h\left(x\right)-\min_{1\le j\le J}\left(-x_{j}^{'}\beta\right)_{+}\right]_{+},\nonumber\\ g_{-,\beta,h}\left(x\right) & :=\left[-h\left(x\right)-\min_{1\le j\le J}\left(x_{j}^{'}\beta\right)_{+}\right]_{+}. \end{align} CGW2025 then constructs the following criterion function \[ W\left(\beta\right):=W_{+}\left(\beta\right)+W_{-}\left(\beta\right),\qquad W_{\pm}\left(\beta\right):=\mathbb{E}\left[g_{\pm,\beta,h_{0}}\left(X_{i}\right)\right]. \] which is shown to ensure $\beta_{0}\in\arg\max_{\beta\in \mathbb{S}^{d-1}}W\left(\beta\right)$ under the MISC condition. CGW2025 then considers two estimation procedures: a two-step procedure in which a first-stage nonparametric estimator of $h_0$ is plugged into (ref) as well as a joint machine learning procedure by formulating the criterion $W$ as a special layer in a deep neural network (DNN) training problem. On a theoretical level, CGW2025 shows that the MISC criterion $W(\beta)$ also features “thin-set” identification with information concentrated on a $(d-1)$-dimensional submanifold that piecewisely takes the form of a hyperplane defined by $\{x_j'\beta_0=0\}$ across $j=1,...,J$. Then the submanifold integral results in this paper are utilized to establish the convergence rate and asymptotic normality of the ReLU-based generalized maximum score estimator in CGW2025.
example[Generalized Partial Means and Nonparametric Regression with Generated Covariates] Suppose that we are interested in the conditional mean of $\psi\left(Y_{i},X_{i}\right)$ given $g\left(X_{i}\right)=c\in\mathbb{R}^{d-m}$, where $\psi$ is some known transformation of $Y_{i}$ and $X_{i}$: \[ \mathbb{E}\left[\rest{\psi\left(Y_{i},X_{i}\right)}g\left(X_{i}\right)=c\right]=\int_{\left\{ g\left(x\right)=c\right\} }h_{0}\left(x\right)w_{0}\left(x\right)d{\cal H}^{m}\left(x\right) \] with $h_{0}\left(x\right):=\mathbb{E}\left[\rest{\psi\left(Y_{i},X_{i}\right)}X_{i}=x\right]$, $w_{0}\left(x\right):=p_{0}\left(\rest xg\left(x\right)=c\right)=\frac{p_{0}\left(x\right)}{\mathcal{J}_g(x)p_{g\left(X\right)}\left(c\right)}$, $p_{0}\left(x\right)$ being the density of $X_{i}$, and $p_{g\left(X\right)}$ is the density of $g(X_i)$. In particular, when $g$ extracts a subvector of $x\in\mathbb{R}^{d}$, say, $g\left(x\right)=\left(x_{1},...,x_{d-m}\right)$ for some $1\leq m<d$, the above specializes to the partial mean functional of the form $\mathbb{E}\left[\rest{Y_{i}}X_{i,1}=c_{1},...,X_{i,d-m}=c_{d-m}\right]$ as studied in \citet*{newey1994kernel}. The weighted boundary average treatment effect considered in cattaneo2025dist is another example with the level set function $g$ given by a known distance function. When $g$ is nonparametrically estimated, the above can be also viewed as the nonparametric regression function with nonparametrically generated covariates as considered in mammen2012nonparametric.\footnote{mammen2012nonparametric develops the asymptotic theory for generated covariates using kernel first-stage estimators, and establishes that, under appropriate conditions, the estimation errors from the first-stage are “smoothed out” partially, which can be viewed as a special case of the rate acceleration result developed in our paper.}
example[Marginal Treatment Effects and Policy Relevant Treatment Effects] \citet*{heckman2005structural,heckman2007econometric} proposes the marginal treatment effect (MTE) as a unifying concept that underlies a wide variety of treatment effects studied in causal inference and program evaluations. It is well-known that MTE can be written as a derivative with respect to the propensity score: \begin{align} MTE\left(x,p\right) & :=\mathbb{E}\left[\rest{Y_{i}\left(1\right)-Y_{i}\left(0\right)}X_{i}=x,p_{0}\left(Z_{i}\right)=p\right],\nonumber \\ & = \frac{\partial}{\partial p}\mathbb{E}\left[\rest{Y_i} X_i= x, p_0(Z_i)=p\right] \end{align} where $p_{0}\left(z\right):=\mathbb{E}\left[\rest{D_{i}=1}Z_{i}=z\right]$ is the propensity score (for the binary treatment $D\in\{0,1\}$) given the instruments $Z_{i}$, $f_{0}\left(z\right|x)$ is the conditional density of $Z_{i}$ given the observed covariates $X_i=x$, and $Y_i$ denotes the observed outcome. It is clear from (ref) that MTE takes the form of the derivative of a submanifold integral, where the submanifold is given by the propensity score level set $p_{0}\left(z\right)=p$. Here, both the $\text{CATE}$ function and the propensity score function can be nonparametrically estimated. As argued in \citet*{carneiro2010evaluating}, in many scenarios it is useful to consider incremental policy reforms and focus on the analysis of marginal policy changes, whose effects are usually concentrated on individuals “at the margin”. The average marginal treatment effect (AMTE) summarizes the mean benefit of the treatment for the subpopulation who is indifferent between participation in the treatment and nonparticipation. Formally, this “at-the-margin” population $\left\{ \left(z,u\right)\in\mathbb{R}^{d_{z}+1}:p_{0}\left(z\right)-u=0\right\} $ is a submanifold defined by the level set of the propensity score function $p$, and the average over this subpopulation can again be represented by a submanifold integral in the form of \[ \text{AMTE}\left(x\right):=\int_{\left\{ p_{0}\left(z\right)=u\right\} }\text{MTE}\left(x,u\right)f_{0}\left(\rest zx\right)d{\cal H}_{\left(z,u\right)}^{d_{z}}. \] In addition, MTE is also used to define a wide range of causal parameters, notably the policy-relevant treatment effect (PRTE) proposed in \citet*{heckman2001policy}. Let $F_{P}$ denote the CDF of $P_{i}:=p_{0}\left(Z_{i}\right)$ and suppose that a proposed policy changes this CDF from $F_{P}$ to $\tilde{F}$. Then the PRTE is defined by \[ \text{PRTE}\left(x,\tilde{F}\right):=\int_{0}^{1}\text{MTE}\left(x,u\right)\omega_{PRTE}\left(u,\tilde{F}\right)du,\quad\omega_{PRTE}\left(u,\tilde{F}\right):=\frac{F_{P}\left(u\right)-\tilde{F}\left(u\right)}{\mathbb{E}_{\tilde{F}}\left(P_{i}\right)-\mathbb{E}_{F_{P}}\left(P_{i}\right)}. \] Since $\text{MTE}$ is defined as a derivative, PRTE takes the form of an average derivative parameter, for which a parametric estimation rate of $n^{-1/2}$ is possible under a certain “vanishing-on-boundary” condition (see Example (ref) for a related discussion). \citet*{carneiro2010evaluating} pointed out that this type of “vanishing-on-boundary” condition may be violated in plausible economic policy scenarios in the analysis of PRTE, which makes it not estimable at a parametric rate $n^{-1/2}$ in general. Our theoretical results can be useful in scenarios where “vanishing-on-boundary” conditions are not imposed in estimation and testing PRTE in various applications.
example[Weighted Average Derivatives] Related to Example (ref) above, suppose that we are interested in the following weighted average derivative of $h_0$: \[ \text{WAD}\left(h_{0}\right):=\int\nabla_x h_{0}\left(x\right)w\left(x\right)dx. \] For example, \citet*[PSS thereafter]{powell1989} considers a density-weighted average derivative, which is a special case of the above with $w\left(x\right):=p^{2}\left(x\right)$, with $p\left(x\right)$ denoting the marginal density of $x$ on its support ${\cal X}$. The standard analysis of $\text{WAD}\left(\hat{h}\right)$ exploits the following integration-by-parts formula (a special case of the Divergence Theorem): \begin{equation} \int_{{\cal X}}\nabla_x h_{0}\left(x\right)w\left(x\right)dx=\int_{\partial{\cal X}}\overrightarrow{{\bf n}}\left(x\right) h_{0}\left(x\right)w\left(x\right)d{\cal H}^{d-1}\left(x\right)-\int_{{\cal X}}h_{0}\left(x\right)\nabla_x w\left(x\right)dx, \end{equation} where the first term is a Hausdorff integral over the submanifold defined by the support boundary $\partial{\cal X}$. The standard analysis of $\text{WAD}\left(h_{0}\right)$, such as in PSS and \citet*{newey1993efficiency}, exploits a “vanishing-on-boundary” assumption, say, $h_{0}\left(x\right)w\left(x\right)=0$ for all $x\in{\cal \partial X}$, so that the first term in (ref) degenerates to $0$. For example, PSS assumes that $p\left(x\right)=0$ on $\partial{\cal X}$ (Assumption 2 in their paper), and hence $\text{WAD}\left(h_{0}\right)=-2\int h_{0}\left(x\right)\nabla_x p\left(x\right)p\left(x\right)dx$. PSS then proceeds to establish that $\sqrt{n} \left(\text{WAD}\left(\hat{h}\right)-\text{WAD}\left(h_{0}\right)\right)=O_p (1)$ and is asymptotic normal. As clear from (ref), the parametric convergence rate of $n^{-1/2}$ relies crucially on the “vanishing on boundary” assumption, without which the leading asymptotic term in $\text{WAD}\left(\hat{h}\right)-\text{WAD}\left(h_{0}\right)$ would be the first term (the submanifold integral) that converges at a slower-than-$n^{-1/2}$ rate, as established in our paper.
example[Structural Functions in NPIV Regression] The previous examples focus on exogenous nonparametric regression functions. Here we present a canonical example of a nonparametric function with endogeneity: the structural functions in the nonparametric instrumental variables (NPIV) model, i.e., the $h_{0}$ function below: \[ Y_{i}=h_{0}\left(X_{i}\right)+\epsilon_{i},\quad\mathbb{E}\left[\rest{\epsilon_{i}}Z_{i}\right]=0. \] Previous work by, for example, \citet*{chenpouzo2015sieve}, \citet*{chen2018optimal}, and \citet*{chen2025adaptive}, provides theoretical results on the estimation and inference on $h_{0}$, as well as various (families of) linear and nonlinear functionals of $h_{0}$, including point evaluation, average derivatives/elasticities, and consumer surplus/welfare functionals. In principle, the $h_{0}$ function in all the submanifold functionals defined in Examples (ref)-(ref) above could be replaced by the structural function in the NPIV model here, and the theory about submanifold integrals developed in the current paper would still be relevant.
example[Nonparametric Quantile] Similarly, the nonparametric function $h_{0}$ need not be restricted as a conditional expectation function; it can be broadly defined as any estimable nonparametric function, for example, nonparametric quantile and density functions. These general nonparametric estimation problems (without endogeneity) have been studied widely in statistics, econometrics, and beyond in various settings. In the presence of endogeneity issues, \citet*{HL07econometrica}, \citet*{chenpouzo2015sieve}, and \citet*{chen2019penalized} have also provided theoretical results for the nonparametric quantile IV regression. The established theoretical results from these previous studies can be combined with the theory of submanifold integrals developed here for the study of new functional parameters of interest.

In summary, we hope that the above examples illustrate the need for estimation and inference on general submanifold integrals of unknown functions, which has not been systematically studied in the existing literature. In the rest of the paper, we shall focus on the case where $h_0$ is a nonparametric regression function, and provide a general analysis of the estimation and inference for linear submanifold integrals (ref), nonlinear submanifold integrals (ref) whose first-order linear approximation takes the form of (ref), as well as integrals on upper contour set of $h_0$ (ref). We present additional results when $h_0$ is a nonparametric density function and a nonparametric instrumental variables regression in the Appendix.

Minimax Rate of Estimation: Lower Bound

We first establish minimax lower bound rates for estimation of linear ($L(h_0)$) and nonlinear ($\Gamma(h_0)$) submanifold integrals when $h_0$ is a nonparametric regression $\mathbb{E}[Y|X]$ or a nonparametric density of $X$ in Section (ref), and establish similar lower bound results when $h_0$ is an NPIV structural function in Subsection (ref).

Throughout this paper we assume that $h_0$ belongs to H\"{o}lder class of functions. We now recall the definition of H\"{o}lder class of functions. Let $\mathcal{X}=\mathcal{X}_{1}\times...\times\mathcal{X}_{d}$ be the Cartesian product of compact intervals $\mathcal{X}_{1},\dots,\mathcal{X}_{d}$, say, ${\cal X}=\left[0,1\right]^{d}$ for simplicity. A real-valued function $h$ on $\mathcal{X}$ is said to satisfy a H\"{o}lder condition with exponent $\gamma\in(0,1]$ if there is a positive number $c$ such that $\left|h(x)-h(y)\right|\leq c\left\| x-y \right\|^{\gamma}$ for all $x,y\in\mathcal{X}$; here $\left\| x-y \right\|=\bigl(\sum_{l=1}^{d}x_{l}^{2}\bigr)^{1/2}$ is the Euclidean norm of $x=\left(x_{1},\dots,x_{d}\right)\in\mathcal{X}$. Given a $d$-tuple $\alpha=\left(\alpha_{1},\dots,\alpha_{d}\right)$ of nonnegative integers, set $\left[\alpha\right]=\alpha_{1}+\dots+\alpha_{d}$ and let $\nabla^{\alpha}$ denote the differential operator defined by \[ \nabla^{\alpha}:=\nabla_x^{\alpha}=\frac{\partial^{\left[\alpha\right]}}{\partial x_{1}^{\alpha_{1}}\dots\partial x_{d}^{\alpha_{d}}}. \] Let $\lfloor s\rfloor$ be a nonnegative integer that is smaller than $s$, and set $s=\lfloor s\rfloor+\gamma$ for some $\gamma\in(0,1]$. A real-valued function $h$ on $\mathcal{X}$ is said to be $s$-$smooth$ if it is $\lfloor s\rfloor$ times continuously differentiable on $\mathcal{X}$ and $\nabla^{\alpha}h$ satisfies a H\"{o}lder condition with exponent $\gamma$ for all $\alpha$ with $\left[\alpha\right]=\lfloor s\rfloor$. Denote the class of all $s$-smooth real-valued functions on $\mathcal{X}$ by $\Lambda^{s}(\mathcal{X})$ (called a H\"{o}lder class), and the space of all $\lfloor s\rfloor$-times continuously differentiable real-valued functions on $\mathcal{X}$ by $C^{\lfloor s\rfloor}(\mathcal{X})$. Define a H\"{o}lder ball with smoothness $s=\lfloor s\rfloor+\gamma$ as \[ \Lambda_{c}^{s}\left(\mathcal{X}\right)=\left\{ h\in C^{\lfloor s\rfloor}(\mathcal{X}):\sup_{[\alpha]\leq\lfloor s\rfloor}\sup_{x\in\mathcal{X}}\left|\nabla^{\alpha}h(x)\right|\leq c,\sup_{[\alpha]=\lfloor s\rfloor}\sup_{x,y\in\mathcal{X},x\neq y}\frac{\left|\nabla^{\alpha}h\left(x\right)-\nabla^{\alpha}h\left(y\right)\right|}{\left\| x-y \right\|^{\gamma}}\leq c\right\} . \]

remIf $h_0\in \Lambda^s$ (for H\"older smooth $s>0$), then $h_0|_{\mathcal{M}}$ is pointwise well-defined and hence the submanifold integrals $L(h_0)$ and $\Gamma(h_0)$ are well-defined. If $h_0$ is only assumed to belong to a Sobolev/Besov space of smoothness $s>0$ but not sup-norm bounded, then one needs to impose a Sobolev/Besov trace condition $s>(d-m)/2$ to ensure that $L(h_0)$ and $\Gamma(h_0)$ are well defined. This is why we assume that $h_0\in \Lambda^s$ for $s>0$ instead of a Sobolev/Besov space of smoothness $s>(d-m)/2$ in this paper.

In the proofs of the minimax lower bound rate results, especially in the proofs of lower bounds for nonlinear submanifold integrals with possibly unknown submanifolds that could depend on $h_0$, we assume that the submanifold function $g$ is H\"{o}lder smooth.

assumption[Smoothness of Submanifolds] $g\in \Lambda^{s_1}({\cal X})$ for some $s_1 \geq \max\{1,s\}$.

We wish to stress that the minimax lower bound rates established in our paper are still valid lower bound rates without imposing the extra smoothness Assumption (ref). Nevertheless, Assumption (ref) is satisfied in typical economics applications with unknown submanifolds.

Case 1: $h_0$ is a Regression or a Density

In Subsection (ref), we first establish a minimax rate lower bound for estimating linear integrals on submanifolds. In Subsection (ref) we then present a general rate lower bound for estimating general nonlinear integrals over submanifolds whose first derivatives take the form of the linear integrals on submanifolds.

For Linear Functional $\theta_0=L(h_0)$

In this subsection we present the minimax lower bound rates for estimation of $\theta_0=L(h_0)$ when $h_{0} (\cdot)$ is the unknown regression function $\mathbb{E}\left[\rest{Y_{i}}X_{i}=\cdot\right]$ (a nonparametric regression) or the unknown density of $X$.

Recall that $L:\Lambda^s \mapsto \mathbb{R}$ is the linear submanifold integral functional

equation[equation omitted — 127 chars of source]

where both the level set function $g$ (consequently the manifold ${\cal M}$) and the weight function $w$ are assumed to be known. We first establish lower bounds for the minimax convergence rates for estimating $\theta_{0}=L\left(h_{0}\right)$ when $h_{0}$ is assumed to belong to the H\"{o}lder class of functions with smoothness $s>0$. The established lower bound rates hold for any possible consistent estimators of $L(h_0)$.

We impose the following standard assumptions on the density of the covariates, the known weight function and the regression error term.

assumption[Density on ${\cal X}$] ${\cal X}=\left[0,1\right]^{d}$, and the density $p_0$ of $X_{i}$ on ${\cal X}$ is uniformly bounded away from zero and infinity.
assumption[Weight Regularity] The weight function $w:\mathcal{X} \mapsto \mathbb{R}$ is uniformly bounded on $\mathcal{M}$ and not identically zero on $\mathcal{M}$.
assumption[Nondegeneracy] $\inf_{x\in{\cal X}}\mathbb{E}\left[\rest{\epsilon_{i}^{2}}X_{i}=x\right]>0$ with $\epsilon_{i}:=Y_{i}-h_{0}\left(X_{i}\right)$.
thm[Lower Bound Rate for $\theta_0=L(h_0)$ When $h_0$ is a Regression] Under Assumptions (ref), (ref), (ref), (ref) and (ref), the rate $r^*_{n}=n^{-\frac{s}{2s+d-m}}$ is the minimax rate lower bound for the estimation of $\theta_{0}\left(P,w\right):=\int_{{\cal M}}h_{0}\left(x\right)w\left(x\right)d{\cal H}^{m}\left(x\right)$, i.e. \[ \liminf_{n\longrightarrow\infty}\inf_{\tilde{\theta}}\sup_{P,w}\mathbb{E}_{P}\left[n^{\frac{2s}{2s+d-m}}\left(\tilde{\theta}-\theta_{0}\left(P,w\right)\right)^{2}\right]\geq c,\quad\text{for some constant }c>0. \] where $P$ is any joint probability distribution of $\left(X_{i},Y_{i}\right)$ that satisfies $h_{0}\left(\cdot\right):=\mathbb{E}_{P}\left[Y_{i}|X_{i}=\cdot\right]\in\Lambda_{c}^{s}\left(\mathcal{X}\right)$ along with Assumptions (ref) and (ref), $w$ is any weight function satisfying Assumption (ref), and $\tilde{\theta}$ is any estimator of $\theta_{0}:=\theta_{0}\left(P,w\right)$.

Note that the lower bound rate $r^*_{n}=n^{-\frac{s}{2s+d-m}}$ reproduces several well-known results in the literature as special cases: When $m=d$, the lower bound rate becomes $r^*_{n}=n^{-\frac{1}{2}}$, reproducing the standard parametric convergence rate for regular (full-dimensional) integral functionals.\footnote{When $m=d$, the Hausdorff measure ${\cal H}^{d}$ coincides with Lebesgue measure in $\mathbb{R}^{d}$.} When $m=0$, the lower bound becomes $r^*_{n}=n^{-\frac{s}{2s+d}}$, reproducing the well-known stone1980optimal minimax optimal rate for point evaluation functionals of the $d$-dimensional nonparametric regression.\footnote{When $m=0$, the Hausdorff measure ${\cal H}^{0}$ becomes a point counting measure.}

Theorem (ref) can be easily adapted to establish a minimax lower bound for $L(h_0)$ for other nonparametric object $h_0$ such as nonparametric density function of $X$ or nonparametric quantile regression of $Y$ on $X$. The following result shows that $r^*_n$ is also a minimax lower bound rate when $h_0$ is a nonparametric density function of $X$.

cor[Lower Bound Rate for $\theta_0=L(h_0)$ When $h_0$ is a Density] Under Assumptions (ref), (ref) and (ref), the rate $r^*_{n}=n^{-\frac{s}{2s+d-m}}$ is the minimax rate lower bound for the estimation of $\theta_{0}\left(P,w\right):=\int_{{\cal M}}h_{0}\left(x\right)w\left(x\right)d{\cal H}^{m}\left(x\right)$, i.e. \[ \liminf_{n\longrightarrow\infty}\inf_{\tilde{\theta}_n}\sup_{P,w}\mathbb{E}_{P}\left[n^{\frac{2s}{2s+d-m}}\left(\tilde{\theta}-\theta_{0}\left(P,w\right)\right)^{2}\right]\geq c,\quad\text{for some constant }c>0. \] where $P$ is any probability distribution of $X_i$ with its density satisfying $h_0\in\Lambda_{c}^{s}\left(\mathcal{X}\right)$ and Assumption (ref), $w$ is any uniformly bounded weight function, and $\tilde{\theta}_n$ is any estimator of $\theta_{0}:=\theta_{0}\left(P,w\right)$.

For Nonlinear Functional $\theta_0=\Gamma(h_0)$

In the previous subsection we focus on the linear integral functionals over a known submanifold, but the key minimax rate $r^*_n = n^{-\frac{s}{2s+d-m}}$ continues to be relevant for a general class of nonlinear submanifold integral functionals.

Specifically, we consider the class of nonlinear integral functionals $\Gamma\left(h_{0}\right)$ whose pathwise derivative takes the form of a linear submanifold integral of form (ref).

assumption[Linearization] (i) $\Gamma:\mathbb{H}\to\mathbb{R}$ is pathwise differentiable with \begin{equation} D\Gamma\left(h_{0}\right)\left[v\right] :=\rest{\frac{\partial}{\partial t} \Gamma(h_0+tv)}_{t=0} =\int_{\left\{ x:g\left(x\right)={\bf 0}\right\} }v\left(x\right)\ol{w}\left(x\right){\cal H}^{m}\left(x\right),\ \forall v\in{\cal V}:=\mathbb{H}-\left\{ h_{0}\right\} , \end{equation} for some uniformly bounded function $\ol{w}$ with $\int_{\{x:g(x)={\bf 0}\}}\ol{w}^2(x)d\mathcal{H}^m(x)\in(0,\infty)$, some level set function $g:{\cal X}\to\mathbb{R}_{d-m}$ satisfying Assumption (ref) uniformly in $h_0\in\Lambda^s_c$, where both $\ol{w}$ and $g$ may depend on the true (and unknown) $h_{0}$; \newline (ii) There is a constant $0<C_r<1$ such that for any $u\in\Lambda^s_c$ with small $\left\| u \right\|_{\infty}$, \[ \left| \Gamma(h_0+u)-\Gamma(h_0)-D\Gamma(h_0)[u] \right|\le C_r\left| D\Gamma(h_0)[u] \right|. \]

Assumption (ref)(ii) is used to prevent the linearization remainder from canceling the first-order separation in the lower bound proof. It is easily satisfied by a Taylor series expansion to second order. It is also satisfied if $\Gamma$ is Fr\'echet differentiable at $h_0$ with respect to the sup-norm, i.e., there exists a function $\omega:(0,\infty)\to(0,\infty)$ with $\omega(t)\to 0$ as $t\downarrow 0$ such that for all $u\in\Lambda^s_c$ with $\left\| u \right\|_\infty\le t$, \[ \left| \Gamma(h_0+u)-\Gamma(h_0)-D\Gamma(h_0)[u] \right| \;\le\; \omega(t)\,\left\| u \right\|_\infty. \]

Assumption (ref) for nonlinear functional $\Gamma(h_0)$ replaces the role of Assumption (ref) for the linear functional $L(h_0)$. We emphasize that the level set function $g$ (and thus the submanifold) is not necessarily known in Assumption (ref). For instance, Examples 3-6 in Section (ref) feature unknown and estimated submanifolds. Also, the upper contour integral $V(h_0)$ of the form (ref) and the surface integral $S(h_0)$ of the form (ref) are examples in which $g=h_0$ (nonparametrically defined).

Under Assumption (ref), the minimax rate lower bound established in Theorem (ref) for linear submanifold integrals continues to be valid for nonlinear submanifold functionals $\Gamma$ as well. This is because the pathwise derivative $D\Gamma(h_0)[\cdot]$ is a linear submanifold integral with minimax convergence rate given by $r^*_n$ in Theorem (ref), and the linearization remainder can only (potentially) slow down the rate. We summarize this observation in the following Corollary.

cor[Lower Bound Rate for $\theta_0=\Gamma(h_0)$ When $h_0$ is a Regression or a Density] Let $\Gamma$ be a potentially nonlinear functional satisfying Assumption (ref), and $h_0$ be a nonparametric regression or a density of $X$. Let Assumptions (ref), (ref), and (ref) (when $h_0$ is a regression) hold. Then, a minimax rate lower bound for estimating $\theta_0:=\Gamma(h_0)$ is given by $r^*_n = n^{-s/(2s + d - m)}$, i.e., \[ \liminf_{n\longrightarrow\infty}\inf_{\tilde{\theta}}\sup_{P}\mathbb{E}_{P}\left[n^{\frac{2s}{2s+d-m}}\left(\tilde{\theta}-\theta_{0}\left(P\right)\right)^{2}\right]\geq c,\quad\text{for some constant }c>0. \]

When $m=d-1$, the lower bound in Corollary (ref) becomes $r^*_{n}=n^{-\frac{s}{2s+1}}$ for nonlinear submanifold integrals. This could be viewed as an extension of horowitz1993et's result on the minimax lower bound rate for smoothed maximum score estimation.

Case 2: When $h_0$ is an NPIV Structural Function

For Linear Functional $\theta_0=L(h_0)$

We now consider the rate lower bound for estimating the linear $\theta_0=L(h_0)$ when $h_0$ is the structural function in the nonparametric instrumental variable (NPIV) model:

align[align omitted — 99 chars of source]

where $X_i \in \mathcal{X} \subset \mathbb{R}^d$ are endogenous regressors and $Z_i \in \mathcal{Z} \subset \mathbb{R}^{d_z}$ are instruments.

Let $T: L^2(\mathcal{X}) \rightarrow L^2(\mathcal{Z})$ be the conditional expectation operator $(Th)(z)=\mathbb{E}[h(X)|Z=z]$. We now introduce some basic conditions on the NPIV model ((ref)).

assumption(i) $X_i$ has compact rectangular support $\mathcal{X} \subset \mathbb{R}^d$ with nonempty interior and the density of $X_i$ is uniformly bounded away from $0$ and $\infty$ on $\mathcal{X}$; (ii) $Z_i$ has compact rectangular support $\mathcal{Z} \subset \mathbb{R}^{d_z}$ and the density of $Z_i$ is uniformly bounded away from $0$ and $\infty$ on $\mathcal{Z}$; (iii) $\inf_{z\in \mathcal{Z}} \mathbb{E}[\epsilon_i^2|Z_i = z] \geq \underline \sigma^2 > 0$ with $\epsilon_i:= Y_i - h_0(X_i)$ satisfying NPIV model (ref); (iv) $h_0 \in \Lambda^s_c$.
assumptionFor any $h_1,h_2\in \Lambda_c^s$, $Th_1 = Th_2 \ \text{in } L^2(\mathcal{Z})$ implies $L(h_1)=L(h_2)$.

To establish a lower bound, we require a link condition that relates smoothness of $T$ to the parameter space for $h_0 \in \Lambda^s_c$. For simplicity we use the condition similar to that used in chen2011rate and chen2018optimal. Let $\{\psi_{\tilde{\jmath},k,\tilde{G}}\}$ denote a tensor-product CDV wavelet basis for $L^2(\mathcal{X})\cap L^{\infty}(\mathcal{X})$ of regularity strictly larger than $s$. See, e.g., Appendix A of chen2018optimal for details on the construction and properties of this basis.

assumption[NPIV Link Condition] There is a positive decreasing function $\nu:(0, \infty)\to(0, \infty)$ such that \begin{equation} \|T h\|_{L^2(\mathcal{Z})}^2 \;\lesssim\; \sum_{\tilde{\jmath},\tilde{G},k} \nu(2^{\tilde{\jmath}})^2\, \bigl\langle h, \psi_{\tilde{\jmath},k,\tilde{G}}\bigr\rangle_{L^2(\mathcal{X})}^2, \qquad \forall h \in \Lambda_c^s(\mathcal{X}). \end{equation}

In the NPIV literature (see HallHorowitz and chen2011rate), the mildly ill-posed case corresponds to $\nu(t) = t^{-\varsigma}$, roughly saying that the operator $T$ converts $s$-smooth functions of $X$ into $(\varsigma+s)$-smooth functions of $Z$. The severely ill-posed case, which corresponds to choosing $\nu(t) = \exp(-\frac{1}{2}t^{\varsigma})$ and roughly says that $T$ maps smooth functions of $X$ into “supersmooth” functions of $Z$.

thm[Lower Bound Rate for $\theta_0=L(h_0)$ when $h_0$ is a NPIV] Let Assumptions (ref), (ref), (ref) and (ref) hold for the NPIV model. We have \begin{align*} \liminf_{n \to \infty} \inf_{\tilde{\theta}_n} \sup_{h \in \Lambda_c^s(\mathcal{X})} \mathbb{E}_{P_{h}} \left[r^{-2}_{NPIV,n} (\tilde{\theta}_n - L(h))^2 \right] \gtrsim C>0 \end{align*} where \begin{equation} r_{NPIV,n} := \begin{cases} n^{-\frac{s}{2(s+\varsigma) + d - m} }& if mildly ill-posed, \\ (\log n)^{-\frac{s}{\varsigma}} & if severely ill-posed, \end{cases} \end{equation} and $\inf_{\tilde{\theta}_n}$ denotes the infimum over all estimators of $L(h)$ based on the sample of size $n$, $\sup_{h \in \Lambda_c^s(\mathcal{X})} \mathbb{E}_{P_{h}}$ denotes the sup over $h \in \Lambda_c^s(\mathcal{X})$ and probability distributions of $(X_i,Z_i,u_i)$ that satisfy Assumptions (ref) and (ref) with fixed $\nu$, and the finite constant $C>0$ does not depend on $n$.

According to Theorem (ref), for severely ill-posed NPIV problems, the minimax lower bound for estimating linear submanifold integral $L(h_0)$ coincides with those in $L^2(X)$ (chen2011rate) and in sup-norm (chen2018optimal) estimation of the NPIV function $h_0$ itself, which is also the same as estimating a NPIV function in $L^2(X)$ estimation of a NPIV function with $(d-m)$ dimensional endogenous regressors. For mildly ill-posed NPIV problems, the minimax lower bound for estimating $L(h_0)$ coincides with those in $L^2(X)$ estimation of a NPIV function with $(d-m)$ dimensional endogenous regressors.

rem[Bounded completeness is not needed for the minimax lower bound for estimating $L(h_0)$] Note that Assumption (ref) is equivalent to \[ L(h)=0\ \text{for all } h\in \ker(T)\cap \Lambda_c^s(\mathcal{X}), \ \text{with~} \ker(T):=\{h\in L^2(\mathcal{X}): Th=0 \text{ in } L^2(\mathcal{Z})\}, \] which is strictly weaker than bounded completeness of $T$. Indeed, bounded completeness implies $\ker(T)=\{0\}$ and hence implies Assumption (ref), but the proof of Theorem (ref) does not use injectivity of $T$: it only constructs two alternatives $h_0,h_1$ such that (a) the induced distributions are close in Kullback--Leibler distance (through small $\|T(h_1-h_0)\|_{L^2(\mathcal{Z})}$), while (b) $|L(h_1)-L(h_0)|$ is separated at the stated rate. Some identification is nevertheless needed to make estimation of $L(h_0)$ nontrivial. If Assumption (ref) fails, then there exists $h^\dagger\in \ker(T)\cap \Lambda_c^s(\mathcal{X})$ with $L(h^\dagger)\neq 0$. For any $h_0\in\Lambda_c^s(\mathcal{X})$ and sufficiently small $t\neq 0$ such that $h_1:=h_0+t h^\dagger\in\Lambda_c^s(\mathcal{X})$, we have $Th_1=Th_0$ and hence the models indexed by $h_0$ and $h_1$ are observationally indistinguishable, but $L(h_1)-L(h_0)=t\,L(h^\dagger)\neq 0$. Hence, for any estimator $\tilde\theta_n$, \[ \sup_{h\in\{h_0,h_1\}} \mathbb{E}_{P_h}\!\left[(\tilde\theta_n-L(h))^2\right] \;\ge\; \frac{(L(h_1)-L(h_0))^2}{4} \;=\; \frac{t^2\,L(h^\dagger)^2}{4}, \] so no consistent estimation of $L(h_0)$ is possible. Thus, Assumption (ref) is a natural minimal identification condition for $L(h_0)$, whereas bounded completeness of $T$ is unnecessary for establishing the lower bound rate in Theorem (ref).

For Nonlinear Functional $\theta_0=\Gamma(h_0)$

cor[Lower Bound Rate for $\theta_0=\Gamma(h_0)$ when $h_0$ is a NPIV] Let $\Gamma$ be a potentially nonlinear functional satisfying Assumption (ref). Let Assumptions (ref), (ref) (for $D\Gamma(h_0)[\cdot]$) and (ref) hold for the NPIV model. We have \begin{align*} \liminf_{n \to \infty} \inf_{\tilde{\theta}_n} \sup_{h \in \Lambda_c^s(\mathcal{X})} \mathbb{E}_{P_{h}} \left[r^{-2}_{NPIV,n} (\tilde{\theta}_n - \Gamma(h))^2 \right] \gtrsim C>0 \end{align*} where $r_{NPIV,n} $ is the same minimax lower bound rate given in Theorem (ref).

We note that Assumption (ref)(ii) and Assumption (ref) (for $D\Gamma(h_0)[\cdot]$) together implies local identification of $\Gamma(h_0)$ only. We could impose a global identification assumption of $\Gamma(h_0)$ as follows:

assumption[Nonlinear analogue of Assumption (ref): global reduced-form identification] For any $h_1,h_2\in\Lambda^s_c$, if $T h_1 = T h_2$ in $L^2(P_Z)$ then $\Gamma(h_1)=\Gamma(h_2)$.

Assumption (ref) states that $\Gamma:\Lambda^s_c \mapsto \mathbb{R}$ is constant on each observational equivalence class $[h]=h+\ker(T):=\{h'\in L^2(\mathcal{X})\cap \Lambda^s_c: Th'=Th \text{ in } L^2(\mathcal{Z})\}$, so $\Gamma(h)$ is point-identified from the reduced-form regression model. This Assumption (ref) is automatically satisfied if the conditional expectation operator $T$ is bounded complete. It is interesting to note that, for the minimax lower bound rate calculation, it suffices to have the local identification of $\Gamma(h_0)$ in a sup-norm neighborhood of $h_0\in\Lambda^s_c$.

Rate-Optimal Estimation

In this section we present simple sieve estimators that achieve the minimax lower bound rate of $r^*_{n}=n^{-\frac{s}{2s+d-m}}$ for $\theta_0 =L(h_0),~\Gamma(h_0),~V(h_0)$ when $h_0 \in\Lambda^s$ is a nonparametric regression function.\footnote{We could also present additional estimators to achieve the minimax lower bound rates for $\theta_0$ when $h_0$ is a nonparametric density of $X$ or a nonparametric instrumental variables regression. We skip those due to the lack of space.} We impose the following regularity condition on a nonparametric regression model in this and the next sections.

assumption[Bounded Variance of Errors] $\sup_{x\in{\cal X}}\mathbb{E}\left[\rest{\epsilon_{i}^{2}}X_{i}=x\right]<\infty$.

For any $h\in\Lambda^s$ we let $\left\| h \right\|_{\infty}:=\sup_{x\in\mathcal{X}} |h(x)|$ denote the sup-norm and $\left\| h \right\|_2:=\sqrt{\mathbb{E}[h(X)^2]}$ denote the $L_2(X)$-norm. Here $L_2(X):=L_2 (P_X )$ denote the space of square integrable functions, which is a Hilbert space under the $L_2(X)$-norm with inner product $\inprod{h_1,h_2}_2\equiv\mathbb{E}\left[h_1\left(X_{i}\right)h_2\left(X_{i}\right)\right]$.

Minimax Rate-Optimal Estimation of $L(h_0)$

We show that the lower bound $r^*_{n}=n^{-\frac{s}{2s+d-m}}$ established in Theorem (ref) for $\theta_0=L(h_0)$ can be attained by the following plug-in sieve estimator:

equation[equation omitted — 157 chars of source]

Hence we conclude that the rate $r^*_n$ is minimax rate optimal, and plug-in estimator with sieve first stage attains this optimal rate.\footnote{In Appendix (ref), we show that this optimal rate is also attained by a plug-in kernel estimator.}

For the sake of concreteness, let ${\mathbb{H}}_{K}$ denote the closed linear span (clsp) in $L_2(X)$-norm of a sieve basis functions $b^{K}(x):=\left(b_{1}^{K}(x),...,b_{K}^{K}(x)\right)^{'}$. By default we use a tensor product basis (see, e.g., chen2007sieve): $b^{K}(x)=\text{vec}\left(\otimes_{\ell=1}^{d}\left(b_{1}^{K}\left(x_{\ell}\right),...,b_{J_{n}}^{K}\left(x_{\ell}\right)\right)\right)$ so $K=J^d$. The sieve dimension is chosen to increase with $n$, with $J\equiv J_{n}\nearrow\infty$ and thus $K \equiv K_{n}\nearrow\infty$ as $n\to\infty$, though we will suppress the subscript $n$ in $J$ and $K$ for notational simplicity. Let $G:=\mathbb{E}\left[b^{K}\left(X_{i}\right)b^{K}\left(X_{i}\right)^{'}\right]$ denote the (population) Gram matrix with its minimum eigenvalue $\lambda_{min} (G)>0$. Let $\ol b^{K}\left(x\right):=G^{-1/2} b^K(x)$ denote the orthonormalized vector of basis functions of $b^K$, then ${\mathbb{H}}_{K}=clsp\left (b^K\right )=clsp \left(\ol b^K \right)$. The sieve LS estimator $\hat{h}$ for $h_0=\mathbb{E}[Y|X]$ is given by

equation[equation omitted — 293 chars of source]

where $\hat{G}:=\frac{1}{n}\sum_{i=1}^{n} b^{K_{n}}\left(X_{i}\right) b^{K_{n}}\left(X_{i}\right)^{'}$ and ${\ol G}_n:=\frac{1}{n}\sum_{i=1}^{n} \ol b^{K_{n}}\left(X_{i}\right) \ol b^{K_{n}}\left(X_{i}\right)^{'}$ denote the empirical Gram matrices. Let \[ h_{0,n}:=\arg\min_{h\in \mathbb{H}_{K_{n}}}\left\| h-h_0 \right\|_2 = \ol b^{K}\left(x\right)^{'}\mathbb{E}[\ol b^{K}\left(X_i\right)h_0(X_i)]= \ol b^{K}\left(x\right)^{'}\mathbb{E}[\ol b^{K}\left(X_i\right)Y_i] \] be the population LS projection of $h_0$ onto the sieve space ${\mathbb{H}}_{K}=clsp \left(\ol b^K \right)$. Following \citet*{chen2014sieve}, we define a sieve Riesz representer $v_{K_n}^{*}\in{\mathbb{H}}_{K}$ as follows:\footnote{Since the minimax lower bound rate $r^*_n$ for estimating $\theta_0=L(h_0)$ is slower than $n^{-1/2}$, there does not exist a Riesz representer for $\theta_0=L(h_0)$ in $L^2(X)$. Following \citet*{chen2014sieve}, a sieve Riesz representer for a linear functional is always well-defined and has a closed-form expression in any finite dimensional linear sieve space ${\cal V}_{K_n}$. The sieve Riesz representer is not unique as it depends on the choice of the linear sieve space, and is used for the variance characterization for $L(\hat{h})$ as well as the sieve influence function. Since we are presenting a linear sieve plug-in estimator $L(\hat{h})$ for $\theta_0=L(h_0)$ in this Section, we simply use the same linear sieve space ${\cal V}_{K_n}={\mathbb{H}}_{K_n}$ for our sieve variance characterization.}

equation[equation omitted — 261 chars of source]

Moreover, it can be solved in closed form as:

align[align omitted — 195 chars of source]

where $L\left[\ol b^{K_{n}}\right]:=\left(L\left[\ol b_1^{K_{n}}\right],..., L\left[\ol b_{K_n}^{K_{n}}\right]\right)^{'}$. Define the sieve variance as

equation[equation omitted — 207 chars of source]

The following lemma establishes the rates at which $\left\| v_{K_{n}}^{*} \right\|_{2}$ and $\left\| v_{K_{n}}^{*} \right\|_{sd}$ grow to infinity. Its proof utilizes the decomposition (ref) in Appendix (ref), which converts the $m$-dimensional Hausdorff integral to a finite sum of Lebesgue integrals on $\mathbb{R}^{m}$.

lem[Growth Rate of Sieve Riesz Representer for $L(h_0)$] Under Assumptions (ref), (ref), (ref) and $\lambda_{min} (G) >0$ for each $K$, we have: \begin{equation} \left\| v_{K_{n}}^{*} \right\|_{2}^{2}=\left\| L\left[\ol b^{K_{n}}\right] \right\|^{2}\equiv\sum_{k=1}^{K_{n}}\left(L\left[\ol b_{k}^{K_{n}}\right]\right)^2\asymp K_{n}^{\frac{d-m}{d}}=J_{n}^{d-m}. \end{equation} Further, under Assumption (ref) and (ref), we have: \[ \left\| v^*_{K_n} \right\|_{sd}^2\asymp\left\| v_{K_{n}}^{*} \right\|_{2}^{2}\asymp K_{n}^{\frac{d-m}{d}}. \]

Lemma (ref) is the core result that demonstrates the rate acceleration provided by the submanifold integral. It shows that the growth rates of the norm of the sieve Riesz representer and of the sieve variance are asymptotically proportional to $J^{(d-m)}$, where the exponent $(d-m)$ is the codimension of the $m$-dimensional manifold, but, importantly, not $d$, the dimension of the ambient space $\mathbb{R}^{d}$ in which the manifold is embedded. In other words, even though the dimensionality of the first-stage estimation of $h_{0}$ is $d$, it is reduced by the integration over the $m$-dimensional submanifold.

For any $h\in \Lambda^s(\mathcal{X})$, we let $P_{K,n}h$ be the empirical LS projection of $h$ onto the sieve space ${\mathbb{H}}_{K}$, which is given by $P_{K,n}h:= b^{K}\left(x\right)^{\prime}\hat{G}^{-1}\frac{1}{n}\sum_{i=1}^{n} b^{K}\left(X_{i}\right)h\left(X_{i}\right)$ with the sup operator norm defined by $\left\| P_{K,n} \right\|_{\infty}:=\sup_{h:0<\norm h_{\infty}<\infty}\frac{\left\| P_{K,n}h \right\|_{\infty}}{\norm h_{\infty}}$. Let $\zeta_{K}:=\sup_{x\in{\cal X}}\left\| b^{K}\left(x\right) \right\|$. We next impose some mild conditions on the linear sieve $b^{K}$ that are satisfied by some commonly used sieve bases such as B-splines and wavelets; see, e.g., chen2015optimal.

assumption[Conditions on Linear Sieves] The linear sieve $b^{K}$ satisfies (i) $\lambda_{min} (G)\geq const. >0$ uniformly in $K$, (ii) $\zeta_{K}=O\left(\sqrt{K}\right)$, (iii) $K\log K/n=o(1)$, (iv) $\left\| P_{K,n} \right\|_{\infty}=O_{p}\left(1\right)$, and (v) $\inf_{h\in {\mathbb{H}}_{K}}\left\| h-h_{0} \right\|_{\infty}\leq K^{-s/d}$ for $h_{0}\in\Lambda^{s}\left(\mathcal{X}\right)$.
thm[Convergence Rate under Sieve First Stage] Let Assumptions (ref), (ref), (ref), (ref), (ref) and (ref)(i)(ii)(iii) hold. Then: (1) \[ \frac{\sqrt{n}L\left(\hat{h}-P_{K_n,n}h_0\right)}{\left\| v_{K_{n}}^{*} \right\|_{sd}} = \frac{1}{\sqrt{n}}\sum_{i=1}^n \frac{v_{K_n}^*(X_i)}{\left\| v_{K_{n}}^{*} \right\|_{sd}}\epsilon_i +o_p(1)=O_p (1),~~~\text{with}~~~\left\| v^*_{K_n} \right\|_{sd}\asymp \sqrt{K_{n}^{\frac{d-m}{d}}}. \] (2) In addition, if Assumption (ref)(iv)(v) holds, then: \[ \hat{\theta}-\theta_{0}\equiv L\left(\hat{h}-h_{0}\right)=O_{p}\left(\sqrt{K_{n}^{\frac{d-m}{d}}/n}+K_{n}^{-s/d}\right). \] Further, if $s>m/2$, then by setting $K_{n}^{*}\asymp n^{\frac{d}{2s+d-m}}$, we obtain \[ \hat{\theta}-\theta_{0}=O_{p}\left(n^{-\frac{s}{2s+d-m}}\right)=O_{p}\left(r^*_{n}\right), \] which attains the lower bound rate $r^*_{n}$ in Theorem (ref) and is thus minimax rate-optimal.

Theorems (ref) and (ref) together demonstrate that the problem of estimating linear integrals over an $m$-dimensional submanifold is akin to the $(d-m)$-dimensional nonparametric regression problem in terms of convergence rates. Intuitively, the integration over the $m$-dimensional submanifold effectively “aggregates out” those $m$ dimensions, leaving only $d-m$ effective dimensions in the nonparametric estimation problem.

In many economic applications as discussed in Section (ref), the level set function $g$ is scalar-valued. Correspondingly the level set submanifold is of dimension $m=d-1$ and co-dimension $d-m = 1$. Here all but one dimensions are “aggregated out” and the minimax optimal convergence rate becomes $L(h_0)$ is $n^{-s/(2s+1)}$, the same rate as a 1-dimensional nonparametric regression problem. This illustrates the (potentially) significant “dimension reduction” achieved through integration.

remIn the proof of Theorem (ref) on the lower bound $r^*_n$ fo $L(h_0)$, it suffices to assume that the density $p_0$ of $X\in \mathcal{X}\subset \mathbb{R}^d$ is known, and also that level set function is smooth $g\in \Lambda^{s_1}({\cal X})$ for some $s_1 \geq \max\{1,s\}$ so that Assumption (ref) is satisfied. Theorem (ref) provides a feasible sieve estimator for $L(h_0)$ that attains the minimax lower bound rate of $r^*_n$ when $p_0$ is unknown and without imposing the extra smoothness condition $g\in \Lambda^{s_1}({\cal X})$. In Appendix Subsection (ref) we present an oracle estimator that uses the extra information of a known density $p_0$ and the extra smooth submanifold condition $g\in \Lambda^{s_1}({\cal X})$. It is interesting that the oracle estimator attains the same lower bound rate of $r^*_n$ as our feasible sieve estimator presented here. This means that in terms of minimax optimal rate of estimating $L(h_0)$, the extra information of knowing/using the marginal density of $X$ as well as the smoothness of submanifold function $g$ does not improve the speed of convergence. This is consistent with the fact that the optimal rate $r^*_n=n^{-s/(2s+d-m)}$ only depends on the smoothness of $h_0$ and the codimension $c:=d-m$ of the submanifold $\mathcal{M}$.

Minimax Rate-Optimal Estimation of $\Gamma(h_0)$

In general, whether the lower bound rate $r_n^*=n^{-s/(2s + d - m)}$ can be attained by any consistent estimator for the nonlinear integral functional $\Gamma(h_0)$ depends on the specific functional form of $\Gamma$. Here we consider a specific but still quite general class of submanifold integrals of nonlinear transformations of $h_0$. We verify Assumption (ref) in this context, and provide lower-level conditions for the attainability of the lower bound rate $r_n^*$, which then implies minimax optimality of the rate $r_n^*$ for the estimation of $\Gamma(h_0)$. Specifically, consider $\Gamma$ defined by

equation[equation omitted — 150 chars of source]

where $\phi\left(t,x\right)$ is a known nonlinear transformation that is $L$-Lipschitz in $t$ (so that its partial derivative $\phi_{1}\left(t,x\right):=\frac{\partial}{\partial t}\phi\left(t,x\right)$ is well-defined almost everywhere and uniformly bounded by the Lipschitz constant $L$ whenever well-defined). We assume that $\phi_{11}\left(t,x\right):=\frac{\partial^2}{\partial t^2}\phi\left(t,x\right)$ is well-defined.

We now provide an umbrella theorem for $\Gamma(h)$, which relates it to the linear submanifold integral we studied earlier, and establishes the attainability of the linear minimax optimal rate $r_n^*$ under appropriate smoothness requirement.

assumption[Regularity Conditions for Submanifold Integral of Nonlinear Transformations] (a) $\phi(\cdot,x)$ is Lipschitz uniformly in $x$, (b) $\phi_{1}\left(t,x\right)$ is Lipschitz, and (c) Assumption (ref) holds.
lem[Pathwise Derivative of $\Gamma$] Let $\theta_0=\Gamma(h_0)$ be defined in (ref) with a known $m$-dimensional submanifold ${\cal M}$ satisfying Assumption (ref). If Assumption (ref)(a) holds, then Assumption (ref) holds with \begin{equation} D\Gamma(h_0)[v] = \int_{\mathcal{M}} v(x)\, \phi_1(h_0(x),x)\, w(x)\, d\mathcal{H}^m(x). \end{equation}
thm[Minimax Rate Optimal Estimation of $\Gamma(h_0)$] Let $\theta_0=\Gamma(h_0)$ be defined in (ref) with a known $m$-dimensional submanifold ${\cal M}$. Let Assumptions (ref), (ref), (ref), (ref), (ref) and (ref) hold. \begin{itemize} • If $s > m/2$, then $\hat{\theta}_{SS}-\theta_0=O_p(r_n^*)$, where $\hat{\theta}_{SS}$ is a Split-Sample estimator: \begin{align} \hat{\theta}_{SS} := \Gamma\left(\bar{h}\right) - \frac{1}{8} \int_{\mathcal{M}} \phi_{11}\left(\bar{h}(x),x\right) \left(\hat{h}_1(x) - \hat{h}_2(x)\right)^2 w(x) d\mathcal{H}^m(x), \end{align} where $\hat{h}_1$ and $\hat{h}_2$ are split-sample sieve estimators with sieve dimension $K_{n}^{*}\asymp n^{\frac{d}{2s+d-m}}$, and $\bar{h} := (\hat{h}_1 + \hat{h}_2)/2$. • If $s > m/2$, then $\hat{\theta}_{LOO}-\theta_0=O_p(r_n^*)$, where $\hat{\theta}_{LOO}$ is a Leave-One-Out estimator: \begin{equation} \hat\theta_{\mathrm{LOO}} := \Gamma(\hat h) -\frac{1}{2n^2}\sum_{i=1}^n\int_{\mathcal{M}} \phi_{11}\!\left(\hat h(x),x\right)\,\hat s_i(x)^2\,\big(\hat\epsilon_i^{(-i)}\big)^2\, w(x)\,d\mathcal{H}^m(x), \end{equation} where $\hat h$ is the sieve LS estimator (ref) with sieve dimension $K_{n}^{*}\asymp n^{\frac{d}{2s+d-m}}$; for $i$-th observation, $\hat s_i(x):= b^{K_n}(x)^\prime \hat G^{-1} b^{K_n}(X_i)$, and $\hat\epsilon_i^{(-i)} := Y_i - \hat h^{(-i)}(X_i)$ is the LOO residual with the LOO sieve LS estimator $\hat h^{(-i)}$ computed using the sample $\{(Y_j,X_j):j\neq i\}_{j=1}^n$. • If $s\geq m$, then $\hat{\theta}-\theta_0=O_p(r_n^*)$, where $\hat{\theta}=\Gamma(\hat{h})$ is a plug-in sieve LS estimator with sieve dimension $K_{n}^{*}\asymp n^{\frac{d}{2s+d-m}}$. \end{itemize} In either (1) or (2) or (3), the rate $r_n^*=n^{-s/(2s+d-m)}$ is minimax-optimal.

Both the SS and the LOO sieve debiased estimators are constructed by subtracting the diagonal parts of the quadratic remainder in the Taylor expansion of $\Gamma(\hat h)$ at $h_0$. Both estimators are shown to attain the minimax lower bound rate $r_n^*$. We note that the split-sample debiased sieve estimator $\hat{\theta}_{SS}$ for $\Gamma(h_0)$ is used to attain an upper bound rate that matches the lower bound rate $r_n^*$ under the mild smoothness condition $s > m/2$. The plug-in sieve estimator $\hat{\theta}=\Gamma\left(\hat{h}\right)$ is extremely easy to compute, but attains the minimax optimal rate only under stronger smoothness condition $s\geq m$.

Theorem (ref) applies to quadratic submanifold functional $Q$ as a special case, where

align[align omitted — 83 chars of source]

and the split-sample debiased estimator for $Q(h_0)$ becomes $\hat{\theta}^Q_{SS}=\int_{{\cal M}}\hat{h}_{1}\left(x\right)\hat{h}_{2}\left(x\right)w\left(x\right)d{\cal H}^{m}\left(x\right)$. The study of $Q(h_0)$ is of conceptual importance given that it represents the simplest type of nonlinear integral functionals: when the standard Taylor expansion of a nonlinear function applies, the quadratic term in the expansion often constitutes the first-order approximation of the form of nonlinearity. For this reason, (full-dimensional) quadratic functionals have long been studied in statistics and econometrics as a leading example of nonlinear integral functionals: for a few examples out of many, see bickel1973some, bickel1988estimating, fan1991estimation, efromovich1996optimal, cai2005nonquadratic for earlier work in statistics on quadratic functionals of nonparametric density and regression functions, as well as chen2018optimal and breunig2019simple for more recent work on quadratic functionals of NPIV structural functions in econometrics. It has also been noted that the estimation of a quadratic functional is tightly related to hypothesis testing; see, for example, cai2005nonquadratic,cai2006adaptive in the statistics literature, and breunig2024adaptive in the adaptive minimax testing for NPIV models.

Minimax Rate-Optimal Estimation of $V(h_0)$

We now turn to another specific type of nonlinear integral functional, defined as an integral on the upper contour set of $h_0$, i.e.

equation[equation omitted — 136 chars of source]

where $h_{0}:\mathbb{R}^{d}\to\mathbb{R}$ is an unknown (scalar-valued) function and $w\left(x\right)$ is known.

Even though here $V$ is a full-dimensional (Lebesgue) integral per se, its pathwise derivative w.r.t. $h$ becomes a submanifold integral over the level set $\{x\in\mathcal{X}:~h_0(x)=0\}$ as we show below. Hence, in this case the “level set function” $g=h_{0}$ is taken to be unknown and needs to be nonparametrically estimated. Notice that such “upper contour integrals” show up in several examples in Section (ref) that feature subpopulations defined through inequalities, such as the welfare/value under optimal treatment assignment. We provide a corresponding theorem for upper contour integrals below.

assumption[Regularity Conditions for $V(h_0)$] (a) $h_{0}$ is continuously differentiable and $\left\| \nabla_{x}h_{0}\left(x\right) \right\|\geq\ul{\epsilon}>0$ on the level set ${\mathcal{M}}_0:=\left\{ x\in{\cal X}:h_{0}\left(x\right)=0\right\}$, (b) $w$, $\nabla_{x} w$, and the second-order (partial) derivatives of $h_{0}$ are all uniformly bounded.
lem[Pathwise Derivative of $V$] Let $\theta_0=V(h_0)$ be defined in (ref). If Assumption (ref)(a) holds, then ${\cal M}_0$ is an $m=(d-1)$ dimensional submanifold satisfying Assumption (ref), and Assumption (ref) holds with \begin{equation} DV(h_0)\left[v\right]=\int_{\mathcal{M}_0}v\left(x\right)\frac{w\left(x\right)}{\left\| \nabla_{x}h_{0}\left(x\right) \right\|}d{\cal H}^{d-1}\left(x\right). \end{equation}

We note that the calculation of the pathwise derivative (ref) is nonstandard, since it involves differentiation w.r.t. the changing boundary of the region of integration.

thm[Minimax Rate Optimal Estimation of $V(h_0)$] Let $\theta_0=V(h_0)$ be defined in (ref). Let Assumptions (ref), (ref), (ref), (ref) and (ref) hold. \begin{itemize} • If $s\geq m/2 +1=(d+1)/2$, then $\hat{\theta}_{SS}-\theta_0=O_p(r_n^*)$, where $\hat{\theta}_{SS}$ is a split-sample debiased estimator \begin{align} \hat{\theta}_{SS} := V\left(\bar{h}\right) - \frac{1}{8} D^2V(\bar{h})\left[\hat{h}_1 - \hat{h}_2, \hat{h}_1 - \hat{h}_2\right], \end{align} where $D^2V$ is the second functional derivative defined in Lemma (ref) in the Appendix, $\hat{h}_1$ and $\hat{h}_2$ are split-sample sieve estimators with sieve dimension $K_{n}^{*}\asymp n^{\frac{d}{2s+1}}$, and $\bar{h} := (\hat{h}_1 + \hat{h}_2)/2$. • If $s\geq m/2 +1=(d+1)/2$, then $\hat{\theta}_{LOO}-\theta_0=O_p(r_n^*)$, where $\hat{\theta}_{LOO}$ is a Leave-One-Out debiased estimator: \begin{equation} \hat\theta_{\mathrm{LOO}} := V(\hat h) -\frac{1}{2n^2}\sum_{i=1}^n D^2V(\hat h)\big[\hat s_i,\hat s_i\big]\cdot \big(\hat\epsilon_i^{(-i)}\big)^2, \end{equation} where $\hat{h}$, $\hat s_i$ and $\hat\epsilon_i^{(-i)}$ are defined in Theorem (ref)(2) with sieve dim $K_{n}^{*}\asymp n^{\frac{d}{2s+1}}$. • If $s\geq m+1=d$, then $\hat{\theta}-\theta_0=O_p(r_n^*)$, where $\hat{\theta}=V(\hat{h})$ is a plug-in sieve estimator with sieve dimension $K_{n}^{*}\asymp n^{\frac{d}{2s+1}}$. \end{itemize} In either (1) or (2) or (3), the rate $r_n^*=n^{-s/(2s+d-m)}=n^{-s/(2s+1)}$ is minimax optimal.

Sieve Inference on Submanifold Integrals

The minimax lower bound results from Section (ref) imply that our functionals of interest, $\theta_0=L(h_0),~\Gamma(h_0),~V(h_0)$ are all irregular functionals of the unknown function $h_0$, where $h_0$ could be a nonparametric regression, a nonparametric density and a NPIV function. Consequently, the well-known $L^2(P_X)$-Riesz representers and the efficient influence functions do not exist for these functionals. Nevertheless, applying the general theories in \citet*{chen2014sieve,chen2014sieveM} for sieve M-estimation and in chenpouzo2015sieve for penalized sieve MD-estimation,\footnote{M-estimation includes nonparametric regression, quantile regression and density estimation as special cases. MD-estimation includes M-estimation, nonparametric IV regression and nonparametric quantile IV as special cases.} the sieve Riesz representers and sieve influence functions are still well-defined for these irregular functionals. The sieve Riesz representer approach enables us to show that the sieve student t statistics for these irregular functionals are still asymptotically standard normal, which can be used to construct simple confidence intervals for $\theta_0$.

For the sake of simplicity, in this section we specialize their general theory to our functionals of interest when $h_0$ is a nonparametric regression. In this section we denote $\theta_0:=\Phi(h):\Lambda^s\mapsto\mathbb{R}$, where $\Phi(h_0)=L(h_0),~\Gamma(h_0),~V(h_0)$. By exploring the submanifold structures, we can characterize the growth rate of the norm of the sieve Riesz representer for the first-order expansions of $\Gamma(h_0)$ and $V(h_0)$, provide a tighter control of nonlinear remainder terms, and establish the asymptotic normality under mild sufficient conditions.

Asymptotic Normality via Sieve Riesz Representers

Let $\mathbb{H}_{K_{n}}:=\left\{ h\in{\cal B}_{K_{n}}:\left\| h-h_{0} \right\|_{\infty}<\ul{\epsilon}\right\} $. Let $\hat{h}_n \in \mathbb{H}_{K_n}$ be any sup-norm consistent sieve estimator of the regression function $h_0\in\Lambda^s(\mathcal{X})$. Let ${\cal V}_{K_n}$ denote a $K_n$-dimensional Hilbert space that contains $h_{0,n}:=\arg\min_{h\in \mathbb{H}_{K_{n}}}\left\| h-h_0 \right\|_2$, and can be equivalently expressed as $clsp\left(\ol \psi^{K_n}\right):=clsp\{\psi_k:k=1,...,K_n\}$, where $\{\psi_k:k\geq 1\}$ is an orthonormal basis for $L_2(P_X)$. Following \citet*{chen2014sieve,chen2014sieveM}, we define the sieve Riesz representer $v_{K_n}^{*}$ and the sieve variance $\left\| v^*_{K_n} \right\|_{sd}^2$ in the same way as in (ref) and (ref).

The following Remark follows from Lemma (ref), which allows us to control the growth rate of the sieve Riesz representer and the sieve variance for the nonlinear submanifold functionals $\Phi=\Gamma$ and $V$.

rem[Growth Rate of Sieve Riesz Representer for $D\Phi(h_0)$] Under Assumptions (ref), (ref), (ref) and $\lambda_{min} \left( R_{K}\right) >0$ for each $K$, we have: \begin{equation} \left\| v_{K_{n}}^{*} \right\|_{2}^{2}=\left\| D\Phi(h_0)\left(\ol \psi^{K_{n}}\right) \right\|^{2}\asymp K_{n}^{\frac{d-m}{d}}=J_{n}^{d-m}. \end{equation} Further, if Assumptions (ref) and (ref) hold, then: $\left\| v_{K_{n}}^{*} \right\|_{sd}^2\asymp \left\| v_{K_{n}}^{*} \right\|_{2}^{2} \asymp K_{n}^{\frac{d-m}{d}}$.

The following is a standard Lindeberg condition for CLT.

assumption[Lindeberg] $\sup_{x\in{\cal X}}\mathbb{E}\left[\rest{\epsilon_{i}^{2}\left\{ \left|\epsilon_{i}\right|>c\right\} }X_{i}=x\right]\to0$ as $c\to\infty$.

In the following, we let $\hat{\theta}=\Phi(\hat{h})$ be the plug-in sieve estimator of $\theta_0=\Phi(h_0)$. We let $\hat{\theta}_{SS}$ and $\hat{\theta}_{LOO}$ respectively denote the Split-Sample sieve estimator and Leave-One-Out sieve estimator of $\theta_0=\Phi (h_0)$ when $\Phi=\Gamma,~V$.

prop[Asymptotic Normality] Let $\hat{h}$ be the sieve LS estimator of $h_0$. Let Assumptions (ref), (ref), (ref), (ref) and (ref) hold. Let the sieve dimension $K_{n}$ satisfy $K_{n}\asymp K^*_n [\log(n)]^{\kappa}$ for some $\kappa \in (0,1)$ and $K_{n}^{*}\asymp n^{\frac{d}{2s+d-m}}$. Then: \[ \frac{\sqrt{n}\left(\tilde{\Phi}-\Phi(h_0)\right)}{\left\| v_{K_{n}}^{*} \right\|_{sd}} = \frac{1}{\sqrt{n}}\sum_{i=1}^n \frac{v_{K_n}^*(X_i)}{\left\| v_{K_{n}}^{*} \right\|_{sd}}\epsilon_i +o_p(1)\overset{d}{\longrightarrow}\mathcal{N}\left(0,1\right)\quad\text{with }\left\| v_{K_{n}}^{*} \right\|_{sd}\asymp\sqrt{K_{n}^{\frac{d-m}{d}}}, \] provided the following conditions hold for the specific functional forms of $\Phi:\Lambda^s \mapsto \mathbb{R}$: \begin{itemize} • For $\Phi=L$: Assumptions (ref) and (ref) hold with $s>m/2$ for $\tilde{\Phi}=L(\hat{h})$. • For $\Phi=\Gamma$: Assumptions (ref) and (ref) hold, $s>m/2$ for $\tilde{\Phi}=\hat{\theta}_{SS},~\hat{\theta}_{LOO}$; and $s>m$ for $\tilde{\Phi}=\Gamma(\hat{h})$. • For $\Phi =V$: Assumption (ref) holds, $s>(d+1)/2$ for $\tilde{\Phi}=\hat{\theta}_{SS},~\hat{\theta}_{LOO}$; and $s>d$ for $\tilde{\Phi}=V(\hat{h})$. \end{itemize}

Inference based on Estimated Sieve Riesz Representers

We note that $v_{K_n}^{*}\left(\cdot\right)$ can be estimated by \[ \hat{v}_{K_{n}}^{*}\left(\cdot\right):=\ol \psi^{K_{n}}\left(\cdot\right)'(\ol {R}_{K_{n}})^{-1}D\Phi\left(\hat{h}\right)\left[\ol \psi^{K_{n}}\right],~~~\text{with}~~\ol {R}_{K_n} :=\frac{1}{n} \sum_i{\ol \psi}^{K_{n}}(X_i){\ol \psi}^{K_{n}}(X_i)'. \] Denote $\sigma^2_{*,K_{n}} :=\left\| v_{K_{n}}^{*} \right\|_{sd}^2$ and $\widehat{\epsilon}_i:=Y_{i}-\hat{h}\left(X_{i}\right)$. We can estimate the sieve variance as\footnote{In small samples,it might be more accurate to estimate the sieve variance using $\widehat{\sigma}^2_{*,K_{n}}=\frac{1}{n}\sum_{i=1}^n \left[ \hat{v}_{K_{n}}^{*}(X_{i})\widehat{\epsilon}_i - \frac{1}{n-1}\sum_{j\neq i,j=1}^n\hat{v}_{K_{n}}^{*}(X_{j})\widehat{\epsilon}_j\right]^{2}$.}

align[align omitted — 322 chars of source]

with $\widehat{\Omega}_{K_n}:= \ol {R}_{K_n}^{-1}\left(\frac{1}{n}\sum_{i=1}^{n}(\widehat{\epsilon}_i)^{2}{\ol \psi}^{K_{n}}(X_i){\ol \psi}^{K_{n}}(X_i)'\right) \ol {R}_{K_n}^{-1}$. The consistency of the sieve variance estimator $\widehat{\sigma}^2_{*,K_{n}}$ can be easily established under the following assumption.

assumption[Higher Error Moments] $\mathbb{E}\left[\left| \epsilon_i \right|^{2+\delta}\right] < \infty$ for some $\delta>2d/(2s-m)$.
prop[Sieve t statistics] Let Assumptions for Proposition (ref) and Assumption (ref) hold. Let $K_{n}\asymp K^*_n [\log(n)]^{\kappa}$ for some $\kappa \in (0,1)$ and $K_{n}^{*}\asymp n^{\frac{d}{2s+d-m}}$. Then: \[ \left| \frac{\widehat{\sigma}_{*,K_{n}}}{\left\| v_{K_{n}}^{*} \right\|_{sd} } - 1 \right| =o_p(1);~~~ \frac{\sqrt{n}\left(\tilde{\Phi}-\theta_0\right)}{\widehat{\sigma}_{*,K_{n}}} = \frac{1}{\sqrt{n}}\sum_{i=1}^n \frac{v_{K_n}^*(X_i)}{\left\| v_{K_{n}}^{*} \right\|_{sd}}\epsilon_i +o_p(1)\overset{d}{\longrightarrow}\mathcal{N}\left(0,1\right), \] where $\tilde{\Phi}$ is the estimator specified in Proposition (ref) (i.e., $L(\hat{h})$, $\hat{\theta}_{SS}$, $\hat{\theta}_{LOO}$, or $\Phi(\hat{h})$, depending on the case and smoothness regime).

As an application of Proposition (ref), we can construct valid $.95$ confidence interval using normal critical value for $\theta_0=\Phi(h_0)$ as follows: \[ \text{CI}:=\left[\tilde{\Phi}-1.96\hat{\sigma}_{\theta}~,~\tilde{\Phi}+1.96\hat{\sigma}_{\theta}\right]~,~~\text{with}~~\hat{\sigma}_{\theta}^{2}=\widehat{\sigma}^2_{*,K_{n}}/n \] provided that the sieve variance is computed using undersmoothed sieve dimension $K_n$, say $K_{n}\asymp K^*_n [\log(n)]^{1/2}$. In practice, we implement an adaptive inference procedure that (i) selects the sieve dimension in a data-driven way (bootstrap--Lepski), with an undersmoothed choice used for valid confidence interval; and (ii) uses a multiplier bootstrap to approximate the distribution of the studentized sieve estimator. See Appendix (ref) for details.

Monte Carlo Simulations

We present two Monte Carlo studies to check finite sample performance of the sieve estimation and inference of linear and nonlinear submanifold integrals. Subsection (ref) considers a linear submanifold integral $\theta_0=L(h_0)$, and Subsection (ref) considers an upper contour integral $\theta_0=V(h_0)$. In both examples, $h_0 (x)=\mathbb{E}[Y|X=x]$ is a nonparametric regression, and is estimated using a tensor-product B-spline sieve LS estimator $\hat{h}$ based on a random sample $\{(Y_i,X_i)\}_{i=1}^n$, drawn from the joint distribution of $(Y,X)$ satisfying: \[ Y_i=h_0(X_i) + \epsilon_{i},~~~\text{with}~~\epsilon_{i}\sim\mathcal{N}\left(0,1\right),~~X_{i1}\sim_{i.i.d.}X_{i2}\sim_{i.i.d.}Uniform\left[-2,2\right], \] and $\epsilon_{i}$ is independent of the covariates $X_i=(X_{i1},X_{i2})'\in \mathcal{X}=[-2,2]^2$ with $d=2$.

In both Monte Carlo studies, we implement multiplier bootstrap Lepski data-driven choice of sieve dimensions $K$ as described in Appendix (ref), with $B_{Lepski}=500$ multiplier bootstrap draws. We report simulation results over $B=1000$ Monte Carlo simulation replications, using performance metrics in terms of Monte Carlo root-mean-squared error (RMSE), bias and standard deviation (SD), as well as coverage rates and coverage length. We report the simulation results for different sample sizes $n=500, 1000, 2000, 4000, 8000$.

Integral on Known Submanifold

We first report simulations results for the estimation of an integral functional of $h_0$ over a known submanifold defined by the unit circle, i.e., $\theta_{0}=L\left(h_{0}\right)$ with $L\left(h\right)=\int_{\mathbb{S}^{1}}h\left(x\right)d{\cal H}^{1}\left(x\right)$ and $h_{0}\left(x\right)=x_{1}^{2}+2\sin\left(x_{1}\right)x_{2}$ with $x=(x_1,x_2)^{\prime} \in{\cal X}=\left[-2,2\right]^{2}$. Under the transformation $x_{1}=\cos\left(\beta\right)$ and $x_{2}=\sin\left(\beta\right)$, we obtain the true value \[ \theta_{0}=L\left(h_{0}\right)=\int_{0}^{2\pi}\left(\text{cos}\left(\beta\right)^{2}+2\sin\left(\cos\left(\beta\right)\right)\sin\left(\beta\right)\right)d\beta=\pi. \] In this exercise, $L(h_0)$ is a linear functional of $h_0$ on a known $m=1$ -dimensional submanifold (i.e., the unit circle $\mathbb{S}^{1}$). Theorem (ref) is directly applicable, and the optimal convergence rate of a spline sieve plug-in estimator $L(\hat{h}) - L(h_0)$ is given by $n^{-s/(2s+1)}$ with the optimal sieve dimension $K^* \asymp n^{2/(2s+1)}$ and $d=2,m=1$.

Given the spline LS estimator $\hat{h}$, the plug-in estimator $\hat{\theta}=L(\hat{h})=\int_{\mathbb{S}^{1}} \hat{h}\left(x\right)d{\cal H}^{1}\left(x\right)$ is numerically computed using the sample average over $M=5000$ Sobol sequence points\footnote{The Sobol sequence sampling, proposed by sobol1967distribution, is a well-known quasi-random Monte Carlo sampling method that generates a deterministic sequence of points, whose distribution asymptotically converges to the uniform distribution, but achieves better finite-sample approximation of the population expectation (integral) by the sample average over $M$ Sobol points.} in the angle space $\left[0,2\pi\right]$. We obtain the plug-in estimator $\hat{\theta}$ and construct the confidence interval $\text{CI}:=\left[\hat{\theta}\pm1.96\hat{\sigma}_{\theta}\right]$ with $\hat{\sigma}_{\theta}^{2}=\widehat{\sigma}^2_{*,K_{n}}/n$ as in (ref), where the pathwise derivative $D\Phi\left(\hat{h}\right)\left[v\right]=L(v)=\int_{\mathbb{S}^{1}}v\left(x\right)d{\cal H}^{1}\left(x\right)$ is computed numerically in the same manner as described above.

Table (ref) reports the RMSE, bias, standard deviation, and the average selected sieve dimension $\bar{\hat{K}}$ for the plug-in estimator using the Lepski sieve dimension choice $\hat{K}$. Note that the rate-optimal $\hat{K}$ is chosen to minimize estimation risk, and does not involve the undersmoothing that is needed for valid confidence interval construction.

table[table omitted — 711 chars of source]

Next, we turn to the inference results. Table (ref) reports the finite-sample performance of the plug-in estimator and the corresponding 95% confidence interval, constructed with bootstrap-Lepski adaptive choice of the undersmoothed sieve dimension $\tilde{K}$ and normal critical values. The average $\bar{\tilde{K}} \approx 51$ is substantially larger than the rate-optimal $\bar{\hat{K}} \approx 17$ in Table (ref), reflecting the undersmoothing needed for valid inference.

table[table omitted — 914 chars of source]

In Table (ref) and subsequent tables, we report the square root of the mean squared error (RMSE), the bias (Bias), the standard deviation (SD)\footnote{The standard deviation is calculated using the standard $\frac{1}{B-1} \sum_{b=1}^B$ formula. For this technical reason, “SD” can be larger than “RMSE”, which is calculated based on the $\frac{1}{B} \sum_{b=1}^B$ formula.} of the estimator, the average lower and upper bounds of the confidence interval (CI_L and CI_U), the average length of the confidence interval (U-L), and the realized coverage probability of the CI (Coverage). Overall, the plug-in estimator and the corresponding CI perform very well under all five sample sizes: the RMSE shrinks (almost at $\sqrt{n}$ rate) as the sample size increases, the bias is of negligible order relative to the standard deviation, and the realized coverage probability is close to the nominal 95% level.

Table (ref) reports the estimator RMSE and the CI coverage depending on: (i) whether the bootstrap-Lepski adaptive choice of sieve dimension or a fixed sieve dimension is used, and (ii) whether normal or bootstrap quantile critical values\footnote{The bootstrap quantiles are also computed based on $B_{Lepski} = 500$ draws.} are used, when constructing the CI. The first column “Adapt+Normal” is identical to that in Table (ref), while the remaining three columns contain results about the three alternative CI constructions. We observe that all four versions produce coverage close to the 95% nominal level.

table[table omitted — 1,058 chars of source]

Integral on Estimated Upper Contour Set

In this subsection, we analyze the estimation of an integral functional of a nonparametric function over the (estimated) unit disk: $\theta_{0}:=V\left(h_{0}\right)$ with

align*[align* omitted — 262 chars of source]

Under the above construction, $h_{0}\left(x\right)\geq0$ if and only if $\norm x\leq1$, and thus $\theta_{0}=\int\mathbb 1\left\{ \norm x\leq1\right\} dx=\pi$. We compute the spline nonparametric regression estimator $\hat{h}\left(x\right)$ for $h_0(x)= \mathbb{E}[Y_i|X_i=x]$, and the plug-in estimator $V\left(\hat{h}\right)=\int\mathbb 1\left\{ \hat{h}\left(x\right)\geq0\right\} dx$ using the sample average over $M=50,000$ Sobol sequence points, and $\text{CI}:=\left[\hat{\theta}\pm1.96\hat{\sigma}_{\theta}\right]$ with the sieve variance $\hat{\sigma}_{\theta}^{2}=\widehat{\sigma}^2_{*,K_{n}}/n$ as in (ref), where the pathwise derivative $DV\left(\hat{h}\right)\left[v\right]$ now takes the following form: \[ DV(h)\left[v\right]=\int_{\left\{x\in{\cal X}:~ h\left(x\right)=0\right\} }\frac{v\left(x\right)}{\left\| \nabla_{x}h\left(x\right) \right\|}d{\cal H}^{d-1}\left(x\right), \] which we numerically approximate via $\hat{D}V\left(h\right)\left[v\right]=\frac{1}{2\epsilon}\int_{\left\{x\in{\cal X}:~ -\epsilon<h\left(x\right)<\epsilon\right\} }v\left(x\right)dx$ based on the mathematical result\footnote{See Theorem 3.13.(iii) of \citet*{evans2015measure}.} that \[ \lim_{\epsilon\searrow0}\frac{1}{2\epsilon}\int_{\left\{x\in{\cal X}:~ -\epsilon<h\left(x\right)<\epsilon\right\} }v\left(x\right)dx=\int_{\left\{x\in{\cal X}:~ h\left(x\right)=0\right\} }\frac{v\left(x\right)}{\left\| \nabla_{x}h\left(x\right) \right\|}d{\cal H}^{d-1}\left(x\right). \] We set $\epsilon=0.01$ in our simulation, and, given that $\left\{x\in{\cal X}:~ -\epsilon<h\left(x\right)<\epsilon\right\} $ may occur infrequently for small $\epsilon$, we use the sample average from $M=50,000$ Sobol sequence points to approximate $\hat{D}V\left(h\right)\left[v\right]$.

Table (ref) reports the rate-optimal estimation results for the upper contour set integral. The average rate-optimal $\bar{\hat{K}}$ grows with $n$ (from approximately 17 at $n=500$ to 36 at $n=8000$), in contrast to the stable $\bar{\hat{K}} \approx 17$ in Simulation 1. The plug-in and LOO-debiased estimators display similar RMSE at the rate-optimal $\hat{K}$ for $n\geq 1000$.

table[table omitted — 1,029 chars of source]
table[table omitted — 1,505 chars of source]

Next, we turn to the inference results. Table (ref) reports the finite-sample performance of the plug-in estimator as well as the Leave-One-Out (LOO) estimator, along with the CIs constructed using bootstrap-Lepski adaptive choice of sieve dimension and normal critical values. The results are again based on $B= 1000$ Monte Carlo replications and (within each replication) $B_{Lepski} = 500$ bootstrap draws when applicable. The results in Table (ref) display a similar pattern as those in Table (ref): Again, the RMSE shrinks (almost at $\sqrt{n}$ rate) as the sample size increases, the bias is of smaller order relative to the standard deviation, and the realized CI coverage is close to the nominal 95% level. One noticeable difference between the two tables, however, lies in that the RMSEs in Table (ref) are substantially smaller than those in Table (ref). Heuristically, this might have been due to the fact that the upper contour set integral being estimated in Table (ref) is itself a full-dimensional integral: even though its asymptotic behavior is theoretically driven by the lower-dimensional submanifold, in finite sample the full-dimensional nature of the integral may have made it overall easier to estimate. Lastly, for $n\geq 1000$, the LOO-debiased estimator does result in (slightly but noticeably) smaller biases relative to the plug-in estimator, yielding comparable or better CI coverage as anticipated.\footnote{At $n=500$ the LOO estimator exhibits a large RMSE outlier, likely reflecting instability the leave-out procedure at this small sample size.}

Table (ref) again reports the RMSE and coverage comparison under the two sieve dimension choice methods and the two critical value choice methods. We observe that all four versions of CI constructions produce coverage close to the 95% nominal level, and the improvement of the LOO estimator over the plug-in one in terms of size control is consistent across all four versions of CI construction procedures.

table[table omitted — 1,543 chars of source]