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.
33,322 characters · 7 sections · 59 citation commands
Simple Adaptive Estimation of Quadratic Functionals in Nonparametric IV Models
{ Keywords: nonparametric instrumental variables, ill-posed inverse problem with an unknown operator, quadratic functional, minimax estimation, leave-one-out, adaptation, Lepski's method.}
\onehalfspacing
Long before the recent popularity of instrumental variables in modern machine learning causal inference, reinforcement learning and biostatistics, the instrumental variables technique has been widely used in economics. For instance, instrumental variables regressions are frequently used to account for omitted variables, mis-measured regressors, endogeneity in simultaneous equations and other complex situations in economic observational data. In economics and other social sciences, as well as in medical research, it is very difficult to estimate causal effects when treatment assignment is not randomized. Instrumental variables are commonly used to provide exogenous variation that is associated with the treatment status, but not with the outcome variable (beyond its direct effect on the treatments).
To avoid mis-specification of parametric functional forms, nonparametric instrumental variables (NPIV) regressions have gained popularity in econometrics and modern causal inference in statistics and machine learning. The simplest NPIV model assumes that a random sample $\{(Y_i,X_i,W_i)\}_{i=1}^n$ is drawn from an unknown joint distribution of $(Y,X,W)$ satisfying
where $h_0$ is an unknown continuous function, $X$ is a $d$-dimensional vector of continuous endogenous regressors in the sense that $\operatorname{\mathbb{E}}[U|X]\neq 0$, $W$ is a vector of conditioning variables (instrumental variables) such that $\operatorname{\mathbb{E}}[U|W] = 0$. The structural function $h_0$ can be identified as a solution to an integral equation of first kind with an unknown operator:
where the conditional density $f_{X|W}$ (and hence the conditional expectation operator $T$) is unknown. Under mild conditions, the conditional density $f_{X|W}$ is continuous and the operator $T$ smoothes out “low regular” (or wiggly) parts of $h_0$. This makes the nonparametric estimation (recovery) of $h_0$ a difficult ill-posed inverse problem with an unknown smoothing operator $T$. See, for example, NP03, HH, \citet*{CFR06handbook}, \citet*{BCK07econometrica}, ChenReiss2011 and DFFR. For a given smoothness of $h_0$, the difficulty of recovering $h_0$ depends on the smoothing property of the conditional expectation operator $T$. The literature distinguishes between the mildly and severely ill-posed regimes, and the optimal convergence rates for nonparametrically estimating $h_0$ are different in the two regimes.
This paper considers adaptive, minimax rate-optimal estimation of a quadratic functional of $h_0$ in the NPIV model ((ref)):
for a known positive, continuous weighting function $\mu$, which is assumed to be uniformly bounded below from zero and from above on some subset of of the support of $X$. Let $\widehat{h}$ be a sieve NPIV estimator of the NPIV function $h_0$ (see e.g., BCK07econometrica). chen2013 and ChenChristensen2017 considered inference on a slightly more general nonlinear functional $g(h_0)$ using plug-in sieve NPIV estimator $g(\widehat{h})$. However, there is no result on any adaptive, minimax rate-optimal estimation of any nonlinear functional $g(h_0)$ of the NPIV function $h_0$ yet. Since a quadratic functional is a leading example of a smooth nonlinear functional in $h_0$, ChenChristensen2017 established the minimax lower bound for estimating a quadratic functional $f(h_0)$ in a NPIV model. They also point out that a plug-in sieve NPIV estimator $f(\widehat{h})$ of the quadratic functional $f(h_0)$ can achieve the lower bound in the severely ill-posed regime, but fails to achieve the lower bound in the mildly ill-posed regime. Moreover, none of the existing work considers adaptive minimax rate-optimal estimation of the quadratic functional $f(h_0)$ in a NPIV model.
In this paper, we first propose a simple leave-one-out sieve NPIV estimator $\widehat{f}_J$ for the quadratic functional $f(h_0)$, and establish an upper bound on its convergence rate. By choosing the sieve dimension $J$ optimally to balance the squared bias and the variance parts, we show that the resulting convergence rate of $\widehat{f}_J - f(h_0)$ coincides with the lower bound of ChenChristensen2017. In this sense the estimator $\widehat{f}_J$ is minimax rate-optimal for $f(h_0)$ regardless whether the NPIV model is severely ill-posed or mildly ill-posed. In particular, for the severely ill-posed case, the optimal convergence rate is of the order $(\log n)^{-\alpha}$, where $\alpha>0$ depends on the smoothness of the NPIV function $h_0$ and the degree of severe ill-posedness. For the mildly ill-posed case, the optimal convergence rate of $\widehat{f}_J - f(h_0)$ exhibits the so-called elbow phenomena: the rate is of the parametric order $n^{-1/2}$ for the regular mildly ill-posed case, and is of the order $n^{-\beta}$ for the irregular mildly ill-posed case, where $\beta\in (0,1/2)$ depends on the smoothness of $h_0$, the dimension of $X$ and the degree of mild ill-posedness.
The minimax optimal estimation rate of $\widehat{f}_J - f(h_0)$ is achieved by the optimal choice of the sieve dimension $J$ (a key tuning parameter) that depends on the unknown smoothness of $h_0$ and the unknown degree of ill-posedness. We next propose a data driven choice $\widehat{J}$ of the sieve dimension based on a modified Lepski method.\footnote{See Lepski90, lepski1997 and \citet*{lepski1997optimal} for detailed descriptions of the original Lepski principle.} The modification is needed to account for the estimation of the unknown degree of ill-posedness. The adaptive, leave-one-out sieve NPIV estimator $\widehat{f}_{\widehat{J}}$ of $f(h_0)$ is shown to attain the minimax optimal rate in the severely ill-posed case and in the regular mildly ill-posed case, but up to a multiplicative $\sqrt{\log n}$ in the irregular mildly ill-posed case. We note that even for adaptive estimation of a quadratic functional of a direct regression in a Gaussian white noise model, efromovich1996 already shown that the extra $\sqrt{\log n}$ factor is the necessary price to pay for adaptation to the unknown smoothness of the regression function.
Previously for the nonparametric estimation of $h_0$ in the NPIV model ((ref)), horowitz2014adaptive considers adaptive estimation of $h_0$ in $L^2$ norm using a model selection procedure. breunig2016 consider adaptive estimation of a linear functional of the NPIV function $h_0$ in a root-mean squared error metric using a combined model selection and Lepski method. These papers obtain adaptive rate of convergence up to a multiplicative factor of $\sqrt{\log(n)}$ (of the minimax optimal rate) in both severely ill-posed and mildly ill-posed cases. \citet*{chen2021} propose adaptive estimation of $h_0$ in $L^\infty$ norm using a modified Lepski method and tight random matrix inequalities to account for the estimated measure of ill-posedness. They show that their data-driven procedure attains the minimax optimal rate in $L^\infty$ norm and is fully adaptive to the unknown smoothness of $h_0$ in both severely and mildly ill-posed regimes. Our data-driven choice of the sieve dimension is closest to that of chen2021, which might explain why we also obtain minimax optimal adaptivity for the quadratic functional $f(h_0)$ in both severely and mildly ill-posed regimes.
While horowitz2014adaptive, breunig2016 and chen2021 use plug-in sieve NPIV estimators in their adaptive estimation of a linear functional of $h_0$, we use a leave-one-out sieve NPIV estimator $\widehat{f}_J$ for the quadratic functional $f(h_0)=\int h_0^2(x) \mu(x) dx$. Recently BC2020 propose a test statistic that is based on a standardized leave-one-out estimator of a quadratic distance for a null hypothesis of $\operatorname{\mathbb{E}}[(h_0(X)-h^R(X))^2\mu(X)]=0$ in a NPIV model (for some parametric, semiparametric or shape restricted $h^R$). They construct an adaptive minimax test using a random exponential scan procedure. We use the unstandardized leave-one-out estimator $\widehat{f}_J$ in our modified Lepski procedure for adaptive minimax estimation of $f(h_0)$ in a NPIV model. It is well-known that adaptive minimax testing and adaptive minimax estimation are related but different (see, e.g., nicklbook). In particular, while both papers apply a tight Bernstein-type inequality for U-statistics (\citet*{houdre2003}) in the proofs, the adaptive optimal rates are different. For instance, the adaptive minimax $L^2$ separation rate of testing in BC2020 is always slower than $n^{-1/2}$, while our adaptive minimax estimation for $f(h_0)$ can achieve the parametric rate of $n^{-1/2}$ for regular mildly ill-posed NPIV models.
Minimax rate-optimal estimation of a quadratic functional in density and direct regression (in Gaussian white noise) settings has a long history in statistics. See, for example, bickel1988estimating, donoho1990minimax, Fan91, efromovich1996, laurent2000, CL06, gine2008, \citet*{collier2017} and the references therein. To the best of our knowledge, there are not many published papers on minimax estimation of a quadratic functional in difficult inverse problems. See butucea2007, butucea2011, che201 and kroll2019rate for deconvolutions and inverse regressions in Gaussian sequence models. Moreover, che201 seems the only published work on adaptive estimation of a quadratic functional in a special deconvolution (with a known operator). Our paper is the first to propose a simple estimator that is adaptive minimax rate-optimal for a quadratic functional in a NPIV model, and also contributes to inverse problems with unknown operators.
The rest of the paper is organized as follows. Section (ref) presents the leave-one-out sieve NPIV estimator of the quadratic functional $f(h_0)$, and derives its optimal convergence rates. Section (ref) first presents a simple data-driven procedure of choosing the sieve dimension using a modified Lepski method. It then establishes the optimal convergence rates of our adaptive estimator of the quadratic functional. Section (ref) provides a brief conclusion and discusses several extensions. All proofs can be found in the Appendices (ref)--(ref).
This section consists of three parts. The first subsection introduces model preliminaries and notation. Subsection (ref) introduces a simple leave-one-out, sieve NPIV estimator of the quadratic functional $f(h_0)$. Subsection (ref) establishes the convergence rate of the proposed estimator, and shows that the convergence rate coincides with the lower bound and hence is optimal.
We first introduce notation that is used throughout the paper. For any random vector $V$ with support $\mathcal V$, we let $L^2(V)=\{\phi:\mathcal V \to \mathbb{R}, \|\phi\|_{L^2(V)}<\infty\}$ with the norm $\|\phi\|_{L^2(V)}=\sqrt{\operatorname{\mathbb{E}}[\phi^2(V)]}$. If $\{a_n\}$ and $\{b_n\}$ are sequences of positive numbers, we use the notation $a_n \lesssim b_n$ if $\limsup_{n\to\infty}a_n/b_n<\infty$ and $a_n\sim b_n$ if $a_n\lesssim b_n$ and $b_n\lesssim a_n$.
We consider a known positive, continuous weighting function $\mu$, which is assumed to be uniformly bounded below from zero and from above on some subset of $\mathcal X$, denoted by $X_\mu$. Denote $L_\mu^2 =\{h:\mathcal X_\mu \to \mathbb{R},\|h\|_\mu<\infty\}$ with the norm $\|h\|_\mu =\sqrt{\int h^2(x)\mu(x)dx}$. We consider basis functions $\{\psi_j\}_{j\geq 1}$ to approximate the NPIV function $h_0$. Its orthonormalized analog with respect to $\|\cdot\|_\mu$ is denoted by $\{\widetilde \psi_j\}_{j\geq 1}$. We assume that the structural function $h_0$ belongs to the Sobolev ellipsoid
Let $T : L^2(X) \mapsto L^2(W)$ denote the conditional expectation operator given by $(T h)(w) = \operatorname{\mathbb{E}}[h(X)|W = w]$. Finally let ${\left\lbrace \psi_1,...,\psi_J\right\rbrace }$ and ${\left\lbrace b_1,...,b_K\right\rbrace }$ be collections of sieve basis functions of dimension $J$ and $K$ for approximating functions in $L^2(X)$ and $L^2(W)$, respectively. We define the sieve measure of ill-posedness which, roughly speaking, measures how much the conditional expectation operator $T$ smoothes out $h$. Following BCK07econometrica the sieve $L_\mu^2$ measure of ill-posedness is
where $\Psi_J = \text{clsp}\{\psi_1,...,\psi_J\} \subset L^2(X)$ denotes the sieve spaces for the endogenous variables. We call a NPIV model (ref)
Let $\{(Y_i,X_i,W_i)\}_{i=1}^n$ denote a random sample from the NPIV model (ref). The sieve NPIV (or series 2SLS) estimator $\widehat h$ of $h_0$ can be written in matrix form as follows (see, e.g., ChenChristensen2017) \[ \widehat h(\cdot) = \psi^J(\cdot)'[\Psi'P_B\Psi]^- \Psi'P_B{\textbf{Y}} = \psi^J(\cdot)'\widehat A B'{\textbf{Y}}/n \] where $P_B= B(B'B)^-B'$ and ${\textbf{Y}} = (Y_1,\ldots,Y_n)'$,
and $\widehat A=n[\Psi'P_B\Psi]^- \Psi' B(B'B)^-$ is an estimator of $A=[S'G_b^{-1}S]^{-1}S'G_b^{-1}$, with $S=\operatorname{\mathbb{E}}[b^K(W_i)\psi^J(X_i)']$ and $G_b=\operatorname{\mathbb{E}}[b^K(W_i)b^K(W_i)']$.
As pointed out by ChenChristensen2017, although one could estimate $f(h_0)$ by the plug-in sieve NPIV estimator $f(\widehat h)$, it fails to achieve the minimax lower bound. We propose a leave-one-out sieve NPIV estimator for the quadratic functional $f(h_0)$ as follows:
where $G_\mu=\int \psi^J(x)\psi^J(x)' \mu(x)dx$. We will show that this simple leave-one-out estimator $\widehat{f_J}$ can achieve the lower bound for estimating $f(h_0)$.
Based on many simulation results in BCK07econometrica and ChenChristensen2017, the crucial regularization parameter in sieve NPIV estimation of $h_0$ is the dimension $J$ of the sieve space used to approximate unknown function $h_0$. In this paper, we simply let $K(J)=c_KJ$ for some constant $c_K\geq 1$. Further, we let $\zeta_{\psi,J}=\sup_x\|G_\mu ^{-1/2}\psi^J(x)\|$ and $\zeta_{b,K}=\sup_w\|G_b^{-1/2}b^K(w)\|$. For instance, $\zeta_{\psi,J} = O(\sqrt J)$ and $\zeta_{b,K} = O( \sqrt K)$ for (tensor-product) polynomial spline, wavelet and cosine bases. Denote $\zeta_J=\max(\zeta_{\psi,J}, \zeta_{b,K})$ for $K=K(J)$. In the rest of the paper we restrict sieve bases to the ones such that $\zeta_J = O(\sqrt J)$.
We first introduce assumptions that are used to derive our rate of convergence of the estimator $\widehat{f_J}$. We denote the sieve Least Squares (LS) projection of $h$ onto $\Psi_J=\text{clsp} \{\psi_1,...,\psi_J\}$ as $\Pi_J h(x)=\psi^J(x)'G_\mu^{-1} \langle \psi^J, h\rangle_\mu$. For $h_0\in \mathcal H_2(p, L)$ we have $\|h_0-\Pi_J h_0\|_{\mu} \leq L J^{-p/d}$ which is used throughout this paper. This implies that $\sqrt{J(\log J)} \|h_0-\Pi_Jh_0\|_\mu=o(1)$ as $J$ goes to infinity (since $p>d/2$).
Below we let $\Pi_K g (w)=b^K(w)'G_b^{-1} \operatorname{\mathbb{E}}[b^K (W) g(W)]$ denote the sieve LS projection of $g\in L^2 (W)$ onto $B_K = \text{clsp} \{b_1,...,b_K\}$.
For a $r\times c$ matrix $M$ with $r \leq c$ and full row rank $r$ we let $M_l^-$ denote its left pseudoinverse, namely $(M'M)^-M'$ where $'$ denotes transpose and $^-$ denotes generalized inverse. Below, $\|\cdot\|$ respectively denotes the vector $\ell_2$ norm when applied to a vector and the operator norm $\|A\|:=\sup_{x:\|x\|=1}\|Ax\|$ when applied to a matrix $A$. Let $(s_1,\dots,s_J)$ denote the singular values, in non-increasing order, of $G_b^{-1/2}SG_\mu^{-1/2}$. In particular $s_J=s_{\min}(G_b^{-1/2}S G_\mu^{-1/2})$.
Discussion of Assumptions: Assumption (ref)(i) ensures identification of the nonlinear functional $f(h_0)$. Assumption (ref) restricts the growth of the sieve dimension $J$. Assumption (ref)(i) is a mild condition on the approximation properties of the basis used for the instrument space and is first imposed in chen2021. In fact, $\|(\Pi_K T-T)h\|_{L^2(W)}= 0$ for all $h \in \Psi_J$ when the basis functions for $B_K$ (with $K\geq J$) and $\Psi_J$ form either a Riesz basis or an eigenfunction basis for the conditional expectation operator. Assumption (ref)(ii) is the usual $L^2$ “stability condition” imposed in the NPIV literature (cf. Assumption 6 in BCK07econometrica). Note that Assumption (ref)(ii) is also automatically satisfied by Riesz bases. Assumption (ref) is a modification of the sieve measure of ill-posedness and was used by EfromovichKoltchinskii2001. Assumption (ref) is also related to the extended link condition in breunig2016 to establish optimal upper bounds in the context of minimax optimal estimation of linear functionals in NPIV models. Finally we note that by definition, $s_J$ satisfies
for all $K=K(J)\geq J>0$. Assumption (ref)(i) further implies that
for some constant $c_\tau>0$. We shall maintain Assumption (ref)(i) and use the equivalence of $s_J$ and $\tau_J^{-1}$ in the paper.
The next result provides an upper bound on the rate of convergence for the estimator $\widehat{f}_J$.
Theorem (ref) presents an upper bound on the convergence rates of $\widehat f_J$ to $f(h_0)$. When the sieve dimension $J$ is chosen optimally, the convergence rate (ref) coincides with the minimax lower bound in ChenChristensen2017 for the mildly ill-posed case, while the convergence rate (ref) coincides with the minimax lower bound in ChenChristensen2017 for the severely ill-posed case. Moreover, within the mildly ill-posed case, depending on the smoothness of $h_0$ relatively to the dimension of $X$ and the degree of mildly ill-posedness $a$, either the first or the second variance term in (ref) dominates, which leads to the so-called elbow phenomenon: the regular case with a parametric rate of $n^{-1/2}$ when $p> a+d/4$; and the irregular case with a nonparametric rate when $p\leq a+d/4$. In particular, Theorem (ref) shows that the simple leave-one-out estimator $\widehat f_J$ is minimax rate optimal provided that the sieve dimension $J$ is chosen optimally.
ChenChristensen2017 actually established lower bound for estimating a quadratic functional of a derivative of $h_0$ in a NPIV model as well. Using Fourier, spline and wavelet bases, we can easily show that our simple leave-one-out, sieve NPIV estimator of the quadratic functional of a derivative of $h_0$ also achieve the lower bound, and hence is minimax rate-optimal. We do not present such a result here since it is a very minor extension of Theorem (ref).
The minimax rate of convergence depends on the optimal choice of sieve dimension $J$, which depends on the unknown smoothness $p$ of the true NPIV function $h_0$ and the unknown degree of ill-posedness. In this section we propose a data-driven choice of the sieve dimension $J$ based on a modified Lepski method; see Lepski90, lepski1997 and lepski1997optimal for early development of this popular method.
In this section we follow chen2021 and let $\Psi_J$ be a tensor-product Cohen-Daubechies-Vial (CDV) wavelet (see, e.g., chapter 4.3.5 of nicklbook) or dyadic B-spline sieve (see, e.g., Appendix A.1 of chen2021) for $\mathcal H_2(p,L)$. Let $\mathcal T$ denote the set of possible sieve dimensions $J$. For example for (order $r$) B-splines, $\mathcal T =\{J=(2^l + r -1)^d:l\in \mathbb{N}\cup\{0\}\}$. Since $\widehat f_J$ is based on a sieve NPIV estimator, we can simply use a random index set $\widehat{\mathcal I}$ that is proposed in chen2021 for their sup-norm rate adaptive sieve NPIV estimation of $h_0$:
where
$\widehat s_J$ is the smallest singular value of $(B'B/n)^{-1/2}(B'\Psi/n)G_\mu^{-1/2}$, and $J^+=\min\{j\in\mathcal T:\,j > J\}$.
We define our data driven choice $\widehat J$ of “optimal” sieve dimension for estimating $f(h_0)$ as follows:
for some constant $c_0>0$ and
where $a\vee b:=\max\{a,b\}$. The random index set $\widehat{\mathcal I}$ is used to compute our data driven choice (ref) since the unknown measure of ill-posedness $\tau_J$ is estimated by $\widehat s_J^{-1}$.
We introduce a non-random index set $\mathcal I=\{J\in\mathcal T:\, J\leq \overline J\}$, where $\overline {J}=\sup \left\{J\in\mathcal T: \tau_{J} J \sqrt{(\log J) / n} \leq \bar{c} \right\}$ for some sufficiently large constant $\bar{c}>0$. Let $\mathcal B=\{h\in L_\mu^2: \|h\|_\infty\leq L\}$ and $\overline p>\underline p\geq 3d/4$. The following assumption strengthens some conditions imposed in the previous section.
The next result establishes an upper bound for the adaptive estimator $\widehat f_{\widehat J}$.
Theorem (ref) shows that our data-driven choice of the key sieve dimension can lead to fully adaptive rate-optimal estimation of $f(h_0)$ for both the severely ill-posed case and the regular mildly ill-posed case, while it has to pay a price of an extra $\sqrt{\log n}$ factor for the irregular mildly ill-posed case (i.e., when $p\leq a+d/4$). We note that when $a =0$ in the mildly ill-posed case, the NPIV model ((ref)) becomes the regression model with $X=W$. Thus our result is in agreement with the theory in efromovich1996, which showed that one must pay a factor of $\sqrt{\log n}$ penalty in adaptive estimation of a quadratic functional in a Gaussian white noise model when $p \leq d/4$.
In adaptive estimation of a nonparametric regression function $\operatorname{\mathbb{E}}[Y|X=\cdot ]=h(\cdot)$, it is known that Lepski method has the tendency of choosing small sieve dimension, and hence may not perform well in empirical work. We wish to point out that due to the ill-posedness of the NPIV model ((ref)), the optimal sieve dimension for estimating $f(h_0)$ is smaller than the optimal sieve dimension for estimating $f(\operatorname{\mathbb{E}}[Y|X=\cdot])$. Therefore, we suspect that our simple adaptive estimator of a quadratic functional of a NPIV function will perform well in finite samples.
In this paper we first show that a simple leave-one-out sieve NPIV estimator of the quadratic functional $f(h_0)$ is minimax rate optimal. We then propose an adaptive leave-one-out sieve NPIV estimator of the $f(h_0)$ based on a modified Lepski method to account for the unknown degree of ill-posedness. We show that the adaptive estimator achieves the minimax optimal rate for the severely ill-posed case and for the regular mildly ill-posed case, while a multiplicative $\sqrt{\log n}$ term is the price to pay for the irregular mildly ill-posed NPIV problem.
Like all existing work using Lepski method, implementation of our data-driven choice relies on a calibration constant. To improve finite sample performance over the original Lepski method, spokoiny2009parameter suggest a propagation approach, \citet*{CCK} and spokoiny2019bootstrap propose bootstrap calibrations in kernel density estimation and in linear regressions with Gaussian errors respectively. chen2021 propose a bootstrap implementation of a modified Lepski method in their minimax adaptive sup-norm estimation in a NPIV model, and show its good performance in finite samples. Their bootstrap implementation can be easily extended to calibrate the constant in our adaptive estimation of the quadratic functional in a NPIV model. We leave this to future refinement.
Our results can be extended in several directions. First, we can relax the Sobolev ball assumption imposed on $h_0$ in the NPIV model. We can let the NPIV function $h_0$ belong to a bump algebra space. The result by collier2017 on minimax estimation of a quadratic functional under sparsity constraints can be useful for this extension. Second, we focus on adaptive estimation of a quadratic functional of the NPIV function $h_0$ in this paper. There are works on minimax-rate estimation and adaptive estimation for more general smooth nonlinear functionals of densities and of nonparametric regressions; see, e.g., birge1995, \citet*{TT2021} and the references therein. We can combine our approach here with those in the literature for extensions to other smooth nonlinear functionals of the NPIV function $h_0$. Such an extension will allow for adaptive minimax estimation of nonlinear policy functionals in economics and modern causal inference.