Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
72,388 characters · 19 sections · 50 citation commands
Inference in a class of optimization problems: Confidence regions and finite sample bounds on errors in coverage probabilities
\def\spacingset#1{ {#1}} \spacingset{1}
\if00 \fi
\if10 {
} \fi
{\it Keywords:} partial identification, normal approximation, sub-Gaussian distribution, finite-sample bounds
\spacingset{1.45}
\doublespacing
This paper presents three methods for carrying out non-asymptotic inference about a function of partially identified structural parameters of an econometric model. The methods apply to models that impose shape restrictions FH:2015, HL:17, a variety of partially identified models Manski:book, Tamer:10 that include discrete games Ciliberto:Tamer:09, and models in which a continuous function is inferred from the average values of variables in a finite number of discrete groups BDM:98,KT16. The specific inference problem consists of finding upper and lower bounds on the partially identified function $f(\psi)$ under the restrictions $g_1(\psi,\mu) \leq 0$ and $g_2(\psi,\mu) = 0$, where $\psi$ is a vector of structural parameters, $\mu$ is a vector of unknown population means of observable random variables, $f$ is a known, real-valued function, and $g_1$ and $g_2$ are known possibly vector-valued functions. The inequality $g_1(\psi,\mu) \leq 0$ holds component-wise.
Most existing methods for inference in our framework are based on asymptotic approximations. They provide correct inference in the limit of an infinite sample size but do not provide information about the accuracy of the asymptotic approximations in finite samples. We provide three methods for obtaining finite-sample lower bounds on the coverage probability of a confidence interval for $f (\psi)$. One method uses asymptotic approximations to obtain a confidence interval. The other two methods do not use asymptotic approximations. All the methods provide information about the accuracy of finite-sample inference.
There are several approaches to carrying out non-asymptotic inference in our framework. Sometimes a statistic with a known finite-sample distribution makes finite-sample inference possible. For example, the \citetalias{CL34} confidence interval for a population probability is obtained by inverting the binomial probability distribution function. Manski:07 used the Clopper-Pearson interval to construct finite-sample confidence sets for counterfactual choice probabilities. Our methods apply to parameters that are not necessarily probabilities. A second existing method consists of using Hoeffding's inequality to obtain a confidence interval. STZ:2018 used this inequality to construct a confidence interval for a partially identified population moment. Hoeffding's inequality requires the underlying random variable to have a known bounded support. Our methods do not require the underlying random variable to have a known or bounded support. Minsker2015 developed a confidence set for a vector of population means by using a method called “median of means.” The bounds provided by this method are looser than the bounds provided by our methods. In addition, Minsker's method depends on certain user-selected tuning parameters. There are no data-based, efficient ways to choose these parameters in applications.
Our first method consists of making a normal approximation to the unknown distribution of the sample average. This method makes certain assumptions about low-order moments of the underlying random variable but does not restrict its distribution in other ways. A variety of results provide finite-sample upper bounds on the errors made by normal approximations. The Berry-Ess\'{e}en inequality for the average of a scalar random variable is a well-known example of such a bound. Bentkus03 provides a bound on the error of a multivariate normal approximation to the distribution of the sample average of a random vector. Other normal approximations for random vectors are given by spokoiny2015; CCK:17; and zhilova2020; among others. Our first method uses a bound on the error of the multivariate normal approximation that is due to Raic2019. \citetalias{Raic2019} bound is a refined and tighter version of the bound of Bentkus03.
The bound of CCK:17 may be tighter than that of Raic2019 when the dimension of $\mu$ exceeds the sample size, but the bound of Raic2019 is tighter when the dimension of $\mu$ is small compared to the sample size, which is the case we treat in this paper. In contrast to conventional asymptotic inference approaches, our first method provides a finite-sample lower bound on the coverage probability of a confidence interval for the partially identified function $f(\psi)$.
The bound provided by our first method is loose in samples of the moderate sizes that occur in most economics applications, though not necessarily in very large samples. This is because it places only weak restrictions on the distribution of the underlying random variable, which may be far from normal. Our second method obtains a tighter bound in moderate size samples by assuming that the distributions of the components of the possibly vector-valued underlying random variable are sub-Gaussian. The sub-Gaussian assumption places stronger restrictions on the thickness of the tails of the relevant distributions than do the assumptions of the bound based on \citetalias{Raic2019} inequality. Our third method tightens the bound obtained with our second method by assuming that if the underlying random variable is vector-valued, then its distribution is sub-Gaussian in a vector sense. This assumption is stronger than the assumption that the components of a random vector are individually sub-Gaussian. The bounds obtained with the second and third methods are identical if the underlying random variable is a scalar.
The bounds provided by all the methods depend on unknown population parameters. This dependence is unavoidable and can be removed only in special cases. The parameters of the bounds of the second and third methods can be estimated, however, which makes it possible to estimate the bounds in applications. We describe how to do this. The resulting estimated bounds are asymptotic. They do not have finite-sample validity but can provide useful, though possibly rough, indications of the magnitudes of the finite-sample bounds. We present the results of Monte Carlo experiments that illustrate the relation between the exact finite-sample bounds and the consistent estimates.
Our work is broadly related to the literature on inference in partially identified models. Tamer:10, canay_shaikh_2017, ho_rosen_2017, and Molinari provide recent surveys. CCT:18 describe a Monte Carlo method for carrying out asymptotic inference for a class of models that includes our framework. BCS:17 and KMS19 develop asymptotic inference methods for subvectors of partially identified parameters in moment inequality models. CCK and BBC construct confidence regions by inverting pointwise tests of a hypothesis about the (sub)vector of parameters that are partially identified by a large number of moment inequalities. The inference problem we treat is different from those in the foregoing papers in that we focus on inference about parameters that are solutions to a class of optimization problems that is different from moment inequality problems. Our methods and results do not apply to moment inequalities. Kline:Tamer describe Bayesian inference in a class of models that includes a special case of the models we treat. Two more closely related papers are HSS:2017 and Shi:Shum:15, who propose a method for asymptotic inference about estimators defined by mathematical programs. However, the class of estimation problems they treat is different from ours and overlaps with ours only under highly restrictive assumptions about both classes. HSS:2017 and Shi:Shum:15 do not provide finite-sample bounds on the errors of their asymptotic approximations.
Our work is also related to the econometrics literature on finite-sample inference. \, STZ:2018 consider finite-sample inference in auction models. Their framework and method are very different from those in this paper. In a different context, CHJ:09 and Rosen:Ura:19 propose finite-sample inference for quantile regression models and for the maximum score estimand, respectively. Their methods and the classes of models they treat are distinct from ours.
The remainder of this paper is organized as follows. Section (ref) describes the inferential problem we treat, our methods for obtaining a confidence interval for $f(\psi)$, and the three methods for obtaining a finite-sample lower bound on the coverage probability of a confidence interval. Section (ref) also describes two empirical studies that illustrate how the inferential problem arises in applications. Section (ref) describes computational procedures for implementing our methods. Section (ref) presents an empirical application of the methods. Section (ref) reports the results of a Monte Carlo investigation of the numerical performance of our methods, and Section (ref) gives concluding comments. The proofs of theorems are presented in online Appendix (ref). Online Appendices (ref)--(ref) provide additional technical information about our methods, a description of Minsker's (2015) method, an additional empirical application, and additional Monte Carlo results.
Section (ref) presents an informal description of inferential problem we address. Section (ref) gives two examples of empirical applications in which the inferential problem arises. Section (ref) provides a formal description of our methods for constructing confidence intervals and bounds on coverage probabilities.
Let $\{ X_i : i=1,\ldots,n \}$ be an independent random sample from the distribution of the random vector $X \in \mathbb{R}^p$ for some finite $p \geq 1$. Define $\mu = \mathbb{E}(X)$ and $\Sigma = \text{cov}(X)$. We assume that both exist. Let $\psi$ be a finite-dimensional parameter and $f (\psi )$ be a real-valued, known function. We assume throughout this section that $f (\psi )$ is only partially identified by the sampling process, though our results also hold if $f (\psi )$ is point identified. We seek a confidence interval for $f (\psi )$, which we define as a data-based interval that contains $f (\psi )$ with probability exceeding a known value. Let $g_1(\psi , \mu)$ and $g_2 (\psi , \mu)$ be possibly vector valued known functions satisfying $g_1(\psi, \mu) \leq 0$ and $g_2(\psi, \mu) = 0$ component-wise. Define
subject to the component-wise constraints:
where $\Psi$ is a compact parameter set. Online Appendix (ref) extends (ref)-(ref) to the case in which $g_1$ and $g_2$ depend on a continuous covariate in addition to $\psi$ and $\mu$.
We are interested in the identification interval $J_{-} \leq f(\psi) \leq J_{+}$. However, this interval cannot be calculated in applications because $\mu$ is unknown. Therefore, we estimate $\mu$ by the sample average $\bar{X} = n^{-1} \sum_{i=1}^n X_i$, and we estimate $J_{+}$ and $J_{-}$ by
subject to the constraints
where $\mathcal{S}$ is a set for which $n^{1/2} (\bar{X} - \mu) \in \mathcal{S}$ with high probability. Since $\mu$ is unknown, we replace it with the variable of optimization $m$ in (ref)--(ref) but require $m$ to satisfy (ref). The resulting confidence interval for $f(\psi)$ is
This is also a confidence interval for the identified set containing $f(\psi)$. Section (ref) provides three different finite-sample lower bounds on the probability that this interval contains $f(\psi)$. That is, Section (ref) provides three finite-sample lower bounds on
The three bounds correspond to increasingly strong assumptions about the distribution of $X$ and are increasingly tight with samples of the moderate sizes found in most economics applications, though not necessarily with very large samples.
The two leading examples of $\mathcal{S}$ in constraint (ref) are a box and an ellipsoid. If $\mathcal{S}$ is a box, let $D$ be a diagonal matrix whose diagonal elements are strictly positive. For example, $D$ might be the diagonal elements of $\Sigma$ if $\Sigma$ is known or the diagonal elements of a consistent estimate, $\widehat{\Sigma}$, if $\Sigma$ is unknown. Choose $\kappa_b(1-\alpha)$ so that the following holds, uniformly in $j = 1,\ldots,p$, with probability $1-\alpha$:
where the subscript $j$ denotes the $j$'th component of a vector or the $(j, j)$ component of a matrix. In this case, (ref) becomes $p$ constraints and can be viewed as a sample analog of $\left| \mathbb{E}(X_j) - \mu_j \right| \leq 0$ with a relaxed constraint on $\mu_j$ for each $j = 1,\ldots,p$. Section (ref) presents methods for choosing $\kappa_b(1-\alpha)$.
If $\mathcal{S}$ is an ellipsoid, let $\Upsilon$ denote a positive definite $p \times p$ matrix, possibly $\Sigma$ or $\widehat{\Sigma}$ if those matrices are non-singular, or the identity matrix. Choose $\kappa_e(1-\alpha)$ so that $$ n (\bar{X} - \mu)' \Upsilon^{-1} (\bar{X} - \mu) \leq \kappa_e(1-\alpha) $$ with probability $1-\alpha$. In this case, (ref) is a single constraint. Section (ref) presents a method for choosing $\kappa_e(1-\alpha)$. When $\Sigma$ is difficult to estimate or is singular, we may use a sphere by choosing a critical value $\kappa_s(1-\alpha)$ such that
with probability $1-\alpha$. In general, the implementation of our methods is simpler if $\mathcal{S}$ is indexed by a scalar critical value $\kappa(1-\alpha)$.
It is straightforward to allow the objective function $f(\psi)$ to depend on $\mu$. For the lower bound $\hat{J}_{-}(\bar{X})$, we introduce an auxiliary variable $t$ that acts as an upper bound on $f(\psi, \mu)$ and solve: $ \min_{\psi,m,t} t $ subject to $f(\psi, m) \leq t$ and (ref). For the upper bound $\hat{J}_{+}(\bar{X})$, we introduce a lower bound $s$ on $f(\psi, \mu)$ and solve: $ \max_{\psi,m,s} s $ subject to $f(\psi, m) \geq s$ and (ref). We focus on the original form (ref)--(ref) in the remainder of this paper because the form with the objective function $f(\psi, \mu)$ can be rewritten in the form (ref)--(ref) by redefining $f$, $g_1$ and $g_2$.
BDM:98 use grouped data to estimate labor supply effects of tax reforms in the United Kingdom. To motivate our setup, we consider a simple model with which BDM:98 explain how to use grouped data to estimate $\beta$ in the following labor supply model with no income effect:
In this model, $h_{it}$ and $w_{it}$, respectively, are hours of work and the post-tax hourly wage rate of individual $i$ in year $t$, and $U_{it}$ is an unobserved random variable that satisfies certain conditions. The parameter $\beta$ is identified by a relation of the form $ \beta = \beta( h_{gt}, lw_{gt}), $ where $h_{gt}$ and $lw_{gt}$ are the mean hours and log wages in year $t$ of individuals in group $g$. There are 8 groups defined by four year-of-birth cohorts and level of education. The data span the period 1978-1992.
A nonparametric version of (ref) is $h_{it} = \xi( w_{it}) + U_{it}$, where $\xi \in \Xi$ is an unknown continuous function and $\Xi$ is a function space. A nonparametric analog of $\beta$ is the weighted average derivative \[ \tilde{\beta} = \int \frac{\partial \xi(u)}{\partial u} w(u) du, \] where $w$ is a non-negative weight function. The average derivative $\tilde{\beta}$ is not identified non-parametrically by the mean values of hours and wages for finitely many groups and time periods. It can be partially identified, however, by imposing a shape restriction such as weak monotonicity on the labor supply function $\xi$. Assume, for example, that $\mathbb{E}[h_{it} - \xi(w_{it}) | g,t] = 0$. {BDM:98 set $\mathbb{E}[h_{it} - \xi(w_{it}) | g,t] = a_g + m_t$, where $a_g$ and $m_t$, respectively, are group and time fixed effects. These are accommodated by our framework but we do not do this in the present discussion.}
The identification interval for $\tilde{\beta}$ is $\tilde{\beta}_{-} \leq \tilde{\beta} \leq \tilde{\beta}_{+}$, where
subject to
The continuous mathematical programming problem (ref)-(ref) can be put into the finite-dimensional framework of (ref)-(ref) by observing that under mild conditions on $\Xi$, $\xi$ can be approximated very accurately by the truncated infinite series
where the $\psi_j$'s are constant parameters, the $\phi_j$'s are basis functions for $\Xi$, and $K$ is a truncation point. In an estimation setting, $K$ can be an increasing function of the sample size, though we do not undertake this extension here. {The approximation error of (ref) can be bounded. Here, however, we assume that $K$ is sufficiently large to make the error negligibly small.} The finite-dimensional analog of (ref)-(ref) is
subject to
$J_{+}$ and $J_{-}$ can be estimated, thereby obtaining $\hat{J}_{+}$ and $\hat{J}_{-}$, by replacing $h_{gt}$ and $w_{gt}$ in (ref)-(ref) with within-group sample averages and adding the constraint (ref).
Ho:Pakes:2014 use the theory of revealed preference to develop an estimator of hospital choices by individuals. HP use data on privately insured births in California. We consider a simplified version of the HP model.
Using the notation of HP, let $p(c, h)$ denote the price an insurer is expected to pay at hospital $h$ for a patient with medical condition $c$. Let $i$ index patients. Then $c_i$ is the medical condition of patient $i$, and $p(c_i, h)$ is the price an insurer is expected to pay at hospital $h$ for patient $i$. Let $l_i$ denote patient $i$'s location, $l_h$ hospital's location, and $d(\cdot,\cdot)$ the distance between the two locations. For hospitals $h \neq h'$, define
That is, $\Delta p(c_i, h, h')$ is the price difference between hospitals $h$ and $h'$ given patient condition $c_i$ and $\Delta d( l_{i}, l_{h}, l_{h'} )$ is the distance difference between hospitals $h$ and $h'$ given patient location $l_i$. Define
where $\psi$ is a scalar parameter that determines price sensitivity relative to distance. Note that the coefficient for distance is normalized to be $-1$. $\psi$ is the key parameter in HP.
Define the four-dimensional vector of instruments based on distance:
Here, the instruments are based on distance measures and constructed to be positive to preserve the signs of the inequalities below in (ref).
Let $S(h,h',s)$ be the set of patients with severity $s$ who chose hospital $h$ but had hospital $h'$ in their choice set. The identifying assumption in HP is that
for all $s, h, h'$ such that $h \neq h'$. We can rewrite (ref) as
where
To see the connection between our general framework and HP's inequality estimator, let $f(\psi) = \psi$, $\mu = (\mu_p, \mu_d)$, and $g_1$ be a collection of inequalities such that
There is no element in $g_2$ (no equality constraints here). Since each element in $\mu$ can be estimated by a suitable sample mean, our general framework includes HP's estimator as a special case.
This section presents our three methods for forming finite-sample lower bounds on $$\mathbb{P} \left[ \hat{J}_{-}(\bar{X}) \leq J_{-} \leq f(\psi) \leq J_{+} \leq \hat{J}_{+}(\bar{X}) \right].$$ The three bounds make assumptions of differing strengths about the distribution of $X$. The bounds are tighter in samples of moderate size with stronger assumptions. All proofs are in Online Appendix (ref). We begin with the following theorem, which applies to all the methods and forms the basis of our approach.
Now define
Then $\mathbb{E}( \bar{Z} ) = 0$. Note that $\Sigma = \text{cov}(Z_i) = \text{cov}(\bar{Z}) = \text{cov}(X)$. We make the following assumption throughout the remainder of the paper.
Suppose for the moment that $\Sigma$ is known. Section (ref) discusses the case in which $\Sigma$ is unknown. If $\Sigma$ is non-singular, let $[ \Sigma^{-1/2} (X_i - \mu) ]_j$ denote the $j$'th component of $\Sigma^{-1/2} (X_i - \mu)$.
This method approximates the distribution of $\bar{Z}$ by a normal distribution. To do this, make the following assumption.
Define the independent random $p$-vectors $W_i \sim N(0, \Sigma)$ $(i=1,\ldots,n)$ and $\bar{W} := n^{-1/2} \sum_{i=1}^n W_i \sim N(0, \Sigma)$. The multivariate generalization of the Lindeberg-L\'{e}vy central limit theorem shows that $\bar{Z}$ is asymptotically distributed as $N(0, \Sigma)$, so the distribution of $\bar{Z}$ can be approximated by that of $W$. The following theorem bounds the error of this approximation.
Theorem (ref) approximates the distribution of $\bar{Z}$ by a multivariate normal distribution and uses a multivariate generalization of the Berry-Ess\'{e}en theorem Raic2019 to bound the approximation error. Theorem (ref) implies that for any $0 < \alpha < 1$,
where $\kappa_{\chi^2_{p}}(1-\alpha)$ is the $(1-\alpha)$ quantile of the chi-square distribution with $p$ degrees of freedom and
The term $\textrm{B}(n,p,\overline{\mu}_3)$ in (ref) is asymptotically negligible but can be large in samples of moderate size because it accommodates “worst case” distributions of $\bar{X}$ that may be far from normal. If $n$ is large enough that $\textrm{B}(n,p,\overline{\mu}_3) < \alpha $, then it follows from (ref) that
It follows from Theorem (ref) that (ref) and (ref) provide lower bounds on the coverage probabilities of confidence intervals for $f(\psi)$ when $\mathcal{S}$ is an ellipsoid.
Table (ref) shows numerical values of the bound $(1-\alpha) - \textrm{B}(n,p,\overline{\mu}_3)$ and critical value $\kappa_{\chi^2_{p}} \big[ 1-\alpha + \textrm{B}(n,p,\overline{\mu}_3) \big]$ for different values of $n$ and $p$ at $\alpha = 0.05$ and $\overline{\mu}_3 = 2$. To have a bound close to $1-\alpha$ and a finite critical value, $n$ must be very large, especially if $p$ is large. This is because (ref) accommodates worst case distributions of $\bar{X}$. Methods 2 and 3, which are discussed next in this section, provide tighter bounds and smaller critical values when $n$ is smaller, though Method 1 can provide a smaller critical value when $n$ is very large and $p$ is small. However, Method 1 is hard to use in applications even when $n$ is large if $\Sigma$ and $\overline{\mu}_3$ are unknown, because the resulting bounds depend on population parameters that are difficult to estimate. This problem is discussed in Section (ref).
Method 2 obtains bounds that are much tighter than those of Method 1 when $n$ is smaller than in Table (ref), and Method 2 does not require $\Sigma$ to be invertible. This is accomplished by assuming that the distributions of the components of $X$ are sub-Gaussian. Specifically, make the following assumption.
Assumption (ref) requires that the distribution of $\tilde{Z}_{ij}$ be thin-tailed. Sub-Gaussian random variables include Gaussian, Rademacher, and bounded random variables as special cases. See, e.g., wainwright2019book. There is a tradeoff between $\Upsilon$ and $\sigma_j^2$. In particular, $\sigma_j^2$ may be larger if $\Upsilon$ is the $p \times p$ identity matrix than if $\Upsilon = \Sigma$ and $\Sigma$ is non-singular.
Define $\sigma^2 := \max_{1 \leq j \leq p} \sigma_j^2$. The following theorem, combined with Theorem (ref), provides a lower bound on a confidence interval for $f(\psi)$ when $\mathcal{S}$ is an ellipsoid.
Theorem (ref) and Method 2 make use of the sub-Gaussianity of $\tilde{Z}_{ij}$, whereas Theorem (ref) and Method 1 allow the tails of the distribution of $X$ to be thicker than sub-Gaussian tails. The critical values of Methods 1 and 2 are compared later in this section after the description of Method 3. Estimation of $\sigma^2$ is discussed in Section (ref).
Method 3 makes the stronger assumption that the distribution of $X$ is sub-Gaussian in a vector sense. Specifically, Method 3 makes the following assumption.
Assumption (ref) is stronger than Assumption (ref), because Assumption (ref) requires the entire vector $X$ to be sub-Gaussian. If $X$ is multivariate normal and $\Upsilon = \Sigma$, Assumption (ref) holds with $\sigma^2 = 1$. In general, however, it is difficult to find simple conditions under which Assumption (ref) is satisfied without assuming that the elements of $\bar{Z}$ are independent of one another.
An application of Theorem 2.1 of HKT:2021 gives the following theorem which, combined with Theorem (ref), provides a lower bound on the coverage probability of a confidence interval for $f(\psi)$ when $\mathcal{S}$ is an ellipsoid.
This theorem and Method 3 yield a smaller critical value and confidence set than Theorem (ref) and Method 2 do, but they require the stronger Assumption 4. Table (ref) shows the critical values of Methods 2 and 3 and chi-square critical values with p degrees of freedom. The critical values are all for $\alpha = 0.05$. None of the critical values depends on $n$. The chi-square critical value achieves an asymptotic coverage probability of 0.95 and yields the smallest confidence set but does not ensure a finite-sample coverage probability of at least 0.95. The sub-Gaussian critical values and resulting confidence sets are larger but ensure finite-sample coverage probabilities of at least 0.95 under the regularity conditions of Theorems (ref) and (ref). Assumption (ref) is easier to satisfy, and its variance proxy is easier to estimate in applications, but the critical value of Method 2 increases more rapidly than the critical value of Method 3 as $p$ gets large.
The critical values of Methods 2 and 3 in Table (ref) can be compared with those of Method 1 in Table (ref). The critical value of Method 1 converges to the chi-square critical value as $n \rightarrow \infty$. The critical values of Methods 2 and 3 do not depend on $n$. Therefore, the critical value of Method 1 is smaller than those of Methods 2 and 3 when $n$ is very large and $p$ is small enough. However, the critical value of Method 1 is infinite with moderate values of $n$, whereas the critical values of Methods 2 and 3 are finite at all values of $n$.
In applications, $\Sigma$ and the sub-Gaussian variance proxy $\sigma^2$ are unknown except in special cases. This section explains how to estimate these quantities and discusses the effect of estimation on the bounds presented in Section (ref).
We begin with the bound of Theorem (ref) and (ref). Let $\widehat{\Sigma}$ be the following estimator of $\Sigma$:
Let $\Sigma_{jk}^{-1}$ denote the $(j,k)$ component of $\Sigma^{-1}$. For each $j,k = 1,\ldots,p$, let
Assumption (ref) is stronger than Assumption (ref). In particular, (ref) implies that the $(X_{ij} - \mu_j)(X_{ik} - \mu_k)$ is sub-exponential. Therefore, $X_{ij}$ is sub-Gaussian because a random variable is sub-Gaussian if and only if its square is sub-exponential vershynin2018high. Also, the product of two sub-Gaussian variables is sub-exponential vershynin2018high. We use Assumption (ref)(iii) to apply Bernstein's inequality to the bound of Method 1 buhlmann2011statistics.
Define the random vector $\widehat{W} \sim N(0, \widehat{\Sigma})$. We approximate the distribution of $W$ by the distribution of $\widehat{W}$ with $\widehat{\Sigma}$ treated as a non-stochastic matrix. Define $\mathbf{P} (\mathcal{S}, \Sigma) := \mathbb{P} (W \in \mathcal{S})$ for $W \sim N(0, \Sigma)$ and
The following lemma gives a finite-sample bound on the error of the approximation.
The condition (ref) is a mild technical condition that can be satisfied easily. The conclusion of Lemma (ref) holds only if $\widehat{\Sigma}$ satisfies certain conditions that are stated in the proof of the lemma in Online Appendix (ref). These conditions are satisfied with probability at least $1-2 e^{-t}$, not with certainty.
Define
Now combine Theorem (ref) and Lemma (ref) to obtain the following theorem.
Theorem (ref) provides a finite-sample upper bound on the error made by approximating $\mathbb{P} \left[ n^{1/2} (\bar{X} - \mu) \in \mathcal{S} \right]$ by $\mathbf{P} (\mathcal{S}, \widehat{\Sigma})$. Combining Theorems (ref) and (ref) yields
Theorem (ref) provides a finite-sample lower bound on $\mathbb{P} \left[ \hat{J}_{-}(\bar{X}) \leq J_{-} \leq f(\psi) \leq J_{+} \leq \hat{J}_{+}(\bar{X}) \right]$ that takes account of random sampling error in $\hat{\Sigma}$. It is not difficult to choose $\mathcal{S}$ so that the right-hand side of (ref) is $1-\alpha$ for any $0 \leq \alpha \leq 1$ if $n$ is large enough and $p$ is small enough to make the term in square brackets on right-hand side of the inequality less than $\alpha$. However, the presence of $\delta_n^*$ greatly decreases the right-hand side of (ref) relative to what it is when $\Sigma$ is known, thereby increasing the size of the confidence region $\mathcal{S}$ for any given value of $\alpha$. In addition, the right-hand side of (ref) depends on several population parameters that are difficult to estimate. Therefore, (ref) is of limited use for applications. The bounds for Methods 2 and 3 with an unknown variance proxy also depend on unknown population parameters, but these are easier to estimate, as is discussed in Section (ref). The dependence of finite-sample bounds on unknown population parameters is unavoidable except in special cases. For example, if $\Upsilon = I_p$ in Theorem (ref) and each element of $X-\mu$ is contained in $[-1, 1]$, then $\sigma_j^2 = 1$ for all $j$, and the first inequality in Theorem (ref) becomes the Hoeffding inequality.
We now consider Method 2, which is the sub-Gaussian case of Assumption (ref). In some special cases, the sub-Gaussian variance proxies, $\sigma_j^2$ are known. For example, if $\Upsilon = I_p$ in Theorem (ref) and each element of $X - \mu$ is contained in $[-1,1]$, then $\sigma_j^2 = 1$ for all $j$. Here, we derive bounds for the case in which $\sigma^2 = \max_{1 \leq j \leq p} \sigma_j^2$ is unknown and must be estimated.
For all $i=1,\ldots,n$ and $j=1,\ldots,p$, let
and the sample variance of $\tilde{X}_{ij}$
The following lemma establishes a finite-sample probability bound on the absolute difference between the sample and population variances of $\tilde{X}_{ij}$. Recall that $\tilde{Z}_{ij} = [ \Upsilon^{-1/2} (X_i - \mu) ]_j$ defined in Assumption (ref).
The next lemma establishes a link between the population variance of $\tilde{Z}_{ij}$ to the variance proxy $\sigma_j^2$.
Now define
which we use to estimate $\sigma^2 = \max_{1 \leq j \leq p} \sigma_j^2$. The following theorem holds for the sub-Gaussian case of Assumption (ref) (Method 2).
A smaller choice of $t$ in (ref) makes the confidence set in (ref) tighter but results in a lower probability bound (ref). The right-hand side of (ref) can be close to zero and depends on the unknown parameters $\sigma_j^2$. Therefore, the bound (ref), like the bound (ref), is of limited use for applications. As is noted in the discussion of (ref), dependence of finite-sample bounds on unknown population quantities is unavoidable except in special cases. Section (ref) describes a practical approach to dealing with this problem that can be used in applications.
A result similar to Theorem (ref) can be obtained for Method 3. We do not undertake this here, however, because the resulting bound analogous to (ref) is loose, depends on unknown population parameters, and is of limited use in applications. Instead, in Section (ref), we describe a way of dealing with unknown population parameters in Methods 2 and 3 that can be used in applications.
Except in special cases, it is not possible to obtain finite-sample inequalities for Methods 2 and 3 that do not depend on unknown population parameters. It is possible, however, to estimate lower bounds on these parameters consistently. It follows from Lemma (ref) that $\tilde{\Sigma}_{jj}$ and $\max_{1 \leq j \leq p} \tilde{\Sigma}_{jj}$, respectively, are consistent estimates of lower bounds on the variance proxies $\sigma_j^2$ and $\sigma^2$ of Method 2. The differences between the lower bounds and the variance proxies are often small. Arguments like those used to prove Lemma (ref) show that the largest eigenvalue of the covariance matrix of $\tilde{Z}_{ij} (j =1,\ldots,p)$ is a lower bound on the variance proxy of Method 3. This can be estimated consistently by the largest eigenvalue of the sample covariance matrix.
Standard methods can be used to obtain asymptotic confidence intervals for the variance proxies of Method 2. These methods do not provide information about differences between true and nominal coverage probabilities in finite samples, but they provide practical indications of the magnitudes of the Method 2 bounds that can be implemented in applications.
Obtaining a useful asymptotic confidence interval for the largest eigenvalue of the covariance matrix of $\tilde{Z}_{ij}$ is difficult. A wide interval whose true coverage probability is likely to be much greater than the nominal probability can be obtained from the Frobenius norm of the difference between the estimated and true covariance matrices.
Table (ref) in Section (ref) compares values of $\hat{J}_{+}$ obtained using estimates and population values of the variance proxies of Methods 2 and 3.
Recall that our general framework is to obtain the bound $[\min_{\psi, m} f(\psi), \max_{\psi, m} f(\psi)]$ subject to $g_1(\psi, m) \leq 0, g_2(\psi, m) = 0, \psi \in \Psi,$ and $n^{1/2} (\bar{X} - m) \in \mathcal{S}$. In many examples, $f(\psi)$ is linear in $\psi$. For example, $\psi$ is the vector of all the parameters in an econometric model and $f(\psi)$ is just one element of $\psi$ or a linear combination of elements of $\psi$.
The restrictions $g_1(\psi, \mu) \leq 0$ include shape restrictions among the elements of $\psi$. Equality restrictions are imposed via $g_2(\psi, \mu) = 0$. The easiest case is that $g_j(\psi, \mu)$ is linear in $(\psi, \mu)$ for each $j=1,2$. In some of examples we consider, $g_j(\psi, \mu)$ is linear in $\psi$, holding $\mu$ fixed, and linear in $\mu$, keeping $\psi$ fixed, but not linear in $(\psi, \mu)$ jointly. This corresponds to the case of bilinear constraints. For example, $g_j(\psi, \mu)$ may depend on the product between one of elements of $\psi$ and one of elements of $\mu$. In practice, $\Psi$ can always be chosen large enough that the constraint $\psi \in \Psi$ is not binding additionally and can be ignored. For example, suppose that $\psi$ is a probability and the constraints in $g_1(\psi, \mu) \leq 0$ and $g_2(\psi, \mu) = 0$ impose a restriction on $\psi$ such as $[a, b]$ for $0 \leq a < b \leq 1$. Then, it is not necessary to impose $\psi \in \Psi = [0,1]$ additionally.
Recall that the leading cases of $\mathcal{S}$ include an ellipsoid and a box. For brevity, we focus on the scenario that normal distributions are used in obtaining $\mathcal{S}$. When $\mathcal{S}$ is a box, the critical value $\kappa_b(1-\alpha)$ can be easily simulated from the $N(0, \widehat{\Sigma})$ and the restriction $n^{1/2}(\bar{X} - \mu) \in \mathcal{S}$ can be written as linear constraints. When $\mathcal{S}$ is an ellipsoid, the critical value $\kappa_e(1-\alpha)$ can be obtained from the $\chi^2( d_\mu )$ distribution, where $d_\mu$ is the dimension of $\mu$. Then, the restriction $n^{1/2}(\bar{X} - \mu) \in \mathcal{S}$ can be written as
This is a convex quadratic constraint in $\mu$.
When some of the constraints $g_1(\psi, \mu) \leq 0$ and $g_2(\psi, \mu) = 0$ are bilinear, the resulting feasible region may not be convex. To deal with the bilinear constraints, we solve optimization problems using mixed integer programming (MIP) with Gurobi in R. By virtue of the developments in MIP solvers and fast computing environments, MIP has become increasingly used in recent applications. For example, bertsimas2016 adopted an MIO approach for obtaining $\ell_0$-constrained estimators in high-dimensional regression models and Reguant:2016 used mixed integer linear programming for computing counterfactual outcomes in game theoretic models.
AE1998 use data from the 1980 and 1990 U.S. census to estimate a model of the relation between the number of weeks per year a woman works and the number of children she has. A simplified but nonparametric version of their model is
where $Y$ is the number of weeks a woman works in a year; $\phi$ is an unknown function; and $D = 0,1$, or $2$ according to whether a woman has 2, 3, or 4 or more children. $D$ is endogenous. $Z$ is a binary instrument for $D$ equal to 1 if the first two children are of the same sex and 0 otherwise. We obtain bounds on the partially identified parameter $\phi(0)-\phi(1)$, which measures the change in the number of weeks a woman works when the number of children she has increases from two to three. We use data consisting of 394,840 observations $(Y,D,Z)$ from the 1980 U.S. census IPUMS.
We assume that $\phi$ is monotone non-increasing and focus on the parameter $f(\psi) = \phi(0)-\phi(1)$, where $\psi = [\phi(0), \phi(1), \phi(2)]'$. The population mean vector is $ \mu = [\textrm{vec} (\Pi), \textrm{vec} (\nu)], $ where
$p_{DZ} (d, z) := \text{Pr}( D = d, Z = z)$, and $I(\cdot)$ is the indicator function. The dependent variable $Y$ is contained in the interval [0,52]. In the analysis described below, we divide $Y$ by 52, so it is contained in the interval $[0,1]$, and denote the empirical analog of $\mu$ by $\bar{X}$. The inequality constraints are
The equality constraints are
The instrumental variable constraints in (ref) are bilinear in the sense that $\Pi \psi$ contains 6 bilinear terms.
We estimate the following 8 identification or 95% confidence intervals for $f(\psi)$.
We do not consider a confidence interval based on (ref) because $\textrm{B}(n,p,\overline{\mu}_3) \approx 2.49 \overline{\mu}_3$ (here, $n = 394,840$ and $p = 7$) and as a result, it is extremely unlikely that $\textrm{B}(n,p,\overline{\mu}_3) < 0.05$, although $\overline{\mu}_3$ is unknown. Computation was carried out on a MacBookPro laptop with an Apple M1 chip and 16 GB of memory.
The results, including computing times, are shown in Table (ref). The intervals are in weeks, not weeks divided by 52. The estimates obtained by methods 1-8 above are labelled in typewriter font: Sample; Box; Chi-Square; three types of Method 2; Method 3; and Minsker; respectively. As expected, the Sample interval is narrower than the Box and \texttt{Chi-Square} intervals, and the \texttt{Chi-Square} interval is narrower than the \texttt{Box} interval. The \texttt{Method 2 ($\Upsilon = \hat{\Sigma}$, $\sigma^2 = 1$)} and \texttt{Method 3} intervals are wider than the \texttt{Chi-Square} interval. The \texttt{Method 2 ($\Upsilon = \hat{\Sigma}$, $\sigma^2 = 1$)} interval is wider than the \texttt{Method 3} interval, which is based on a stronger assumption about the distribution of the observed random variables. The foregoing methods are motivated by consistency of the estimates of the population parameters they depend on. They do not take account of random sampling error in the estimates. The \texttt{Method 2 ($\Upsilon = I_7$, $\sigma^2 = 1$)} and \texttt{Minsker} intervals ensure finite-sample coverage probabilities of 0.95 but are wider than the other intervals. The \texttt{Method 2 ($\Upsilon = I_7$, $\sigma^2 = \hat{\sigma}^2 (\log n)$)} interval provides a tighter upper bound of 26.73 by estimating the variance proxy instead of relying on the bounded support of $X-\mu$. Only the \texttt{Sample} estimate yields an informative lower bound on $\phi(0) - \phi(1)$. Excluding the \texttt{Sample} bounds, which do not take account of random sampling error in the estimate of $\mu$, the upper bounds indicate that an increase in the number of children a woman has from 2 to 3 reduces her annual employment at most by approximately 11-37 weeks, depending on the estimation method. All of the computing times are small.
This section presents the results of Monte Carlo experiments that investigate the widths and coverage probabilities of nominal 95% confidence intervals for $f (\psi)$ that are obtained by using the methods described in Section (ref). For reasons explained in Section (ref), we concentrate on the upper confidence limit $\hat{J}_{+}$ and do not investigate $\hat{J}_{-}$. The experiments are designed to mimic the empirical application of Section (ref).
Data were generated by simulation from model (ref) as follows. Using the notation of Section (ref), set
The values of $\phi$ are similar to those used in FH:2015, who considered larger support of $D$, and the values of $\Pi$ are the same as the sample probabilities of $\Pi$ in the empirical example of Section (ref). Simulated values of $(D,Z)$ were drawn from $\Pi$. The outcome variable $Y$ was generated from
The parameter of interest is $f(\psi) = 52 [ \phi(0) - \phi(1)] = 5$, which is not point-identified. We impose that each element of $(\psi, \mu)$ is non-negative and impose the monotonicity constraint (ref). The population value of $[J_{-}, J_{+}]$ is $[0.381, 5.438]$. Thus, the data are not very informative about $J_{-}$. Consequently, we focus on $\hat{J}_{+}$, the estimate of $J_{+}$, in the experiments and do not investigate $J_{-}$.
We carried out experiments with $\mathcal{S}$ an ellipsoid and $\mathcal{S}$ a box. We report the averages and empirical coverage probabilities of $\hat{J}_{+}$ obtained by using 8 methods described in Section (ref). The nominal coverage probability was 95% and there were 500 Monte Carlo replications per experiment. In experiments in which $\mathcal{S}$ is a box, $\kappa_b(1-\alpha)$ was computed using $10^6$ random draws from $N(0,\widehat{\Sigma})$.
Figure (ref) shows only 100 (out of 500) realized values for each different estimated bound by sample size $n \in \{ 10^4, 10^5, 10^6 \}$. Each symbol corresponds to one Monte Carlo realization, and the population value $J_{+}$ is shown as a horizontal line. In the figure, both $x$ and $y$ axes are shown in the log scale with base 10. Table (ref) summarizes the results of the Monte Carlo experiments. The quantities in parentheses in the table are the values of $\hat{J}_{+}$ for the Method 2 ($\Upsilon = \hat{\Sigma}$, $\sigma^2 = 1$) and Method 3 methods obtained with the population values of the sub-Gaussian variance proxies, not estimates. {In the examples treated in these experiments, inference based on the true variance proxies and on consistent estimates of the proxies would be almost identical.}
The Sample bounds converge to the true upper bound as $n$ increases; however, they are not suitable for inference since sampling errors are ignored. All bounds get smaller as $n$ gets large, but the Chi-Square, Method 2 ($\Upsilon = \hat{\Sigma}$, $\sigma^2 = 1$), and Method 3 methods provide much tighter bounds than the other methods. The values of $\hat{J}_{+}$ for the Method 2 ($\Upsilon = \hat{\Sigma}$, $\sigma^2 = 1$) and Method 3 methods obtained with estimated and true values of the sub-Gaussian variance proxies are nearly equal. In every case except the \texttt{Sample} bounds, the estimated bounds are larger than $J_{+}$, resulting in the 100% empirical coverage probabilities. Overall, the simulation results are consistent with those of the empirical application in Section (ref).
Table (ref) shows that the bounds on $J_{+}$ vary greatly, depending on the method used and assumptions made about population parameters. In an application, a researcher might compute the bounds using all the methods and assumptions that are relevant to the application. The researcher can then decide which bounds to use or report, depending on his/her beliefs about the accuracy of the asymptotic approximations s/he is making and how much risk of inaccuracy s/he is willing to accept.
To check sensitivity to sub-Gaussian assumptions, we carried out additional Monte Carlo experiments by replacing $V$ in (ref) with $V \sim {(\chi^2(5) - 5)}/{\sqrt{10}}$. In other words, we use a standardized chi-square random variable, which is not sub-Gaussian but sub-exponential. The results of the additional Monte Carlo experiments are similar to those reported here with the standard normal $V$. See Online Appendix (ref) for details.
This paper has described a method for carrying out inference on partially identified parameters that are solutions to a class of optimization problems. The parameters arise, for example, in applications in which grouped data are used for estimation of a model's structural parameters. Inference consists of obtaining confidence intervals for the partially identified parameters. The paper has presented three methods for obtaining finite-sample lower bounds on the coverage probabilities of the confidence intervals. The methods correspond to three sets of assumptions of increasing strength about the underlying random variable. With the moderate sample sizes found in most economics applications, the bounds become tighter as the assumptions strengthen. The paper has also described a computational algorithm for implementing the methods. The results of Monte Carlo experiments and an empirical example illustrate the methods' usefulness. The paper has focused on the case in which a vector of first moments is the only unknown population parameter in the population version of the optimization problem. It might be useful to extend our formulation to the case in which the population optimization problem contains other population parameters. It is an open question how to obtain finite-sample lower bounds on coverage probabilities in such a case.