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.
54,040 characters · 19 sections · 79 citation commands
Adaptive estimation for some nonparametric instrumental variable models
\allowdisplaybreaks
\parindent0ptAbstract: The problem of endogeneity in statistics and econometrics is often handled by introducing instrumental variables (IV) which fulfill the mean independence assumption, i.e. the unobservable is mean independent of the instruments. When full independence of IV's and the unobservable is assumed, nonparametric IV regression models and nonparametric demand models lead to nonlinear integral equations with unknown integral kernels. We prove convergence rates for the mean integrated square error of the iteratively regularized Newton method applied to these problems. Compared to related results we derive stronger convergence results that rely on weaker nonlinearity restrictions. We demonstrate in numerical simulations for a nonparametric IV regression that the method produces better results than the standard model.
\vskip .2in {\sl MSC: AMS 2010 subject classification.} primary 62G08, secondary 62G20\\ {\sl Keywords and phrases:} Nonparametric regression, instrumental variables, nonlinear inverse problems, regularization
\onehalfspacing \parindent8pt
Dependence of an unobservable error term and covariates is a frequent problem in statistical and econometrical modeling known as endogeneity. An efficient way to deal with endogeneity is to use instrumental variables (IV) in the estimation. These are additional variables which can assumed to be independent or mean independent of the unobservable. In the context of nonparametric estimation the IV approach usually leads to ill-posed problems with an unknown operator that needs to be estimated. The solution $\varphi$ of the nonparametric IV problem can be characterized by a possibly nonlinear operator equation
In some regression models $\psi = 0$, in others $\psi$ is a function that has to be estimated from observations by some estimator $\widehat{\psi}$. The operator $\mathcal{F}: \mathbb{X} \to \mathbb{Y}$ is an integral operator between some Banach or Hilbert spaces $\mathbb{X}$ and $\mathbb{Y}$ which is unknown in applications. Only an estimator $\widehat{\mathcal{F}}$ is available. The inverse of the operators $\mathcal{F}$ or $\widehat{\mathcal{F}}$ is usually not continuous. Even with an arbitrarily small variance in $\widehat\psi$ and $\widehat\mathcal{F}$ we usually have $\mathop{\rm {\mathbb V}ar}\nolimits(\|\widehat{\mathcal{F}}^{-1}\widehat\psi\|_\mathbb{X}) = \infty$ and the straightforward estimator $\widehat\varphi = \widehat{\mathcal{F}}^{-1}\widehat\psi$ is typically inconsistent. We discuss specific examples for nonparametric IV models and the related operators together with the respective literature in Section (ref).
In this paper we describe and analyze a consistent estimator for this type of problem, when $\mathcal{F}$ is an operator between Hilbert spaces. The estimator is based on the iteratively regularized Gau{\ss}-Newton method (IRGNM) with iterated Tikhonov regularization defined below in (ref). Details about the method will be given in Section (ref).
This method was suggested by Baku:92. Important monographs on this topic are BK:04b and KNS:08. These contributions consider only problems with known operators and deterministic right hand side in equation (ref). The use of IRGNM for nonparametric IV problems was proposed and analyzed by DFHJM:14. They derived rates for convergence in probability with a priori parameter choice using variational methods.
The novelty of this paper is that we prove significantly faster convergence rates for the mean integrated squared error (MISE) rather than convergence in probability under a different set of assumptions. In addition, we propose adaptive estimation with Lepski\u\i's principle and prove rates for this case. Furthermore, we assume a significantly weaker nonlinearity condition for the operator $\mathcal{F}$ which has a clear interpretation and is reasonable for most applications while the nonlinearity condition in DFHJM:14 is difficult to interpret and to check. We also prove faster rates of convergence when the regression function is smooth enough. Our proofs do not use variational methods. Instead we rely on spectral methods as in BauHohMun:09. We also use a modification of Hoeffding's inequality from McD:89.
The paper is organized as follows. We discuss in Section (ref) some IV models which fit into the framework of this paper and explain the estimator. The estimator is introduced in Section (ref). Section (ref) contains convergence rate theorems. Finally, we present some numerical simulations in Section (ref). All proofs are in the Appendix.
In our general framework a function $\varphi^{\dagger}$ is characterized by the possibly nonlinear operator equation
i.e. $\varphi^{\dagger}$ is the true solution. Here $\mathcal{F} :B_{2R}(\varphi_0) \subseteq \mathbb{X} \rightarrow \mathbb{Y}$ is an operator between Hilbert spaces with norms $\|\cdot\|_\mathbb{X}$ and $\|\cdot\|_\mathbb{Y}$ respectively. Tyical examples for $\mathbb{X}$ and $\mathbb{Y}$ are $L^2$ and $L^2$ based Sobolev spaces $H^i$ for $i=1,2,\ldots$. A ball $B_{2R}(\varphi_0)$ with radius $2R$ around an initial guess $\varphi_0$ is contained in the domain of $\mathcal{F}$. In practice, large values of $R$ are possible. The operator equation is allowed to be ill-posed, i.e. $\mathcal{F}^{-1}$ may not be continuous. Furthermore, the operator $\mathcal{F}$ is not known in applications. Only a series of estimators $\widehat{\mathcal{F}}_n :B_{2R}(\varphi_0) \subseteq \mathbb{X} \rightarrow \mathbb{Y}$ are available where $n$ denotes the sample size. We assume that $\varphi^{\dagger}$ is a unique solution to (ref) in $B_{2R}$, i.e. the problem is locally identified. In the following we discuss econometric examples for this setup.
\paragraph{Mean independence} The simplest nonparametric IV regression model has a separable error term and a mean independence condition
Here and in all following models $Y$ and $U$ are univariate random variables, while $X$ and $Z$ can be multivariate and their dimensions do not have to coincide. The regressor $X$, the instrument $Z$ and the response $Y$ are observed, while the error term $U$ is unobservable.
This model was proposed by NewPow:03 and florens:03. It was further studied and applied in HalHor:05, BluCheKri:07, CheRei:11, DFFR:11 FJV:11, Hor:11, JVV:11, GS:12_reg, ChePou:12, Horowitz14, ChenChrist:15, BreJoh:15, ChenChrist:18, as well as Babii:20 among others. For an overview see Hor:14_survey.
We can write (ref) equivalently as $\mathbb{E}[\varphi(X)|Z] = \mathbb{E}[Y|Z]$ and if the conditional densities $f_{X|Z}$ and $f_{Y|Z}$ exist, as
We define the linear integral operator $(\mathcal{F}_{ce}\varphi)(z) := \int f_{X|Z}(x|z) \varphi(x) dx$ with integral kernel $f_{X|Z}(x|z)$ and the function $\psi(z) := \int y f_{Y|Z}(y|z)dy$. Model (ref) can be given in operator form $(\mathcal{F}_{ce}\varphi)(z) = \psi(z)$. The integral kernel $f_{X|Z}$ and thereby $\mathcal{F}_{ce}$ as well as the function $\psi$ are unknown and have to be estimated from a sample of $Y,X,Z$. An Density estimators $\widehat f_{X|Z}$ and $\widehat f_{Y|Z}$ give estimators $\widehat\mathcal{F}_{ce}$ and $\widehat \psi$ in a natural way. While the main focus of this paper is on nonlinear operator equations, we use model (ref) as a benchmark for the IRGNM applied to model (ref) below.
The model identifies the regression function $\varphi$ if and only if $\mathcal{F}_{ce}$ is injective. This property is called completeness, see Hault:11, HF:15, Andrews:17, and BF:20.
\paragraph{Full independence} In many applications the error term can be assumed to be independent of the instrument. Hence, mean independence of the instrument can be replace by full independence as proposed in DFHJM:14
Since the new assumptions $U \protect\mathpalette{\protect\independenT}{\perp} Z$ and $\mathbb{E}[U] = 0$ imply $\mathbb{E}[U|Z] = 0$ but not vice versa model (ref) makes stronger assumptions than model (ref). Consequently, whenever (ref) identifies the solution so does (ref). Furthermore, there are cases in which (ref) can identify a solution, while (ref) fails. This is for example the case with discrete instruments and continuous regressors as discussed in DFHJM:14, Torgo, HF:15, CFF:19, and Loh:19.
We can translate model (ref) into an operator equation by defining the operator
When $Y,X,Z$ have a joint density $f_{YXZ}$, taking the derivative with respect to $u$ yields the alternative operator
Model (ref) is equivalent to the operator equations $\widetilde\mathcal{F}_{ind}(\varphi) = 0$ or $\mathcal{F}_{ind}(\varphi) = 0$. Note that the operators are nonlinear due to the first line of (ref) or (ref). Furthermore, the operators are not known and have to be estimated. A density estimator $\widehat f_{YXZ}$ gives a straight forward estimator $\widehat \mathcal{F}_{ind}$.
For any $\varphi$ that sets the first line of the operator $\mathcal{F}_{ind}$ to $0$ also $c+ \varphi$ with $c\in\mathbb{R}$ sets it to $0$. In addition, for any $\varphi$, the second line of $\mathcal{F}_{ind}$ is set to $0$ by $\varphi -\mathbb{E}[Y-\varphi(X)]$. Hence, for any solution $\varphi$ of the first line of the operator we have $\mathcal{F}_{ind}\big(\varphi-\mathbb{E}[Y-\varphi(X)]\big)=0$. The nonlinear inverse problem is to find a $\varphi$ that solves the first line of $\mathcal{F}_{ind}$. The second line is a parametric problem that can be estimated with the parametric rate. When we discuss this example below we will only consider the first line of the operator as this is dominating the convergence rate.
Let us denote the integral kernel of the first line of the operator (ref) and its estimator by $k_{ind}(y,x,z) := f_{YXZ}(y,x,z) - f_{YX}(y,x) f_Z(z)$ and $\widehat k_{ind}(y,x,z) := \widehat f_{YXZ}(y,x,z) - \widehat f_{YX}(y,x) \widehat f_Z(z)$ respectively. Then the first component of the operator reads $(\mathcal{F}_{ind}(\varphi))(u,z) = \int k_{ind}(u+\varphi(x),x,z)dx$.
\paragraph{Nonparametric IV quantile regression} Another model that leads to a different nonlinear operator equation is nonparametric IV quantile regression proposed by HorLee:07. For $q \in [0,1]$ the $q$-th quantile regression function $\varphi_q$ is characterized by
If the joint density $f_{YXZ}$ exists, the model is equivalent to an operator equation $\mathcal{F}_q(\varphi_q)=0$ with
where $F_{YXZ}(y,x,z):=\int_{-\infty}^y f_{YXZ}(\tilde{y},x,z)\,d\tilde{y}$. Different estimation procedures for this model were proposed and analyzed in HorLee:07, ChePou:12, GS:12_quant DFHJM:14, and Breunig:15. Local identification properties of this and related models are discussed in CCLN:14.
We can write $(\mathcal{F}_q(\varphi))(z) = \int k_q(\varphi(x),x,z)dx$ with integral kernel \[ k_q(y,x,z) := F_{YXZ}(y,x,z) - q f_{XZ}(x,z). \] Replacing $q f_Z(z)$ by $\int q f_{XZ}(x,z)dx$ is impractical in applications but makes it easier to discuss properties of (ref) in this paper. $\mathcal{F}_q$ and $k_q$ are unknown and have to be estimated. If we plug-in a density estimator $\widehat f_{YXZ}$, we get straight forward estimators $\widehat k_q$ and $\widehat{\mathcal{F}}_q$.
\paragraph{Non-separable model} A related example that falls in our framework is nonparametric IV regression with unseparable error, wich was proposed in CheImbNew:07. See also CheHan:05. The model is
It was pointed out in HorLee:07, and CheImbNew:07 that this model is already contained in model (ref). Let $F_U$ be the cumulative distribution function of $U$. Normalize $\widetilde U:= F_U(U)$ and $\widetilde\phi(x, \tilde u) := \phi(x,F_U^{-1}(\tilde u))$. Then $\widetilde U$ is uniformly distributed on $[0,1]$. The value of $\widetilde U$ corresponds to a quantile in model (ref). This reduces (ref) to model (ref) with $\varphi_q(x) = \widetilde\phi(x,q)$.
\paragraph{Further examples} We briefly comment on further econometric models that fall into the framwork of this paper. A problem that has a similar mathematical structure as IV regression appears in some nonparametric demand models for differentiated products. It was considered with mean independence assumption similar to (ref) in BH:11, BH:14 and with full independence similar to (ref) in DHK:14. Some models for games of incomplete information lead to a nonlinear inverse problem with deterministic operator, see for example FS:10. Nonlinear inverse problems with deterministic operators also occur in functional linear quantile regression (without instrumental variables) as in Kato:12. The estimator in this paper can be applied to these type of problems. However, the error analysis would be different since there is no randomness in the operator. Also related are nonparametric ARCH($\infty$) models which can be treated as linear inverse problem, see LM:05. Further linear inverse problems in econometrics are discussed in CFR:07.
Remember that $\varphi^{\dagger}$ denotes the true solution and let $\varphi_0$ be an initial guess. Our method is based on linearizing $\mathcal{F}$ which motivates the following assumption.
The iteratively regularized Gau\ss-Newton method with iterated Tikhonov regularization consists of two nested iterations. The outer iteration is a Newton method. It starts at $\varphi_0$ and produces in the $j$-th step the estimate $\widehat{\varphi}_{j+1}$. In the $j$-th step the operator is linearized as $\widehat{\mathcal{F}}(\varphi) \approx \widehat{\mathcal{F}}_n'[\widehat{\varphi}_j](\varphi - \widehat{\varphi}_j) + \widehat{\mathcal{F}}_n(\widehat{\varphi}_j)$. A regular Newton method would invert the linear operator $\widehat{\mathcal{F}}_n'[\widehat{\varphi}_j]$ to compute the next step. Due to the ill-posedness, this would be unstable and we use a regularized inverse instead. The regularized inverse is computed by $m$-times iterated Tikhonov regularization which is the inner iteration of the method. In the following scheme the Newton iteration is indexed by $j$ and the Tikhonov iteration by $i$, and $\alpha_j > 0$ is a regularization parameter
As usual for Newton methods, convergence can fail if the initial guess $\widehat\varphi_0$ is too far from the true solution $\varphi^\dag$. In practice and in simulations the method proves to be quite robust to the choice of $\widehat\varphi_0$. If no a priori information about $\varphi^\dag$ is available, $\widehat{\varphi}_0 = 0$ is usually a good choice.
With a small $\alpha_j$ the method has a large variance due to the ill-posedness of $\mathcal{F}$. While a larger $\alpha_j$ controls the variance but adds some bias. We choose $\alpha_0$ large enough to stabilize the problem and let $\alpha_j$ decay in every Newton step by
to reduce the bias. A second parameter that has to be chosen is the number of inner iterations $m$. A large $m$ is of advantage for very smooth $\varphi^{\dagger}$. We will address the choice of $\alpha_0$ and $m$ in Section (ref) and Assumption (ref). The Newton iteration needs to be stopped at an appropriate iteration step. The size of the regularization parameter is linked to the number of steps. Hence, the number of steps corresponds to a bias variance trade-off. We will investigate parameter choice with a priori knowledge in Section (ref) and fully data driven in Section (ref).
We introduce the following notations for shorter formulas \[ T_\dag := \mathcal{F}'[\varphi^{\dagger}] \qquad \widehat T_{n,j} := \widehat{\mathcal{F}}_n'[\widehat{\varphi}_j] \qquad \widehat T_{n \dag} := \widehat{\mathcal{F}}'_n[\varphi^{\dagger}]. \] An alternative formulation of the method can be obtained by using the functional calculus. Let $\widehat T_{n,j}^*$ denote the adjoint operator of $\widehat T_{n,j}$ and set
Then (ref) is equivalent to
The convergence theory is presented in four steps. We start by introducing assumptions for the general operator equation (ref) as well as for the IV regression models (ref) and (ref). Afterwards, we state convergence rate result for the MISE with a priori choice of the stopping parameter $j$. Then we compare the result to HorLee:07. The last step is a theorem with data driven choice of $j$ by Lepski\u\i's principle.
As usual for nonparametric methods a smoothness assumption has to be imposed on the true solution $\varphi^{\dagger}$ to get convergence rates. In our setup with an ill-posed operator equation (ref) it is necessary to link the smoothness of $\varphi^{\dagger}$ to the smoothing properties of the operator $\mathcal{F}$. An efficient and popular way to formulate this is a source condition. The following definition uses the functional calculus.
When $T_\dag$ is a linear integral operator with kernel $\frac{\partial}{\partial y}k(\varphi^{\dagger}(x),x,z)$ as in Example (ref), this definition can be interpreted in the following way. We assume for simplicity that $T_\dag$ is compact which is for example the case if $\frac{\partial}{\partial y}k(\varphi^{\dagger}(x),x,z)$ is continuous. It was shown in reade841 and reade842 that the singular values of such an operator decay at least polynomially if $\frac{\partial}{\partial y}k(\varphi^{\dagger}(x),x,z)$ belongs to a Sobolev space, and exponentially if $\frac{\partial}{\partial y}k(\varphi^{\dagger}(x),x,z)$ is analytic.
Let $(\sigma_t,\, u_t,\, v_t)$ be a singular system for $T_\dag$. The source condition (ref) implies for $e_0 = \varphi_0 - \varphi^{\dagger}$ \[ \omega = \sum_{t\in \mathbb{N}}\frac{\langle e_0,\, v_t \rangle}{\Lambda(\sigma_t^2)} u_t \in \mathbb{X} \qquad \mbox{and thereby} \qquad \sum_{t=1}^\infty \left(\frac{\langle e_0,\, v_t \rangle}{\Lambda(\sigma_t^2)}\right)^2 < \infty. \] Hence, a $\omega$ fulfilling (ref) only exists if $\Lambda$ compensates the decay of the singular values in a way that $\langle e_0,\, v_t \rangle\Lambda(\sigma_t^2)^{-1}$ is square summable. The decay of singular values describes the smoothing properties of the $T_\dag$ with respect to the singular vectors. While the decay of $\langle e_0,\, v_t \rangle$ describes the smoothness of $e_0$ with respect to the singular vectors. Thus, the rate of decay for $\Lambda(x)$ when $x \searrow 0$ compares these two degrees of smoothness. For the examples above the source condition compares the smoothness of $f_{YXZ}$ with the smoothness of the regression function $\varphi^{\dagger}$.
When $\sigma_t$ and $\langle e_0,\, v_t \rangle$ both decay polynomially or both decay exponentially, i.e. $\sigma_t \lesssim \exp(-c_\sigma t)$ and $\langle e_0,\, v_t \rangle \lesssim \exp(-c_{e_0} t)$ with some constants $c_\sigma$ and $c_{e_0}$, the source condition is fulfilled with $\Lambda(x) = x^\mu$. Where $\mu >0$ is a sufficiently small constant. A source condition with polynomial $\Lambda$ is called a H\"older source condition, which is a concept that goes back to lavrentev62 and morozov68. For exponential decay of $\sigma_t$ but only polynomial decay of $\langle e_0,\, v_t \rangle$ the source condition holds when the operator is rescaled to $\|T_\dag\| < 1$ and $\Lambda(x) = (-\ln(x))^{-p}$ with some $0<p$. This choice of $\Lambda$ was proposed by mair94 and hohage:97 and is called logarithmic source condition.
Despite the word “condition” in the name “source condition” it is rather a relation that selects an index function. Corollary 2 in mh08 shows that for any compact injective operator $T_\dag$ and any $e_0$ exists an index function $\Lambda$ such that a source condition is fulfilled.
In this paper we focus on H\"older source conditions with $\mu > 1/2$. Notice that this implies $e_0 \in \text{Range}(\mathcal{F}'[\varphi^{\dagger}]^*)$. The case of $\mu \le 1/2$ and logarithmic source conditions was analyzed in DFHJM:14. We make the formal assumption:
Closely related to the smoothness of the true solution is the choice of the parameters $\alpha_0$ and $m$ for the IRGNM. In the following assumption $\|T\|_{\mathcal{L}(\mathbb{X},\mathbb{X})} := \sup_{\varphi} \{\|T\varphi\|_\mathbb{X}~| ~\|\varphi\|_\mathbb{X}=1\}$ denotes the usual operator norm for linear operators.
Both parameters need to be large enough but there is not much harm in choosing them larger than necessary. Any $\alpha_0$ and $m$ fulfilling Assumption (ref) will lead to comparable estimates. However, increasing $\alpha_0$ will lead to a few more Newton steps. The lower bound of $\alpha_0$ depends on the derivative of the estimated operator and is thereby random but not unknown.
As usual for nonparametric methods the rate of convergence increases if the true solution $\varphi^{\dagger}$ is smoother, i.e. if $\mu$ is larger. But this increase is only realized if $m \ge \mu$. However, $m$ does not act as a regularization parameter. Since the inner iteration is numerically cheap it is save to chose a larger value for $m$ without having a significant disadvantage.
The non-linearity of $\mathcal{F}$ needs to be restricted for the algorithm to work. We use a Lipschitz condition on the derivative for this purpose.
The special structure of the of $\mathcal{F}_{ind}$ and $\mathcal{F}_q$ allows us to replace Assumption (ref) for the IV regression examples by the following alternative. Lemma (ref) in the appendix shows that Assumption (ref) implies Assumption (ref).
Common nonparametric density estimators are strongly consistent. Assumption (ref) implies for the operators (ref) and (ref) \[ \sup_{y, x, z} \left| \frac{\partial^2}{\partial y^2} f_{YXZ}(y, x, z) \right| < \infty \qquad \mbox{or} \qquad \sup_{y, x, z} \left| \frac{\partial}{\partial y} f_{YXZ}(y, x, z) \right| < \infty \] respectively.
The estimation error in the operator and its derivative needs to be bounded by exponential inequalities.
The following lemma shows that Assumption (ref) holds for the IV regression applications (ref) and (ref) under mild conditions when $\mathbb{Y}$ is a $L^2$ space and $\|\cdot\|_D$ is the Hilbert-Schmidt norm. The Hilbert-Schmidt norm bounds the operator norm from above and is denoted by $\|\cdot\|_{HS}$. For linear integral operators it coincides with the $L^2$ norm of the integral kernel.
Our first convergence rate theorem assumes that $\mu$ in Assumption (ref) is known, i.e. the smoothness of the true solution is known. Adaptive estimation will be discussed in the next section.
In the special cases of the IV regression examples in $L^2$ spaces, the convergence rate can be given more explicitly.
We can compare our result to the rates for nonparametric quantile regression in HorLee:07. They use nonlinear Tikhonov regularization for the operator $\mathcal{F}_q$ and proved optimal rates under assumptions which are more restrictive than ours. In contrast to our rates, their rates do not depend on the derivative $\widehat T_{n \dag}$. The main challenge for nonlinear Tikhonov regularization \[ \widehat{\varphi} = \mathop{\mathrm{argmin}}_\varphi\|\mathcal{F}_q \varphi \|_\mathbb{Y}^2 + \alpha\|\varphi\|_\mathbb{X}^2 \] is to find the minimizer of the nonlinear functional $\|\mathcal{F}_q \varphi \|_\mathbb{Y}^2 + \alpha\|\varphi\|_\mathbb{X}^2$ which usually has multiple local minima. In HorLee:07 it is assumed that this minimum is known exactly which is unrealistic in practice. A convergence analysis which takes the performance of a minimization algorithm into account would typically lead to a different rate which also depends on some derivative depending on the particular minimization algorithm. The IRGNM does not have this problem since we only have to solve a linear least squares problem in every Newton step. For a fair comparison of the convergence rates, we will assume that the term $\left(\mathbb{E}\big[\|\widehat{\mathcal{F}}_n(\varphi^{\dagger})\|_\mathbb{Y}^2\big]\right)^{\frac{2\mu}{2\mu + 1}}$ dominates our rate. For sake of simplicity we assume $d_z=1$. Hence, our rate in the case of $\mathbb{X}=L^2$ is \[ \mathbb{E}\big[\|\widehat{\varphi}_{J^*} -\varphi^{\dagger}\|_{L^2}^2\big] = \mathcal{O} \left((n^{-1} h^{-2} + h^{2r} )^{\frac{2\mu}{2\mu + 1}}\right), \] while the rate in HorLee:07 is $\mathcal{O}\left(n^{-(2\beta-1)/(2\beta+a)}\right)$ with $a$ and $\beta$ defined an their paper as $\alpha$ and $\beta$. HorLee:07 restricts the bandwidth choice to \[ h = C_h n^{-\gamma} \quad \text{ with }\quad \frac{2\beta+a -1)}{2r(2\beta+a} <\gamma< \frac{a}{2(2\beta+a)}. \] Under their assumptions on $a$ and $\beta$, the density estimator achieves the rate \[ n^{-1} h^{-2} + h^{2r} = \mathcal{O}\left(n^{-\frac{2\beta +a-1}{2\beta-1}}\right). \] Note that this will only coincide with the optimal rate $n^{-\frac{2r}{2r+2}}$ in special cases. Futhermore, in their notation the source condition in our Assumption (ref) holds for any $\mu < (\beta-\frac{1}{2})/a$. Hence, $2\mu/(2\mu+1) < (2\beta-1)/(2\beta +a -1)$, where $\mu$ can be chosen such that the left hand side is arbitrarily close to the right hand side. Therefore, under their assumptions our rate is arbitrarily close to \[ \left(n^{-\frac{2\beta +a-1}{2\beta-1}}\right)^\frac{2\beta-1}{2\beta +a -1} = n^{-\frac{2\beta-1}{2\beta+a}}, \] which is also the optimal rate in Theorem 2 and 3 in HorLee:07.
Our assumptions are less restrictive which leads to faster rates of convergence in many cases compared to HorLee:07. Firstly, we have no restrictions on $h$ which means we can chose the optimal bandwidth $h = \mathcal{O}\left(n^{-\frac{1}{2r+2}}\right)$ and achieve \[ \mathbb{E}\big[\|\widehat{\varphi}_{J^*} -\varphi^{\dagger}\|_{L^2}^2\big] = \mathcal{O} \left(n^{-\frac{2r}{2r+2}\frac{2\mu}{2\mu + 1}}\right) \] This has also the advantage that we can choose $h$ by standard data driven bandwidth selectors like cross-validation. It is pointed out in HorLee:07 that a bandwidth selector for their assumptions does not yet exist.
Secondly, we have no upper bound on $\mu$ while their method requires $\mu \le 1$. Hence, our rate can be significantly better for smooth $\varphi^{\dagger}$. This is not a contradiction to their optimality result, as their result only holds under their more restrictive assumptions.
There is no uniform optimality theory for nonlinear inverse problems. Our setup is too general to settle optimality for the rates above. However, this can be achieved in special cases with further assumptions. We briefly discuss optimality for some special cases in this section. First of all, the rate in (ref) is known to be optimal in the case of a linear operator $\mathcal{F}$ with separable noise. See for example Tautenhahn:98.
For the nonparametric instrumental quantile regression, optimal minimax rates for convergence in probability were derived in HorLee:07. If in Corollary (ref) the influence of $\delta^{noi}_n$ and $\sigma^{noi}_n$ dominates $\delta^{der}_n$ and $\sigma^{der}_n$, the same rate is proved for the risk of the IRGNM. This result for the risk is even stronger than convergence in probability. Nevertheless, when $\delta^{der}_n$ and $\sigma^{der}_n$ dominate the convergence, the rate in Corollary (ref) can be slower than the one in HorLee:07. As discussed in Section (ref) the analysis of the nonlinear Tikhonov regularization does not consider the influence of an optimization algorithm that solves the minimization problem in (ref). If the properties of such an algorithm were added to the convergence analysis, the derivative of $\mathcal{F}_q$ will likely show up in the convergence rate in a similar way as in Corollary (ref).
Furthermore, we like to compare our results for the nonparametric instrumental regression with full independence to those with mean independents. When the influence of $\delta^{noi}_n$ and $\sigma^{noi}_n$ dominates $\delta^{der}_n$ and $\sigma^{der}_n$ and a source condition with the same $\mu$ holds for both models, the rates coincide. These rates are known to be optimal for the nonparametric instrumental regression with mean independents as shown in HalHor:05 and CheRei:11. However, the source condition for both models are linked only implicitly. The $\mu$ and thereby the rates can differ in actual examples. The advantage of model (ref) together with the IRGNM is not to improve upon the rates, but to be applicable in cases where model (ref) fails to identify the solution. \fi
The parameter choice (ref) in Theorem (ref) and Corollary (ref) depends on the unknown $\mu$ , which is unfeasible in practice. We present in this section convergence rates for a data driven choice of the stopping parameter $J$ by Lepski{\u\i}'s principle. This is a popular parameter choice rule in the context of statistical inverse problems, see Tsybakov:00, BauHoh:05, Mathe:06, BauHohMun:09, and HW:16. In the context of nonparametric IV, Lepski{\u\i}'s principle was used for adaptive estimation in ChenChrist:15. The following theorem gives convergence rates of the MISE with a Lepski{\u\i} type parameter choice. We lose a logarithmic factor compared to Theorem (ref). The constant $C_d$ and $\gamma_{nl}$ used in the theorem are specified in the appendix in formula (ref) and in Lemma (ref) respectively.
We evaluate the small sample behavior of the estimator based on the IRGNM with simulated data. As a test problem we use a nonparametric IV regression consistent with models (ref) and (ref), with univariate covariate $X$ and instrument $Z$. This setup allows to compare our estimator with an estimator which solves (ref) with iterated Tikhonov regularization.
The test problem described in the next section is solved on the domain \[ \mathrm{supp}\,(Y) \times \mathrm{supp}\,(X) \times \mathrm{supp}\,(Z) = [-1/2 , 1/2] \times [0 , 1] \times [0 , 1] \] discretized by an equidistant grid with $100 \times 100 \times 100$ nodes. The joint density is estimated on this grid by a standard adaptive density estimator. Trying different density estimators, we found that both method are quite robust with respect to the density estimate. They tolerate some undersmoothing of the density as long as the stopping index $J$ and the regularization parameter $\alpha$ are chosen properly. In the simulations below the same density estimate is used for both the IRGNM and the iterated Tikhonov regularization which allows for a fair comparison of the methods. The initial guess for both methods is the constant function with the value $\mathbb{E}[Y]$, and the penalty functional for both methods is the squared $H^1$ norm.
The Fr\'{e}chet derivative is implemented as in Example (ref). The partial derivative of the density and the derivatives for the $H^1$ norm are computed by the central differencing scheme. Operators and norms are evaluated using numerical integration. We tried rectangle rule, trapezoid rule and Simpson's rule but could not find a significant difference in the output of the estimator.
The least squares problems in each step of the iterated Tikhonov regularization and in the inner iteration of the IRGNM are computed by QR decomposition. Note that only one QR decomposition is needed in every Newton step. We tried different numbers of iterations $m$ in the inner iteration of the IRGNM and the iterated Tikhonov regularization for the test example below. No significant difference in the results was observed, which indicates that $\mu$ is not large.
The regularization parameters for the example below are $\alpha_0 = 1$ and $\alpha_{n+1} = 0.9 \alpha_n$. Lepski{\u\i}'s principle is used to find the stopping parameter of the Newton iteration. For the alternative estimator using model (ref) the regularization parameter $\alpha$ has to be chosen instead, which is done by Lepski{\u\i}'s principle as well. The iterated Tikhonov regularization is computed for a large number of different $\alpha$. Then one of these approximation is chosen by Lepski{\u\i}'s principle. Hence, both methods are fully data driven.
The regressor of the test example is generated by some function $g$ and a random variable $V$ such that $X = g(Z) + V$ and $V \protect\mathpalette{\protect\independenT}{\perp} Z$. In addition, an exact solution $\varphi^{\dagger}$ and an error term $U$ depending on $V$ but not on $Z$ are chosen. Then $Y$ is defined as $Y := \varphi^{\dagger}(X) + U$. With this construction both models (ref) and (ref) identify the true solution. The functions and probability densities that were used for the test example are
The densities of $V$ and $U$ are constructed with Gaussians in a way that the expectation of $U$ depends on $v$. Figure (ref) shows the exact solution (blue) compared to the solution a nonparametric regression without instrumental variables would yield asymptotically (green).
Both methods were tested on samples of $500$ and $1000$ observations. For each of the two sample sizes 1000 samples were generated and the joint density $f_{YXZ}$ was estimated by a kernel method.
Figures (ref) and (ref) show histograms for the $L^2$ error of the reconstructions for both methods and different sample sizes. The values are normed by the initial error, so that on this scale the initial error becomes $1$.
In Figure (ref) we compare the errors of both methods for the sample size $n=500$. Both methods produce acceptable results. The variance as well as the number of outliers observed for the method with independent instrument are significantly smaller than the variance or number of outliers of the method with the conditional mean assumption. The latter method produces a considerable number of outliers with the same or even larger errors than the initial guess. This cannot be observed for the IRGNM. In addition, the mean error of the IRGNM is smaller.
Similar histograms for sample size $n=1000$ are displayed in Figures (ref). Both methods perform well. The advantages of the IRGNM with less outliers, smaller variance and smaller mean error can be observed again. The following table provides the mean and some quantiles of the errors normed by the initial error.
We close this section with examples of median reconstructions for both sample sizes displayed in Figure (ref). They illustrate the advantage of the regression model with independent instruments solved with the IRGNM.
These results suggest that both methods give consistent estimators for the nonparametric instrumental regression with clear advantages for the regression model with independent instruments (ref) and the IRGNM.
The author would like to thank Thorsten Hohage and Johannes Schmidt-Hieber for interesting and fruitful discussions on this topic. He also would like to thank two anonymous referees for valuable comments that improved the paper.