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.
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.
Debiased Machine Learning of Set-Identified Linear Models
abstractThis paper provides estimation and inference methods for an identified set's boundary (i.e., support function) where the selection among a very large number of covariates is based on modern regularized tools. I characterize the boundary using a semiparametric moment equation. Combining Neyman-orthogonality and sample splitting ideas, I construct a root-$N$ consistent, uniformly asymptotically Gaussian estimator of the boundary and propose a multiplier bootstrap procedure to conduct inference. I apply this result to the Partially Linear Model, the Partially Linear IV Model and the Average Partial Derivative with an interval-valued outcome.
Introduction and Motivation.
Interval-valued outcomes are ubiquitous in economic research. Examples of such outcomes include bidders' valuation in English auctions (HaileTamer), income and wages (Trostel, Gafarov), house prices (GRT, BerSasaki), and county-level employment rates (Dorn). An outcome is interval-valued if the actual outcome $Y$ is missing, but there exist an observable lower bound $Y_L$ and an upper bound $Y_U$ so that
align[align omitted — 62 chars of source]
When an outcome is interval-valued, the parameter of interest is a set, where each point corresponds to a possible random variable $Y$ in the band (ref).
The main contribution of this paper is to provide an estimator of the identified set's boundary, where the selection among high-dimensional controls is based on modern machine learning/regularized methods. The paper focuses on identified sets whose boundary can be represented by a moment equation as in BM, BMM. In this paper, the equation depends on an identified functional nuisance parameter, for example, a conditional mean function. A naive approach would be to plug-in a machine learning estimate of the nuisance parameter into the moment equation and solve for the boundary. However, modern regularized methods (machine learning techniques) have bias converging slower than the parametric rate, which cannot be made small by classic techniques (e.g., undersmoothing). As a result, plugging such estimates into the moment equation produces a biased, suboptimal estimate of the boundary itself.
To overcome the transmission of the bias into the second stage, I adjust the moment equation to make it insensitive or, formally, Neyman-orthogonal, to the biased estimation of the nuisance parameter. While orthogonality has been extensively studied in the point-identified case, set identification presents several challenges. The first one is the non-smoothness of $x \rightarrow \min (x,0)$ function at $x=0$, often occurred in censored LAD (e.g., Powell). The second one is to establish uniformly valid inference over the boundary in addition the pointwise one.
When the identified set is multi-dimensional, its boundary consists of continuum points. As a result, economists are interested in uniform inference in addition to the pointwise inference. Establishing uniform inference is not trivial. To control the speed at which an empirical sample average concentrates around the population mean, I invoke maximal inequalities of CCKAS instead of Markov inequality that is typically sufficient in the point-identified case. I propose multiplier bootstrap algorithm to conduct inference. By virtue of orthogonality, only the moment function (the second stage), not the nuisance parameter estimate (the first stage), needs to be resampled in simulation. As a result, this multiplier bootstrap is faster to compute than the weighted bootstrap of CCMS, which is based on a non-orthogonal moment and involves resampling of both stages. I demonstrate the proposed approach in a simulation exercise and give a brief empirical illustration.
\paragraph{Literature review. Set identification} This paper bridges the gap between three literatures: set-identified models, debiased/orthogonal machine learning, and non-smooth models. Set identification is a vast area of research, encompassing a wide variety of approaches: linear and quadratic programming, random set theory, support function, and moment inequalities
(Manski90, ManskiPepper, Manski:2002, HaileTamer, CHT, BM, Molinari2008, CilibertoTamer, LeeBound, Stoye, AndrewsShiECMA, BMM2, CCMS, BMM3, CherRigStoker, BMM, CLR, FanPark, KaidoWhite, KaidoSantos, KaidoWhite2, Pakesetal, ShiShum, Kaido:2016, Kasy2016, KlineTamer, AndrewsShi, CanayBugniShi, Kaido, ChenTamerChristensen, GafarovMeierOlea, Shi, Gafarov, KaidoMolinariStoye, SyrgkanisTamer, Torgovitsky,MolinariStoye, BerSasaki, Honore, AndrewsRothPakes, kallus2020localized, FanTao, FanTao2, MolinariMolchanovPeng, HsiehShiShum, DongHsiehShum), see e.g. Tamer:2010 or Molinari:2018 for a review. This paper generalizes the BMM's model by allowing its components to depend on a functional nuisance parameter, covering e.g., CCMS and Kaido as special cases.
\paragraph {Orthogonality. }Next, this paper contributes to a large body of work on debiased inference for parameters following regularization or model selection (Neyman:1959, Neyman:1979, HardleStoker1989, NeweyStoker, Newey1994, Robins, robinson:88, ZhangZhang, JM, chernozhukov2016double, LRSP, Program, sasaki2018estimation, sasaki2020unconditional, Sasaki, chiang2019multiway, ning2020doubly, chernozhukov2021debiased, chernozhukov2021automatic, CherSem, NSS, singh2020debiased, Colangelo, Lieli, ZimLech). A basic idea is to make the moment condition insensitive, or, formally, Neyman-orthogonal, to the biased estimation of the nuisance parameter. For a semiparametric GMM setting, the work by AckerbergChen derives an orthogonal moment condition whose nuisance functions are identified by conditional moment restriction, such as conditional mean and conditional quantiles. Combining Neyman-orthogonality and sample splitting, LRSP and chernozhukov2016double derive a root-$N$ consistent and asymptotically normal estimator for a single target parameter. This has idea has been extended for many functional parameters in $Z$-estimation framework, in the context of distribution regression (BelCherWei) and quantile regression (sasaki2020unconditional). Next, the paper is related to literature on the non-smooth estimating equations (Powell, Powell2, PowellStockStoker, kaplan_sun_2017, franguridi2021conditional). Finally, the paper contributed to a small, but growing literature on machine learning for bounds and partially identified models (kallus2019assessing, jeong2020robust, SemSupp2, Bonvini_2021).
\paragraph{Structure of the paper.} The paper is organized as follows. Section (ref) demonstrates main points for the partially linear model of robinson:88. Section (ref) states theoretical results. Section (ref) applies the results to models with an interval-valued outcome. Section (ref) presents finite-sample evidence. Section (ref) contains an empirical illustration. Section (ref) contains the proofs of main results.
Set-Up.
General Framework
I focus on parameters that are linear in an unobserved scalar outcome $Y$. The identified set takes the form
align[align omitted — 125 chars of source]
where the random vector $V(\eta_0) \in \mathrm{R}^{d}$ depends on a nuisance function $\eta_0$. Examples of the nuisance functions
$$
\eta_0 = \eta_0(X)
$$
include the propensity score, the conditional density, and the regression function, among others. The matrix $\Sigma \in \mathrm{R}^{d \times d}$ is identified by a moment equation
align[align omitted — 61 chars of source]
and is assumed invertible. The key innovation of this framework is to allow the vector $V(\eta)$, the matrix function $A(W, \eta)$, and the bounds $Y_L, Y_U$ to depend on a functional nuisance parameter $\eta$, covering the models in BMM, CCMS, Kaido, and many others as special cases.
Examples
example[Partially Linear Model]
Consider the partially linear model of robinson:88
\begin{align}
Y &= D' \beta_0 + f_0(X) + U, \quad \mathbb{E} [U\mid D,X] = 0,
\end{align}
where $D \in \mathrm{R}^d$ is a treatment (policy) variable, $\beta_0$ is a causal (structural) parameter, $$X = (1,X_1, X_2, \dots, X_{p_X}) $$ is a vector of covariates whose dimension may be large relative to the sample size (e.g., $p_X \gg N$), and $f_0(\cdot)$ is an integrable function. The parameter $\beta_0$ can be represented as the minimizer of the least squares criterion function
\begin{align}
\beta_0& =\arg \min_{b \in \mathrm{R}^d, f \in L_2(P)} \mathbb{E} (Y - D' b - f(X))^2.
\end{align}
As in Lemma (ref), $\beta_0$ coincides with the minimizer of a shorter criterion function
\begin{align}
\beta_0& =\arg \min_{b \in \mathrm{R}^d} \mathbb{E} (Y - (D - \eta_0(X))' b)^2.
\end{align}
Thus, the identified set $\mathcal{B}$ is a special case of model (ref)-(ref) with $V(\eta)$ and $A(W,\eta)$ defined as follows. The treatment regression function is
\begin{align}
\eta_0(X) = \mathbb{E}[ D \mid X],
\end{align}
the treatment residual is
\begin{align}
V (\eta) = D-\eta(X),
\end{align}
the matrix function is
\begin{align}
A(W,\eta) = (D - \eta(X))(D - \eta(X))'.
\end{align}
example[Partially Linear IV Model]
Consider the following partially linear IV model
\begin{align}
Y &= D' \beta_0 + f_0(X) + U, \quad \mathbb{E}[ U \mid Z, X]=0 \\
D &= m_0(X) + E, \quad \mathbb{E} [ E \mid X ] =0, \\
Z &= \eta_0(X) + V, \quad \mathbb{E}[ V \mid X ] =0,
\end{align}
where $Z \in \mathrm{R}^d$ is the instrument for $D$. The parameter $\beta_0$ can be characterized by a moment equation
\begin{align*}
\mathbb{E} (Z - \eta_0(X)) ( Y - (D - m_0(X))' \beta_0) = 0,
\end{align*}
which gives a closed-form expression for $\beta_0$
\begin{align*}
\beta_0 = (\mathbb{E} (Z - \eta_0(X)) (D - m_0(X))' )^{-1} \mathbb{E} (Z - \eta_0(X)) Y.
\end{align*}
The model is a special case of (ref)-(ref) with
\begin{align}
V (\eta) &= Z - \eta(X) \\
A (W, \eta,m) &= (Z - \eta(X)) (D - m(X))',
\end{align}
where $V = V(\eta_0)$ in (ref). If $Z = D$, the model (ref)-(ref) coincides with (ref)-(ref).
example[Average Partial Derivative]
An important parameter in economics is the average partial derivative. This parameter shows the average effect of a small change in a variable of interest $D$ on the outcome $Y$ conditional on the covariates $X$. To describe this change, define the conditional expectation function of an outcome $Y$ given the variable $D$ and exogenous variable $X$ as $$\mu(D,X):= \mathbb{E} [Y|D,X]$$ and its partial derivative with respect to $D$ as $\nabla_D\mu(D,X):=\nabla_D\mu(d,X)|_{d=D}$. Then, the average partial derivative is defined as
\begin{align}
\beta &= \mathbb{E} \nabla_D\mu(D,X).
\end{align}
For example, when $Y$ is the logarithm of consumption, $D$ is the logarithm of price, and $X$ is the vector of other demand attributes, the average partial derivative stands for the average price elasticity.
Assume that the variable $D$ has bounded support $\mathcal{D} \subset \mathcal{R}^d$. Furthermore, the conditional density $f ( D \mid X)$ has positive density on this support a.s. in $X$. HardleStoker1989 have shown that the average partial derivative can be represented as $\beta = \mathbb{E} V Y,$
where $V = - \nabla_D \log f(D|X) = - \frac{\nabla_D f(D|X)}{f(D|X)}$ is the negative partial derivative of the logarithm of the density $f(D|X)$. When $Y$ is interval-valued, the identified set $\mathcal{B}$ for $\beta$ is a special case of (ref)-(ref) with $A(W, \eta) = \Sigma = I_d$, the nuisance function $\eta_0(D, X) = \frac{\nabla_D f(D|X)}{f(D|X)}$ and the vector $V ( \eta) = -\eta$. Kaido studies a special case of this problem without covariates.
Single treatment
In this section, I derive an orthogonal moment equation for the upper bound $\beta_U$ on the causal parameter $\beta_0$ in Example (ref). Because the bias due to sign mistake in $Y^{\text{best}}(\eta) \neq Y^{\text{best}}(\eta_0)$ proves to be second-order, I derive the orthogonal moment treating $Y^{\text{best}}(\eta_0)$ as observed. Next, I formally control the bias due to the sign mistake.
The following two subsections elaborate on these points.
\paragraph{Moment equation for $\beta_U$. }Consider Example (ref) with a single treatment (i.e., $d=1$). The identified set $\mathcal{B}$ becomes a closed interval $[\beta_L, \beta_U]$. The upper bound $\beta_U$ is
align[align omitted — 172 chars of source]
To maximize the numerator of (ref), take $Y=Y_U$ for positive values of $D - \eta_0(X)$ and $Y=Y_L$ otherwise. Define the best-case outcome
align[align omitted — 147 chars of source]
Plugging (ref) into (ref) gives the moment function for $\beta_U$
align[align omitted — 124 chars of source]
\paragraph{ Establishing orthogonality. } In what follows, $\eta$ corresponds to an instance of the nuisance parameter whose true value is $\eta_0$. Consider an infeasible moment function
align[align omitted — 152 chars of source]
where the best-case outcome $Y^{\text{best}}(\eta_0) $ is treated as if it was observed, and $ \pmb{\eta}$ indicates the instance of the nuisance parameter treated as unknown. The moment equation (ref) is not orthogonal to the perturbations of $\widehat{\eta}(X) - \eta_0(X)$
align*[align* omitted — 168 chars of source]
Therefore, the bias of the estimation error $\widehat{\eta}(X)-\eta_0(X)$ translates into the moment (ref). To overcome the transmission of this bias, robinson:88 proposes an orthogonal moment equation
align*[align* omitted — 168 chars of source]
where the true value of $\gamma_U(x)$ is $$\gamma_{U,0}(X) = \mathbb{E}[ Y^{\text{best}}(\eta_0) \mid X]. $$ The orthogonality condition is
align[align omitted — 265 chars of source]
Thus, the bias of the estimation error, $\widehat{\eta}(X) - \eta_0(X)$, does not translate into the moment (ref). The orthogonality condition with respect to $\gamma_{U,0}$ can be verified in a similar way.
\paragraph{ Bias due to sign mistake. } I now discuss whether the mistake in the best-case outcome $Y^{\text{best}} (\eta) \neq Y^{\text{best}} (\eta_0)$, which has been so far ignored, has any effect on the support function estimate. Consider the difference between the feasible moment $g (W, \beta_U, \{ \pmb{\eta}, \pmb{\gamma}_U \})$ and its infeasible analog
align*[align* omitted — 351 chars of source]
Define the first-order bias $B_1( \eta, \eta_0)$ as
align*[align* omitted — 118 chars of source]
and the second-order one
align*[align* omitted — 127 chars of source]
Below, I describe the conditions under which $B_1( \eta, \eta_0)$ and $ B_2( \eta, \eta_0)$ are negligible, (i.e., $o(N^{-1/2})$). If they hold, the feasible moment equation $g (W, \beta_U, \{ \pmb{\eta}, \pmb{\gamma}_U \})$ is insensitive to the biased estimation of $\eta$. As a result, the support function estimator based on $g (W, \beta_U, \{ \pmb{\eta}, \pmb{\gamma}_U \})$ is asymptotically unbiased under plausible conditions.
Consider the best-case outcome $Y^{\text{best}} (\eta)$
$$
Y^{\text{best}} (\eta) = Y_L + (Y_U - Y_L) 1 \{ D- \eta(X) \geq 0 \}.
$$
The sign mistake $Y^{\text{best}} (\eta) \neq Y^{\text{best}} (\eta_0)$ occurs on the events
align[align omitted — 190 chars of source]
On these events, the residual cannot exceed estimation error in absolute value
align[align omitted — 167 chars of source]
If the width $Y_U - Y_L$ is bounded by $M_{UL}$, the estimation error of $Y^{\text{best}}(\eta) $ is bounded as
align[align omitted — 240 chars of source]
Suppose the conditional density $ h_{V (\eta_0) \mid X} (t,X) $ of $V(\eta_0)$ is bounded by $M_h$ a.s.. Invoking (ref) gives
align[align omitted — 303 chars of source]
As a result, the bias terms $B_1( \eta, \eta_0)$ and $B_2( \eta, \eta_0)$ shrink at the quadratic rate. Combining (ref) and (ref) gives a feasible moment function
align[align omitted — 188 chars of source]
Multi-dimensional case
In this section, I derive an orthogonal moment for the support function, starting from a non-orthogonal one due to BMM, BM.
\paragraph{Moment Equation for Support Function. } As shown in BM, the identified set $\mathcal{B}$ in (ref) is a compact and convex set. Thus, it can be described by its projections onto a unit sphere
align[align omitted — 99 chars of source]
For any direction $q \in \mathcal{S}^{d-1}$, define the support function as the upper bound on $q' \beta_0$
align[align omitted — 78 chars of source]
As proposed in BM and BMM, define the projected weighting vector
align[align omitted — 52 chars of source]
the best-case outcome $Y(p,\eta)$
align[align omitted — 85 chars of source]
and the projection parameter $p(q)$
align[align omitted — 48 chars of source]
Then, the moment equation for $\sigma(q)$ is
align[align omitted — 103 chars of source]
\paragraph{Orthogonal Moment for Support Function. } The moment equation (ref) is sensitive to the biased estimation of the nuisance parameter $\eta$. To avoid the transmission of this bias into the second stage, I construct another moment function $g(W, p, \xi),$ where $\xi(p)$ is the nuisance parameter, in the following steps.
enumerate• Starting from an infeasible, smooth moment
\begin{align}
m_0(W, p, \pmb{\eta}) := z (p,\pmb{\eta}) Y (p,\eta_0)
\end{align}
derive an infeasible orthogonal moment $g_0(W, p, \xi(p)) $ obeying (ref) for each $p$.
• Invoke Lemma (ref) to bound the bias
\begin{align}
\sup_{p \in \mathcal{P}} \mathbb{E} | z(p,\eta) \left(Y (p,\eta)- Y (p,\eta_0) \right) | = O( \mathbb{E} \| \eta(X) - \eta_0(X) \|^2),
\end{align}
where $p$ belongs to a compact bounded set $\mathcal{P}$ defined below. When $\Sigma = I_d$, $\mathcal{P} = \mathcal{S}^{d-1}$.
• Combine (ref) and (ref) to obtain the feasible orthogonal moment
\begin{align}
g(W, p, \xi(p)) = g_0(W, p, \xi(p)) + z(p,\eta) (Y (p,\eta)- Y (p,\eta_0)).
\end{align}
example*[Example (ref), cont.]
Consider Example (ref). The projected weighting vector (ref) is
\begin{align*}
z(p, \eta) = p' (Z - \eta(X)).
\end{align*}
The infeasible orthogonal moment $g_0(W, p, \xi(p)) $ is
\begin{align}
g_0(W, p, \xi(p)) = z(p,\eta) (Y (p,\eta_0) - \gamma (p, X)),
\end{align}
where
\begin{align}
\gamma_0(p,x) &= \mathbb{E} [ Y (p,\eta_0) \mid X] \\
&= \mathbb{E} [ Y_L \mid X=x] + \mathbb{E} [ (Y_U - Y_L) 1\{ p' V(\eta_0) > 0 \} \mid X=x] \nonumber
\end{align}
The nuisance function $\xi_0(p) = (\eta(\cdot), \gamma(p,\cdot))$. Invoking (ref) gives a feasible orthogonal moment
\begin{align}
g(W, p, \xi(p)) = z(p,\eta) (Y (p,\eta) - \gamma (p, X)).
\end{align}
Corollary (ref) establishes the asymptotic theory for the support function estimator based on (ref).
example*[Example (ref), cont.]
Consider Example (ref). The projected weighting vector is
\begin{align*}
z(q, \eta) = - q' \partial_D \log f(D \mid X).
\end{align*}
The infeasible orthogonal moment $g_0(W, p, \xi(p)) $ is
\begin{align}
g_0(W,q,\xi(q)) = z(q, \eta) Y(q, \eta_0) + q' \partial_D \log f(D \mid X) \mu(q,D,X) + q' \nabla_D \mu(q,D,X),
\end{align}
where
\begin{align*}
\mu_0(q,D,X) &= \mathbb{E} [ Y_L \mid D, X] + \mathbb{E} [ (Y_U - Y_L) \mid D,X ] 1\{ -q' \partial_D \log f(D \mid X) > 0 \} \\
&= \gamma_{L,0}(D, X) + \gamma_{UL,0}(D, X) 1\{ -q' \partial_D \log f(D \mid X) > 0 \}.
\end{align*}
Invoking (ref) gives a feasible moment equation
\begin{align}
g(W, p, \xi(p)) &= g_0(W,q,\xi(q)) + z(q,\eta) (Y (q,\eta) -Y (q,\eta_0)) \\
&= z(q,\eta) Y(q, \eta) +q' \partial_D \log f(D|X) \mu(q,D,X) + q' \nabla_D \mu(D,X) \nonumber.
\end{align}
In absence of the conditioning covariates, (ref) coincides with the efficient score in Kaido. Corollary (ref) establishes the asymptotic theory for the support function estimator based on (ref).
Overview of Main Results
The Support Function Estimator $\widehat{\sigma}(q)$ has two stages. In the first stage, I construct an estimate $\widehat{\xi}$ of the nuisance parameter $\xi_0$ using some regularized estimator. In the second stage, I compute the estimated values $(\widehat{\xi}_i)_{i=1}^N$ and the support function estimate. I use different samples in the first and the second stage in the form of cross-fitting. The number $K$ of cross-fit partitions is assumed to be fixed/finite relative to $N$.
definition[Cross-Fitting]
\begin{compactenum}
• For a random sample of size $N$, denote a $K$-fold random partition of the sample indices $[N]=\{1,2,...,N\}$ by $(J_k)_{k=1}^K$, where $K$ is the number of partitions and the sample size of each fold is $n = N/K$. For each $k \in [K] = \{1,2,...,K\}$ define $J_k^c = \{1,2,...,N\} \setminus J_k$.
• For each $k \in [K]$, construct estimates $\widehat{\xi}_k = \widehat{\xi}( W_{i \in J_k^c})$ and $\widehat{\eta}_k = \widehat{\eta}( W_{i \in J_k^c})$ of the nuisance parameters $\xi_0$ and $\eta_0$ using only the data $\{ W_{j}: j \in J_k^c \}$. For any observation $i \in J_k$, define $\widehat{\xi}_i = \widehat{\xi}_k (W_i)$ and $\widehat{\eta}_i = \widehat{\eta}_k (W_{i})$.
\end{compactenum}
Definition (ref) introduces cross-fitting. Cross-fitting plays an essential role in modern debiased inference in semi-parametric models; see, e.g., bch:2010,zheng:laan,chernozhukov2016double for recent examples and hasminskii:debiased and schick1986asymptotically for early, classical uses of simpler sample-splitting methods for debiased inference. Cross-fitting is proposed for the cases where the nuisance parameter $\xi_0(\cdot)$ does not depend on $p \in \mathcal{P}$, such as Example (ref) and Examples (ref)--(ref) under Assumption (ref).
definition[Support Function Estimator]
Let $\widehat{\xi}$ and $\widehat{\eta}$ be the estimates of $\xi_0$ and $\eta_0$. Define
\begin{align}
\widehat{\Sigma} &= \dfrac{1}{N} \sum_{i=1}^N A(W_i, \widehat{\eta}_i), \quad \widehat{p}(q) = (\widehat{\Sigma}^{-1})' q \\
\widehat{\sigma}(q)&= \dfrac{1}{N} \sum_{i=1}^N g(W_i, \widehat{p}(q) , \widehat{\xi}_i( \widehat{p}(q) )).
\end{align}
definition[Multiplier Bootstrap]
Let $(e_i)_{i=1}^N: e_i $ are i.i.d. truncated exponential random variables $\text{Tr Exp}(1)$ on $[0, \bar{M}]$ independent of the data. Define the bootstrap analog of $ \widehat{\sigma}(q)$ as
\begin{align}
\widetilde{\Sigma} &= \dfrac{1}{N} \sum_{i=1}^N \dfrac{e_i}{\bar{e}} A(W_i, \widehat{\eta}_i), \quad \widetilde{p}(q) = (\widetilde{\Sigma}^{-1})' q \\
\widetilde{\sigma}(q)&= \dfrac{1}{N} \sum_{i=1}^N \dfrac{e_i}{\bar{e}} g(W_i, \widetilde{p}(q) , \widehat{\xi}_i( \widetilde{p}(q) )).
\end{align}
Under mild conditions on $\xi$, the Support Function Estimator delivers a high-quality estimate $\widehat{\sigma}(q)$ of the support function $\sigma(q)$ with the following properties
enumerate• With probability (w.p.) $\rightarrow 1$, the estimator converges uniformly over the unit sphere $\mathcal{S}^{d-1}$
\begin{align}
\sup_{ q \in \mathcal{S}^{d-1} } | \widehat{\sigma} (q) - \sigma(q) | = O_P(1/\sqrt{N}) = o_P(1).
\end{align}
• The estimator $\widehat{\sigma}(q)$ is asymptotically Gaussian
\begin{align}
S_N(q) :=\sqrt{N} (\widehat{\sigma}(q) - \sigma(q)) = \mathbb{G}_N(q)+ o_P(1) \quad uniformly in \mathcal{S}^{d-1},
\end{align}
where the empirical process $\mathbb{G}_N(q)$ is approximated by a Gaussian process $\mathbb{G}(q)$, which is a tight $P$-Brownian bridge in $ \ell^{\infty}(\mathcal{S}^{d-1})$.
Define the bootstrap statistic
align*[align* omitted — 88 chars of source]
\paragraph{Pointwise asymptotics. } The sharp identified set for $q' \beta_0$ is $[ - \sigma(-q), \sigma(q)]$. Its $(1-\tau)$-pointwise confidence region (CR) is
align*[align* omitted — 171 chars of source]
where the critical values $\widehat{C}_{\tau/2} (q)$ and $\widehat{C}_{1-\tau/2} (q) $ are the $\tau/2$ and $1-\tau/2$ quantiles of the bootstrapped statistic $|\widetilde{S}_N(q)|$. Plugging $$q= \mathbf{e}_k = (0,0,\dots, 0, \underbrace{1}_{k}, 0, \dots, 0) \in \mathrm{R}^d, \quad k = 1,2, \dots, d$$ gives the CR for the projection of the identified set $\mathcal{B}$.
\paragraph{Uniform asymptotics. } The $(1-\tau)$-uniform confidence region (CR) for $[ - \sigma(-q), \sigma(q)]$ is
align[align omitted — 195 chars of source]
where $\widehat{C}^{*}_{\tau/2}$ and $\widehat{C}^{*}_{1-\tau/2}$ are the quantiles of the bootstrapped statistic $\sup_{q \in \mathcal{S}^{d-1}} |\widetilde{S}_N(q)|$. Likewise, for any function $f(\cdot)$ and a critical value $\widehat{c}_N = c_N + o_P(1)$ and $c_N = O_P(1)$,
align*[align* omitted — 130 chars of source]
where ${\mathrm{P}}^{e} (\cdot)$ is the probability conditional on the data.
Theoretical Results.
\setcounter{example}{0}
\paragraph{Notation.} I use the empirical process notation. For a generic function $f$ and a generic sample $(W_i)_{i=1}^N$, denote the empirical sample average by $$ {\mathbb{E}_{N}} f(W_i) := \dfrac{1}{N} \sum_{i=1}^N f(W_i) $$ and the
scaled, demeaned sample average by $$\mathbb{G}_N f(W_i) := 1/\sqrt{N} \sum_{i=1}^N [f(W_i) - \int f(w) d P(w)].$$
For two sequences of random variables $\{ a_N, b_N, N \geq 1\}: a_N \lesssim_{P} b_N$ means $a_N = O_{P} (b_N)$. For two sequences of numbers $\{a_N, b_N, N \geq 1\}$, $a_N \lesssim b_N$ means $a_N = O (b_N)$. Let $a \wedge b = \min \{ a, b\}, a \vee b = \max \{ a, b\} $. The $\ell_2$ norm of a vector is denoted by $\| \cdot \|$, the $\ell_1$ norm is denoted by $\| \cdot \|_1$, the $\ell_{\infty}$ norm is denoted by $\| \cdot \|_{\infty}$, and $\ell_0$ norm is denoted by $\| \cdot \|_{0}$. For a matrix $Q$, let $\|Q\|$ be the maximal eigenvalue of $Q$ and $\| Q \|_F$ be the Frobenius norm of $Q$. For a random vector $W$, let $\| W \|_{P,c}:= (\int |W|^c d P)^{1/c}$. The random sample $(W_i)_{i=1}^N$ is a sequence of independent copies of a random element $W$ taking values in a measurable space $(\mathcal{W}, \mathcal{A}_{\mathcal{W}})$ according to a probability law $P$. The $\| f \|_{P_N, 2}$ is the empirical $\ell_2$-norm, denoted as $\| f \|_{P_N, 2}:=(N^{-1} \sum_{i=1}^N f^2(W_i))^{1/2}$. Define the projection set
align[align omitted — 185 chars of source]
and let $C_P:= 2 \max \operatorname{eig} (\Sigma^{-1}) $. Let $\mathcal{F}_c$ be the space of continuous functions obeying two conditions: (1) $f(Z)$ has a continuous distribution when $Z$ is a tight Gaussian process with non-degenerate covariance function and (b) $f(\xi_N + c) -f(\xi_N) = o(1)$ for any $c=o(1)$ and any $\| \xi_N \|= O_P(1)$. Note that $[p_X]:=\{1,2,\dots, p_X\}$.
Assumptions
Assumption (ref) is a standard identification condition. It ensures that the eigenvalues of $\Sigma$ are bounded from above and below.
assumption[Identification]
There exist constants $\lambda_{\min} >0$ and $\lambda_{\max} < \infty $ that bound the eigenvalues of $\Sigma$ in (ref) from above and below $0 < \lambda_{\min} \leq \min \operatorname{eig} (\Sigma) \leq \max \operatorname{eig}(\Sigma) \leq \lambda_{\max}$.
Assumption (ref) ensures that the support function is differentiable on the unit sphere. It requires the distribution of the weighting vector $V(\eta_0)$ to be sufficiently smooth. For example, if the vector $V(\eta_0)$ has a symmetric (i.e., spherical) distribution around zero, the normalized vector $\| V(\eta_0) \|^{-1} V(\eta_0)$ is uniformly distributed on the surface of the unit sphere $\mathcal{S}^{d-1}$, and Assumption (ref) holds. Assumption (ref) is a common regularity condition in set-identified models (e.g., Condition C.1 in CCMS) and censored median regression (e.g., Assumption R.2 in Powell).
assumption[Smooth boundary]
Let $d \geq 2$. There exists a finite constant $\mathcal{C}_V$ such that
\begin{align}
\sup_{q \in \mathrm{S}^{d-1}} {\mathrm{P}} \left( \dfrac{ |q' \Sigma^{-1/2} V (\eta_0) | }{ \| \Sigma^{-1/2} V (\eta_0) \| } \leq \delta \right) \leq \mathcal{C}_V \delta.
\end{align}
example[Gaussian Noise]
Consider Example (ref) with $$V(\eta_0) = Z - \eta_0(X) \sim N(0, \Sigma)$$ independent of $X$. Then, $\Sigma^{-1/2} V(\eta_0) \sim N(0, I_d)$ is standard Gaussian vector and
$\Sigma^{-1/2} V(\eta_0)/ \| \Sigma^{-1/2} V (\eta_0) \|$ is uniformly distributed over the unit sphere $\mathcal{S}^{d-1}$. For any $q \in \mathcal{S}^{d-1}$, $|q'\Sigma^{-1/2} V(\eta_0)|/ \| \Sigma^{-1/2} V (\eta_0) \|$ is uniformly distributed on $[0,1]$ (Pitman). As a result, (ref) holds with $\mathcal{C}_V = 1$ conditional on $X$ uniformly in $X$.
Assumption (ref) is a sufficient condition for the smoothness of the boundary. If it holds, the moment equation (ref) is differentiable in $p$. The gradient
align[align omitted — 71 chars of source]
is a uniformly continuous function of $p$ (see Lemma (ref) in Online Appendix). As a result, there exists a uniform Gaussian approximation for the support function estimator.
remarkConsider Example (ref). If $D$ and $X$ consist of discrete variables only, the distribution of $V (\eta_0)$ cannot be continuous, and Assumption (ref) fails. Discrete distributions imply flat surfaces on the identified set, which may not be compatible with uniform Gaussian approximation. In this case,
CCMS suggests adding a small amount of continuously distributed noise and work with a slightly expanded set with smooth boundary, while Gafarov provides an alternative approach.
The moment functions $A(W, \eta)$ in (ref) and $g(W, p, \xi(p))$ depend on the nuisance parameters $\eta_0$ and $\xi_0$, respectively. Definition (ref) introduces a sequence of nuisance realization sets $\Xi_N \subseteq \Xi$ and $\mathcal{T}_N \subseteq \mathcal{T}$ that contain $\xi_0$ and $\eta_0$, as well as their estimators $\widehat{\xi}$ and $\widehat{\eta}$, with probability $1-\phi_N$. As the sample size increases, the sets $\Xi_N$ and $\mathcal{T}_N$ shrink. The shrinkage speed is measured by the rates below.
definition[Moment Rates]
Let $\{ \Xi_N, N \geq 1\}$ and $\{ \mathcal{T}_N, N \geq 1\}$ be sequences of subsets of $\Xi$ and $\mathcal{T}$, respectively, obeying the following conditions. (1) The true values $\xi_0$ and $\eta_0$ belong to $\Xi_N$ and $\mathcal{T}_N$ for all $N \geq 1$. There exists a sequence of numbers $\phi_N = o(1)$ such that the first-stage estimators $\widehat{\xi}$ of $\xi_0$ and $\widehat{\eta}$ of $\eta_0$ belong to $\Xi_N$ and $\mathcal{T}_N$ with probability at least $1-\phi_N$. Define the sequences $\mu_N, r_N'', r_N', A_N, \delta_N$ as
\begin{align}
\sup_{\xi \in \Xi_N} \sup_{p \in \mathcal{P}} | \mathbb{E} [ g(W, p, \xi(p)) - g (W,p, \xi_0(p)) ] | = \mu_N \\
\sup_{\eta \in \mathcal{T}_N} \| \mathbb{E} [ A (W, \eta) - A(W, \eta_0) ] \| = A_N \\
\sup_{\xi \in \Xi_N} \sup_{p \in \mathcal{P}} (\mathbb{E} (g(W, p, \xi(p)) - g (W,p, \xi_0(p)) )^2 )^{1/2} = r_N” \\
\sup_{p, p_0 \in \mathcal{P} \| p - p_0 \| \lesssim \tau_N} (\mathbb{E} ( g(W, p, \xi_0(p)) - g(W,p_0, \xi_0(p_0)) )^2 )^{1/2} = r_N' \\
\sup_{\eta \in \mathcal{T}_N} (\mathbb{E} \| A(W, \eta) - A(W, \eta_0) \|^2)^{1/2} = \delta_N,
\end{align}
where $\tau_N:= N^{-1/2} \log N$ for the case when $A(W,\eta_0)$ is being estimated and $\tau_N:= r_N':=0$ when $A(W, \eta_0) = \Sigma$ is known.
assumption[Regularity Conditions]
(1) There exist absolute constants $c'>2$ and $\bar B_A< \infty$ such that the following matrix norms are bounded
$$\sup_{\eta \in \mathcal{T}_N} \mathbb{E} \| A(W, \eta) \|^{c'} \leq \bar B_A, \quad \mathbb{E} \| A(W, \eta_0) \|_F^{2} \leq \bar B_A$$ (2) There exists a sequence $v_N = o(N^{-1/4})$ such that $$ \| \mathbb{E}_N A(W_i, \eta_0) - \Sigma \| = O_P (v_N) = o_P(1). $$
Assumption (ref) requires the rates of Definition (ref) to decay sufficiently fast. The matrix moment function $A(W, \eta)$ is assumed to be already orthogonal with respect to $\eta$. This is the case for covariance matrices in Examples (ref) and (ref).
assumption[Convergence Rates]
(1) For $c_3>2$ in Assumption (ref), suppose the numerator moment rates $\mu_N, r_N'', r_N'$ obey the following bounds: $\mu_N = o(N^{-1/2})$ and $$(r_N'' + r_N' ) \log^{1/2} (1/(r_N'' + r_N')) + N^{-1/2+1/c_3} \log N = o(1).$$
(2) The matrix rates $A_N$ and $\delta_N$ obey the following bounds: $A_N = o(N^{-1/2})$ and $\delta_N = o(1)$.
Assumption (ref) bounds the complexity of the function class $$\mathcal{G}_{\xi} = \{ g (W,p, \xi(p)), p \in \mathcal{P}\}.$$
assumption[Complexity Conditions]
(1) There exists a measurable envelope function $G_{\xi} = G_{\xi}(W)$ that almost surely bounds all elements in the class $$ \sup_{p \in \mathcal{P}} | g(W,p,\xi(p)) | \leq G_{\xi}(W)
\quad \text{a.s.}$$
There exists $c>2$ such that $\|G_{\xi}\|_{P,c} := \left(\int_{w \in \mathcal{W}} (G_{\xi}(w))^c \right)^{1/c} < \infty$. (2) There exist constants $a,v$ that do not depend on $N$ such that the uniform covering entropy of the function class $\mathcal{G}_{\xi}$ is bounded
\begin{align}
\log \sup_{Q} N(\epsilon \| G_{\xi} \|_{Q,2}, \mathcal{G}_{\xi} , \| \cdot \|_{Q,2}) \leq v \log (a/\epsilon), \quad for all 0 < \epsilon \leq 1.
\end{align}
Results
Define the influence function $h_g(W,q)$ as
align[align omitted — 77 chars of source]
and the influence function for the matrix estimation
align[align omitted — 106 chars of source]
where $G(p)$ is the gradient defined in (ref). Finally, define
align[align omitted — 57 chars of source]
theorem[Limit Theory for the Support Function Process]
Suppose Assumptions (ref)-(ref) hold. Then, the support function process $S_N(q) = \sqrt{N} (\widehat{\sigma} (q) - \sigma(q))$ is asymptotically linear uniformly on $ \mathcal{S}^{d-1}$
$$S_N(q) = \mathbb{G}_N [h(W,q)] + o_{P} (1)\text{ uniformly on } \mathcal{S}^{d-1},$$
where $h(W, q)$ is as in (ref). Furthermore, the process $S_N(q)$ admits the following approximation
$$
S_N(q) =_d \mathbb{G}[h(q) ] + o_P(1) \quad \text{ in } \ell^{\infty} (\mathcal{S}^{d-1}),$$
where the process $ \mathbb{G}[h(q)] $ is a tight $P$-Brownian bridge in $ \ell^{\infty} (\mathcal{S}^{d-1})$ with a non-degenerate covariance function
\begin{align*}
\Omega(q_1,q_2 ) = \mathbb{E} [h(W,q_1)h(W,q_2 )] - \mathbb{E}[h(W,q_1)] \mathbb{E}[h(W,q_2 )], \quad q_1, q_2 \in \mathcal{S}^{d-1}.
\end{align*}
Theorem (ref) is my first main result. It says that the Support Function Estimator is asymptotically equivalent to a tight Gaussian process with a non-degenerate covariance function. Previous work (e.g., CCMS, Kaido) has derived similar Gaussian approximations for the estimators based on classic nonparametric methods. Combing Neyman-orthogonality and sample splitting, Theorem (ref) allows to accommodate both classic nonparametric and modern regularized/machine learning estimators. Corollary (ref) states that the inference properties of support function estimator.
corollary[Limit Inference on Support Function Process]
Suppose Assumptions (ref)--(ref) hold. For any $f \in \mathcal{F}_c$, $\widehat{c}_N = c_N + o_P(1)$ and $c_N = O_P(1)$,
\begin{align*}
{\mathrm{P}} ( f(S_N) \leq \widehat{c}_N) - {\mathrm{P}} ( f( \mathbb{G}[h(q)]) \leq \widehat{c}_N) \rightarrow 0.
\end{align*}
If $c_N (1-\tau)$ is the $(1-\tau)$-quantile of $ f( \mathbb{G}[h(q)]) $ and $\widehat{c}_N(1-\tau) = c_N(1-\tau) + o_P(1)$ is any consistent estimate of this quantile, then
\begin{align*}
{\mathrm{P}} ( f(S_N) \leq \widehat{c}_N(1-\tau)) \rightarrow 1-\tau.
\end{align*}
theorem[Limit Theory for the Bootstrap Support Function Process]
Under conditions of Theorem (ref), the bootstrap support function process $\widetilde{S}_N(q) = \sqrt{N} (\widetilde{\sigma} (q) - \widehat{\sigma}(q))$ is asymptotically linear uniformly on $ \mathcal{S}^{d-1}$
\begin{align*}
\widetilde{S}_N(q) &= \mathbb{G}_N [ (e-1) h (W,q)]+ o_{P}(1) uniformly on \mathcal{S}^{d-1}.
\end{align*}
Furthermore, the bootstrap support function process admits an approximation conditional on the data:
$$
\widetilde{S}_N(q) = \widetilde{\mathbb{G}[h(q)]} + o_{P^e}(1) \quad \text{ in } \ell^{\infty} (\mathcal{S}^{d-1}), \text{ in probability } P,$$
where $\widetilde{\mathbb{G}[h(q)]}$ is a tight $P$-Brownian bridge in $\ell^{\infty} (\mathcal{S}^{d-1})$ with the same distribution as the process $\mathbb{G}[h(q)]$ defined in Theorem (ref), and independent of $ \mathbb{G}[h(q)] $.
Theorem (ref) is my second main result. It establishes the validity of multiplier bootstrap for uniform inference on the support function. In contrast to the weighted bootstrap of CCMS, the multiplier bootstrap does not require re-estimating the first-stage nuisance parameter in each bootstrap repetition. Instead, the nuisance parameter is estimated on an auxiliary sample once and plugged into the bootstrap sampling procedure.
corollary[Limit Inference on Bootstrap Support Function Process]
For any $\widehat{c}_N = c_N + o_P(1)$ and $c_N = O_P(1)$,
\begin{align*}
{\mathrm{P}} ( f(S_N) \leq \widehat{c}_N) - {\mathrm{P}}^e ( f(\widetilde{S}_N) \leq \widehat{c}_N) \rightarrow_P 0.
\end{align*}
If $\widehat{c}_N(1-\tau)$ is the $(1-\tau)$-quantile of $ f(\widetilde{S}_N) $ under $P^e$, then
\begin{align*}
{\mathrm{P}} (f(\widetilde{S}_N) \leq \widehat{c}_N(1-\tau)) \rightarrow_P 1-\tau.
\end{align*}
Partially Linear IV Model
\setcounter{assumption}{0}
\setcounter{lemma}{0}
\setcounter{corollary}{0}
The First Stage
In this section, I give examples of the first-stage estimators as well as the nuisance realization sets. I focus on the partially linear IV model of Example (ref). This discussion automatically covers Example (ref) that is a special case of Example (ref) with $D=Z$.
Definition (ref) introduces sequences of nuisance realization sets $\{ \mathcal{T}_N, N \geq 1\}$ and $\{ \mathcal{M}_N, N \geq 1\}$ that contain the true value of $\eta_0$ and $m_0$ and their estimates $\widehat \eta$ and $\widehat m$ with probability $1-o(1)$. As the sample size $N$ increases, the sets $\mathcal{T}_N$ and $\mathcal{M}_N$
shrink. The shrinkage speed is measured by mean square rates $\eta_N$ and $m_N$, respectively.
definition[Mean Square Rates]
Define the mean square rate for the expectation functions $\eta_0(x)$ and $m_0(x)$
\begin{align*}
\sup_{\eta \in \mathcal{T}_N} \left(\mathbb{E} \|\eta (X)- \eta_0(X)\|^2\right)^{1/2} = \eta_N \\
\sup_{m \in \mathcal{M}_N} \left(\mathbb{E} \| m (X)-m_0(X)\|^2\right)^{1/2} = m_N \end{align*}
Definition (ref) introduces mean square convergence rates for the expectation functions $\eta_0$ and $\gamma_0$. The bounds on $\eta_N$ are established for
a wide variety of regularized estimators, including partition estimators (CattaneoFarrell, CattaneoFarrellFeng), $\ell_2$-boosting (Luo), deep neural networks (Schmidt, Farrell), linear and nonlinear sieve estimators Chen2007, penalized sieve estimators Chen2011, random forest in small (WagerWalther) and high (syrgkanis2020estimation) dimensions with sparsity structure. For the sake of brevity, the examples of nuisance realization sets $ \mathcal{M}_N$ and $\mathcal{T}_N$ are omitted in this text, but they can be found in SemCher, Appendix B or CGST, Section 5.
\paragraph{First-Stage Fitted Values: Basic Case.} In this paragraph, I consider the case when the nuisance parameter $\xi_0$ does not depend on $p$.
assumption[Independent and Symmetric Residual]
The following conditions hold. (1) The interval width $Y_U - Y_L$ is independent of $V(\eta_0)$ conditional on $X$:
\begin{align}
(Y_U - Y_L) \quad \rotatebox[origin=c]{90}{$\models$} \quad V(\eta_0) \mid X.
\end{align}
(2) The vector $V(\eta_0)$ is independent of $X$, for any $p \in \mathcal{P}$. (3) $V(\eta_0)$ has a spherical distribution, which implies ${\mathrm{P}} (p' V(\eta_0) >0) = 1/2$ for any $p \in \mathcal{P}$.
Suppose Assumption (ref) holds. Define
$$
\gamma_{L,0}(X):= \mathbb{E}[ Y_L \mid X], \quad \gamma_{UL,0} (X) = \mathbb{E}[ (Y_U - Y_L) \mid X].
$$
Then, the Riesz representer function takes the form
align[align omitted — 90 chars of source]
and the first-stage fitted values take the form
align[align omitted — 147 chars of source]
where $(\widehat{\gamma}_L (X_i), \widehat{\gamma}_{UL}(X_i))_{i=1}^N$ are the cross-fit first-stage estimates of Definition (ref).
\paragraph{First-Stage Fitted Values: High-Dimensional Sparse Case.} In this paragraph, I sketch a possible estimator of the Riesz representer function without relying on Assumption (ref). Define
$$
\rho_0(p,X) = \mathbb{E}[ \mathcal{Y} (p,\eta_0) \mid X] = \mathbb{E}[ (Y_U - Y_L) 1\{ p' V(\eta_0) > 0 \} \mid X].
$$
Assume that the expectation function is approximated as
align[align omitted — 89 chars of source]
where $Z: \mathcal{X} \rightarrow \mathrm{R}^{p_X}$ is a set of $p_X$ measurable basis functions of the covariates $X$, $\nu_0(p)$ is a $p_X$-dimensional vector, and $R_p(\eta_0,X)$ is an approximation error, and $\Lambda: \mathrm{R} \rightarrow \mathrm{R}$ is the linear link function (Example (ref)) and logistic link function (Example (ref)).
example[Linear Lasso Estimator]
Let $$\mathcal{Y}_i(p, \widehat \eta):= (Y_{U,i} - Y_{L,i}) 1 \{ p' V_i(\widehat \eta) >0 \}, \quad i=1,2,\dots, N$$ be an estimate of $\mathcal{Y}_i(p, \eta_0)$. Given the penalty level
$$
\lambda_{Z} = \dfrac{1.1}{N^{1/2}} \Phi^{-1} \left(1 - \dfrac{\bar{\gamma}}{ 2 p_X N^d} \right), \quad \bar{\gamma} = .1/ \log N
$$
and the diagonal matrix of penalty loadings $\widehat \Psi_{p}$, define
$$
\widehat \nu (p) \in \arg \min_{\nu \in \mathrm{R}^{p_X}} N^{-1} \sum_{i=1}^N ( \mathcal{Y}_i(p, \widehat \eta)- Z(X_i)' \nu )^2 + \lambda_{Z} \| \widehat \Psi_{p} \nu \|_1
$$
and the fitted value as
$$
\widehat \rho(p,X_i):= Z(X_i)' \widehat \nu (p).
$$
Given the cross-fit estimate $ \widehat \gamma_L (\cdot)$ of the $\gamma_{L,0}(\cdot)$, define
\begin{align}
\widehat \gamma(p, X_i) = \widehat \gamma_L (X_i) +Z(X_i)' \widehat \nu (p).
\end{align}
An important special case occurs when the brackets have constant width
$$
Y_U - Y_L = \Delta \quad \text{a.s.},
$$
and the nuisance parameter reduces to the conditional probability
$$
\mathbb{E} [ \mathcal{M} (p, \eta_0) \mid X]:={\mathrm{P}} (p' V(\eta_0) > 0 \mid X).
$$
Suppose this probability can be approximated as
align[align omitted — 125 chars of source]
where $ \Lambda (t) = \exp t/ (\exp t+1)$ is the logistic link function.
example[Logistic Lasso Estimator]
Let $$L(y,t):= -1 (1\{ y=1 \} \log \Lambda(t) + 1\{ y=0 \} \log (1-\Lambda(t)) )$$ be the logistic loss function and let $\mathcal{M}_i(p, \widehat \eta) = 1\{ p' V_i (\widehat \eta)>0\}$ be the estimated outcome.
Given the penalty level $\lambda_Z$ and the diagonal matrix of penalty loadings $\widehat \Psi_{p}$, define
$$
\widehat \nu (p) \in \arg \min_{\nu \in \mathrm{R}^{p_X}} N^{-1} \sum_{i=1}^N L (\mathcal{M}_i(p, \widehat \eta) , Z(X_i)' \nu ) + \lambda_{Z} \| \widehat \Psi_{p} \nu \|_1
$$
The final estimator of fitted values is
\begin{align}
\widehat \gamma(p, X_i) &= \widehat \gamma_L (X_i) + \widehat \rho(p,X_i)=\widehat \gamma_L (X_i) + \Delta \Lambda( Z(X_i)' \widehat \nu(p)).
\end{align}
Lemma (ref) provides the first-stage convergence rates for linear Lasso estimator. In contrast to regular linear Lasso model, the outcome $\mathcal{Y}(p, \eta_0)$ depends on the nuisance parameter $\eta_0$ and therefore has to be estimated. I show that the linear Lasso estimator remains valid as long as the first-stage mean square rate decays sufficiently fast.
lemma[Validity of Linear Lasso with Estimated Outcome]
The following conditions hold for $N$ large enough and a sequence $\zeta_N=o(1)$ and $\Lambda(t)=t$. (i) The model (ref) is approximately sparse with $s_{\nu} =s_{\nu}(N)$
\begin{align*}
\sup_{ p \in \mathcal{P}} \| \nu_0(p) \|_0 \leq s_N
\end{align*}
and $\log (p_X \vee N) \leq \zeta_N N^{1/3}$. (ii) Heteroscedasticity. There exists a constant $c_{\zeta} >0$ so that $0<c_{\zeta} \leq \mathbb{E}[ \zeta^2(p, \eta_0) \mid X ] \text{ a.s. }$. (iii) Lipschitz property of $\nu_0(p)$. For some finite constant $C_L, \sup_{p_1, p_2 \in \mathcal{P}} \| \nu_0(p_1) - \nu_0(p_2) \|_1 \leq C_L \| p_1 - p_2 \|$. (iv) The approximation error decays fast:
$$ \sup_{p \in \mathcal{P}} (\mathbb{E} R^2_0(p,X))^{1/2} = o (\bar{\sigma}_N), \quad \bar{\sigma}_N:=\sqrt{s_N \log (p_X \vee N)/N}.$$
(vi) First Stage. (a) The first-stage mean square rate is fast enough $\eta_N = o (\bar{\sigma}_N)$ and (b) There exists $\eta^{\infty}_N = o(1)$ so that
$\sup_{\eta \in \mathcal{T}_N} \sup_{x \in \mathcal{X}} \| \eta(x) - \eta_0(x) \| \leq \eta^{\infty}_N = o(1)$. Then, under additional conditions on $Z(X)$ in Assumption (ref),
the estimate $\widehat \nu (p)$ of Example (ref) is uniformly sparse, that is $\sup_{p \in \mathcal{P}} \| \widehat \eta(p) \|_0 \leq C_X s_N$,
and the following performance bounds hold:
\begin{align}
\sup_{ p \in \mathcal{P}} \| Z(X)' (\widehat \nu (p) - \eta_0(p)) \|_{P_N, 2} &\leq C_X \sqrt{\dfrac{s_N \log p_X }{N}} \\
\sup_{ p \in \mathcal{P}} \| \widehat \nu (p) - \eta_0(p) \|_1 &\leq C_X \sqrt{\dfrac{s^2_N \log p_X }{N}}.
\end{align}
lemma[Validity of Logistic Lasso with Estimated Outcome]
Suppose the conditions of Lemma (ref) hold for (ref) with logistic link function $\Lambda(t)$. Then, under Assumption (ref),
the estimate $\widehat \nu (p)$ of Example (ref) is uniformly sparse, that is $\sup_{p \in \mathcal{P}} \| \widehat \eta(p) \|_0 \leq \widetilde{C} s_N$, and the bounds (ref)--(ref) hold.
The Second Stage
In this section, I describe the Support Function Estimator for the Partially Linear IV Model. The estimator for Example (ref) is obtained by replacing Steps 1 and 2 and 4 of Algorithm (ref) by their analogs in Algorithm (ref). As an input, Algorithm (ref) takes a direction $q \in \mathcal{S}^{d-1}$ and the first-stage fitted values $(\widehat{\eta}(X_i), \widehat{m}(X_i), \widehat{\gamma}(p,X_i))_{i=1}^N$.
algorithm[algorithm omitted — 1,132 chars of source]
algorithm[algorithm omitted — 716 chars of source]
Results
In this section, I verify Assumptions (ref)-- (ref) for Example (ref). Then, I establish asymptotic Gaussian approximation for the Support Function Estimator.
assumption[Bounded Width $Y_U - Y_L$]
The width $Y_U - Y_L$ is bounded by a finite constant $M_{UL}$ a.s., that is
$$
Y_U - Y_L \leq M_{UL} \text{ a.s. }
$$
assumption[Regularity Conditions]
(1) For every $p \in \mathcal{P}$, the residual vector $V(\eta_0) $ has a conditional density $h_{\{p' V(\eta_0) \mid X=x\}}(\cdot,x)$ that is bounded uniformly over $X$ by $M_h$. (2) The residual vector norms $\| V(\eta_0) \| =\| Z-\eta_0(X)\|$ and $ \| E(m_0) \|= \|D-m_0(X)\|$ are subGaussian random variables conditional on $X$. (3) The random variable $Y_L$ has finite conditional second moment $\sup_{x \in \mathcal{X}} \mathbb{E}[ Y^2_L \mid X=x] \leq \bar C_L$ and $\| Y_L \|_{P,4} \leq \bar C_L$ for some finite $ \bar C_L< \infty$. Finally, for some $c'>2$, $\| \| V(\eta_0) \| (Y_L - \gamma_{L,0} (X))\|_{P,c'}$ and $\| V (\eta_0) \|_{P,c'}$ are finite.
Lemma (ref) shows that the moment equation (ref) incurs only a second-order bias due to the sign mistake of $V(\eta_0)$ as long as $V(\eta_0)$ is continuously distributed.
lemma[First-Order Bias]
Under Assumptions (ref) and (ref) (1), the first-order bias shrinks at the quadratic speed
\begin{align*}
\sup_{p \in \mathcal{P}} | \mathbb{E} [ z(p, \eta_0) ( Y(p, \eta) - Y(p, \eta_0) )] | \leq 2 C_P^2 M_{UL} M_h \mathbb{E} \| \eta(X) - \eta_0(X) \|^2.
\end{align*}
Likewise, the second-order bias shrinks at the quadratic speed
\begin{align*}
\sup_{p \in \mathcal{P}} | \mathbb{E} [ (z(p, \eta) - z(p, \eta_0)) ( Y(p, \eta) - Y(p, \eta_0) )] | \leq 2 C_P^2 M_{UL} M_h \mathbb{E} \| \eta(X) - \eta_0(X) \|^2.
\end{align*}
lemma[Verification of Assumption (ref)]
Suppose Assumption (ref)(2) holds. For any sequence $\ell_N \rightarrow \infty$, Assumption (ref) is satisfied with $v_N = \sqrt{ \ell_N /N}$, in particular, one can take $\ell_N = \log N$ to ensure that $v_N = o(N^{-1/2} \log N)$.
Lemma (ref) verifies the Assumption (ref) for Example (ref). Define the mean square rates for the expectation functions
align*[align* omitted — 270 chars of source]
lemma[Verification of Assumptions (ref)--(ref)]
Let $\eta_0(X)$ be as in (ref), $V(\eta)$ be as in (ref), the matrix function $A(W, \eta, m)$ be as in (ref) and the orthogonal moment function be as in (ref). Suppose Assumptions (ref) and (ref) and (ref) hold. Furthermore, suppose (1) the elements of $\Gamma_{L,N}$ are bounded by finite constant $M_{\gamma}$ in the sup-norm: $\sup_{ \gamma_L \in \Gamma_{L,N} } \sup_{x \in \mathcal{X}} | \gamma_L(x) | \leq M_{\gamma}$ and $\sup_{ \gamma_{UL} \in \Gamma_{UL,N} } \sup_{x \in \mathcal{X}} | \gamma_{UL}(x) | \leq M_{UL\gamma}$ (2) the elements of $\mathcal{M}_N$ and $\mathcal{T}_N$ are bounded by finite constant $M_{\eta}< \infty$ in the sup-norm: $$\sup_{ m \in \mathcal{M}_N } \sup_{x \in \mathcal{X}} \| m(x) \| \leq M_{\eta} , \quad \sup_{ \eta \in \mathcal{T}_N } \sup_{x \in \mathcal{X}} \| \eta(x) \| \leq M_{\eta}. $$ Then, Assumption (ref)(1) holds. Furthermore, the bias rates in Definition (ref) can be bounded as follows with $\gamma^2_N := 2(\gamma^2_{L,N} + \gamma^2_{UL, N})$
\begin{align*}
\mu_N &\leq 4 C^2_P M_{UL} M_h \eta^2_N + C_P \eta_N \cdot \gamma_N \\
A_N &\leq \eta_N \cdot m_N.
\end{align*}
In addition, for $\tau_N := N^{-1/2} \log N$, the sequences $r_N'', r_N', \delta_N$ obey
$$
r_N'' = O (\eta_N + \gamma_N), \quad r_N' = O ( N^{-1/4} \log^{1/2} N), \quad \delta_N = O (\eta_N + m_N).
$$
Thus, if $\eta_N = o(N^{-1/4})$ and $ (m_N + \gamma_N) \cdot \eta_N = o(N^{-1/2})$ and $m_N = o(1)$ and $\gamma_N = o(1)$, Assumption (ref) holds.
Combining the statements in Lemmas (ref)--(ref), I obtain the following corollary.
corollary[Asymptotic Theory for Partially Linear IV Model with Interval-Valued Outcome]
Suppose Assumptions (ref), (ref) and (ref) and the conditions of Lemma (ref) hold with the fast enough first-stage rates $\eta_N, m_N$ and $\gamma_N:= 2(\gamma_{L,N} + \gamma_{UL, N})$ such that Assumption (ref) holds. Then, Theorems (ref) and (ref) and Corollaries (ref) and (ref) hold for the Support Function Estimator of the Algorithm (ref) with the first-stage fitted values (ref) and with the influence function $h(W,q)$ equal to (ref).
corollary[Asymptotic Theory for Partially Linear IV Model with Interval-Valued Outcome]
Suppose Assumptions (ref), (ref) and the conditions of Lemma (ref) and Lemma (ref) hold with the fast enough first-stage rates $\eta_N, m_N$ and $\gamma_N:= \gamma_{L,N}$ such that Assumption (ref) holds. Then, Theorems (ref) and (ref) and Corollaries (ref) and (ref) hold for the Support Function Estimator of the Algorithm (ref) with the influence function $h(W,q)$ equal to (ref), where the first-stage fitted values are as in (ref) (assuming Lemma (ref) holds with $\bar{\sigma}_N = o(N^{-1/4})$ ) or as in (ref) (assuming Lemma (ref) holds with $\bar{\sigma}_N = o(N^{-1/4})$).
remarkConsider Example (ref). Note that the short least squares regression (ref) uses only $d$ out of (infinitely) many restrictions implied by the exogeneity restriction (ref). Therefore, the identified set
$$
\mathcal{B}_1:=\bigg\{ b \in \mathrm{R}^d: \exists f_b \in \mathcal{L} \text{ and } Y \in [Y_L, Y_U]: \quad \mathbb{E}[ Y - D b - f_b(X) \mid D, X] =0 \bigg\}
$$
is a convex subset of $\mathcal{B}$, but may not coincide with $\mathcal{B}$.
Average Partial Derivative
In this section, I present the identification, estimation and inference results for Example (ref).
\setcounter{assumption}{0}
\setcounter{lemma}{0}
\setcounter{corollary}{0}
\setcounter{example}{0}
\setcounter{definition}{0}
Assumption (ref) states the sufficient conditions for compactness and convexity of the identified set $\mathcal{B}$.
assumption[Regularity Conditions for Average Partial Derivative]
The following conditions hold. (1) There exists a compact and convex set $\mathcal{D}$ with nonempty interior containing the support of $D$, such that $f_0 (d \mid X) = 0$ on the boundary of $\mathcal{D}$ a.s. in $X$. (2) The random variables $D$ and $f_0(D \mid X)$ and $\nabla f_0 (D \mid X)$ have the density conditional on $X$. (3) The random vector $$\eta_0(D,X) = \nabla f_0 (D \mid X)/ f_0 (D \mid X)$$ is $L_{P,2}$-integrable, that is, there exists a finite constant $C_{\text{APD}}< \infty$ such that
$$
\| \| V(\eta_0) \| \|_{P,2} := \| \| \nabla f_0 (D \mid X)/ f_0 (D \mid X) \| \|_{P,2} \leq C_{\text{APD}}.
$$
(4) The functions $\gamma_{L,0}(D,X)$ and $\gamma_{UL,0}(D, X)$ are bounded by some finite constant $C_{\gamma}< \infty$ on the support of $D$ and $X$. (5) Furthermore, $\gamma_{L,0}(D,X)$ and $\gamma_{UL,0}(D, X)$ are continuously differentiable w.r.t to $D$ with bounded derivatives a.s. in $D$ and $X$. (5) The random variables $\| V(\eta_0) \|$ and $
\| V(\eta_0) \| (Y_L - \gamma_{L,0}(D,X))
$ are $L_{P,2}$-integrable and $L_{P,c'}$-integrable for $c'>2$.
lemmaSuppose Assumption (ref) (1)-(4) hold. Then, (a) the identified set $\mathcal{B}$ for $\beta_0$ is compact and convex and (b) the support function of $\mathcal{B}$ is given by (ref).
In addition, if Assumption (ref) holds, which implies that $\mathcal{B}$ is strictly convex.
Lemma (ref) extends the Theorem 2.1 of Kaido to allow the density $D$ to be conditioned on $X$.
assumption[Margin Condition]
(1) There exists an absolute constant $\bar C_f <\infty$ so that, in some neighborhood $(0, \bar t)$ of zero,
\begin{align*}
\sup_{q \in \mathcal{S}^{d-1}} {\mathrm{P}} \left( |q' \nabla_D f_0(D \mid X)/ f_0 (D \mid X) | \leq \delta \right) \leq \bar C_f \delta , \quad \delta \in (0, \bar{t}).
\end{align*}
The margin condition is commonly used in classification analysis (MammenTsybakov, Tsybakov) and empirical welfare maximization (KitagawaTetenov, MbakopTabord) and bounds SemSupp2. This paper is the first one to introduce it for the study of average partial derivative with an interval-valued variable.
definition[Worst-Case Rates]
Let $f_0(D \mid X)$ and $\nabla_D f_0(D \mid X)$ be the conditional density and its derivative. Let $\{ F_N, \quad N \geq 1 \}$ and $\{ F^1_N, \quad N \geq 1 \}$ be sequences of realization sets of the estimates of $f_0 (D \mid X)$ and $\nabla_D f_0 (D \mid X)$, respectively. Assume that these sequences shrink at the following worst-case rates $f^{\infty}_N$ and $f^{\infty}_{1,N}$
\begin{align*}
\sup_{ f \in F_N } \sup_{d, x} | f (d \mid x) - f_0( d \mid x) | \leq f^{\infty}_N \\
\sup_{ f \in F^1_N } \sup_{d, x} \| \nabla_{d} f (d \mid x) - \nabla_d f_0( d \mid x) \| \leq f^{\infty}_{1,N}
\end{align*}
Define
$$
\eta^{\infty}_N:= f^{\infty}_{1,N} + f^{\infty}_N.
$$
Furthermore, assume that the elements of $F_N$ are bounded by some $\bar B_f< \infty$:
$$
\sup_{ f \in F_N } \sup_{d, x} \{ | f^{-1} (d \mid x ) |, | f (d \mid x ) |\} \leq \bar B_f.
$$
and
$$
\sup_{ f \in F^1_N } \sup_{d, x} \| \nabla_{d} f (d \mid x) \| \leq \bar B_f.
$$
The Algorithm
algorithm[algorithm omitted — 1,082 chars of source]
Results
Lemma (ref) shows that the moment equation (ref) incurs only a second-order bias due to the sign mistake of $V(\eta_0)$ as long as $V(\eta_0)$ is continuously distributed. In contrast to the setup of Lemma (ref), the estimation error is not orthogonal to the space of $V(\eta_0)$. As a result, the first-order bias is bounded by invoking the margin assumption and the $\ell_{\infty}$-rate.
lemma[First-Order Bias]
Suppose Assumptions (ref) and (ref) holds. Then, the first-order bias shrinks at quadratic speed, that is, for some constant $\bar C$ large enough,
\begin{align*}
\sup_{ \eta \in \mathcal{T}_N } | \mathbb{E} [z(q, \eta_0) (Y(q, \eta) - Y(q, \eta_0)) ] | = O ((\eta^{\infty}_N)^2).
\end{align*}
Likewise, the second-order bias shrinks at quadratic speed
\begin{align*}
\sup_{ \eta \in \mathcal{T}_N } |\mathbb{E} [ (z(q, \eta) - z(q, \eta_0)) (Y(q, \eta) - Y(q, \eta_0)) ] | = O ((\eta^{\infty}_N)^2).
\end{align*}
definition[Mean Square Rates]
Define the mean square rates for the expectation functions
\begin{align*}
\sup_{\gamma_L \in \Gamma_{L,N}} \left(\mathbb{E} (\gamma_{L} (D,X) - \gamma_{L,0}(D,X))^2 \right)^{1/2} =: \gamma_{L,N} \\
\sup_{\gamma_{UL} \in \Gamma_{UL,N}} \left(\mathbb{E} (\gamma_{UL} (D,X) - \gamma_{UL,0}(D,X))^2 \right)^{1/2} =: \gamma_{UL,N}
\end{align*}
and for their derivatives
\begin{align*}
\sup_{\gamma_L \in \Gamma_{L,N}} \left(\mathbb{E} \| \nabla_D \gamma_{L} (X) -\nabla_D \gamma_{L,0}(X) \|^2 \right)^{1/2} =: \gamma^1_{L,N} \\
\sup_{\gamma_{UL} \in \Gamma_{UL,N}} \left(\mathbb{E} \| \nabla_D \gamma_{UL} (X) -\nabla_D \gamma_{UL,0}(X) \|^2 \right)^{1/2} =: \gamma^1_{UL,N}
\end{align*}
and let $$
\gamma_N := \gamma_{L,N} + \gamma_{UL,N} + \gamma^1_{L,N} +\gamma^1_{UL,N}.
$$
corollary[Asymptotic Theory for Average Partial Derivative with an Interval-Valued Outcome]
Suppose Assumptions (ref) and (ref) and (ref) hold. In addition, suppose $| \gamma_{L} (D,X) | + \| \nabla_D \gamma_{L} (D,X) \| + | \gamma_{UL} (D,X) | + \| \nabla_D \gamma_{UL} (D,X) \|< M_{\gamma} \text{ a.s.}$. Then, the sequences $\mu_N$ and $r_N'$ can be bounded as follows
$$
\mu_N = O ( (\eta^{\infty}_N)^{3/2}+ \eta^{\infty}_N \cdot \gamma_N).
$$
and $r_N' = O (\gamma_N + (\eta^{\infty}_N)^{1/2} ) $. In particular, if $\eta^{\infty}_N = o (N^{-1/3})$ and $\eta^{\infty}_N \cdot \gamma_N = o (N^{-1/2})$ and $\gamma_N = o(1)$, Assumption (ref) holds. Assumption (ref) holds automatically since $A(W, \eta_0) = \Sigma = I_d$ is a known matrix. Finally, Assumption (ref) holds. Then, Theorems (ref) and (ref) hold for the Support Function Estimator of Algorithm (ref) and the influence function equal to
\begin{align*}
h(W,q) = g(W,q,\xi(q)) - \mathbb{E} [ g(W,q,\xi(q)) ],
\end{align*}
where $g(W,q,\xi(q))$ is given in ((ref)).
Simulation Study
In this section, I compare the performance of the classic (series-based non-orthogonal) and the proposed (lasso-based orthogonal) approaches in a moderate-dimensional sparse design. The first approach is to plug a least squares series first-stage estimator into the non-orthogonal moment equation (ref). The second one is to plug an $\ell_1$-regularized least squares series into the orthogonal moment equation (ref). Regularization helps to leverage the sparsity assumption and to reduce risk.
Consider the partially linear model of Example (ref) with $d=2$ treatments. The outcome equation is generated as
align[align omitted — 97 chars of source]
where $D=(D_1, D_2)$ is the treatment vector, $X \in \mathrm{R}^{p_X}$ is the $p_X$-covariate vector with $p_X=50$, and $U \sim N(0, \sigma_U^2)$ is a normal shock independent of $D$ and $X$ with $\sigma_U =1$. For each treatment $m=1,2$, its reduced form is
align[align omitted — 83 chars of source]
where the coefficients $$\theta_0 = \alpha_1 = \alpha_2 = (1, 1/2^2, \dots, 1/j^2, \dots, 1/p_X^2). $$ The parameters in (ref) and (ref) are chosen as
align*[align* omitted — 79 chars of source]
The covariates $X$ are generated from $N(0, \Omega)$, where $\Omega$ is a Toeplitz matrix with correlation coefficient $\rho = 0.5$. That is, for every $(i,j) \in \{1,2,\dots, p_X\}^2$, $\Omega_{ij} = \rho^{|i-j|}$. The vector $V$ is independent of $(X,U)$ and is drawn from the bivariate normal distribution
align[align omitted — 129 chars of source]
The outcome $Y$ is not included into the data. Instead, the support of $Y$ is partitioned into the bins $ \cup_{s=1}^{S} [b_s, b_{s+1})$ of width $\Delta$:
$$
b_{s+1} = b_s + \Delta, \quad s=1, 2, \dots, S-1,
$$
where $b_0 = -\infty$ and $b_{S+1} = \infty$. The observed bounds $Y_L$ and $Y_U$ are taken to be
align*[align* omitted — 94 chars of source]
Thus, the observed data vector $W=(X, D, Y_L, Y_U)$ but does not contain $Y$.
I now derive the true (population) support function. Plugging $Y_U - Y_L = \Delta$ into (ref) gives
$$
\sigma(q) = q' \Sigma^{-1} \mathbb{E} V(\eta_0) Y(q, \eta_0) = q' \Sigma^{-1} \mathbb{E} V Y_L + \Delta \mathbb{E} \max (q' \Sigma^{-1} V,0).
$$
Invoking (ref) gives $q' \Sigma^{-1} V \sim N(0, q' \Sigma^{-1} \Sigma \Sigma^{-1} q) = \sqrt{q' \Sigma^{-1} q} N(0, 1)$, which implies
$$
\mathbb{E} \max (q' V, 0) = \sqrt{q' \Sigma^{-1} q/2\pi}.
$$
Thus, the support function is
align[align omitted — 146 chars of source]
The classic and the proposed approaches are implemented via Algorithm (ref) with different sets of the first-stage fitted values. For each treatment $m \in \{1,2\}$, the first-stage regression parameter $\alpha$ is estimated as
align[align omitted — 163 chars of source]
where $\lambda_{D_m} = 0$ in the series-based case and $\lambda_{D_m} > 0$ in the lasso-based case. For the lasso estimator, the penalty parameter $\lambda_{D_m}$ is chosen according to Algorithm 1 in Program (i.e., the default value of \url{rlasso} package). In both cases, the treatment fitted values are
align*[align* omitted — 78 chars of source]
In the orthogonal case (the proposed approach), the best-case outcome fitted values are
align*[align* omitted — 82 chars of source]
where $\gamma_L$ is estimated by Lasso regression of $Y_L$ on $X$ similarly to (ref). Since the classic case does not require partialling out,
the fitted values $\widehat{\gamma}_U(X)$ are set to zero.
The estimator's performance is summarized in terms of its risk and coverage. The total risk of the estimator is defined as
align[align omitted — 104 chars of source]
It coincides with the Hausdorff distance $d(\widehat{\mathcal{B}},\mathcal{B})$ between the estimated ($\widehat{\mathcal{B}}$) and the true ($\mathcal{B}$) sets. Furthermore, the outer and the inner risks are defined as
align[align omitted — 201 chars of source]
The rejection frequency is the share of rejected simulation draws
align[align omitted — 96 chars of source]
where $c^{*}_{1- \alpha} $ is the $(1-\alpha)$-quantile of $R_H^b$ of the bootstrap process
align[align omitted — 122 chars of source]
Table (ref) compares the classic and the proposed estimators in terms of risk and coverage. Across the board, the proposed estimator has a smaller total risk than the classic one, by a factor of ranging from $2.0$ for $\Delta=1$ to $3.0$ for $\Delta=3$. The classic estimator has higher risk due to the excessive noise of series estimators $\widehat{\alpha}_1$ and $\widehat{\alpha}_2$ in a regime with $p_X =50$ covariates.
I investigate the importance of the first-stage regularization and second-stage orthogonalization, applying one at a time. Table (ref) compares the non-orthogonal (Columns (1)-(4)) and the orthogonal (Columns (5)-(8)) estimators based on the true treatment first stage. Across the board, the orthogonal estimator has smaller total risk than the non-orthogonal one, by a factor of $2.5$ on average. The variance reduction occurs due to partialling out the relevant controls from the outcome $Y^{\text{best}}(\eta_0)$. In contrast, the risks of series-based non-orthogonal (Table (ref), Column (3)) and orthogonal (Table (ref), Column (7)) estimators are close to each other. In particular, the high risk of the series-based estimator cannot be improved by orthogonalization.
Next, I compare the ortho (Table (ref), Columns 5--8) and non-ortho (Table (ref), Columns 1--4) estimators based on the lasso-based first stage.
Across the board, the risk of the ortho version is substantially smaller than the non-ortho one, by a factor ranging from $2$ to $5.5$. The non-ortho version has a higher risk, because the estimates of non-zero coefficients in $\alpha_m$ are shrunk to zero. While the shrinkage bias could be reduced by invoking the post-lasso instead of the lasso estimator, post-single-selection inference may not be robust to moderate deviations from zero.
Empirical illustration
This section demonstrates the proposed approach by estimating the gender wage gap with a bracketed wage variable. First, I show that a frequent empirical practice -- midpoint regression -- gives biased results. Instead of this approach, I propose reporting identified set (i.e., the lower and the upper bound) for the parameter of interest and demonstrate how to estimate the set.
The sample for the analysis comes from the U.S. March Supplement of the Current Population Survey (CPS) in 2015, as studied in MulliganRubinstein and CCMS. The selected sample consists of white non-hispanic individuals, aged 25 to 64 years, working more than 35 hours per week during at least 50 weeks of the year. The resulting sample comprises $32,523$ workers, including $18,137$ men and $14,386$ women. The object of interest is the gender wage gap -- the average difference in log wages between men and women after controlling for the observed characteristics. The data set is augmented by the lower bound $Y_L$ and the upper bound $Y_U$ defined as
align*[align* omitted — 94 chars of source]
where $b_0=1$ and $$ b_{s+1} - b_s = \Delta \quad s=1,2,\dots, S.$$
I assume that the log wage variable follows the partially linear regression of Example (ref)
$$
Y = D' \beta_0 + X' \gamma_0+ U, \quad \mathbb{E}[ U \mid X, D] =0.
$$
Here, the outcome variable $Y$ is the logarithm of the hourly wage rate, the treatment/policy variable $D$ is an indicator for female gender, and the vector of $p_X=260$ controls includes demographic indicators, region, and experience indicators, as well as their interactions. The coefficients $\beta_0$ and $\gamma_0$ are the target and the nuisance parameters, respectively. The parameter $\beta_0$ is identified as a minimizer of the least squares loss function
$$
\beta_0:=\arg \min_{b \in \mathrm{R}} \mathbb{E} (Y - (D- \eta_0(X))b)^2.
$$
When $Y$ is unobserved, a frequent approach is to replace $Y$ by the bracket midpoint
$$Y_M:=(Y_L+Y_U)/2.$$
The “mid-point” regression parameter is taken to be
$$
\beta_{M}:=\arg \min_{b \in \mathrm{R}} \mathbb{E} (Y_M - (D- \eta_0(X))b)^2,
$$
which can differ from $\beta_0$. I report the estimate and the $95 \%$ CI of $\beta_0$ and $\beta_M$. The first-stage treatment and the outcome expectation functions are estimated via logistic and linear Lasso regression of Program implemented in the hdm $R$ package, respectively. In addition, I also consider the random forest estimator as implemented in the ranger $R$ package. The final estimator is taken to be the Double Machine Learning estimator of chernozhukov2016double with $K=2$-fold cross-fitting. The Lasso-based first-stage fitted values (Columns (1)-(2)) and the random-forest-based (Columns (3)-(4)).
Instead of reporting $\beta_M$ -- which is a biased measure of $\beta_0$ -- I propose reporting the lower and the upper bound on $\beta_0$. The moment equation for $\beta_U$ is given in (ref) and its estimate is defined in the Algorithm (ref) with the fitted values described below. The symmetry-based specification (SYM) imposes Assumption (ref), and the first-stage fitted values are taken to be (ref). The sparsity-based specification (SPRS) imposes the model (ref), and the first-stage are taken to be as in Example (ref), (ref). The treatment expectation function $\eta_0(\cdot)$ is estimated the same as in the point-identified specifications.
Table (ref) summarizes the findings for the bracket width $\Delta \in \{1,2,3\}$. The “ground-truth” gender wage gap ranges between $18 \%$ (Columns (3)-(4)) and $20 \%$ (Columns (1)-(2)). I call this estimate “ground-truth” because this is the estimate to be reported if the wage $Y$ was observed. When the bracket width $\Delta=1$ is small, the midpoint estimate $\beta_M$ is close to the estimate of $\beta_0$. When $\Delta=2$, the midpoint estimate $\beta_M$ ranges between $25\%$ (RF) and $28 \%$ (Lasso). Furthermore, the 95 $\%$ CI for $\beta_M$ does not contain the “ground-truth“ estimate of $\beta_0$ for either RF or Lasso first-stage method. When $\Delta=3$, the magnitude of the bias remains substantial, which speaks against midpoint regression for large values of bracket width.
To estimate $[\beta_L, \beta_U]$, I consider four specifications: SYM-Lasso, SYM-RF, SPRS-Lasso and SPRS-RF. Across the board, the bounds contain $\beta_0$. Furthermore, the specifications yield close results despite being based on very different starting assumptions. That said, the bounds $[\beta_L, \beta_U]$ come out wide. As discussed in Remark (ref), the bounds $[\beta_L, \beta_U]$ may not be sharp for $\beta_0$, since they are utilize only $d$ (out of infinitely many) moment restriction implied by (ref). The derivation of sharp bounds for $\beta_0$ is left for the future work.
table[table omitted — 2,132 chars of source]
table[table omitted — 2,450 chars of source]
Proofs
\paragraph{Empirical process notation. } Let $\widehat{\eta}_k$ and $\widehat{\xi}_k, \quad k=1,2,\dots, K$ be as in Definition (ref). Define an event
align*[align* omitted — 135 chars of source]
By union bound, this event holds with probability approaching one $${\mathrm{P}} (\mathcal{E}_N ) \geq 1- K \epsilon_N = 1-o(1).$$
For a given partition $k$ in $\{1,2, \dots, K\}$, define the partition-specific averages
align*[align* omitted — 182 chars of source]
Define the function $\psi_0(p)$
align[align omitted — 90 chars of source]
and observe that plugging $p_0(q)$ into $\psi_0(p)$ gives the support function $ \sigma(q) $ at $q$:
$$ \psi_0(p_0(q))= \psi_0(\Sigma^{-1}q)=\sigma(q).$$
For $i \in J_k$ and $k=1,2,\dots, K$, define the partition-specific conditional expectation
align[align omitted — 145 chars of source]
and its weighted sample analog
align[align omitted — 128 chars of source]
Finally, define the matrix error terms
align*[align* omitted — 203 chars of source]
and the weighted matrix error
align*[align* omitted — 214 chars of source]
\paragraph{Empirical process remainder terms.} Define the remainder term
align*[align* omitted — 176 chars of source]
the bias term
align*[align* omitted — 99 chars of source]
the second-order remainder term
align*[align* omitted — 78 chars of source]
and the bootstrap term
align*[align* omitted — 182 chars of source]
Thus, the support function process $S_N(q)$ of Theorem (ref) can be decomposed as
align*[align* omitted — 438 chars of source]
\paragraph{Misclassification events.} For $p, p_0 \in \mathcal{P}$, define the events $\mathcal{E}_{+}(p), \mathcal{E}_{-}(p), \mathcal{E}_{-}(p,p_0), \mathcal{E}_{+}(p,p_0) $
align[align omitted — 199 chars of source]
and
align[align omitted — 207 chars of source]
Proofs of Main Results.
proof[Proof of Lemma (ref)]
Define
\begin{align*}
z(p, \eta):&= p' V(\eta) \\
B_1 (W, \eta, p) :&= p' V ( \eta_0) ( Y(p, \eta) - Y(p, \eta_0) ) \\
B_2 (W, \eta, p) :&= p' ( V (\eta) - V ( \eta_0)) ( Y(p, \eta) - Y(p, \eta_0) ).
\end{align*}
Observe that
\begin{align*}
z(p, \eta) ( Y(p, \eta) - Y(p, \eta_0)) = B_1 (W, \eta, p) + B_2 (W, \eta, p).
\end{align*}
The mistake in $Y(p, \eta)$ can only occur if $p' V ( \eta_0)$ is small enough
\begin{align}
\bigg\{ Y(p, \eta) \neq Y(p, \eta_0) \bigg\} &\Leftrightarrow \bigg\{ \mathcal{E}_{+}(p) or \mathcal{E}_{-}(p)\bigg\} \nonumber \\
&\Rightarrow \bigg\{ 0 < | p' V ( \eta_0) | < | p' (V ( \eta) - V (\eta_0))| \bigg\} \nonumber \\
&\Rightarrow \bigg\{ 0 < | p' V ( \eta_0) | < C_P \| V ( \eta) - V (\eta_0) \| \bigg\}=: \mathcal{E}_{+-}(p) .
\end{align}
Recall that
\begin{align}
Y(p, \eta) -Y(p, \eta_0)=(Y_U - Y_L) 1\{ \mathcal{E}_{+}(p) \cup \mathcal{E}_{-}(p)\}
\end{align}
Invoking $| \mathbb{E} X | \leq \mathbb{E} | X| $ and (ref) gives
\begin{align}
| \mathbb{E} B_1 (W, \eta, p) |&= | \mathbb{E} p' V ( \eta_0) (Y(p, \eta) -Y(p, \eta_0) ) | \nonumber \\
&\leq \mathbb{E} | p' V ( \eta_0) | (Y_U - Y_L) 1\{ \mathcal{E}_{+}(p) \cup \mathcal{E}_{-}(p)\} \nonumber \\
&\leq \mathbb{E} | p' V ( \eta_0) | (Y_U - Y_L) 1\{ \mathcal{E}_{+-}(p) \}
\end{align}
Invoking definition of $\mathcal{E}_{+-}(p)$ in (ref) and $Y_U - Y_L \leq M_{UL} \text{ a.s. }$ gives
\begin{align}
&\mathbb{E} | p' V ( \eta_0) | (Y_U - Y_L) 1\{ \mathcal{E}_{+-}(p) \} \nonumber \\
&\leq C_P M_{UL} \mathbb{E} \| V ( \eta) - V (\eta_0) \| 1\{ \mathcal{E}_{+-}(p) \}.
\end{align}
The second-order bias term is bounded as
\begin{align}
| \mathbb{E} B_2 (W, \eta, p) | \leq C_P M_{UL} \mathbb{E} \| V ( \eta) - V (\eta_0)\| 1\{ \mathcal{E}_{+}(p) \cup \mathcal{E}_{-}(p)\}.
\end{align}
Invoking Assumption (ref) gives
\begin{align}
&\mathbb{E} \| \eta_0(X) - \eta(X) \| 1\{ \mathcal{E}_{+}(p) \cup \mathcal{E}_{-}(p)\} \nonumber \\
&=\mathbb{E}_{X} \| \eta(X) - \eta_0(X) \| \int_{- C_P \| \eta(X) - \eta_0(X) \| }^{C_P \| \eta(X) - \eta_0(X) \| }h_{p' V(\eta_0) \mid X} (t, X) dt \nonumber \\
&\leq 2 C_P M_h \mathbb{E}_{X} \| \eta(X) - \eta_0(X) \|^2 .
\end{align}
Combining the bounds gives
\begin{align}
| \mathbb{E} [B_1 (W, \eta, p) + B_2 (W, \eta, p)] | &\leq |\mathbb{E} B_1 (W, \eta, p) | + | \mathbb{E} B_2 (W, \eta, p) | \nonumber \\
&\leq 4 C_P^2 M_{UL} M_{h} \mathbb{E} \| \eta(X) - \eta_0(X) \|^2 .
\end{align}
proof[Proof of Theorem (ref)]
Step 1. This step is required only if $\Sigma$ is unknown, such as in Example (ref). As shown in Lemma (ref), for $v=1$ (regular case) and $v=e$ (bootstrap case),
\begin{align}
(\widehat{\Sigma}^v (\widehat{\eta}))^{-1} - \Sigma^{-1} = -\Sigma^{-1} ( \widehat{\Sigma}^v (\eta_0) - \Sigma ) \Sigma^{-1} + M^v,
\end{align}
where the remainder matrix $M^v$ obeys $\| M^v \| = o_P (N^{-1/2})$ for both cases. Take $v=1$. Post-multiplying the LHS above by $q$ gives
\begin{align}
\sqrt{N} (\widehat{p}(q) - p_0(q)) &= - \Sigma^{-1} ( \widehat{\Sigma} (\eta_0) - \Sigma ) \Sigma^{-1}q + o_P(1) \nonumber
\end{align}
Likewise, taking $v=e$ gives
\begin{align*}
\sqrt{N} (\widetilde{p}(q) - p_0(q)) &= ((\widehat{\Sigma}^e (\widehat{\eta}))^{-1} - \Sigma^{-1})'q \\
&= - \Sigma^{-1} ( \widetilde{\Sigma}^e (\eta_0) - \Sigma ) \Sigma^{-1}q + o_P(1).
\end{align*}
For some $N$ large enough, $ \widehat{p}(q) \in \mathcal{P} \quad \forall q \in \mathcal{S}^{d-1}$ and $ \widetilde{p}(q) \in \mathcal{P} \quad \forall q \in \mathcal{S}^{d-1}$ with probability $1-o(1)$.
Step 2. Let $k=1,2,\dots, K$ denote the partition index. We bound $R_{1,k} (\widehat{p}(q))$ and $R^e_{1,k} (\widehat{p}(q))$. Define the function class
\begin{align*}
\mathcal{F}_{2k}^v := \{ v\cdot(g(\cdot,p,\widehat{\xi}_k(p)) - g(\cdot,p_0,\xi_0(p_0))), \quad p , p_0 \in \mathcal{P}, \quad \| p - p_0 \| \leq \tau_N \},
\end{align*}
where $v=1$ (regular case) and $v=e$ (bootstrap case). The class $\mathcal{F}_{2k}$ is obtained as $$\mathcal{F}_{2k} \subset \mathcal{G}_{\widehat{\xi}_k} - \mathcal{G}_{\xi_0},$$ where $\mathcal{G}_{\xi_0}$ and $\mathcal{G}_{\widehat{\xi}_k} $ are defined in Assumption (ref). On the event $\mathcal{E}_N$,
\begin{align*}
&\sup_{p \in \mathcal{P} } |g(W_i,p,\widehat{\xi}_k(p)) - g(W_i,p_0,\xi_0(p_0)) | \\
&\leq \sup_{p \in \mathcal{P}} | g(W_i,p,\widehat{\xi}_k(p)) | + \sup_{p \in \mathcal{P}} | g(W_i,p,\xi_0(p)) | \\
&\leq G_{\widehat{\xi}_k} + G_{\xi_0},
\end{align*}
and $G_{\widehat{\xi}_k \xi_0} := G_{\widehat{\xi}_k} + G_{\xi_0}$ is a measurable envelope for the class $\mathcal{F}_{2k}$. Note that $\| G_{\widehat{\xi}_k \xi_0} \|_{P,c} \leq \| G_{\xi_0} \|_{P,c} + \| G_{\widehat{\xi}_k } \|_{P,c} \leq 2 C_1$. The uniform covering entropy of the function class $ \mathcal{G}_{\widehat{\xi}_k} - \mathcal{G}_{\xi_0}$ is bounded as
\begin{align*}
& \log \sup_{Q} N(\epsilon \| G_{\xi} + G_{\xi_0} \|_{Q,2}, \mathcal{G}_{\widehat \xi} - \mathcal{G}_{\xi_0} , \| \cdot \|_{Q,2}) \\
&\leq \log \sup_{Q} N(\epsilon/2 \| G_{\xi} \|_{Q,2}, \mathcal{G}_{\xi} , \| \cdot \|_{Q,2}) + \log \sup_{Q} N(\epsilon/2 \| G_{\xi_0} \|_{Q,2}, \mathcal{G}_{\xi_0} , \| \cdot \|_{Q,2}) \\
&\leq 2v \log (2a/\epsilon)
\end{align*}
by the proof of Theorem 3 in andrews:1994b and Assumption (ref). The class $\mathcal{F}_{2k}^e $ is obtained by multiplication of $\mathcal{F}_{2k}$ by an integrable random variable independent of the data, and therefore retains $P$-Donsker and uniform covering properties of the class $\mathcal{F}_{2}$. In particular,
$G^v_{\widehat{\xi}_k \xi_0} := |v| (G_{\widehat{\xi}_k} + G_{\xi_0})$ is a valid envelope for $ \mathcal{F}_{2k}^v $. Next, for $v=1$ and $v=e$,
\begin{align*}
&\sup_{\xi \in \Xi_N} \mathbb{E} v^2 (g(W,p,\xi(p)) - g(W,p_0,\xi_0(p_0)) )^2 \\
&\leq 2 \bigg( \sup_{\xi \in \Xi_N} \sup_{p \in \mathcal{P}} \mathbb{E} v^2 ( g(W,p,\xi(p)) - g(W,p,\xi_0(p)))^2 \\
&+ \sup_{p_0, p \in \mathcal{P}, \quad \| p - p_0 \| \leq \tau_N} \mathbb{E} v^2 (g(W,p,\xi_0(p)) - g(W,p_0,\xi_0(p_0)) )^2 \bigg) \\
&\leq 2((r_N”)^2 + (r_N')^2).
\end{align*}
Invoking Lemma (ref) conditional on $(W_i)_{i \in J_k}$ and taking $\sigma = 2 (r_N' + r_N'')$ gives
\begin{align*}
\sup_{f \in \mathcal{F}_2^v} | {\mathbb{G}_{n,k}} [f] | &= \sup_{ p_0, p \in \mathcal{P}} | \widehat{\psi}_k (p, \widehat{\xi}_k) - \widehat{\psi}_k(p_0, \xi_0) - (\psi(p, \widehat{\xi}_k) - \psi(p_0, \xi_0)) | \\
&\lesssim_P (r_N”+r_N' )\log^{1/2} (1/(r_N”+r_N' )) + N^{-1/2+1/c} \log N \\
&= o_P(1),
\end{align*}
where the last equality follows from Assumption (ref). By Lemma (ref), $\sup_{f \in \mathcal{F}_2^v} | {\mathbb{G}_{n,k}} [f] | = o_P(1)$ holds unconditionally. as shown in Step 1, wp $1-o(1)$, $\widehat{p}(q) \in \mathcal{P} \quad \forall q$, and
\begin{align*}
\sup_{q \in \mathcal{S}^{d-1}}| R_{1,k} (\widehat{p}(q)) | \lesssim_P \sup_{p \in \mathcal{P}}| R_{1,k} (p) | \lesssim_P \sup_{f \in \mathcal{F}_2} | {\mathbb{G}_{n,k}} [f] | = o_P(1).
\end{align*}
Step 3. Bound on $R_{2,k} (\widehat{p}(q))$. On the event $\mathcal{E}_N$, for any $k=1,2,\dots, K$
\begin{align*}
\sup_{q \in \mathcal{S}^{d-1}}| R_{2,k} (\widehat{p}(q)) | &\leq \sup_{p \in \mathcal{P}} \sqrt{N} | \psi(p, \widehat{\xi}_k) - \psi_0(p) | \\
&\leq \sup_{ \xi \in \Xi_N} \sup_{p \in \mathcal{P}}\sqrt{N} | \psi(p, \xi) - \psi_0(p) | \lesssim \sqrt{N} \mu_N.
\end{align*}
Combining the bounds over a finite set $k=1,2,\dots, K$ gives
\begin{align*}
\dfrac{1}{K} \sum_{k=1}^K \bigg[R_{1,k} (\widehat{p}(q)) + R_{2,k} (\widehat{p}(q)) \bigg] = o_P(1).
\end{align*}
Step 4. Bound on $R(\widehat{p}(q), p_0)$.
By Assumption (ref), Lemma (ref) and Step 1, on the event $\mathcal{E}_N$,
\begin{align*}
\sup_{q \in \mathcal{S}^{d-1}} | R(\widehat{p}(q), p_0(q)) | = o (\sqrt{N} \| (\widehat{\Sigma} (\widehat{\eta}))^{-1} - \Sigma^{-1} \|) = o_P(1).
\end{align*}
Step 5. Conclusion. For the influence function $h(W,q)$ in (ref), the function class
\begin{align*}
\mathcal{H} &= \bigg\{ h(\cdot, q), \quad q \in \mathcal{S}^{d-1} \bigg\} \subseteq \mathcal{G}_{\xi_0} +\mathcal{H}_A
\end{align*}
is included into the sum $\mathcal{G}_{\xi_0} + \mathcal{H}_A $. The function classes $\mathcal{G}_{\xi_0}$ and $\mathcal{H}_A $ are Donsker classes with square integrable envelopes (by Assumption (ref) and Lemma (ref), respectively). The statement of the lemma follows from the Skorohod-Dudley-Wichura construction, as in Skorohod, Dudley and Wichura.
The bootstrap support function process can be decomposed as
align*[align* omitted — 132 chars of source]
The first summand can be decomposed as
align*[align* omitted — 381 chars of source]
The weighted moment can be decomposed as
align*[align* omitted — 363 chars of source]
proof[Proof of Theorem (ref)]
By Comment B.1 in CCMS, if the bootstrap random element converges in probability $P$ unconditionally (i.e, $Z_N=o_P(1)$), then $Z_N = o_{P^e} (1)$
in $L^1(P)$ sense and hence in probability $P$, where $P^e$ denotes the probability measure conditional on the data.
Step 1. As shown in Step 1 of the proof of Theorem (ref), $\widetilde{p}(q) \in \mathcal{P} \quad \forall q$. The bound on $R_{1,k}^e (\widetilde{p}(q))$ is established in the proof of Theorem (ref), Step 2. By construction,
\begin{align*}
\sqrt{N}(\psi(p, \widehat{\xi}_k) - \psi(p, \xi_0)) = R_{2,k} (p) = \mathbb{E} [ v \cdot (g(W, p, \widehat{\xi}_k) - g(W, p, \xi_0)) \mid (W_i)_{i \in J_k} ].
\end{align*}
Thus, the bounds on $R_{2,k}(p)$ and $R(p, p_0)$ are established in Steps 3 and 4 of Theorem (ref).
Step 2. As shown in Step 5 of Theorem (ref), the function class $\mathcal{H} \subseteq \mathcal{G}_{\xi_0} + \mathcal{H}_A$ is a Donsker class with
square-integrable envelopes. Then by the Donsker theorem for exchangeable bootstraps, weak convergence holds conditional on the
data,
\begin{align*}
\mathbb{G}_N [ (e-1) h(W,q) ] /\bar{e} \Rightarrow \widetilde{\mathbb{G}[h(q)]} under P^e in probability P,
\end{align*}
where $ \widetilde{\mathbb{G}[h(q)]} $ is a P-Brownian bridge independent of $G[h(q)]$ with the same distribution as $G[h(q)]$.
abstractThis appendix contains supplementary statements for the paper “Debiased Machine Learning of Set-Identified Linear Models” by Vira Semenova. Appendix A contains useful technical statements. Appendix B contains additional proofs. Appendix C contains supplementary tables.