The exact contents of citations.db main_text.text for this paper — one flattened LaTeX string, title through conclusion, appendix excluded, unmodified except for removing email addresses. This is what our citation measures are computed over.
105,192 characters
Simple subvector inference on sharp identified set in affine models
\begin{frontmatter}
\title{Simple subvector inference on sharp identified set in affine models \protect\thanksref{T1}}
\runtitle{Simple subvector inference on sharp identified set}
\thankstext{T1}{This version is July 12, 2024. Earlier versions of the paper was circulated with the title ``Inference on Scalar Parameters in Set-Identified Affine Models'' and ``Inference in high-dimensional Set-Identified Affine Models''. The work is based on a chapter in my PhD dissertation \citep{gafarov2017essays}. First-draft date: November 10, 2015.}
\begin{aug}
\author{\fnms{Bulat} \snm{Gafarov}\thanksref{t1}\ead[label=e1]{[email removed]}}
\thankstext{t1}{I am extremely grateful to Joris Pinkse and Patrik Guggenberger for their very helpful and detailed comments. I would like to thank (in alphabetical order) Donald Andrews, Andres Aradillas-Lopez, Christian Bontemps, Ivan Canay, Peng Ding, Graham Elliott, Zheng Fang, Dalia Ghanem, Joachim Fryberger, Ronald Gallant, Michael Gechter, Marc Henry, Keisuke Hirano, Sung Jae Jun, Nail Kashaev, Francesca Molinari, Demian Pouzo, Adam Rosen, Thomas Russell, Andres Santos, Xiaoxia Shi, Jing Tao, and Alexander Torgovitsky
for their comments and suggestions.}
\runauthor{B. Gafarov}
\address{University of California, Davis, Department of Agricultural and Resource Economics \\
}
\end{aug}
\begin{abstract}
This paper studies a regularized support function estimator for bounds on components of the parameter vector in the case in which the identified set is a polygon.
The proposed regularized estimator has three important properties: (i) it has a uniform asymptotic Gaussian limit in the presence of flat faces in the absence of redundant (or overidentifying) constraints (or vice versa); (ii) the bias from regularization does not enter the first-order limiting distribution; (iii) the estimator remains consistent for sharp (non-enlarged) identified set for the individual components even in the non-regualar case.
These properties are used to construct \emph{uniformly valid} confidence sets for an element $\theta_{1}$ of a parameter vector $\theta\in\mathbb{R}^{d}$ that is partially identified by affine moment equality and inequality conditions.
The proposed confidence sets can be computed as a solution to a small number of linear and convex quadratic programs, leading to a substantial decrease in computation time and guarantees a global optimum.
As a result, the method provides a uniformly valid inference in applications in which the dimension of the parameter space, $d$, and the number of inequalities, $k$, were previously computationally unfeasible ($d,k=100$).
The proposed approach can be extended to construct confidence sets for intersection bounds, to construct joint polygon-shaped confidence sets for multiple components of $\theta$, and to find the set of solutions to a linear program.
Inference for coefficients in the linear IV regression model with an interval outcome is used as an illustrative example.
\end{abstract}
\begin{keyword}
\kwd{affine-moment inequalities}
\kwd{asymptotic linear representation}
\kwd{higher-order analysis}
\kwd{delta method}
\kwd{interval data}
\kwd{intersection bounds}
\kwd{partial identification}
\kwd{ regularization}
\kwd{ strong approximation}
\kwd{ stochastic programming}
\kwd{subvector inference}
\kwd{ uniform inference}
\end{keyword}
\end{frontmatter}
\newpage
\setcounter{page}{1}
\section{Introduction}\label{sec:Introduction}
Strong econometric assumptions can lead to poor estimates.
Sometimes, moment inequalities can provide alternative estimates under weaker assumptions.
Linear models with interval-valued outcome data are a good example.\footnote{Other examples of affine-moment inequalities include monotone instrumental variables (\citet{manski2000monotone}, \citet{Freyberger201541}) and models with missing data (\citet{manski2003partial}).}
It is common practice to replace income-bracket data with the corresponding midpoints when estimating the returns to schooling (\cite{trostel2002estimates}).
However, the conventional approach is applicable only under strong assumptions about the distribution of the residual term.\footnote{Another common approach is to assume Gaussian distribution for the residuals and apply the maximum likelihood method (\citet{stewart1983least}). }
The affine moment inequality approach to interval-valued data proposed by \cite{manski2002inference} can set-identify the return to schooling without such strong assumptions.
Multiple methods can be used to construct confidence sets (CS) for parameters defined by moment inequalities.
The pioneering procedures of \cite{chernozhukov2007estimation} and \cite{andrews2010inference} (AS) and their subsequent refinements by
\cite{bugni2014inference} (BCS) and \cite{kaido2015inference} (KMS) are powerful statistical methods that solve this inference problem in the small-dimensional case.
Some applications, such as panel or semiparametric regression models with interval-measured outcome variables, have a high-dimension parameter space, which poses a computational challenge for the existing procedures.\footnote{For example, \cite{trostel2002estimates} consider a panel regression with more than 60 variables that include country fixed effects, time effects, exogenous demographic control variables, and their interactions.}
I propose a novel regularized support function estimator for the lower and upper extremes of the identified set for an element $\theta_{1}$ of an unknown parameter vector $\typVector{\theta}\in\R^{d}$ in models defined by affine moment equalities and inequalities.
In the example of returns to schooling, $\theta_{1}$ corresponds to the returns to schooling and $\typVector{\theta}\in\R^{d}$ to the full vector of the regression coefficients.
The novel estimator has a closed-form asymptotic Gaussian distribution, which I use to construct uniformly valid confidence bounds and confidence intervals for $\theta_{1}$.
The proposed set has valid asymptotic coverage probability uniformly over a class of data-generating processes (DGP).
Uniformity in DGP is a desirable property, as it results in better coverage-probability control in small samples compared to point-wise analogs in \emph{nonregular statistical models} such as the affine moment inequality model.
The regularized support function proposed in this paper is a solution to a convex quadratic program that minimizes the sum of $\theta_{1}$ and a penalty $\mu_{n}\left\Vert \typVector{\theta}\right\Vert ^{2}$ with $\mu_{n}\to0$, subject to the sample moment restrictions.
If the set of optima for $\mu_n=0$ is not a singleton, this additional convex term selects the optimum with the minimal norm as $n$ increases.
The standard errors are computed using the sample variance of the weighted moment conditions at the unique optima.
To correct the asymptotic bias resulting from the regularization exactly, I suggest using the argmin of the regularized program with a larger tuning parameter $\kappa_{n}\to0$.
If $\kappa_{n}/\mu_{n}\to\infty$ as $n\to\infty,$ then the bias correction does not affect the asymptotic distribution of the estimator.
To achieve a uniformly valid confidence interval (CI), I replace the exact correction with an upper bound on the maximum of $\mu_{n}\left\Vert \typVector{\theta}\right\Vert ^{2}$ over the argmin set of the nonregularized program.
The proposed CIs have several attractive statistical and computational properties that make them viable in high-dimensional affine moment inequality models.
First, the estimator of the regularized support function has an asymptotically linear (Bahadur-Kiefer) representation that provides an easy-to-compute asymptotic standard error.
Consequently, this paper is the first to propose a closed-form estimator of the bounds on $\theta_{1}$ convex moment inequality models with asymptotic Gaussian distribution in the non-regular case (in the absence of strict convexity).
In contrast, the estimator of the ordinary support function estimator used in the existing literature (\citet{beresteanu2008asymptotic}, \citet{kaido2014asymptotically}, \citeauthor{Freyberger201541} (2015, FH), \citet{gafarov2018delta}, among others) can have a non-Gaussian asymptotic distribution, which complicates uniform inference.
To establish the uniform remainder bound in the Bahadur-Kiefer expansion, I developed a novel second-order directional envelope theorem, which is a theoretical result of independent interest.
Second, the proposed approach requires only a fraction of the computational time of the existing uniform procedures if $\typVector{\theta}$ has a large number of dimensions.
The computational cost is low since it involves only four quadratic programs, it does not require any resampling, and it depends on covariance of the moment conditions at two points.
In the asymptotic analysis, I only consider the case of a fixed dimension of $\theta$ and a fixed number of inequalities for analytical simplicity.
The goal of the present study is to focus on the computational difficulties resulting from the high-dimensional moment inequalities.
(There are recent papers that are concerned with the impact of the growing number of inequalities on the statistical properties of the inference procedures; see, for example, \cite{belloni2018subvector}.)
The computation time for my procedure increases slowly in the dimension of $\theta\in\R^{d}$ and takes only 0.1 second for $d=20$ and $k=40$ moment inequalities and 2.5 seconds for $d=100$ and $k=200$. As a result, the proposed method can address parameter $\theta$ with a large dimension and a large number of moment conditions.
In contrast, the existing uniform-inference methods for moment inequalities proposed by AS, KMS, and BCS are based on costly nonconvex optimization.
\label{reply:R1.1} In fact, despite the fact that the moment inequalities are convex (linear) the standardized moment conditions are not convex which can potentially result in creation of multiple local optima which complicate computation for these statistical procedures (see Appendix Sections \ref{sec:Convergnce rate} and \ref{sec:globalOptimum}).
I provide an example of an affine moment inequality model by showing that the number of local optimal solutions in existing uniform procedures (AS, BCS, and KMS) can grow exponentially with dimension $d$.
As a result, the procedures take more computational time and can produce misleadingly short CIs if the optimization routine fails to find the global optimum.
It takes 100 seconds to compute the CI of AS in an affine model with $d=20$ and $k=40$ moment inequalities which is 1000 times longer than the newly proposed method.
The simulation evidence suggests that {the computational speed gains} come without a substantial loss of statistical power.
In the non-regular cases, the proposed uniform CIs have length properties that are not worse than those of the existing uniform methods. (The proposed \emph{uniform} CI has a length within simulation error from the projection CI of AS in the Monte Carlo (MC) design considered in this paper.) In the regular cases, the novel confidence bounds attain the efficiency bound of the usual support function estimators (as shown in \cite{kaido2014asymptotically}).
\label{reply:R2.1.1}
The proposed idea of regularizing support functions for linear moment inequality inference gained a subsequent development in a recent work \cite{cho2023simple}.\footnote{The first draft of \cite{cho2023simple} was circulated on 27 May 2019 on Arxiv depository two years after the initial distribution of the present paper as a PhD dissertation chapter in \cite{gafarov2017essays}.}
The authors argue that if one is satisfied with an inference on an \emph{enlarged} identified set, one can simply regularize the moment inequalities by adding random noise to their coefficients under milder regularity conditions.
This operation restores a Gaussian limit of the perturbed estimated support function, allowing conventional bootstrap inference on the \emph{enlarged} identified set.
Since the inference is done on the enlarged identified set, one does not need to impose any constraint qualification conditions on the moment inequalities --- they are satisfied automatically with probability 1 after adding the random noise to the coefficients.
Unfortunately, this approach results in confidence sets with zero power against all local alternatives corresponding to other nested enlargements of the identified set.
In practical terms, it means that the confidence sets can be very conservative if too much noise was added or fail to control size in small samples if the added noise was insufficient (there is no theory that would determine the minimal level of required noise for a given sample size).
In contrast, the theory provided in the present paper explicitly studies the impact of the choice of tuning parameters on the size of the regularized identified set.
Such an analysis requires imposing constraint qualification conditions on the moment inequalities (Assumption \ref{ass:ULICQ}).
Under these conditions, Theorem \ref{thm:bounds} shows that for sufficiently small value of the tuning parameters $\mu_n$ and $\kappa_n$ both lower and upper bounds based on the regularized value function coincide exactly with the nonregularized support function in the regular case when the set of primal solutions is a singleton.
More generally, under the maintained assumptions, the regularized estimators remain consistent for the bounds on the sharp (non-enlarged) identified set, thus providing non-trivial local power against the relevant alternatives.
Identified sets defined by affine inequalities appear in various economic applications, in particular those dealing with discrete variables and shape restrictions. Linear models with interval outcome, which were originally studied in \cite{manski2002inference} and \cite{haile2003inference}, are just one example of affine inequalities. Other examples include bounds on marginal effects in dynamic discrete choice panel models (\cite{honore2006bounds}, \cite{torgovitsky2016nonparametric, torgovitsky2018partial}), bounds on average treatment effects (\cite{kasy2016partial}, \cite{Laffer2018}, \cite{russell2017sharp}), nonparametric IV models with shape restrictions (\cite{manski2000monotone}, \cite{Freyberger201541}), errors in variables (\cite{MOLINARI200881}), intersection bounds (\cite{honore2006bounds2}), revealed preference restrictions (\cite{kline2016bounding}), and game-theoretic models (\cite{syrgkanis2017inference}).
This paper focuses on the case with finitely many affine unconditional inequalities that are non-overidentifed and that result in a non-empty identified set.
This, of course, rules out many interesting economic applications that involve overidentifying moment inequalities \citep[for example,][]{shi2018estimating}, or where moment inequalities are nonlinear \citep[for example,][]{pakes2015moment}.
I also contribute to the growing literature on inference on non-differentiable functions and regularized estimators.
Bounds on components of a parameter characterized by linear moment conditions considered in this paper are an example of a nondifferentiable (\emph{nonregular}) function of a parameter (the expectation of the data) that has an asymptotically normal estimator (the sample mean).
Distributions of nondifferentiable functions of the normal estimator are hard to approximate using standard methods.
(See Section~\ref{subsec:asymptoticDistribution} for a detailed discussion.)
I propose differentiable lower and upper bounds that converge to the nondifferentiable function of interest as the sample size grows.
Since the bounds are regular parameters themselves, the standard delta method and bootstrap can be used to conduct one-sided or two-sided inference on the bounds.
In the regular case, these bounds collapse and coincide with the original parameter of interest, which results in a $\sqrt{n}$-consistent and asymptotically normal estimator.
In the non-regular case, the bounds converge at a slower rate and result in a locally biased estimator, which is acceptable for valid one-sided inference.
Another interesting statistical problem that appears in many applications is inference for extrema of finitely many means of random variables.
It is known as \emph{ intersection-bounds problem } (\cite{hall2010bootstrap} and \cite{chernozhukov2013intersection}) and can be framed as a value of a linear program.
The regularized support function estimator can also be used for uniform delta-method CSs in this setting (see Appendix Section~\ref{subsec:overidentification}).
The approach considered here is expected to have statistical properties similar to \cite{chernozhukov2013intersection}, but has the advantage of closed-form standard errors and critical values, which correspond to the standard normal distribution.
The paper is structured as follows. Section~\ref{sec:Setup} describes the setup, gives examples, and summarizes the available literature.
Section~\ref{sec:results} provides the main results, applies them to uniform inference on projections, and discusses extensions (overidentified inequality models, joint CSs, characterization of argmin sets).
Section~\ref{sec:Monte-Carlo} provides the results of the Monte Carlo experiments.
Section~\ref{sec:Conclusion} concludes.
The following notational conventions are used.
$\bydef$ denotes definitions.
$\expect\left[\cdot\right]$ denotes expectation with respect to a probability distribution $\measTrue$.
Uppercase English letters denote random variables (scalar, vector, or matrix valued), and lowercase letters denote the corresponding realizations, for example, $\Data_i$ and $\data_i$.
$\measEmp$ denotes the empirical distribution and $f\left(0^+\right)$ denotes $\lim_{x\downarrow0}f\left(x\right)$.
The vector $\ej\bydef(0,...1,...0)'$ is the $j$-th coordinate vector, where 1 occurs at position $j$. $\ej$ is the projector on the $j$-th coordinate.
The symbol $\jSet$ denotes a finite set of indices $\jSet\bydef\{i_1,...,i_\ell\}\subset\numbers$ and $\jMat\bydef(e_{i_1},...,e_{i_\ell})^\prime$ as a coordinate projection matrix in the corresponding Euclidean space.
$\left|\jSet\right|$ denotes the cardinality of the set $\jSet$. u.h.c. stands for upper-hemicontinuous correspondence.
\section{Setup, motivating examples, and related literature}\label{sec:Setup}
\subsection{\label{subsec:GeneralSetup} Support function and projections of identified sets}
\label{R3.5.c}
Consider a support function for a polygon $\setID\subset\R^{d}$ (a set defined by a system of linear equalities and inequalities) that depends on a data generating process parametrized by a measure $\measTrue$ evaluated at a direction $e_1\bydef(1,0,\dots)^\prime
\in\R^d$,\footnote{See, for example, \citet[][ chapter 13]{rockafellar1970convex}}
\begin{equation}
\begin{array}{ccc}
\funLValO\bydef\displaystyle\min_{ {\theta}\in\setID}\eOne {\theta}.
\end{array} \label{prog:minP}
\end{equation}
The set $\setID$, also referred to as an \emph{identified set} for a parameter vector $ {\theta}\in\R^{d}$, consists of all solutions to the following system of affine moment equalities/inequalities:
\begin{equation}
\begin{cases}
\expect\moment{W}{\theta} =0, & j\in\setEq,\\
\expect\moment{W}{\theta} \leq0, & j\in\setIneq.
\end{cases}\label{eq:GenericMomentConditions}
\end{equation}
Here, $\moment{W}{\theta} \bydef \sum_{\ell=1}^d \Data_{j\ell} \theta_\ell -\Data_{j(d+1)}$.
I consider a setup with finitely many unconditional moment functions (that is, $\left|\setEq\cup\setIneq\right|=k$, $\card{\setEq}=p$, $0\leq p \leq d$, $k<\infty$).
The random matrix corresponding to an individual observation $\Data$ has probability measure $\measTrue$ with sample space $\R^{k\times\left(d+1\right)}$,
\begin{align*}
W &= \left(\begin{array}{c}
W'_1\\
\vdots \\
W'_k
\end{array} \right) ,\quad W_j = (W_{j1},\dots,W_{j(d+1))})', \quad j=1,\dots,k.
\end{align*}
The econometrician observes an i.i.d. sample $\left\{ \data_{i}\in\R^{k\times\left(d+1\right)}|\,i=1,...,n\right\} $ of random matrix $W$.
There is a straightforward way to extend the analysis to the case of dependent data as long as a CLT for averages of $\data_{i}$ remains valid.
This setup reduces to a finite-dimensional parametric statistical model once the support function evaluated at $e_1$ is represented as a linear program,
\begin{align}
\funLValO&=\min_{ {\theta}\in\R^{d} } \eOne {\theta} \label{eq:primalLP} \\
\text{s.t. } & e_j^\prime \Ap \theta= e_j^\prime\bp , \quad j\in\setEq,\nonumber\\
& e_j^\prime \Ap \theta\leq e_j^\prime\bp , \quad j\in\setIneq.\nonumber
\end{align}
Here, the coefficients on the left-hand side, $\Ap$, and the right-hand side, $\bp$, taken together constitute matrix $\expect \Data \bydef (\Ap|\bp)$. It means $W$’s expectation is a $k\times (d + 1)$ matrix whose first $d$
columns is $A_P$ and the last column is the $k$-dimensional vector $b_p$.
Following optimization theory terminology, I occasionally refer to the support function at $e_1$ as \emph{value} of program \eqref{eq:primalLP} and to the corresponding argmin set as the set of \emph{optimal solutions}.
The focus on the first coordinate of $\theta$ as an objective of \eqref{prog:minP} is without loss of generality.
As discussed in Section~\ref{subsec:subvectors}, the support function evaluated at any unit direction vector $a\in\R^d$ can be represented in the form \eqref{prog:minP} after some redefinition of the parameter space (matrices $ \Ap$ and $\bp $ can be functions of the unit vector $a$).
To ensure boundedness of the support function, I assume that system (\ref{eq:GenericMomentConditions}) includes inequalities that make the identified set compact,
\begin{equation}
-\infty<-\underline{c}_{\ell}\leq\theta_{\ell}\leq \bar{c}_{\ell}<\infty \quad \text{for } \ell=1,...,d, \label{eq:Box}
\end{equation}
for some constants $\underline{c},\bar{c}\in\R^d_+$.
Moreover, I impose the following assumption:
\begin{assumption}
\label{assu:Non--empty} $\Theta\left(\measTrue\right)$ is non-empty
for the probability measure $\measTrue$.
\end{assumption}
This assumption implies that $\funLValO<+\infty$.
It is valid whenever the model corresponding to \eqref{eq:GenericMomentConditions} is correctly specified.
Under Assumption \ref{assu:Non--empty}, support functions evaluated at $\pm e_j$ characterize projections of $\setID$ on individual coordinates $j$ of $\theta$.
In particular, the marginal (projected) identified set for the coordinate $\theta_{1}$ can be represented as an interval $\setMarginal=\left[\funLValO,\funUValO\right]$,
where the bounds are support functions evaluated at directions $\eOneVec$ and $-\eOneVec$.
Indeed, the upper bound can also be written as a (minus) support function at $e_1$,
\begin{equation}
\funUValO\bydef\max_{ {\theta}\in\setID} \{ \eOne {\theta}\}=-\min_{ {\theta}\in\setID} \{-\eOne {\theta}\}.
\end{equation}
The analysis for the upper bound is analogous to the one for the lower bound, so from here on I focus on the lower bound.
\subsection{\label{subsec:examples} Motivating example: Linear IV model with interval-valued outcome }
The polygon-shaped identified sets considered in this paper appear in many econometric models that feature discrete data, as mentioned in the introduction.
Many of these applications are concerned with instrumental variables.
I use the linear IV model with interval outcome \citep[for example,][]{manski2002inference,chernozhukov2007estimation} to illustrate ideas throughout the paper. Other applications of the proposed method include monotone IV \citep{manski2000monotone} and nonparametric IV \citep{Freyberger201541}.
Consider a linear model for a random vector $(Y,X,Z)$ satisfying
$$\expect\left[Y- {\theta}^{\prime} X| Z\right] =0,$$
where $\theta\in\R^d$ is the vector of regression coefficients.
The true outcome $Y$ is unobserved; only its a.s. bounds $\left[\underline{Y},\overline{Y}\right]$ are observed.
Suppose that the vector of instrumental variables $Z$ has finite support $\left\{ z_{1},\dots, z_{K}\right\} \subset\R^{d}$.
In this case, the model implies a polygon-shaped identified set $\setID$ for $\theta$ defined by a set of moment inequalities,
\begin{equation}
\begin{cases}
\expect\left[\underline{Y}\Ch{ Z= z_{j}} \right]\leq {\theta}^{\prime} \expect \left[ X\Ch{ Z= z_{j}}\right] , & j=1,...,K,\\
\expect\left[\overline{Y}\Ch{ Z= z_{j-K}} \right]\geq {\theta}^{\prime} \expect \left[X\Ch{ Z= z_{j-K}} \right], & j=K+1,...,2K.
\end{cases}\label{eq:intervalOutcome}
\end{equation}
These inequalities can be represented in the standard form (\ref{eq:GenericMomentConditions})
with $p=0$, $k=2K$, and the following observation matrix:
\begin{align*}
\Data_{j\ell}\bydef\begin{cases}
- X_{\ell}\Ch{ Z= z_{j}} , & \text{ for }j=1,\dots,K, \ell =1,\dots,d,\\
X_{\ell}\Ch{ Z= z_{j-K} } , & \text{ for }j=K+1,\dots,2K, \ell =1,\dots,d,
\end{cases}\\
\Data_{j(d+1)}\bydef\begin{cases}
\underline{Y}\Ch{ Z = z_{j} } , & \text{ for }j=1,...,K,\\
-\overline{Y}\Ch{ Z = z_{j-K} } , & \text{ for }j=K+1,...,2K.
\end{cases}
\end{align*}
If it is known a priori that for some support points $j$ $$\expect\left[\underline{Y}\Ch{ Z= z_{j}} \right] = \expect\left[\overline{Y}\Ch{ Z= z_{j-K}} \right],$$
then one can replace the corresponding pair of inequalities with a single equality,
\begin{equation}
\expect\left[\frac{1}{2}(\underline{Y}+\overline{Y})\Ch{ Z= z_{j}} \right]= {\theta}^{\prime} \expect X\Ch{ Z= z_{j}} .
\end{equation}
In this case, $p$ (the number of equality restrictions) is equal to the number of such support points.
One can further incorporate additional shape restrictions information such as signs of components of $ {\theta}$ in the form of linear inequalities to narrow the identified set.
This example appears in the context of the estimation of return to schooling using survey data.
\citet{trostel2002estimates} study economic returns to schooling
for 28 countries using data from the International Social Survey Programme from 1985 to 1995. They estimate the conventional \citet{mincer1974schooling}
model of earnings (the human capital earnings function), which has $Y$, the log of hourly wages, satisfying
\begin{equation}
\expect\left[Y - {\theta}^{\prime} X| Z\right]=0 ,\label{eq:schooling}
\end{equation}
where the first regressor, $ X_{1}$, is the years of schooling; the other components of $X$ play a role of additional controls.
The component $\theta_{1}$ is then interpreted as the return to schooling.
It is equal to the percentage change in wages due to an additional year of schooling.
To correct for the endogeneity bias in $\theta_{1}$ resulting from omitting the latent ability variable, one can use an instrument vector $Z$ that correlates with $X_{1}$ (for example, an indicator of whether a good school is in proximity or an indicator of the quarter of birth).
Exact measurements of $Y$ are not available for some countries (including the US); only hourly-income- bracket data $\left[\underline{Y},\overline{Y}\right]$ are
available.
Instrumental variables $Z$ (coinciding with the control variable) considered in \citet{trostel2002estimates} take discrete values.
The variables include annual fixed effects, union status, marital status, age and age squared, and country-year dummies (in the case of the aggregate equation).
With inclusion of the country and time effects, $d$ can be larger than 60.
Because of the large number of support points for $Z$, the corresponding system~\eqref{eq:GenericMomentConditions} would also have a large number of linear moment inequalities $k$.
The conventional technique to estimate IV regression with interval-outcome data is to replace the interval data with the corresponding midpoints and estimate the model using the OLS method.
This technique is valid only under the unreasonably strong condition
\begin{align}
\expect\left[\left(Y-\frac{1}{2}(\underline{Y}+\overline{Y})\right) Z\right] & =0.\label{eq:midpoint}
\end{align}
If Equation \eqref{eq:midpoint} is violated, then the OLS estimator with midpoints is inconsistent for the true parameter $\theta$.\footnote{The OLS is, however, consistent for the best linear predictor of the midpoint that may not have the desirable economic interpretation; see, for example, \citet{shi2020uniform}.}
Without assuming that \eqref{eq:midpoint} is satisfied, support function estimators (defined below) provide consistent bounds on marginal identified sets of the true parameter $\theta$ \citep[see][]{beresteanu2008asymptotic,bontemps2012set}.
When constructing valid CSs on projections of identified sets, strictly speaking, one needs only a lower bound on the support function at $e_1$ in \eqref{prog:minP}.
Other applications may instead require estimators of an upper bound on the (minimum) value of an optimization problem.
The upper bound can be used, for example, to bound the set of optimal solutions or to construct a consistent \emph{inner} estimator of a convex identified set \citep[the inner set estimator has important applications in inference; see, for example, ][]{bugni2017inference}.
\subsection{\label{subsec:asymptoticDistribution} Review of existing results on the asymptotic distribution of support function estimators for polygon-shaped identified sets }
Following \cite{beresteanu2008asymptotic} the parameter $\funLValO$ (the support function at $e_1$) can be estimated using a sample analog,
\begin{align}
\hat{\underline{v}}_n &=\min_{ {\theta}\in\R^{d} } \eOne {\theta} \\
\text{s.t. } & \begin{cases}
\sMean\moment{W_i}{\theta} =0, & j\in\setEq,\\
\sMean\moment{W_i}{\theta} \leq0, & j\in\setIneq,
\end{cases}
\end{align}
where the observations $W_i$ are independent copies of the random matrix $W$.
The asymptotic distribution of this estimator has been extensively studied in the case of the strictly convex identified set.
The strict convexity implies that the support function at $e_1$, $\hat{\underline{v}}_n$ , is differentiable in the sample mean $\sMean W_i$ (under additional regularity conditions discussed below).
As a consequence, its sample analog has an asymptotic Gaussian distribution, admits of bootstrap inference, and attains semiparametric efficiency \citep{kaido2014asymptotically}.
\label{reply:R1.3.1}In contrast to the strictly convex case, the Gaussian limit is not guaranteed anymore in the general (non-deterministic) linear moment inequality model (\cite {kaido2014asymptotically} only allow for deterministic linear inequalities under particular regularity conditions).
The stochastic programming approach \citep{shapiro1991asymptotic} enables an asymptotic analysis of the support function estimators for a fixed direction in the general convex case, which includes the linear moment inequality model.
In this section, I briefly review the two main ideas in this approach, Lagrangian duality of convex programs and the delta method for directionally differentiable functions.
I conclude with an overview of the alternative inference approaches and their relation to stochastic programming methods.
The Lagrangian duality theory provides additional, often more convenient, formulations of convex programs.
Namely, the (primal) program in \eqref{prog:minP} has the same value as its Lagrangian dual formulation,
\begin{align}
\funLValO=\max_{ {\lambda}\in\R^{p}\times\R_{+}^{k-p}} & \left\{ - {\lambda}^{\prime}\bp\right\} \label{eq:DualProgram0},\\
\text{s.t. } & {\lambda}^{\prime}\Ap=-\eOne,\nonumber
\end{align}
given Assumption~\ref{assu:Non--empty}.
Iff the dual program has a bounded set of solutions $\argminL(\measTrue)$, the support function for a given direction has bounded directional derivatives in $\Ap,\bp$ (defined precisely below).
Proposition 5.45 in \cite{bonnans2013perturbation} provides a necessary and sufficient condition for $\argminL(\measTrue)$ to be bounded:
\begin{condition}[Slater's]
\label{con:slater} There exist $ \theta \in \setID$ s.t. $\expect\moment{W}{\theta} < 0$ for all $ j\in\setIneq$.
\end{condition}
Condition~\ref{con:slater}
enables another, Lagrangian min/max, representation of \eqref{prog:minP}, which is particularly convenient for computing derivatives of the support function in a given direction for the delta-method-based inference procedures.
Suppose that $\Lambda\subset{\R}^{p}\times\R_{+}^{k-p}$ is some compact set that contains $\argminL(\measTrue)$.
The support function for a given direction $ e_1$ can then be represented as \citep[see, for example,][p. 437]{bonnans2013perturbation}
\begin{equation}
\funLValO=\min_{ {\theta}\in \Theta }\max_{ {\lambda}\in\Lambda} \{\eOne {\theta} +{\lambda}^{\prime}(\Ap\theta-\bp)\}. \label{prog:minmax}
\end{equation}
The directional derivative for the value of this min/max program is provided by a corresponding envelope theorem \citep[for example,][Theorem 7.28]{shapiro2014lectures}.
Namely, the envelope theorem gives a derivative of any perturbed version of the min/max program,
\begin{equation}
\funLValO[P,t]\bydef\min_{ {\theta}\in\Theta}\max_{ {\lambda}\in\Lambda} \{\eOne {\theta} +{\lambda}^{\prime}((\Ap +t h_{A,t})\theta-(\bp+t h_{b,t}))\}, \label{prog:minmaxPerturbed}
\end{equation}
with respect to a scalar $t$ for any uniformly converging sequence of directions $h_t\bydef(h_{A,t},h_{b,t})\to(h_{A},h_{b})$.
The derivative takes the form
\begin{equation}
\lim_{t\to0} \frac{\funLValO[P,t]-\funLValO[P]}{t}=\min_{ {\theta}\in \argminT(\measTrue)}\max_{ {\lambda}\in\argminL(\measTrue)} \{ {\lambda}^{\prime}( h_{A} \theta- h_{b})\}, \label{prog:minmaxDerivative}
\end{equation}
where $\argminT(\measTrue)$ and $\argminL(\measTrue)$ are sets of primal and dual optima for \eqref{prog:minP}.
Both sets can be nonsingletons as illustrated below.
\begin{example}[Bivariate linear IV with an interval outcome]\label{exa:runningExample}
Suppose that $\theta\in\R^2$, $X=Z$, and regressors $Z_{1}$ and $Z_{2}$ take values in $\left\{ 0,1\right\} $ with $\expGeneric Z_{1}=\expGeneric Z_{2}=\frac{1}{2}$.
As in \eqref{eq:intervalOutcome}, the identified set for $\theta$ can be characterized by eight inequality constraints,
\begin{equation}
\expGeneric\left[\underline{Y}\psi_z(Z)\right]\leq \expGeneric\left[Z_{1}\psi_z(Z)\right] \theta_{1}+\expGeneric\left[Z_{2}\psi_z(Z)\right]\theta_{2}\leq\expGeneric\left[\bar{Y}\psi_z(Z)\right], \label{eq:fullsystem}
\end{equation}
where indicator functions $\psi_z(Z)=\Ch{Z=z}$ correspond to all combinations of $z\in \{0,1\}^2$.
Suppose, for illustrative purposes, we are interested only in the identified set $\setID$ defined by the following subsystem of four inequalities:
\begin{equation}
\expGeneric\left[\underline{Y}Z_{1}\right]\leq\frac{1}{2}\theta_{1}+\theta_{2}\expGeneric\left[Z_{1}Z_{2}\right]\leq\expGeneric\left[\bar{Y}_{i}Z_{1}\right],\\
\expGeneric\left[\underline{Y}\left(1-Z_{1}\right)\right]\leq \theta_{2}\expGeneric\left[\left(1-Z_{1}\right)Z_{2}\right]\leq\expGeneric\left[\bar{Y}\left(1-Z_{1}\right)\right].\label{eq:subsystem}
\end{equation}
In order to represent the identified set on a diagram, suppose further that the a.s. bounds on the outcome variable $Y$ satisfy $\expect[\bar{Y}|Z_1=i]=- \expect[\underline{Y}|Z_1=i] = \frac{1}{2}\Delta_{i}\geq0 $ for $i\in\{0,1\}$ (that is, $\Delta_{i}$ is the average length of the outcome interval depending on $Z_1$).
\begin{figure}[H]
\begin{centering}
\begin{tabular}{ccc}
\includegraphics[scale=0.5,trim={2cm 0 1cm 0 }]{Figures/Example/RhoNegative.pdf}& \includegraphics[scale=0.5,trim={2cm 0 1cm 0 }]{Figures/Example/RhoZero.pdf}& \includegraphics[scale=0.5,trim={2cm 0 1cm 0 }]{Figures/Example/RhoPositive.pdf}\tabularnewline
\end{tabular}
\par\end{centering}
\caption{\label{fig:Identified-sets-Example2}The identified sets in Example~\ref{exa:runningExample} for various values of $\rho$.}
\end{figure}
The shape of the full identified set $\setID$ depends on the value of $\rho\bydef \expGeneric\left(Z_{1}Z_{2}\right)$.
The corresponding marginal identified set for $\theta_1$ can be written in explicit form, $\setMarginal =\left[-\Delta_{1}-2\abs{\rho}\Delta_{0},\Delta_{1}+2\abs{\rho}\Delta_{0}\right]$.
If $\rho=0$, $\argminT_2$---the second coordinate of the set of primal optima of the program in \eqref{prog:minP}---is not uniquely defined.
The case of a non-singleton set of dual optima occurs when the number of binding inequalities at the primal optimum is larger than the dimension of the parameter space.
Suppose that we further restrict the parameter space by imposing $\theta_2=0$.
Consider the following two moment inequalities from system \eqref{eq:fullsystem}:
\begin{equation}
- \frac{1}{2} \theta_{1} \leq -\expGeneric\left[\underline{Y}Z_{1}\right],\\
- \frac{1}{2} \theta_{1} \leq -\expGeneric\left[\underline{Y} \right].\label{eq:subsystemIntersectionBounds}
\end{equation}
The support function at $e_1$ corresponding to \eqref{eq:subsystemIntersectionBounds} takes the explicit form
\begin{equation}
\funLValO = \argminT_1 = 2 \min \{\expGeneric \underline{Y}Z_{1},\expGeneric\underline{Y}\}.
\end{equation}
If $\expGeneric\underline{Y}=\expGeneric\underline{Y}Z_{1}$, then both inequalities in \eqref{eq:subsystemIntersectionBounds} are binding at the optimum, resulting in a non-singleton dual-optimum set (compare with the parameter-on-the-boundary problem and the intersection-bounds problem considered in \citet{andrews2001testing} and \citet{chernozhukov2013intersection}, respectively).
Indeed, the dual formulation of the support function at $e_1$ takes form
\begin{align}
\funLValO=\max_{ {\lambda}\in\R_{+}^{2}} & \left\{ \expGeneric\left[\underline{Y} \right] \lambda_1+ \expGeneric\left[\underline{Y}Z_1 \right] \lambda_2 \right\} = \expGeneric\left[\underline{Y} \right] \max_{ {\lambda}\in\R_{+}^{2}} \left\{ \lambda_1+ \lambda_2 \right\} \label{eq:DualProgramIntersectionBounds},\\
\text{s.t. } & -\frac{1}{2}\lambda_1- \frac{1}{2}\lambda_2 = -1.\nonumber
\end{align}
This program has the constant value of the objective function on the entire optimization domain, which coincides with its argmin set,
\begin{equation}
\argminL =\{\lambda\in\R_+^2| \lambda_1+ \lambda_2 =2\}. \QEDB
\end{equation}
\end{example}
The derivative \eqref{prog:minmaxDerivative} depends on particular optimal primal and dual solutions that are selected by a given perturbation $(h_{A},h_{b})$, unless both sets of solutions are singletons.
I refer to parameter values $\Ap$ and $\bp$ resulting in nonsingleton (primal or dual) solutions as \emph{nonregular}.
\cite{shapiro1991asymptotic} proposes a generalization of the delta method for the directionally differentiable functions of asymptotic normal estimators.
Under Condition~\ref{con:slater}, the sample analog of $\funLValO$ has a solution (and thus is well defined) with probability approaching 1.
Theorems 3.4 and 3.5 in \cite{shapiro1991asymptotic} provide an asymptotic distribution of the sample support function in a given direction.
It takes form
\begin{equation}
\frac{1}{\sqrt{n}}(\hat{\underline{v}}_n-\funLValO[P]) \wTo \min_{ {\theta}\in \argminT(\measTrue)}\max_{ {\lambda}\in\argminL(\measTrue)} \{ {\lambda}^{\prime}( \eG_{A}(P) \theta- \eG_{b}(P))\}, \label{eq:asyDistributionShapiro}
\end{equation}
where $\eG(P)\bydef(\eG_{A}(P),\eG_{b}(P))$ is the limiting zero-mean Gaussian process for $\eG_n(P)\bydef\frac{1}{\sqrt{n}} \sum_{i=1}^n(\data_i-\expect \Data)$ indexed by $\measTrue$.
In general, the limit \eqref{eq:asyDistributionShapiro} is non-Gaussian since either of the two sets $\argminT$ and $\argminL$ can be non-singleton. (It does reduce to a Gaussian distribution when both sets are singletons.)
Several robust methods have been designed for inference on components of $\theta$ defined by moment inequalities (nonlinear, in general).
The dominant approach has been to test individual values of $\theta$ or their subvectors and then to invert the tests \citep[for example,][]{chernozhukov2007estimation,andrews2001testing,andrews2019inference,chernozhukov2019inference, cox2019simple,kaido2019confidence}.
The main advantage of the test-inversion approach is that fixing $\theta$ simplifies the asymptotic distribution of the test statistics and allows for valid inference under weak assumptions (in particular, one can allow for an empty true identified set as in \citet{andrews2019spuriousinference}).
Although having attractive statistical properties, the test-inversion methods can become computationally intractable in settings with high-dimensional parameters $\theta$ since they are based on grid search and resampling methods.\footnote{See further discussion of the computational issues in Appendix \ref{sec:Discussion}.} In contrast, delta-method inference using the support function estimators for a fixed direction proposed in this paper remains a computationally tractable (frequentist) option in high-dimensional settings such as mentioned in Section \ref{subsec:examples} since it fully uses the linear programming structure of the problem.
Inference using the sample support function for a fixed direction in the non-regular case faces three challenges: (i) nondifferentiability makes the standard bootstrap inconsistent \citep[see][]{Fang2018inference}; (ii) discontinuity in the directional derivative \eqref{prog:minmaxDerivative} results in poor (uniform) approximation of \eqref{eq:asyDistributionShapiro} by the numerical bootstrap \citep[see][]{dumbgen1993nondifferentiable,Hong2018numerical};
(iii) the estimator is necessarily (asymptotically) biased \citep[see][]{hirano2012impossibility}.
The regularized support function estimator proposed in the next section addresses these challenges in a robust and computationally tractable way.
\section{Main results}\label{sec:results}
\subsection{\label{subsec:Regularization} Bounds based on the regularized primal program }
The main source of inference complications, the lack of smoothness in the support function, resolves itself when the sets of primal, $\argminT,$ and dual, $\argminL,$ optima are singletons (see equations \eqref{prog:minmaxDerivative} and \eqref{eq:asyDistributionShapiro}). Consequently, my proposal is to consider a regularized support function at $e_1$ that has a unique solution and approximates the original support function from a known direction, either from above or from below.
\footnote{\label{rep1.4}If both primal and dual solutions of the regularized program are unique, all the directional derivatives of the regularized support function at $e_1$ given by \eqref{prog:minmaxDerivative} coincide and are given by formulas
\begin{equation}
\frac{\partial \minV\mnP }{\partial (\Ap)_{ij}}=\argminL_i\mnP\argminT_j \mnP\text{ and } \frac{\partial \minV \mnP}{\partial (\bp)_i}=-\argminL_i\mnP,
\end{equation}
where $\argminT \mnP$ and $\argminL\mnP$ are, correspondingly, the primal and dual solutions of the regularized program with a regularization parameter $\mn$. }
Since the direction of such a regularization bias is known, the corresponding estimators of the regularized support function at $e_1$ admit standard one-sided delta-method and bootstrap inference based on the normal limiting distribution.
Without regularization, both primal and dual solutions of \eqref{prog:minmax} can be non-singletons.
To keep the problem analytically tractable, I focus on the case in which the dual set is a singleton by assumption, and I propose a regularization for the primal problem.
In this way, I can explicitly consider only the bias from the primal regularization and develop the necessary bias correction.
Appendix Section~\ref{subsec:overidentification} discusses the complementary case of dual regularization.
The following \emph{regularized primal program} is strictly convex for any $\mu>0$ and hence has a unique primal solution\label{reply:R1.4} $\argminT (\mu, P)$ and approximates the optimal value of program~(\ref{prog:minP}):
\begin{align}
\funLVal[\mu ] & =\underset{ {\theta}\in\setID}{\mbox{min}}\left\{\eOne {\theta}+\mu \ltwo{ {\theta}}^{2}\right\}.\label{prog:regProg}
\end{align}
The corresponding dual program has a unique solution $\argminL(\mu,P)$ iff the set of constraints satisfies the linear independence constraint qualification (LICQ; see \citet[p.178]{shapiro1991asymptotic} and \cite{ wachsmuth2013licq}),
\begin{condition}[LICQ]
\label{con:LICQ} The matrix of gradients of binding constraints has a full rank for any $ {\theta}\in\setID$.
\end{condition}
\begin{rem}
There is a purely computational reason to ensure that LICQ holds.
If it is violated, Newton-type algorithms, which typically guarantee a (fast) quadratic rate of convergence to a stationary point, have a linear rate of convergence or do not converge at all (see, for example, \cite{golishnikov2006newton}).
\end{rem}
The set of DGPs that satisfy LICQ is not closed since a limit of a sequence of linearly independent matrices can be a reduced rank matrix.
This means that depending on the size of the smallest singular value of the gradients of the set of binding constraints at any point may be arbitrary small while still satisfying LICQ.
As a result, the dual solutions $\argminL(\mu,P)$ may be arbitrary large, implying large derivatives of $\minV\mnP $
and may require very high sample sizes for reasonable precision of the delta-method inference (see equation \eqref{eq:boundsOnLinearApproximation} in Lemma \ref{lem:boundOnVgrowth} in Appendix Section \ref{app:smoothness} for details).
In the next section, I provide sufficient conditions for LICQ that explicitly ensure that the class of DGP under consideration is closed, so that we can uniformly control the quality of the delta-method inference within this class.
Then I illustrate these conditions in the context of Example~\ref{exa:runningExample}.
\subsubsection{Testable sufficient conditions for uniqueness of dual solutions}\label{subsec:UniqueDual}
LICQ is often considered hard to verify \citep[see][]{kaido2019constraint}.
In fact, a direct test of this assumption would face two problems: (i) the set $\setID$ is unknown but can only be estimated with an error; (ii) the set of binding inequalities at each $\theta$ is also unknown.
Because of their multiple-testing nature, both problems would complicate inference beyond the practical level.
Moreover, LICQ does not provide an explicit bound on the dual variables required for uniform validity analysis.
Some authors, including \cite{hsieh2017inference}, make a high-level assumption about boundedness of the dual variables.
Others (KMS and \cite{andrews2019inference}) focus on test inversion and thus only require restriction on non-degenerate covariance matrix of the moment functions.
Both KMS and \cite{andrews2019inference} compute critical values of a test at a particular point $\theta$ which they then invert. (\cite{kaido2019confidence} additionally propose a method for computationally efficient interpolation of confidence sets based on test inversion.)
KMS, for example, only need a bounded dual variable for each of the bootstrap draws for their purposes.
The local linear program that is used in KMS bootstrap has a bounded dual variables almost surely under the aforementioned covariance constraints since its inequality constraints have random bootstrapped coefficients.
In contrast, I need stronger assumptions to be able to estimate the nuisance parameter, the argmin set $\argminT(\mnP)$ and the dual solutions $\argminL(\mnP)$, which allows me to avoid the test inversion stage and achive compuational gains.
I propose a sufficient condition for LICQ in the form of bounds on the values of two auxiliary optimization programs.
The values of these programs, in turn, explicitly characterize an upper bound on the dual variables, allowing for a uniform asymptotic analysis.
Specifically, within this class of DGP satisfying this assumption, we can uniformly control the quality of the delta-method inference by establishing an explicit uniform bound on the relevant second-order directional derivatives.
To introduce the sufficient conditions for LICQ, I use the following notation.
For any $\jSet\subset\setIneq$, let the matrix $\activeMat=(e_{i_1},...,e_{i_\ell})^\prime$ correspond to $ \activeSet\bydef \setEq \cup \jSet = \{i_1,...,i_\ell\} $, a set of active constraints.
Let $\eta_1(\cdot)$ be the smallest left singular value function of a matrix; that is, $\eta_1(A) \bydef \sqrt{ \min_{u}( u^\prime A A^\prime u/ u^\prime u ) }$ .
The sufficient conditions can now be formulated as the following two assumptions on every submatrix $\activeMat (\Ap|\bp)$ that include all $p$ equality constraints and either $d-p$ or $1+(d-p)$ inequality constraints.
\begin{assumption}\label{ass:ULICQ} Measure $\measTrue$ satisfies two conditions:
\begin{enumerate}[label={\textbf{\Alph*.}},
ref={\theassumption.\Alph*}]
\item \label{assu:RankCondition} For any combination $\jSet$ consisting of all $p$ equality constraints and $d-p$ inequality constraints, the corresponding submatrices of coefficients $\activeMat (\Ap|\bp) $ of the full set of constraints $ (\Ap|\bp)$ have singular values that are uniformly bounded from below by a positive number $\eta(\measTrue)$.
\item \label{assu:NoOveridentification} Any combination $\jSet$ consisting of all $p$ equality constraints and $d-p+1$ inequality constraints cannot be simultaneously satisfied as equality at any point $\theta\in\setID$.
\end{enumerate}
\end{assumption}
Assumption \ref{ass:ULICQ} can be summarized using two characteristics:
\begin{align}
\eta(\measTrue)&\bydef\min_{\text{s.t. }\begin{matrix}
\jSet \subset\setIneq \\
\card\jSet = d-p
\end{matrix}}\eta_1\left(\activeMat (\Ap|\bp) \right) >0,\label{eq:rankD}\\
s(\measTrue)&\bydef\min_{\text{s.t. } \begin{matrix}
\jSet\subset \setIneq \\
\card\jSet = d-p+1\\
{\theta}\in \setID
\end{matrix}}\ltwo{\activeMat (\Ap\theta-\bp)} >0.\label{eq:GenericMomentConditions-1}
\end{align}
These numbers measure how close a given DGP $\measTrue$ to a violation of LICQ condition in population and determine how many observations are required to meet LICQ and have non-empty feasible set for a sample analog of program \ref{prog:regProg} with a given probability (see Lemma \ref{lem:wellDefined} in Appendix).
Both characteristics $\eta(\measTrue)$ and $ s(\measTrue)$ can be consistently estimated using a plug-in approach.
As long as a sample analog of $\setID$ is nonempty, sample analogs of $\eta(\measTrue)$ and $ s(\measTrue)$ will be generically positive.
In principle, a critical value can be obtained for a formal test of hypothesis $\eta(\measTrue)\geq \underline{\eta}$ and $s(\measTrue)\geq \underline{s}$ for any given pair of numbers $\underline{\eta},\underline{s}$ along the lines of \cite{cragg1997inferring}.\footnote{See also Appendix Remark \ref{rem:contastLICQ2} for an alternative representation of Assumption \ref{ass:ULICQ} as a single characteristic minimal bound on a different subset of submatrices of $ (\Ap|\bp)$.}
I leave this formal test for future research.
Assumption~\ref{ass:ULICQ} rules out more than $d$ binding inequality constraints at any point $\theta\in\setID$.
There are some special empirical applications resulting in a singleton identified set $\setID$ in which such \emph{overidentifying} inequality constraints can appear, which can result in multiplicity of dual solutions.\footnote{See, for example, \cite{gafarov2014identification} and \citet{shi2018estimating} who consider cases of infinitely many inequality conditions.}
This multiplicity can be eliminated using a regularization of the dual program \eqref{eq:DualProgram0}; that is, Assumption~\ref{assu:NoOveridentification} can be relaxed within the framework proposed in this paper, but such an extension is left for future research.
I briefly discuss this proposal in Appendix Section \ref{subsec:overidentification}.
It is instructive to see what Assumption~\ref{ass:ULICQ} implies for our running example.
\begin{example}[ continues=exa:runningExample ]\label{exa:runningExample5}
In this setup, there are four moment inequality conditions defined in \eqref{eq:subsystem}. They correspond to the following matrix of coefficients:
\begin{equation}
(\Ap|\bp) = \left(\begin{array}{cc|c}
\frac{1}{2} & \expGeneric\left[Z_{1}Z_{2}\right] & \expGeneric\left[\bar{Y} Z_{1}\right]\\
-\frac{1}{2} & -\expGeneric\left[Z_{1}Z_{2}\right] & -\expGeneric\left[\underline{Y}Z_{1}\right]\\
0 &\expGeneric\left[\left(1-Z_{1}\right)Z_{2}\right]& \expGeneric\left[\bar{Y}\left(1-Z_{1}\right)\right]\\
0 & -\expGeneric\left[\left(1-Z_{1}\right)Z_{2}\right]& -\expGeneric\left[\underline{Y}\left(1-Z_{1}\right)\right]
\end{array}\right)
\end{equation}
Here $p=0$, so to check Assumption~\ref{assu:RankCondition} one need to consider all subsets with two rows out of four, six combinations in total.
For example, the submatrix with rows $\jSet=\{1,2\}$ takes form
\begin{equation}
\activeMat (\Ap|\bp) = \left(\begin{array}{cc|c}
\frac{1}{2} & \expGeneric\left[Z_{1}Z_{2}\right] & \expGeneric\left[\bar{Y}Z_{1}\right]\\
-\frac{1}{2} & -\expGeneric\left[Z_{1}Z_{2}\right] & -\expGeneric\left[\underline{Y}Z_{1}\right]
\end{array}\right).
\end{equation}
This matrix has a full row rank iff $\expGeneric\left[\bar{Y}Z_{1}\right]\neq \expGeneric\left[\underline{Y}Z_{1}\right]$ or $\Delta_1>0$. As result, we just verified that for $\jSet=\{1,2\}$ we have $\eta_1\left(\activeMat (\Ap|\bp) \right) >0$.
Similarly, the matrix corresponding to rows $\jSet=\{3,4\}$ has full rank iff both $\Delta_0>0$ and $\expGeneric\left[\left(1-Z_{1}\right)Z_{2}\right]\neq0$.
Under those conditions, submatrices with pairs of rows $\jSet\in\big\{\{1,3\},\{1,4\},\{2,3\}\{2,4\}\big\} $ are also all full rank, and thus their left singular values are all positive.
To summarize, Assumption~\ref{assu:RankCondition} is satisfied iff
\begin{align}
&\Delta_0>0,\Delta_1>0, \label{eq:non-trivial-intervals}\\
&\expGeneric\left[\left(1-Z_{1}\right)Z_{2}\right]\neq0.\label{eq:noMulticollinearity}
\end{align}
What do these conditions mean in practical terms?
Inequalities~\eqref{eq:non-trivial-intervals} imply that the upper and lower bounds on $Y$ are different from each other on average, conditional on $Z_1$, while inequality~\eqref{eq:noMulticollinearity} implies that the instrument $Z_2$ (with values in $\{0,1\}$) is not perfectly correlated with $Z_1$.
In fact, in the case where the bounds $\underline{Y}$ and $\bar{Y}$ coincide with probability 1, one can replace the corresponding pair of moment inequalities with a single equality.
The degenerate case $\expGeneric\left[\left(1-Z_{1}\right)Z_{2}\right]=0 $ is analogous to the multicollinearity problem in the usual linear regression setup.
Inequalities~\eqref{eq:non-trivial-intervals} imply that Assumption~\ref{assu:NoOveridentification} is satisfied.
To illustrate a violation of Assumption~\ref{assu:NoOveridentification}, suppose that we add one more moment condition corresponding to instrument $Z_2$, $$\expGeneric\left[\underline{Y}Z_{2}\right]\leq \theta_{1}\expGeneric\left[Z_{1}Z_{2}\right]+ \frac{1}{2}\theta_{2}.$$
Assumption~\ref{assu:NoOveridentification} would be violated if, for example, this additional inequality was binding at the corner points of the original identified set, $\theta=(-\Delta_{1}\mp 2\rho\Delta_{0}, \pm \Delta_{0})$.
It is worth contrasting Assumption~\ref{ass:ULICQ} with LICQ:
the latter only restricts the gradients of the submatrices of $\Ap$ at any point where the corresponding constraints are active.
It turns out that instead of checking the active constraints at any given point, one can restrict the singular values of submatrices $\activeMat (\Ap|\bp)$.
In this simple example, we can manually verify that LICQ holds under Assumption~\ref{ass:ULICQ} (the general proof is given in Lemma \ref{lem:LICQ} in the Appendix Section \ref{subsec:LICQdiscussion}).
First, let us consider pairs of constraints with collinear gradients that potentially could violate LICQ.
For example, the submatrix of $\Ap$ corresponding to rows $\jSet=\{1,2\}$,
\begin{equation}
\activeMat \Ap = \begin{pmatrix}
\frac{1}{2} & \expGeneric\left[Z_{1}Z_{2}\right] \\
-\frac{1}{2} & -\expGeneric\left[Z_{1}Z_{2}\right]
\end{pmatrix},
\end{equation}
is always a reduced rank matrix.
Nevertheless, the corresponding two inequality constraints cannot be binding simultaneously by the Rouch\'e-Capelli theorem (also known as Kronecker–Capelli theorem) since $ \activeMat (\Ap|\bp)$ for them has rank equal to 2 which is larger than 1, the rank of $\activeMat \Ap$.
Correspondingly, these two constraints do not violate LICQ at any point (since they do not intersect) despite having collinear gradients.
Second, for constraint pairs like $\jSet=\{1,3\}$ we have $ \rank (\activeMat \Ap )= \rank ( \activeMat (\Ap|\bp)) = 2$, and hence the corresponding corner point (intersection of these constraints) exist as a unique solution to $(\activeMat \Ap) \theta = \activeMat \bp$. That corner point does not violate LICQ since $ \rank(\activeMat \Ap )=d=2$.
By this logic, we can prove that all the points in this example have at most 2 binding constraints, all of which have linearly independent gradients.
Third, the faces of the polygon are always defined by a subset of the inequalities defining their corners.
Hence the inequalities defining the faces are also linearly independent (a matrix of full rank has all submatrices of full rank) and any point on a face does not violate LICQ.
Finally, we do not need to check the LICQ in the internal points of $\setID$ since there are no binding constraints at those points by definition.
Thus, in this example, we just verified that LICQ indeed holds under Assumption \ref{ass:ULICQ}. $\QEDB$
\end{example}
Appendix Section \ref{subsec:LICQdiscussion} contains further technical details about Assumption \ref{ass:ULICQ}.
\subsubsection{Tighter bounds on the support function at a given direction} \label{subsec:tighterBounds}
Clearly, for any positive $\mu$, the value of the regularized program \eqref{prog:regProg} $\minV(\mu\measTrue)$ is larger than $\minV (\measTrue)$ by at least $\mu\ltwo{\argminT(\mu,\measTrue)}^{2}$.
A tempting approach would be to correct for this regularization bias by subtracting $\mu\ltwo{\argminT(\mu,\measTrue)}^{2}$.
Unfortunately, this correction makes the corrected value function $\minV(\mu,\measTrue)-\mu \ltwo{\argminT(\mu,\measTrue)}^{2}$ non-differentiable in the paramteters $A_P,b_P$ if the non-regularized value function $\minV (\measTrue)$ itself is non-differentiable.
So to preserve differentiability of the bias-corrected value function, which is crucial for delta-method inference, one cannot use argmin of the actual regualrized programm for bias corrections.
Instead, I suggest to tighten the bounds on $\minV (\measTrue)$ using one the following two corrections,
\begin{align}
\minVin(\mu,\kappa,\measTrue) &\bydef \minV(\mu,\measTrue)-\mu \ltwo{\argminT(\kappa,\measTrue)}^{2},\label{eq:innerBound} \\
\minVout(\mu,\measTrue) &\bydef\minV(\mu,\measTrue)-\mu \ltwo{\theta^*}^{2},\label{eq:outerBound}
\end{align}
where $\argminT(\kappa,\measTrue)$ is an argmin of the regularized program \eqref{prog:regProg} with a larger tuning parameter $\kappa$ instead of $\mu$, $\theta^*$ is any point in $ \argminT(\measTrue)$.
As $\mu$ and $\kappa$ shrink to zero, these bounds continuously shrink to $\minV (\measTrue)$ above and below, respectively. The expressions in \eqref{eq:innerOuterBounds} provide valid conservative bounds from above and below for $\minV(\measTrue)$ that are useful for uniform one-sided inference (coverage probability in this case can be higher than nominal level). For the remainder of the paper, the focus will be on the confidence bounds that cover $\minV (\measTrue)$ from below. So $\minVout(\mu,\measTrue) $ will play the major role. The uses of $\minVin(\mu,\kappa,\measTrue)$ for uniform inference are discussed in Appendix Section \ref{subsec:overidentification}
\begin{thm}
\label{thm:bounds}For any DGP parameterized by $\measTrue$ satisfying Assumption~\ref{assu:Non--empty} and any $\kappa\geq \mu \geq 0$, the following bounds hold:
\begin{equation}
\minVout(\mu,\measTrue) \leq \minV(\measTrue) \leq \minVin(\mu,\kappa,\measTrue). \label{eq:innerOuterBounds}
\end{equation}
If, in addition, $\measTrue$ satisfies Assumption ~\ref{ass:ULICQ}--\ref{assu:NoOveridentification}, then there exist $\bar{\mu}(\measTrue)>0$ such that $\minVin(\mu,\kappa,\measTrue)=\minV(\measTrue)$
for any fixed $\mu<\kappa<\bar{\mu}(\measTrue)$. Furthermore, if $\argminT(\measTrue)$ is a singleton, then $\minVout(\mu,\measTrue)=\minV(\measTrue)$ for any $\mu<\bar{\mu}(\measTrue)$.
\end{thm}
\begin{proof}
See Appendix Section~\ref{sec:Proof-of-Theorem1}.
\end{proof}
Although for a given DGP $P$, characterized by $(A_P,b_P)$, the cutoff $\bar{\mu}(\measTrue)$ is well defined, it can change discontinuously as a result of a small change in $(A_P,b_P)$. It makes a consistent estimation of $\bar{\mu}(\measTrue)$ a challenging task.
So instead of trying to estimate $\bar{\mu}(\measTrue)$, I suggest using two appropriately chosen shrinking sequences $\mu_n$ and $\kappa_n$ instead of fixed tuning parameters $\mu$ and $\kappa$.
The gap between $\minV(\measTrue)$ and $\minVout(\mu_n,\measTrue) $ is at least $\mu_n (\ltwo{\theta^*}^{2}-\ltwo{\argminT(\mu_n,\measTrue)}^2)\geq0$; for $ \minVin(\mu_n,\kappa_n,\measTrue) $ the gap is at most $\mu_n (\ltwo{\argminT(\mu_n,\measTrue)}^{2}-\ltwo{\argminT(\kappa_n,\measTrue)}^2)\leq0$.
In the nonregular case, that is where $\argminT(\measTrue)$ is non-singleton, the gap for $\minVout(\mu_n,\measTrue) $ is shrinking to $0$ at rate $\mu_n$ which still results in a consistent estimators of $\minV (\measTrue)$ at rate $\mu_n$ that can be slightly slower than the regular parametric rate $1/\sqrt{n}$ (for example, $ {\sqrt{\ln{n}}}/{\sqrt{n}}$).
As a result, the corresponding confidence sets cover the projections of the sharp identified set, not an enlargement of it (in contrast to other studies, including \cite{cho2023simple}, who do not use consistent estimators of the sharp bounds on the identified set and study an enlarged identified set instead).
The gap for $\minVin(\mu_n,\kappa_n,\measTrue)$ becomes exactly equal to zero for some sufficiently large $n$ in both regular and nonregular cases, resulting in exact corresponding point-wise one-side confidence sets.
\subsection{\label{subsec:inference} Large-sample results for the estimated regularized support function}
Consider the analog estimator of $\funLVal[\mu_{n}]$ for some sequence $\mu_{n}$,
\begin{equation}
\funLValHat[\mu_{n}]\bydef\underset{ {\theta}\in\setID[\measEmp]}{\mbox{min}}\left\{\eOne {\theta}+\mu_{n}\ltwo{ {\theta}}^{2}\right\} \label{prog:minSSAregularized}.
\end{equation}
As we shall see in this subsection, this estimator admits an asymptotic linear (Bahadur-Kiefer) expansion with an explicit approximation error, which determines the precision of asymptotic normal and bootstrap inference.
This allows me to establish the validity of the corresponding inference methods \emph{uniformly} over a reasonable class of DGP.
Uniform asymptotic validity is crucial for small-sample performance of inference methods in the presence of possible discontinuous changes of the asymptotic distribution of estimators (for example, $\funLValHat[\mu_{n}]$) with respect to DGP parameters (for example, $\measTrue$).
\subsubsection{A class of DGPs under consideration}
I consider the class of all measures $\Measures=\Measures(\underline{\eta},\underline{s},\varepsilon,\bar{M})$ that satisfy Assumptions~\ref{assu:Non--empty}-\ref{assu:Moments} with some positive constants $\underline{\eta}$, $\underline{s}$, $\varepsilon$, and $\bar{M}$.
\begin{assumption}
\label{assu:Moments}There exist $\varepsilon>0$ and $\bar{M}<\infty$ such that
\begin{align}
\expect\left\Vert \Data\right\Vert & ^{2+\varepsilon} <\bar{M}.
\end{align}
\end{assumption}
To summarize, every $\measTrue\in \Measures$ satisfies $\setID\neq \emptyset$,
$\eta(\measTrue)>\underline{\eta}$, $s(\measTrue)>\underline{s}$, and $\expect\left\Vert \Data\right\Vert^{2+\varepsilon} < \bar{M}^{2+\varepsilon}$. \label{reply:R3.5.c.precompact}
The class $\Measures$ includes DGPs (parameterized by measure $\measTrue$) with multiplicity of primal solutions of \eqref{prog:minP} but assumes uniqueness of dual solutions.
This property is necessary to show that program~\eqref{prog:minSSAregularized} has a nonempty sample argmin and a unique vector of Lagrange multipliers in large enough samples with probability approaching 1 uniformly in $\measTrue\in\Measures$ as $n\to\infty$ (see Lemma~\ref{lem:wellDefined} in Appendix Section \ref{sec:Proof-of-Theorem2}).
As discussed in Section~\ref{subsec:UniqueDual}, one can study a case with multiple dual solutions and a unique primal solution in a similar way (see Appendix Section~\ref{subsec:overidentification}).
\subsubsection{A higher-order envelope theorem and a strong approximation of the regularized support function}
A Bahadur-Kiefer expansion of estimator defined in program~\eqref{prog:minSSAregularized} can be established using the envelope theorem \eqref{prog:minmaxDerivative}.
As discussed in Section~\ref{subsec:asymptoticDistribution}, this theorem typically holds for the sample support function at $e_1$ .
However, this theorem does not provide bounds on the higher-order directional derivatives.
So, without additional assumptions, the directional delta method of \cite{shapiro1991asymptotic} does not provide means to evaluate the error of the asymptotic approximation of the sample support function at $e_1$ by \eqref{eq:asyDistributionShapiro}.
Since the limiting distribution changes discontinuously with the DGP $\measTrue$, poor performance of asymptotic methods based on the limit \eqref{eq:asyDistributionShapiro} should be expected when the parameter $\measTrue$ is close to a discontinuity point.
In other words, inference procedures that estimate \eqref{eq:asyDistributionShapiro} directly (for example, subsampling or directional bootstrap) will inevitably be only point-wise valid and can have poor finite-sample performance in such cases. (It is an implication of the impossibility theorem of \citet{hirano2012impossibility} for a functional that is directionally differentiable.)
In contrast, the regularized support function at $e_1$ admits a stronger version of the envelope theorem that provides a bound on the error of the linear approximation uniformly over $\measTrue\in\Measures$.
I developed this novel bound using a \emph{second-order directional Taylor expansion} of a system of generalized inequalities that define the optimal solutions to the regularized program (see Lemmas~\ref{lem:directionalDerivative} and \ref{lem:boundOnVgrowth} in Appendix Section \ref{app:smoothness}).
Using this result, the asymptotic linear (Bahadur-Kiefer) representation of the value function in program~\eqref{prog:minSSAregularized} follows almost immediately (see Lemma~\ref{lem:Bahadur} in Appendix Section \ref{sec:Proof-of-Theorem2}):
\begin{equation}
\sqrt{n} ( \funLValHat[\mu_{n}]- \minV\mnP ) = \frac{1}{\sqrt{n}}\sum^{n}_{i=1}\argminL\mnP^\prime\momentVec{\data_i}{\argminT\mnP} +O_\Measures(\frac{1}{\mn \sqrt{n}}). \label{eq:Bahadur}
\end{equation}
The first term on the right-hand side of this representation is a scaled sample average of a zero-mean random variable, which admits a uniform Gaussian approximation. (Indeed, the binding constraints have zero mean at $\argminT$, while the nonbinding constraints are multiplied by zero dual variables $\argminL$.)
The residual term $O_\Measures(1)$ denotes a uniformly tight sequence of a random process indexed by $\measTrue\in\Measures$. Analogously, I denote any random sequence as $ o_\Measures(1)$ if
$\lim_{n\to\infty}\sup_{\measTrue \in \mathcal{P}}\measTrue\left( \ltwo{\zeta_n(\measTrue) }\geq \epsilon \right) =0$.
The Bahadur-Kiefer representation \eqref{eq:Bahadur} is important for three reasons.
First, it suggests a coupling of $\funLValHat[\mu_{n}]$ with a Gaussian process (through the Yurinsky theorem).
Second, it implies uniform validity of the optimization-free score bootstrap, which is particularly computationally convenient.
Third, it suggests an analog estimator of the asymptotic variance of the regularized support-function estimator,
\begin{equation}
\funLStdHat[\mu_{n}][][2]=\sMean{(\funLLagrangeHat[\mu_{n}][][\prime]\momentVec{\data_i}{\funLArgHat[\mu_{n}][]})^2}, \label{eq:standardErrorFormula}
\end{equation}
where $\funLArgHat[\mu_{n}]$ and $\funLLagrangeHat[\mu_{n}]$ are, respectively, the optimum and the vector of Lagrange multipliers of (\ref{prog:minSSAregularized}), which are provided by common constraint-optimization software packages.
These implications are summarized in the following theorem.
\begin{thm}
\label{thm:consistency}Consider any sequence $\mu_{n}$ such that
$\mu_{n}\to0$ and $\mu_{n}\sqrt{n}\to\infty$. Then with probability approaching 1 uniformly in $\measTrue\in \Measures$,
\begin{align}
\lim_{n\to\infty} \sup_{\measTrue \in \mathcal{P}}&\pi(\sqrt{n}(\funLValHat[\mu_{n}]-\funLVal[\mu_{n}][\measTrue]), N\left(0,\funLStd[\mu_{n}][\measTrue][2]\right)) = 0,\label{eq:coupling}\\
\funLArgHat[\mn] &=\funLArg[\mn] + O_\Measures(\frac{1}{\mu_n\sqrt{n}}),\label{eq:thetaConsistency}\\
\funLLagrangeHat[\mn]&=\funLLagrange[\mn] + O_\Measures(\frac{1}{\sqrt{n}}),\label{eq:lambdaConsistency}\\
\funLStdHat[\mn]& = \funLStd[\mn] + o_\Measures(1). \label{eq:convergenceSpeedStd}
\end{align}
The function $\pi(\cdot,\cdot)$ is the Levy-Prohorov metric, which metricizes the weak topology of probability measures (see \cite{van1996weak}).
\end{thm}
\begin{proof}
The strong approximation result \eqref{eq:coupling} is based on the uniform bound on the higher-order directional derivatives (implied by Assumptions \ref{assu:Non--empty},\ref{ass:ULICQ}) and the generalization of the \cite{yurinskii1978error} coupling proposed in \citet[][Proposition A.5.2 on p. 457]{van1996weak} ( it is implied by \ref{assu:Moments} for i.i.d. data).
See Appendix Section~\ref{sec:Proof-of-Theorem2} for details.
\end{proof}
The coupling with a Gaussian random process \eqref{eq:coupling} is a \emph{strong approximation} result, which can be understood using a geometric interpretation.
The distance between the difference $ \sqrt{n} ( \funLValHat[\mu_{n}]- \minV\mnP )$ and some sequence of zero-mean Gaussian r.v.s $ N\left(0,\funLStd[\mu_{n}][\measTrue][2]\right)$
converges to zero with the same \emph{uniform} rate for all DGP $\measTrue \in \mathcal{P}$.
This property is stronger than conventional CLT-type results, since it does not require the existence of a limiting distribution for the approximating Gaussian variables.
\subsection{Uniformly valid inference\label{subsec:pointwise-ci}}
Theorems~\ref{thm:bounds} and \ref{thm:consistency} can be used to construct uniformly valid confidence bands (one-sided CS) for $\minV(\measTrue)$.
The corresponding generic algorithm takes the following form: \vspace{0.5cm}
\begin{enumerate}
\item[Step 1.] Compute the regularized sample support function $\funLValHat[\mu_{n}]$ at $e_1$ defined in \eqref{prog:minSSAregularized}.
\item[Step 2.] Compute the standard error using \eqref{eq:standardErrorFormula}.
\item[Step 3.] Compute a bias adjustment using sample analogs of either $\mu_{n} \ltwo{\argminT(\kappa_{n},\measTrue)}^{2} $ (for an upper confidence band) or the norm of some point in the argmin set $\mu_{n} \ltwo{\theta^*}^{2}$ (for a lower confidence band).
\end{enumerate}
\vspace{0.5cm}
The following subsections explain specific implementations of this algorithm for uniformly valid CSs for a projection of the identified set on a single coordinate or multiple coordinates and for the argmin set of a linear program with estimated coefficients.
\subsubsection{Application to confidence sets for a scalar projection of the identified set } \label{sec:DriftingDGP}
Let's revisit one of the primary objects of interest in the moment-inequality models: CSs on projections of the identified set, $\setMarginal\bydef\left[\funLValO,\funUValO\right]$ .
Theorems~\ref{thm:bounds} and \ref{thm:consistency} suggest that the outer-bound estimator for the minimal value,
\begin{equation}
\minV^{out}(\mn,\measEmp) \bydef\funLValHat[\mu_{n}]-\mu_{n}\ltwo{\theta^*\left(\measEmp\right)}^{2},\label{eq:biasCorrectedEstimator}
\end{equation}
is asymptotically unbiased (in the regular case) or biased downward (in the nonregular case).
(Some downward bias is acceptable since the purpose of CSs is to cover the lowest point $\funLValO$ from below.)
To achieve this property, the estimator $\theta^*\left(\measEmp\right)$ should have a norm that is not smaller than the minimal norm in $ \argminT(\measTrue)$ with probability approaching 1.
We can use an estimator $\theta^*\left(\measEmp\right)$ that converges to the point with the coordinates
\begin{equation}
\theta^*_i(\measTrue) \bydef \max \{ {\theta^{+}_i(\measTrue)}, {\theta^{-}_i(\measTrue)}\},
\end{equation}
where
\begin{equation}
{\theta}^{\pm}_i(\measTrue)\bydef\abs{ \underset{ {\theta}\in\setID[\measTrue], {\theta_1}\leq\funLValO[\measTrue]+\mn} \min\left\{\pm {\theta_i}\right\} }.\label{eq:defThetaStar}
\end{equation}
By definition, $\ltwo{\theta^*}\geq \ltwo{\theta}$ for any $\theta \in \argminT(\measTrue)$. (See the proof of consistency in Lemma~\ref{lem:argminBound} in Appendix A.)
This bound on $\ltwo{\theta^*}$ has two attractive properties. First, it can be (uniformly) consistently estimated using only $2k$ linear programs, which can be easily computed even in models with a very large-dimension $k$ and a large number of inequalities using interior-point numerical optimization methods.
Second, in the regular case, in which $\argminT(\measTrue)$ is a singleton, this estimator is asymptotically unbiased since it converges to $ \ltwo{\theta^*} = \ltwo{\argminT(\measTrue)} $.
So for any such fixed $\measTrue$, by Theorem~\ref{thm:bounds} we have $\minV^{out}(\mu_{n},\measTrue)=\minV(\measTrue)$ for sufficiently small $\mu_n$; that is, it is possible to construct CIs with correct (nonconservative) coverage in the regular case.
Using the bias-corrected estimator $\minV^{out}(\mn,\measEmp)$ and its analog for the upper bound, $\maxV^{out}(\mn, \measEmp)$, I construct the following delta-method CSs:
\begin{equation}
\begin{cases}
\mbox{CB}_{\alpha,n,\Measures} & =\left[\minV^{out}(\mn,\measEmp)-z_{1-\alpha}n^{-1/2} \hat{ { \underline{\sigma}}}^{reg}_n,\infty\right),\\
\mbox{CI}^{\theta_1}_{\alpha,n,\Measures} & =\left[\minV^{out}(\mn,\kappa_n,\measEmp)-z_{1-\alpha}n^{-1/2}\hat{ { \underline{\sigma}}}^{reg}_n; \maxV^{out}(\mn,\kappa_n,\measEmp)+z_{1-\alpha}n^{-1/2}\hat{ { \bar{\sigma}}}^{reg}_n\right],\\
\mbox{CI}_{\alpha,n,\Measures}^{\mathcal{S}} & =\mbox{CI}^{\theta_1}_{\alpha/2,n,\Measures}.
\end{cases}\label{eq:CIsPoint}
\end{equation}
Here, $z_{1-\alpha}$ is $1-\alpha$ quantiles of the standard Gaussian distribution and
\begin{align}
\hat{ {\underline \sigma}}^{reg}_n &\bydef \max\{\funLStdHat[\mu_{n}][][],\sigma_0\},\\
\hat{ { \bar\sigma}}^{reg}_n &\bydef \max\{\funUStdHat[\mu_{n}][][],\sigma_0\}, \label{eq:sigma0reg}
\end{align}
for some small positive number $\sigma_0$.
$\mbox{CB}_{\alpha,n,\Measures}$ is a one-sided CB for $\funLValO$ (and for $\theta_1$ as well), $\mbox{CI}^{\theta_1}_{\alpha,n,\Measures}$ is a two-sided CI that covers any $\theta_1$ in the identified set, and $\mbox{CI}_{\alpha,n,\Measures}^{\mathcal{S}}$ is a two-sided CI (based on the Bonferroni inequality) that covers the entire maginal identified set $\setMarginal$.
The tuning parameter $\sigma_0$ is introduced to guarantee the nominal coverage in the cases in which only nonstochastic constraints are binding at the optimal solutions corresponding to the support functions; otherwise, this degeneracy can lead to superconsistent support-function estimators with a non-Gaussian limiting distribution.
One can set $\sigma_0=0$ if this concern is not appropriate in a particular application.
As before, let $\Measures$ contain all measures $\measTrue$ that satisfy Assumptions~\ref{assu:Non--empty}--\ref{assu:Moments} with some uniform positive constants $\underline{\eta}$, $\underline{s}$, $\varepsilon$, and $\bar{M}$.
\begin{thm}
\label{thm:TheoremUniform}
Suppose that $0<\alpha<1/2$, $\mu_{n}\to0$, and $\mu_{n}\sqrt{n}\to\infty$.
Then the following results hold:
\begin{eqnarray*}
\liminf_{n\to\infty}\inf_{\measTrue\in\Measures} \measTrue\left(\setMarginal\subset\mbox{CB}_{\alpha,n,\Measures} \right)=\liminf_{n\to\infty}\inf_{\measTrue\in\Measures}\inf_{ {\theta}\in\setID}\measTrue\left(\theta_{1}\in\mbox{CB}_{\alpha,n,\Measures} \right) & \geq & 1-\alpha, \\
\liminf_{n\to\infty}\inf_{\measTrue\in\Measures}\measTrue\left(\setMarginal\subset\mbox{CI}_{\alpha,n,\Measures}^{\mathcal{S}}\right)\geq1-\alpha,\text{ }\liminf_{n\to\infty}\inf_{\measTrue\in\Measures}\inf_{ {\theta}\in\setID}\measTrue\left(\theta_{1}\in\mbox{CI}_{\alpha,n,\Measures}^{\mathcal{S}}\right) & \geq & 1-\alpha.
\end{eqnarray*}
\end{thm}
\begin{proof}
See Appendix \ref{sec:Proof-of-Theorem2}.
\end{proof}
Note that the worst-case (that is, the smallest) asymptotic coverage probability of $\mbox{CB}_{\alpha,n,\Measures}$ is exactly equal to $1-\alpha$ for a fixed regular DGP such that $\argminT(\measTrue)$ is a singleton and $\lim_{n\to\infty}\funLStd[\mu_{n}][\measTrue][2]>\sigma_0\geq 0 $.
For a nonregular DGP (or sequence of DGPs), the asymptotic coverage can only be larger than $1-\alpha$.
How conservative is $\mbox{CB}_{\alpha,n,\Measures}$ in the nonregular case?
We can evaluate it by comparing its average length with that of a confidence bound that has exact point-wise asymptotic coverage probability in a Monte Carlo study (see Section \ref{sec:Monte-Carlo}).
Theorem \ref{thm:coveragePW} in Appendix Section \ref{subsec:pointwise-results} provides an alternative bias correction that remains non-conservative in the non-regular case and results in shorter point-wise confidence bounds.
This approach is based on using an estimator of the inner bound $\minVin(\mu,\kappa,\measTrue)$.
I conclude this section with a brief discussion of the tuning parameters.
The theory of optimal choice of tuning parameter is beyond the scope of this paper. The following considerations, however, can provide some guidance for the optimal choice.
Theorem~\ref{thm:bounds} suggests that the tuning parameters should be smaller than $ \bar{\mu}(\measTrue)$ to avoid the bias in the first-order asymptotic distribution.
This choice is infeasible since $\bar{\mu}(\measTrue)$ is unknown, so one has to let $\kappa_n$ and $\mn$ go to zero.
The optimal rates of $\kappa_n$ and $\mn$ should balance the higher-order variance and the worst-case bias.
A specific choice of the tuning parameters is discussed in Section \ref{sec:Monte-Carlo}.
\subsubsection{Joint confidence sets for general subvectors \label{subsec:subvectors}}
It is trivial to extend the analysis to
$$\underline{v}(\measTrue;a)=\displaystyle\min_{ {\theta}\in\setID} a^\prime {\theta}$$
for any $a\in R^d$ with $\ltwo{a}=1$.
Indeed, Assumptions~\ref{assu:Non--empty}--\ref{assu:Moments}
are invariant with respect to orthogonal transformations of the coordinates; that is, they are satisfied for the following program (with $\tilde{\theta}=U^\prime \theta$, $\tilde{\Ap} = \Ap U$, and $a^\prime= \eOne U $ for any orthogonal matrix $U$):
\begin{align}
\underline{v}(\measTrue;a) &=\min_{ {\theta}\in\R^{d} } \eOne \tilde{\theta} \\
\text{s.t. } & \begin{cases}
e_j^\prime\tilde{\Ap} \tilde{\theta} =e_j^\prime\bp, & j\in\setEq,\\
e_j^\prime\tilde{\Ap} \tilde{\theta} \leq e_j^\prime\bp, & j\in\setIneq.
\end{cases}
\end{align}
We can think of $\tilde{\Ap}$ as a coefficient matrix under a different measure $\tilde{\measTrue}$, $A_{\tilde{\measTrue}}$.
The set of measures $\Measures$ from Section~\ref{sec:DriftingDGP} includes $\tilde{\measTrue}$ corresponding to all orthogonal transformations of $\Ap$.
The identified set $\setID$ is convex, so any projection of it can be characterized using support functions for a corresponding direction. One can construct a joint CS for $\setID$ as follows. For any set of directions $\mathcal{A}\subset R^d$, take
$$ \text{CS}^\mathcal{A}_{\alpha,n} = \{\theta| a\in\mathcal{A} , a^\prime \theta \leq -\minV^{out} (\mn,\measEmp;-a) + c_{1-\alpha}n^{-1/2} \max\{\underline{\sigma}(\mn,\measEmp;-a),\sigma_0\}\},$$
where $c_{1-\alpha}$ is $1-\alpha$ quantiles of the maximum of the corresponding asymptotic Gaussian variables (also known as $\sup t$ statistics) that can be estimated using multiplier bootstrap enabled by the asymptotic linear representation,
\eqref{eq:Bahadur}.\footnote{
I leave full analysis of mulitplier bootstrap procedure in this setup for future work; see additional discussion in Appendix Section \ref{sec:multiplier boostrap}.}
One can also use the Bonferroni inequality-based standard Gaussian critical value $c_{1-\alpha}= z_{1-\alpha/|\mathcal{A}|}$.
Choosing the set of directions appropriately $\mathcal{A}$, we can construct joint CSs for projections of $\setID$ on any subvectors $\theta$.
If $\mathcal{A}$ has finitely many elements, $ \text{CS}^\mathcal{A}_{\alpha,n}$ is a polygon.
So we can plot it directly without performing test inversion as in the one-dimensional case.
The confidence set $\mbox{CI}_{\alpha,n,\Measures}^{\mathcal{S}}$ is a particular case of $ \text{CS}^\mathcal{A}_{\alpha,n}$ corresponding to $\mathcal{A}=\{e_1,-e_1\}$ and the Bonferroni estimate of $c_{1-\alpha}$.
It seems natural to construct a joint CS for $\theta$ based on directions that correspond to the normal vectors of the moment conditions.
For simplicity, assume that $p=0$.
The original system \eqref{eq:GenericMomentConditions} may have some inequalities that are slack for any point $\theta\in\setID$.
We can characterize the identified set $\setID$ as the solution to a tight system of inequalities
\begin{equation}
e_j^\prime \Ap \theta \leq \underline{b}_j, j\in\setIneq, \label{eq:tightIDSet}
\end{equation}
where
\begin{align}
\underline{b}_j &= \max_{ {(\vartheta,\theta)}\in\R^{d+1} }\vartheta \label{eq:tightBoundsProg} \\
\text{s.t. } & \begin{cases}
\vartheta &= e_j ^\prime \Ap \theta, \\
e_\ell^\prime{\Ap} {\theta} & \leq e_\ell^\prime\bp, \ell\in\setIneq.
\end{cases}
\end{align}
Every inequality in system \eqref{eq:tightIDSet} is active at least at one point in $\setID$ (any point in the argmax of \eqref{eq:tightBoundsProg}). Programs \eqref{eq:tightBoundsProg} meet Assumptions~\ref{assu:Non--empty}--\ref{assu:Moments}, and therefore the outer estimators $ \hat{\underline{b}}_j^\text{out}$ are half-median unbiased with the corresponding standard error estimators $\hat{\sigma}_j$. Then the following polyhedron CS will cover any point $\theta\in\setID$ with asymptotic probability of at least $1-\alpha$ uniformly over $\measTrue\in\Measures$,
$$ \text{CS}^\mathcal{N}_{\alpha,n} =\big \{\theta \big|
e_j^\prime \hat{A}_P \theta \leq \hat{\underline{b}}_j^\text{out} + z_{1-\frac{\alpha}{k}} n^{-1/2} \max\{\hat{\sigma}_j,\sigma_0\} , j\in\setIneq \big\},$$
where $z_{1-\frac{\alpha}{k}}$ is the critical value of the standard normal r.v. and $k$ is the number of (in)equality restrictions.
The generalization to the case $p\neq0$ is straightforward.
\section{Monte Carlo experiments}\label{sec:Monte-Carlo}
\subsection{Overview}
In this Monte Carlo study our goal is to evaluate how the confidence bound's length and the corresponding coverage probability depends on a number of key factors: (i) how close gradients of the relevant moment inequality to the non-regular case (i.e. being collinear with the vector $e_1$) ; (ii) how tuning parameter choice affects the conservativeness of confidence bounds; (iii) how the dimension of the problem and the number of inequalities affect length, coverage and computational time for the proposed methods.
The last exercise in the list also involves a comparison with the AS projection implemented using KMS EM computational algorithm.
\subsection{Proximity to non-regular case}
To illustrate the advantages of uniform coverage compared to point coverage, we can study a simple design with $k=4$ moment inequalities where the angle $\omega$ between the gradient of one of the faces of the identified set and vector $e_1$ varies continuously.
The expectation of the moment inequalities can be parametrized as follows
\begin{align*}
\expect W &= \expect (\Ap|\bp)=\left(\begin{array}{cc|c}
-\cos(\NEangle) & -\sin\left(\NEangle \right)& \cos(\NEangle) + \sin\left(\NEangle \right) \\
\cos(\NEangle) & \sin(\NEangle)& \cos(\NEangle) + \sin\left(\NEangle \right) \\
0 & -1 & 1 \\
0 & 1 & 1
\end{array} \right) .
\end{align*}
The shape of this set is a parallelogram analogous to the one in Figure~\ref{fig:Identified-sets-Example2} in with $\rho\geq0$.
The parameter $\omega\in[0^{\circ},36^{\circ}]$ defines the angle between the normal vectors of the rear sides of the parallelogram and the horizontal axis.
The value $\omega=0^{\circ}$ corresponds to a square-shaped identified set, also referred to as a \emph{non-regular case} because the sides are orthogonal to the gradient of the objective function.
All other values are termed \emph{regular}.
In the vicinity of $\omega=0^{\circ}$ pointwise valid $\mbox{CB}_{\alpha,n}$ may have a coverage probability below the nominal level because it lacks uniform validity.
The expectation $\expect W$ is parameterized to guarantee $\argminT_1=\argminT_2=-1$ and $\argmaxT_1=\argmaxT_2=1$ for all values of $\omega$.
The components of $W_i$ are independent Gaussian random variables with variance $s_2^{2}=0.01$.
For each value of $\omega$, I compute the frequency of coverage and the excess average length over the identified set for $\mbox{CB}_{\alpha,n}$ and $\mbox{CB}_{\alpha,n,\Measures}$ based on the sample sizes $n\in\{100,10000\}$.
The number of MC simulations is $1000$ for every combination of $n$ and $\omega$.
The focus is on the nominal coverage probability $\alpha=0.95$.
As the main choice of tuning parameters, I use $\mu_n = \hat\mu_0 \sqrt{n^{-1}\ln\ln n}$ and $\kappa_n = \hat\mu_0 \sqrt{n^{-1}\ln n}$ , where $$\hat\mu_0 = \sqrt{\frac{1}{n}\sum_{i=1}^n \big( \argminL^\prime(0,\measEmp)\data_i e_1 - \frac{1}{n}\sum_{i=1}^n \argminL^\prime(0,\measEmp)\data_i e_1 \big)^2
}.$$
While a theory of optimal choice of $\hat\mu_0$ is beyond the scope of this paper, this particular choice has some advantages: (i) $\hat\mu_0$ depends only on the behavior of the relevant moment inequalities as selected by non-zero components of $ \argminL^\prime(0,\measEmp)$; (ii) it does not depend on the dimension of the parameters $\theta$, since only the first column of $w_i$ is involved; (iii) it has the same scale as the standard deviation of the relevant components of $w_i$ which justifies our perturbation analysis (the impact of the regularization term $\mu_n \ltwo{\theta}^2$ is larger than the sample variation, $\sim \hat\mu_0 n^{-1/2}$).
As a result, this choice resulted in a good alignment of the theoretical predictions with the simulations, as can be seen below.
To appreciate the importance of uniformly valid confidence bounds, we start our study with point-wise valid confidence bounds $\mbox{CB}_{\alpha,n}$.
Figure \ref{fig:pwFreq_opt1} panel (a) shows the corresponding coverage frequency for a small sample size $n=100$ and a large sample with $n=10,000$.
For values of $\omega$ that are far from $0$, both sample sizes have observed a frequency of covering the correct bound of the identified set for $\theta_1$ close to the nominal coverage probability of $\alpha=0.95$.
The same is true for the nonregular case with $\omega = 0$.
However, in the vicinity of the nonregular case $\omega \in (0,4^\circ]$, large sample sizes are necessary to achieve a satisfactory coverage frequency.
In contrast, the coverage frequency for the corresponding sample sizes for the uniformly valid confidence bounds $\mbox{CB}_{\alpha,n,\Measures}$ given on Figure \ref{fig:pwFreq_opt1} panel (b) is uniformly at least as large as $\alpha$ for all values of $\omega$.
One can see that the only point on Figure \ref{fig:pwFreq_opt1} panel (b) with coverage is significantly higher than $\alpha=0.95$ for $n=10,000$ is $\omega = 0$, that is, for most designs (or equivalently, directions of the support function), the uniformly valid confidence bounds $\mbox{CB}_{\alpha,n}$ have exact coverage.
The difference in coverage frequencies of $\mbox{CB}_{\alpha,n}$ and $\mbox{CB}_{\alpha,n,\Measures}$ can be understood by considering the behavior of the corresponding average excess length of the bounds, i.e. the Monte Carlo average of the difference between the corresponding confidence bounds and the true bound of the identified set $\argminT_1=-1$.
Figure \ref{fig:pwFreq_opt1} panel (c) compares the lengths of the two bounds for the small sample case $n=100$ where the point-wise inference becomes unreliable.
One can see that the confidence bounds should have an excess length approximately equal to $0.025$ to match the coverage frequency with the nominal probability $\alpha$ in the proximity of $\omega =0$.
The length of $\mbox{CB}_{\alpha,n,\Measures}$ is larger than necessary for $\omega=0$, while $\mbox{CB}_{\alpha,n}$ has the correct length for $\omega=0$, but is too short in the neighborhood $ (0,4^\circ]$.
For values $\omega\geq 5$ both confidence bounds have approximately the same length, which also corresponds to the correct coverage probability predicted by Theorems \ref{thm:bounds} and Appendix Theorem~\ref{thm:coveragePW}.
\subsection{Sensitivity to choice of $\mu_n$ and $\kappa_n$} \label{subsec:sensitivityMC}
First, consider the effect of the tuning parameter $\mu_n$ on the coverage probability for the uniformly valid confidence bounds $\mbox{CB}_{\alpha,n,\Measures}$.
Figure \ref{fig:lengthCompare_opt1} panel (a) compares the coverage frequency of $\mbox{CB}_{\alpha,n,\Measures}$ for the sample size $n=100$ as a function of $\omega$ for two choices: (i) $\hat\mu_0 \sqrt{n^{-1}\ln\ln n}$ (baseline) and (ii) $\hat\mu_0 \sqrt{n^{-1}\ln n} $ (large). The baseline option has slightly less conservative coverage in the neighborhood of the nonregular design ($\omega=0$) and similar performance for other values of $\omega$.
As a result, the baseline option is used for all the other simulations.
Unlike its uniform counterpart, the point-wise valid confidence interval $\mbox{CB}_{\alpha,n}$ depends on additional tuning parameter sequence $\kappa_n$.
Theorem \ref{thm:coveragePW} claims that for one-sided bounds $\mbox{CB}_{\alpha,n}$, the asymptotic coverage is exactly equal to $\alpha$ as long as $\kappa_n$ shrinks to zero slower than $\mu_n$.
This allows values of $\kappa_n$ to be smaller than $\mu_n$ and still result in (point-wise) valid inference.
Figure \ref{fig:lengthCompare_opt1} panel (b) compares three choices: (i) $\kappa_n= \hat\mu_0 \sqrt{n^{-1}\ln n} > \mu_n$; (ii) $ \kappa_n=\mu_n$; (iii) $ \kappa_n=0<\mu_n$.
From the practitioner's perspective, options (ii) and (iii) are attractive since they only require the specification of one tuning parameter $\mu_n$.
Both additional choices (ii) and (iii) result in valid point-wise coverage in large samples (see Appendix Figures
\ref{fig:pwFreq_opt2} and \ref{fig:pwFreq_opt3}).
For $n=100$, $\omega>4^\circ$ all three options have nearly indistinguishable coverage frequency.
The behavior is notably different between the three for $ \omega \in [0,4^\circ]$.
Choice (ii) still results in a lower probability of coverage than the required probability, but the problematic neighborhood $[0,2^\circ]$ is smaller than for choice (i).
The choice $\mu_n=\kappa_n$ for $n=10,000$ is particularly in good alignment with the nominal coverage $\alpha=0.95$.
The choice (iii) essentially results in $\mbox{CB}_{\alpha,n}$ being the same length as $\mbox{CB}_{\alpha,n,\Measures}$ (see Appendix Figure \ref{fig:lengthCompare_opt3}).
The similar performance in simulation suggests that $\mbox{CB}_{\alpha,n}$ with $\kappa_n=0$ may have uniformly valid coverage like $\mbox{CB}_{\alpha,n,\Measures}$.
However, a formal study of uniform validity for $\mbox{CB}_{\alpha,n}$ with $\kappa_n=0$ is more difficult to conduct than for $\mbox{CB}_{\alpha,n,\Measures}$.
\begin{figure}[H]
\begin{centering}
\caption{\label{fig:pwFreq_opt1} Coverage frequency (a,b) and average excess length (c) in the $2$-dimensional design for $\mbox{CB}_{\alpha,n}$ and $\mbox{CB}_{\alpha,n,\Measures}$ as function of $\omega$ in the $2$-dimensional design.}
\begin{tabular}{cc}
\includegraphics[scale=0.9]{Figures/tunning_case_1/MCresults2d_PW.pdf}& \includegraphics[scale=0.9]{Figures/tunning_case_1/MCresults2d_Uniform.pdf} \tabularnewline
\end{tabular}
\begin{tabular}{c}
\includegraphics[scale=0.9]{Figures/tunning_case_1/MCresults2d_lengthPWvsUniform.pdf}\tabularnewline
\end{tabular}
\par\end{centering}
Note: The dotted lines correspond to the asymptotic uniform 95\% confidence interval for the parameter $p=0.95$ of the Bernoulli random variable based on a random sample of $1000$ simulations based on the Bonferroni correction for 14 hypothesis tests.
Values of $\omega$ close to zero in panel (a) result in nonnegligible under-coverage. As the sample size grows, the problematic area shrinks.
Values of $\omega$ close to zero in panel (b) result in nonnegligible conservative coverage.
\end{figure}
\begin{figure}[H]
\begin{centering}
\caption{\label{fig:lengthCompare_opt1} Sensitivity of coverage frequency to tuning parameter choices.}
\begin{tabular}{cc}
\includegraphics[scale=0.9]{Figures/tunning_compare/MCresults2d_compareMU.pdf}& \includegraphics[scale=0.9]{Figures/tunning_compare/MCresults2d_compareKappa.pdf}\tabularnewline
\end{tabular}
\par\end{centering}
Note: The dotted lines correspond to the asymptotic uniform 95\% confidence interval for the parameter $p=0.95$ of Bernoulli random variable based on a random sample of $1000$ simulations based on Bonferroni correction for 14 hypothesis tests.
\end{figure}
\subsection{Effect of high dimensions}\label{subsec:highDMC}
In this section, we study the effect of the dimension of $\theta$ on the coverage probability and the length of the confidence bounds.
The design of the constraints matrix $\expect W=\expect (\Ap|\bp)$ is given by
\begin{align*}
\Ap=\left(\begin{array}{c c}
1 & -a\\
0 & -\I[d-1]\\
0 & \I[d-1]
\end{array}\right) \text{ and } \bp=\left(\begin{array}{ c}
0\\
\iota_{2}\\
\iota_{d-3}/\sqrt{d-3}\\
\iota_{2}\\
\iota_{d-3}/\sqrt{d-3}
\end{array}\right).
\end{align*}
The first line corresponds to an equality constraint that depends on a unit vector with direction $a\in\R^{d-1}$.
This equality constraint is deterministic, that is, $\operatorname{Var}W_{1j}=0$ for all $j=1,\dots,d+1$. Other constraints are inequality constraints with coefficients that are i.i.d Gaussian r.v. with $\operatorname{Var}W_{ij}=0.01$ for all $j=1,\dots,d+1$ and $i=2,\dots,(1+2d)$.
The identified set for the first coordinate $\theta_1$ in this design corresponds to values of a support function of a $d-1$ dimensional rectangular box $[-1,1]^2\times[-1/\sqrt{d-3},1/\sqrt{d-3}]^{d-3}$ in direction $a$ (and $-a$).
As before, the focus is on the lower bound, since by design the upper bound is symmetric.
When $d=3$ and $a=(1,0)^\prime$, this design reduces to the two-dimensional design considered in the previous subsection with $\omega=0$.
We consider two possible values for $a$: (a) regular case, $a=\iota_{d-1}/\sqrt{d-1}$; (b) nonregular case, $a=(1,0,\dots,0)'\in\R^{d-1}$.
In the regular case $\argminT = (-\iota_2,-\iota_{d-3}/\sqrt{d-3})$. In the nonregular case, the argmin set $\argminT$ is a convex hull of $2^{d-1}$ corner points with coordinates $(-1,\pm1,\pm1/\sqrt{d-3})$.
Every such corner point has the same distance to $0$, equal to $\sqrt{3}$, for $d>2$ regardless of the dimension $d$.
This normalization was chosen to make performance comparisons as dimension grows while keeping the diameter of the identified set fixed. In this way, we isolate the effect of the number of dimensions and constraint from the effect of the diameter of the identified on the dual variables.
By design, the asymptotic variance of the regularized estimators will not change with the dimension.
We used $n=1000$ data points and $N=1000$ Monte Carlo simulations. We use the first option for the tuning sequences, namely $\mu_n = \hat\mu_0 \sqrt{n^{-1}\ln\ln n}$ and $\kappa_n = \hat\mu_0 \sqrt{n^{-1}\ln n}$.
First, consider the coverage probability of the confidence bounds $\mbox{CB}_{\alpha,n}$ and $\mbox{CB}_{\alpha,n,\Measures}$. Panel (a) in Figure \ref{fig:ndFreq} shows that in the regular case both the uniform and the point-wise valid confidence set have a coverage frequency that is within the simulation error bound from the nominal coverage probability of $\alpha=0.95$.
At the same time, the coverage of the two procedures becomes noticeably different from each other for $d>40$.
It suggests that the size of the neighborhood where the uniform inference becomes important gets wider with $d$.
The nonregular design provided on panel (b) in Figure \ref{fig:ndFreq} shows that uniformly valid confidence bounds have a coverage frequency nearly equal to $100\%$.
The point-wise bounds $\mbox{CB}_{\alpha,n}$ are less conservative for all $d$, but still reach the $100\%$ coverage frequency for $d\geq21$.
It suggests that the sample size required for exact asymptotic coverage for $\mbox{CB}_{\alpha,n}$ gets larger as dimension grows.
The probability of coverage above the nominal level $\alpha=0.95$ in the nonregular case with high dimension $d$ reflects the fact that as the number of corner points increases exponentially with $d$.
As a result, the nonregularized support function gets increasingly biased as it has asymptotic distribution of a minimum of $2^{d-1}$ Gaussian r.v. as evident from \eqref{eq:asyDistributionShapiro}.
Fortunately, this extreme conservatism only appears in the worst possible case.
For a randomly chosen direction of a support function $a\in\R^{d-1}$, the performance is expected to be closer to that of panel (a) of Figure \ref{fig:ndFreq}, at least for sufficiently large sample sizes $n$.
As a benchmark, I use $\mbox{CB}_{\alpha,n,AS}$---one-sided confidence bounds for $\argminT_1$ based on \cite{andrews2010inference} with Bonferroni critical values implemented using the fast E-A-M algorithm of \cite{kaido2015inference}.
This choice of benchmark is one of the fastest available uniformly valid procedure in the literature. The two alternative approaches, \cite{bugni2014inference} and \cite{kaido2015inference}, can provide uniformly valid CSs with a potentially shorter average length than $\mbox{CB}_{\alpha,n,AS}$.
However, both are expected to take considerably more time to compute because they add time-consuming profiling or calibration steps.
AS results are only available in the nonregular case because of the software restrictions.\footnote{The code is available at https://molinari.economics.cornell.edu/programs.html.
This code is highly optimized for the case of inference on bounds on individual coordinates of the identified set. The code does not allow making confidence sets for arbitrary directions $a$.}
Nevertheless, AS is not adaptive and is expected to have a similar length regardless of whether we compute the bounds in regular or nonregular directions $a$.
The coverage probability equal to the confidence sets to $100\%$ also applies to $\mbox{CB}_{\alpha,n,AS}$.
Figure \ref{fig:ndLen} shows that for $d\geq9$ for $n=1000$, the uniformly valid $\mbox{CB}_{\alpha,n,\Measures}$ is less conservative than $\mbox{CB}_{\alpha,n,AS}$.
\subsection{Computational time}
The main advantage of $\mbox{CB}_{\alpha,n,\Measures}$ compared to $\mbox{CB}_{\alpha,n,AS}$ is the computational gain that is achieved because of two factors: (i) using only linear and convex quadratic programs; (ii) no need to use multi-start procedures for global optimization of non-convex programs.
Other alternative procedures like KMS and BCS would have additional burden of computing simulation-based critical values and are expected to be much slower.
Note that the computational cost does not depend on whether the design is regular or not, so Figure \ref{fig:time} compares the average computation time as a function of dimension $d$ on a modern multi-core laptop.\footnote{All simulations are using 2023 Macbook Pro with M2 Max CPU (12 cores) and 96 GB of RAM. The EAM algorithm for $\mbox{CB}_{\alpha,n,AS}$ takes advantage of all 12 cores. }
The computational time for the delta-method for regularized support functions grows very slowly with dimension $d$; the average time for $d=100$ is about 2 seconds.
In contrast, $\mbox{CB}_{\alpha,n,AS}$ is already nearly 1000 times slower for $d=21$.
\begin{figure}[H]
\caption{\label{fig:ndFreq}Coverage frequency in the $d$-dimensional design as function of $d$ in (a) regular and (b) non-regular cases.}
\begin{centering}
\begin{tabular}{cc}
\includegraphics[scale=0.9]{Figures/ndBoxTurned/MCresultsNdLD_Coverage.pdf}& \includegraphics[scale=0.9]{Figures/ndBox/MCresultsNdLD_Coverage.pdf}\tabularnewline
\end{tabular}
\par\end{centering}
Note: The dotted lines correspond to the asymptotic uniform 95\% confidence interval for the parameter $p=0.95$ of Bernoulli random variable based on a random sample of $1000$ simulations based on Bonferroni correction for 11 hypothesis tests.
AS coverage frequency is computed using $100$ simulations and for $d\leq21$, while the delta-method confidence bounds are based on $1000$ simulations. In all cases each simulation is based on $n=1000$ observations.
\end{figure}
\vspace{-0.5cm}
\begin{figure}[H]
\caption{\label{fig:ndLen}Length in the $d$-dimensional design as function of $d$ in (a) the regular and (b) non-regular cases.}
\begin{centering}
\begin{tabular}{cc}
\includegraphics[scale=0.9]{Figures/ndBoxTurned/MCresultsNd_length.pdf}& \includegraphics[scale=0.9]{Figures/ndBox/MCresultsNdLD_length.pdf}\tabularnewline
\end{tabular}
\par\end{centering}
Note: AS length is computed using $100$ simulations and for $d\leq21$, while the delta-method average length are based on $1000$ simulations. In all cases each simulation is based on $n=1000$ observations.
\end{figure}
\vspace{-0.5cm}
\begin{figure}[H]
\begin{centering}
\caption{\label{fig:time}Computational time in the $d$-dimensional design as function of $d$.}
\includegraphics[scale=0.8]{Figures/ndBox/MCresultsNdLD_time.pdf}
\par\end{centering}
Note: The scale is logarithmic. AS average time is computed using $100$ simulations and for $d\leq21$ only, while the delta-method average times are based on $1000$ simulations. In all cases, each simulation is based on $n=1000$ observations.
\end{figure}
\section{Conclusion}\label{sec:Conclusion}
This paper demonstrated that the regularization approach provides a fast way to construct point-wise and uniform CSs for a $\theta_{1}$
that is comparable to or shorter than those of the existing literature.
Monte Carlo simulations showed that the proposed CSs have good finite sample coverage properties.
The computational benefits of the new approach are particularly prominent if the dimension of $\theta$ is large.
The regularization framework can be extended in a number of ways to allow for overidentification and joint inference.
The proposed approach is attractive in applications such as a linear model with an interval-valued outcome variable and a large number of regressors and in problems with parameters represented as intersection bounds.
\label{reply:R1.2a} Focus on the affine inequalities simplifies the large sample analysis. However, there are many interesting applications that characterize the identified set for the structural parameters using non-linear moment conditions, in particular, among the structural models of the industrial organization \citep{pakes2015moment}. Analysis of such non-linear moment inequalities goes beyond the scope of the current study, but would be a promising direction of future research.
\bibliographystyle{ecta}
\bibliography{main}
\newpage{}