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.
111,446 characters · 18 sections · 132 citation commands
High Dimensional Binary Choice Model with Unknown Heteroskedasticity or Instrumental Variables
\onehalfspacing
\onehalfspacing
In this paper, we study the estimation of semiparametric binary choice models in high-dimensional settings. We adopt a standard model specification in which the 0-1 valued dependent variable $Y$ is modeled as \[ Y=\boldsymbol{1}(\boldsymbol{X}_{1}^{s^{\ast}\prime}\boldsymbol{\beta}^{\ast}+\varepsilon>0)=\boldsymbol{1}\left(X_{1}\beta_{1}^{\ast}+X_{2}\beta_{2}^{\ast}+...+X_{s^{\ast}}\beta_{s^{\ast}}^{\ast}+\varepsilon>0\right), \] where $\boldsymbol{X}_{1}^{s^{*}}\equiv(X_{1},X_{2},...,X_{s^{\ast}})'$ represent the explanatory variables (signals), $\boldsymbol{\beta}^{*}\equiv(\beta_{1}^{\ast},\beta_{2}^{\ast},...,\beta_{s^{\ast}}^{\ast})'$ are the unknown true coefficients to estimate, $\varepsilon$ represents an unobserved disturbance (error term), and $\boldsymbol{1}\left(\cdot\right)$ is the indicator function that equals one if $\cdot$ is true and 0 otherwise. In the context of high-dimensional settings, we allow for the possibility that either the number of candidate explanatory variables $\boldsymbol{X}_{1}^{p_{n}}\equiv\left(X_{1},X_{2},...,X_{p_{n}}\right)'$, denoted by $p_{n}$, diverges as the sample size $n$ increases, or the number of available instrumental variables (IV) $\boldsymbol{Z}_{1}^{p_{n}}\equiv\left(Z_{1},Z_{2},...,Z_{p_{n}}\right)'$, also denoted by $p_{n}$, diverges depending on the specific application.
Estimating binary choice models is far more challenging than estimating linear models due to the non-linearity of the indicator function. Researchers often resort to imposing distributional assumptions on the error term $\varepsilon$, such as assuming it follows a normal or logistic distribution, to enable feasible estimation. However, the likelihood of $\varepsilon$ precisely conforming to these specific distributions is very low. Complicating things further, heteroskedasticity in $\varepsilon$, common in economic data, can invalidate these distributional assumptions. Furthermore, when endogeneity issues are in play, which is also frequently encountered in empirical studies, these distributional assumptions provide no help in identifying $\boldsymbol{\beta}^{*}$, even if a sufficient number of valid IVs are available in $\boldsymbol{Z}_{1}^{p_{n}}$. One may address the endogeneity issue using the control function method detailed in standard textbooks like Wooldridge2010. The key idea is to include the residuals $\boldsymbol{u}$ obtained from regressing $\boldsymbol{X}_{1}^{s^{*}}$ on the valid IVs in $\boldsymbol{Z}_{1}^{p_{n}}$ in the regression, with the hope of controlling the endogeneity through additional control variables $\boldsymbol{u}$. However, it is essential to acknowledge an important limitation of the control function approach--it demands not just any set of valid IVs in $\boldsymbol{Z}_{1}^{p_{n}}$ but rather the exact right set of valid IVs. For example, missing any valid IVs of $\boldsymbol{Z}_{1}^{p_{n}}$ can lead to a violation of the requisite assumptions for this approach. We refer interested readers to Section 5 of LewbelDongYang for more details.
Dealing with many candidate explanatory variables or IVs adds another layer of difficulty. The processes of recovering the model's sparsity structure (e.g., Tibshirani1996), choosing the strongest instruments (e.g., Belloni_et_al2012), selecting the appropriate moment conditions (e.g., Liao2013 and ChengLiao2015), and conducting estimation and inferences in high-dimensional GMM (e.g., BelloniEtal2018, GautierTsybakov, and GoldEtal2020) become even more challenging in the context of the binary choice model. To our knowledge, existing literature on high-dimensional binary choice models often assumes independence between the error term $\varepsilon$ and the explanatory variables $\boldsymbol{X}$. Classic works like FanLi2001 and more recent studies such as chetverikov2023 and KhanEtal2023 rely on this assumption. However, these papers exclude the possibility of heteroskedasticity in the error term. Furthermore, we are not aware of any research on high-dimensional binary choice models that allow for endogeneity and instrumental variables in the usual sense.\footnote{Here, the usual sense means the following: A valid IV $Z$ for an endogenous variable $X$ ($\textrm{Cov}\left(X,\varepsilon\right)\neq0$) in linear models satisfies both the relevance and the exogeneity conditions, that is, $\textrm{Cov}\left(X,Z\right)\neq0$ and $\textrm{Cov}\left(\varepsilon,Z\right)=0$.} Another strand of literature considers binary prediction, possibly in high-dimensional settings. These works do not assume distributional assumptions and consider a framework similar to that in Manski1975. Important contributions in this direction include ElliottLieli2013, ChenLee2018, and babiiEtal2021. An interesting observation in babiiEtal2021 is that logistic regression can provide accurate predictions even if the independence between $\boldsymbol{X}$ and $\varepsilon$ fails to hold or $\varepsilon$ does not obey logistic distribution.
In this paper, we propose a feasible special regressor method to address these challenges. Our approach is semiparametric, as it does not rely on any distributional assumptions regarding $\varepsilon$ and allows for general heteroskedasticity. When certain explanatory variables are endogenous, our method only requires the availability of IVs that are valid in the conventional sense. The special regressor method, initially proposed in Lewbel2000, is based on the assumption that there exists a special regressor $V$ in the model
where the coefficient before $V$ is normalized to $1$. $V$ is assumed to be a continuous regressor with support larger than that of $-(\boldsymbol{X}_{1}^{s^{\ast}\prime}\boldsymbol{\beta}^{\ast}+\varepsilon)$ and satisfy $V\perp\left(\varepsilon,\boldsymbol{X}_{1}^{s^{\ast}}\right)|\boldsymbol{Z}_{1}^{p_{n}}$.\footnote{When all candidate explanatory variables are exogenous, $\boldsymbol{Z}_{1}^{p_{n}}=\boldsymbol{X}_{1}^{p_{n}}$ and the exclusion restriction on $V$ can be expressed as $V\perp\varepsilon|\boldsymbol{X}_{1}^{p_{n}}$.} Then, under the classic exogeneity condition, $\mathbb{E}\left(\boldsymbol{Z}_{1}^{p_{n}}\varepsilon\right)=0$, Lewbel2000 showed that
where $\tilde{Y}$ is defined as
Using the moment conditions in equation ((ref)), we can construct the regular GMM estimator to estimate $\boldsymbol{\beta}^{*}$. The main advantage of the special regressor approach is its ability to estimate the coefficients in a manner akin to linear models.
A drawback of the special regressor approach is its requirement to estimate the conditional density $f\left(V|\boldsymbol{Z}_{1}^{p_{n}}\right)$. Even when dealing with a moderate number of elements (e.g., $p_{n}=10$) in $\boldsymbol{Z}_{1}^{p_{n}}$, nonparametrically estimating $f\left(V|\boldsymbol{Z}_{1}^{p_{n}}\right)$ is challenging due to the curse of dimensionality. As a result, the special regressor estimation becomes infeasible if we allow the number of elements in $\boldsymbol{Z}_{1}^{p_{n}}$ to diverge. This paper is motivated by the observation that often only a few elements in $\boldsymbol{Z}_{1}^{p_{n}}$ are “relevant” to $V$. Then, the estimator becomes feasible, provided we can distinguish these relevant $Z$s from the irrelevant ones with high probability.
Interestingly, the foundational study by HallRacineLi revealed a surprising prevalence of “irrelevant” components in estimating conditional densities. Applying the special regressor approach, DongLewbel2015 adopted “negative age” as the special regressor to study inter-state migration. Their sample consists of 23 to 59 years old individuals who had completed education and were not retired in 1990. Of all the regressors, “negative age” is expected to be independent of factors such as “education” (since all individuals had completed their education), “gender”, “race”, and potentially “government benefit”. This can be further seen from the application of XueYangZhou, which studied the migration intention of rural residents in China using the special regressor method. The special regressor adopted in their paper is the “average daily precipitation”. Since precipitation is strongly exogenous and unlikely to be affected by individual characteristics, they assumed $V$ to be independent of all other regressors and only calculate $\hat{f}\left(V\right)$ for ((ref)). Notably, this study did not formally test the validity of this assumption. Our paper proposes a data-driven procedure to address this challenge rigorously.
Dimension reduction for conditional density estimation poses a significant challenge. The cross-validation (CV) approach, as proposed in HallRacineLi, is computationally demanding, making it difficult to apply in high-dimensional scenarios. Efromovich2010 required tensor products of basis over each dimension; this makes his approach infeasible in high-dimensional settings. The approach in ZhangJMP2023 suffers from the same issue because it needs to compute conditional expectations at each dimension. In this paper, we present an innovative dimension reduction technique inspired by the “sure independence screening” method introduced in FanLv2008 for linear models. To mitigate the high computational burden associated with high-dimensional linear models, FanLv2008 suggested regressing the dependent variable on each explanatory variable individually and retaining only those with the strongest correlation to the dependent variable. However, using correlation to measure dependence is inappropriate for our goals, as it captures only linear relationships. Consequently, we adopt the concept of “distance covariance” (DC), first detailed and studied in Szekely_et_al. As a metric of general dependence, DC is zero if and only if the two random vectors under examination are independent.
Leveraging this desirable property, we propose a DC-based screening procedure. When estimating $f\left(V|\boldsymbol{Z}_{1}^{p_{n}}\right)$, the primary strategy is to employ the screening method to reduce the dimension of the conditioning set. Subsequently, the CV procedure from HallRacineLi is applied post-screening to refine the screening outcomes further. Integrating these two methods bypasses the often complicated task of determining optimal tuning parameters. This procedure also reduces the necessity for “perfect” variable selection during the screening phase, offering a more dependable and practical approach to estimate $f\left(V|\boldsymbol{Z}_{1}^{p_{n}}\right)$ in the context of high-dimensional data.\footnote{This two-step selection approach can be extended to other nonparametric estimation scenarios, making it of independent interest.} We extensively explore the theoretical attributes of our proposed selection method and the subsequent conditional density estimation. Our findings demonstrate that our selection technique can identify the components in $\boldsymbol{Z}_{1}^{p_{n}}$ that are relevant to $V$ with a high probability. Furthermore, the post-selection conditional density estimator achieves the “oracle” rate of convergence.
With the feasible estimator for $f\left(V|\boldsymbol{Z}_{1}^{p_{n}}\right)$ in hand, we proceed to address the variable selection and moment (IV) selection problems for high-dimensional binary choice models. Specifically, for the case where there are many exogenous explanatory variables $\boldsymbol{X}_{1}^{p_{n}}$, we introduce two estimators obtained from combining linear regression with either weighted or unweighted $L_{1}$ penalties. These estimators can estimate $\boldsymbol{\beta}^{\ast}$ while selecting true signals within $\boldsymbol{X}_{1}^{p_{n}}$ simultaneously. When there are many candidate IVs $\boldsymbol{Z}_{1}^{p_{n}}$, a mixture of valid and invalid, we propose a GMM with $L_{1}$ penalty procedure, which can simultaneously estimate $\boldsymbol{\beta}^{\ast}$ and detect invalid IVs within $\boldsymbol{Z}_{1}^{p_{n}}$. We establish the asymptotic properties of these approaches, providing a solid theoretical foundation for their application. Furthermore, we demonstrate the practical effectiveness of these methods in finite samples using both simulated and real-world data.
The rest of this paper is organized as follows. In Section (ref), we propose a novel dimension reduction method for the conditional density estimation. With the feasible conditional density estimator, we address the classic variable and moment (IV) selection problems in the context of special regressor estimation in Section (ref). We investigate the small sample properties of our methods through Monte Carlo experiments in Section (ref), where we also provide a practical guide for choosing tuning parameters and kernel functions to implement our proposed procedures. In Section (ref), we further illustrate our approaches by examining migration decisions of rural residents in China. Section (ref) concludes this paper.
Proofs and tables are presented in the online supplementary appendix. Specifically, Appendix (ref) provides a gentle introduction to the DC and CV procedures. The proofs of all the main theorems are collected in Appendices (ref)--(ref), while the proofs of all technical lemmas are deferred to Appendix (ref). An investigation of variable selections with unweighted $L_{1}$ penalty is presented in Appendix (ref). Lastly, Appendix (ref) contains tables for simulation and application results.
For ease of reference, we list the notations maintained throughout this paper here.
Notation. All vectors are column vectors. Unless stated otherwise, we follow the convention of using capital letters for random variables and their corresponding lowercase letters for sample realizations. We use $\rho$ to denote the eigenvalues of a matrix. For instance, $\rho_{\min}$ and $\rho_{\max}$ denote a matrix's minimum and maximum eigenvalues, respectively. The notation $\left\Vert \boldsymbol{x}\right\Vert $ represents the Euclidean norm of a vector $\boldsymbol{x}$. For a matrix $\mathbf{A},$ we define $\left\Vert \mathbf{A}\right\Vert \equiv\sqrt{\text{trace}\left(\mathbf{AA}^{\prime}\right)}$ and $\left\Vert \mathbf{A}\right\Vert _{\infty}\equiv\max_{jl}\mathbf{A}_{jl}$, where $\mathbf{A}_{jl}$ is the $\left(j,l\right)$-th element of $\mathbf{A}$. The symbol $\left\Vert \cdot\right\Vert _{0}$ represents the $L_{0}$-norm, which counts the total number of nonzero elements in a vector. $\left\Vert \cdot\right\Vert _{1}$ denotes the $L_{1}$-norm. For deterministic series $\left\{ a_{n}\right\} _{n=1}^{\infty}$ and $\left\{ b_{n}\right\} _{n=1}^{\infty}$, the notation $a_{n}\propto b_{n}$ means that $0<C_{1}\leq\lim\inf_{n\rightarrow\infty}\left\vert a_{n}/b_{n}\right\vert \leq\lim\sup_{n\rightarrow\infty}\left\vert a_{n}/b_{n}\right\vert \leq C_{2}<\infty$ for some constants $C_{1}$ and $C_{2}$, $a_{n}\lesssim b_{n}$ means $\lim\sup_{n\rightarrow\infty}\left\vert a_{n}/b_{n}\right\vert \leq C<\infty$ for some constant $C$, $a_{n}\gtrsim b_{n}$ if $b_{n}\lesssim a_{n},$ and $a_{n}\ll b_{n}$ if $a_{n}=o\left(b_{n}\right),$ and $a_{n}\gg b_{n}$ if $b_{n}\ll a_{n}$. The term $\mathcal{A}^{c}$ stands for the complement of the set $\mathcal{A}$, and $\left\vert \mathcal{A}\right\vert $ is the number of elements in $\mathcal{A}$. As $n$ tends to infinity, the notations $\overset{P}{\rightarrow}$ and $\overset{d}{\rightarrow}$ indicate convergence in probability and distribution, respectively. The symbol $C$ denotes various positive constants, which may change from one instance to the next.
As previously mentioned, an essential prerequisite for implementing the special regressor estimator is the ability to estimate $f\left(V|\boldsymbol{Z}_{1}^{p_{n}}\right)$, which is used to construct $\tilde{Y}$ as defined in equation ((ref)). To achieve this, we presuppose a certain sparsity in the conditional set of $f\left(V|\boldsymbol{Z}_{1}^{p_{n}}\right)$. Specifically, we assume that only a subset of the $Z$s are included in the conditional set and refer to them as “relevant to $V$”. This assumption is detailed in Assumption (ref). To measure the dependence between $V$ and the $Z$s, we use the distance covariance (DC), first introduced and examined in Szekely_et_al. Section (ref) provides a brief overview of this concept. In Section (ref), we describe our screening process and outline its theoretical properties. We propose a practical dimension reduction technique in Section (ref) by combining our screening method with the approach from HallRacineLi. We further extend the theory from HallRacineLi in Section (ref), allowing for applying the dimension reduction approach to a broader range of scenarios. Section (ref) summarizes the practical procedure for easier reference.
We emphasize that this novel approach to dimension reduction for conditional density is the main innovation of this paper. This method has the potential to be adapted to various scenarios involving high-dimensional conditional density estimation, including the studies discussed in Sections (ref) and (ref).
We measure the dependence between $V$ and another regressor or IV $Z_{l}$ using DC introduced and studied in Szekely_et_al. We briefly review this statistic in this section and refer interested readers to Appendix (ref) and Szekely_et_al for more details.
Let $\phi_{V}\left(t\right)$, $\phi_{Z_{l}}\left(s\right)$ and $\phi_{V,Z_{l}}\left(t,s\right)$ denote the characteristic functions of $V,Z_{l},$ and $\left(V,Z_{l}\right),$ respectively. Specifically, \[ \phi_{V,Z_{l}}\left(t,s\right)=\mathbb{E}\left(e^{\mathrm{i}tV+\mathrm{i}sZ_{l}}\right)=\int\int e^{\mathrm{i}tv+\mathrm{i}sz_{l}}f\left(v,z_{l}\right)dvdz_{l}, \] where $\mathrm{i}$ denotes the imaginary unit, and $\phi_{V}\left(t\right)\textrm{ and }\phi_{Z_{l}}\left(s\right)$ are similarly defined. The DC between univariate random variables $V$ and $Z_{l}$ is defined as \[ \mathcal{V}^{2}\left(V,Z_{l}\right)=\int_{\mathbb{R}^{2}}\left\Vert \phi_{V,Z_{l}}\left(t,s\right)-\phi_{V}\left(t\right)\phi_{Z_{l}}\left(s\right)\right\Vert \omega\left(t,s\right)dtds, \] where $\left\Vert \phi\right\Vert ^{2}=\phi\bar{\phi}$ for a complex-valued function $\phi$, with $\bar{\phi}$ being the conjugate of $\phi$, and \[ \omega\left(t,s\right)=\frac{1}{\pi^{2}t^{2}s^{2}}. \] One nice property of DC is that $\mathcal{V}\left(V,Z_{l}\right)=0$ if and only if $V$ and $Z_{l}$ are independent. While other positive weighting functions can ensure this property as well, this particular choice of $\omega(t,s)$ (also recommended in Szekely_et_al) is notable for yielding a very simple sample analog, as follows.
To estimate $\mathcal{V}^{2}\left(V,Z_{l}\right)$, the sample DC is \[ \mathcal{V}_{n}^{2}\left(V,Z_{l}\right)=S_{n1}\left(V,Z_{l}\right)+S_{n2}\left(V,Z_{l}\right)-2S_{n3}\left(V,Z_{l}\right), \] where
and $(v_{j},z_{lj}),j=1,2,...,n,$ are i.i.d. realizations of $(V,Z_{l})$. Szekely_et_al showed that \[ \mathcal{V}_{n}^{2}\left(V,Z_{l}\right)\rightarrow\mathcal{V}^{2}\left(V,Z_{l}\right)\text{ almost surely.} \] and proposed the test-statistic for null hypothesis $H_{0}:V\perp Z_{l}$ defined as $n\mathcal{V}_{n}^{2}(V,Z_{l})/S_{n2}(V,Z_{l})$, where the denominator makes the test statistic scale-free (e.g., increasing $V$ or $Z_{l}$ by ten times does not alter the test statistic).
As will be detailed in Section (ref), we propose a DC screening procedure in this paper based on the following statistic:
Estimating conditional density in high-dimensional settings is a challenging problem. The well-known paper, HallRacineLi, proposed a data-driven CV procedure to automatically reduce the dimension of conditional density estimation. However, it is infeasible to handle the high dimensional case due to the high computation cost and the curse of dimensionality. To tackle this problem, we use the screening technique initially proposed by FanLv2008 for ultra-high dimensional linear models. The idea of the screening in the linear model case is to calculate the correlation between the dependent variable and one independent variable at a time and keep the independent variables with the highest correlation with the dependent variable for the model, e.g., 5% of the candidate independent variables.
In the context of nonparametric conditional density estimation, the work most closely related to this paper is LiZhongZhu, who extended the screening method to explore more general nonlinear relationships between the dependent variable and independent variables, leveraging the concept of distance correlation.\footnote{The relationship between distance correlation and distance covariance is analogous to that of ordinary correlation and covariance. In LiZhongZhu, the empirical distance correlation is defined as $\mathcal{V}_{n}(V,Z_{l})/\sqrt{\mathcal{V}_{n}(V,V)\cdot\mathcal{V}_{n}(Z_{l},Z_{l})}$.} Using a similar idea, we measure the dependence between $V$ and another regressor or IV $Z_{l}$ one at a time, employing the statistic defined in ((ref)). It is worth noting that LiZhongZhu imposed restrictive sub-Gaussian tail assumption. Our work differs from theirs by allowing for much heavier-tailed distributions, as stated in the second part of Assumption (ref), presented below. This flexibility significantly broadens the applicability of the DC-based screening approach, particularly in economics, where thin-tail conditions like sub-Gaussian distributions can be rather restrictive.
Before presenting the key assumptions and theoretical results, we introduce some technical terms. Let $(\tilde{V},\tilde{Z}_{l})$ denote an independent copy of $(V,Z_{l})$. We define the quantity \[ S_{2}\left(V,Z_{l}\right)=\mathbb{E}(\vert V-\tilde{V}\vert)\mathbb{E}(\vert Z_{l}-\tilde{Z}_{l}\vert). \] Additionally, we introduce the term $\kappa_{n}$, which represents a slowly diverging sequence with the rate of divergence no faster than $\log n$.
{0pt}
{1em}
Equation ((ref)) is essential for the screening procedure. Equation ((ref)) is to ensure we can apply the method in HallRacineLi. Both equations ((ref)) and ((ref)) can be implied by condition (19) in HallRacineLi. Alternatively, we can express the sparsity of the conditional density of $V$ given $\boldsymbol{Z}_{1}^{p_{n}}$ solely through conditional independence, specifically using equation ((ref)) only. However, as discussed intensively in HallRacineLi, the definition of conditional independence per se can be ambiguous.\footnote{That is the following: $V\perp Z_{1}|Z_{2},$ $V\perp Z_{2}|Z_{1},$ but $V\not\perp$ $\left(Z_{1},Z_{2}\right)$ may hold at the same time.} We refer readers to their paper's example of linear combinations of standard normal random variables (on page 1016). To avoid such ambiguity, HallRacineLi focus on conventional independence, which can imply both equations ((ref)) and ((ref)).
The DC screening procedure proposed below can only guarantee to identify $(Z_{p^{*}+1},Z_{p^{*}+2}...,Z_{p_{n}})$ that satisfy both ((ref)) and ((ref)) in Assumption (ref)(2) as irrelevant. There are situations where certain covariates are dependent on $V$, but only indirectly through other covariates. For example, $V$ is a function of $Z_{1}$, and $Z_{2}$ is correlated with $Z_{1}$ but independent of all other determinants of $V$. This violates Assumption (ref)(2), and as a result, the DC screening is likely to identify $Z_{2}$ as relevant. To address this issue, we generalize the findings in HallRacineLi in Section (ref), showing that applying their CV-based bandwidth selection algorithm post-screening can eliminate the influence of such $Z_{2}$ on $V$ when estimating the conditional density of $V$ on $(Z_{1},Z_{2})$.
Assumption (ref)(3) imposes mild moment conditions for $V$ and $\boldsymbol{Z}_{1}^{p_{n}}$. It is apparent that higher dimension $p_{n}$ requires more restrictive moment conditions, i.e., larger $\delta$. We restrict $\delta\geq1$ so that the selection error can be negligible when establishing Theorem (ref). The false discovery rate (FDR) defined in equation ((ref)) below can also be better controlled for larger $\delta$ due to sharper error bounds. We allow $V$ and $Z_{l}$ to have as low as $\left(2+\delta\right)$-th finite moment, but jointly $\mathbb{E}(\left\vert VZ_{l}\right\vert ^{2+\delta})$ should be finite. To accommodate ultra-high dimensional data (i.e., $p_{n}\propto\exp(n^{C})$ for some positive $C$), we do need to impose thin tail restrictions, such as sub-Gaussian, on $V$ and $\boldsymbol{Z}_{1}^{p_{n}}$ (see, e.g., Assumption C1 in LiZhongZhu). However, ultra-high dimensional data are rare in economics. In this sense, we consider our assumption quite general. Assumption (ref)(4), which distinguishes relevant variables from irrelevant ones, is also made in LiZhongZhu
Let $\varsigma_{n}$ denote the threshold for the test statistic $\mathcal{\hat{T}}_{nl}$ defined in equation ((ref)). We determine $Z_{l}$ to be relevant to $V$ if \[ \widehat{\mathcal{J}}_{l}\equiv\mathbf{1}(\mathcal{\hat{T}}_{nl}\geq\varsigma_{n})=1. \] We define the set of indices of the truly relevant variables and the set of indices of selected relevant variables, respectively, as
Following the multiple test literature (see, e.g., BenjaminiHochberg), we adopt the true positive rate (TPR) and the FDR defined as
to measure the performance of the screening procedure proposed above. The following results are on the Type-I and Type-II errors, which are useful for showing the properties of TPR$_{n}$ and FDR$_{n}$ summarized in Theorem (ref).
{0pt}
{1em}
Theorem (ref) shows that with some suitably chosen $\varsigma_{n}$, the proposed test statistic $\mathcal{\hat{T}}_{nl}$ can separate the relevant variables from irrelevant variables with high probability. However, there is some power loss up to $\sqrt{\kappa_{n}\log n}$. The error bound is directly linked to $\delta$. The larger $\delta$ or more restrictive moment conditions are, the sharper the bound becomes. There is some middle ground where our procedure is silent. Built on Theorem (ref), the theorem below shows the properties of our screening procedure in terms of TPR$_{n}$ and FDR$_{n}$.
{0pt}
{1em}
In addition to the properties of TPR$_{n}$ and FDR$_{n}$, Theorem (ref) also establishes the “oracle” property that our screening procedure can select the truly relevant variables with probability approaching one. However, the theoretical rate of $\varsigma_{n}$ only provides a limited guide in practice. For this reason, we propose a practical procedure for choosing the threshold in the following subsection, utilizing the CV approach proposed in HallRacineLi. It is worth noting that this procedure does not rely on the perfect selection.
In this section, we first discuss how to choose $\varsigma_{n}$ based on the results stated in Theorem (ref). Then, we propose a more practical way of handling variable selection, which utilizes the method by HallRacineLi. After briefly reviewing their method, we present our procedure and derive the convergence rate of the post-selection conditional density estimator. This result is summarized in Theorem (ref).
Our procedure involves keeping as many as $Z_{l}$ that the method of HallRacineLi can handle from the initial screening stage, after which we apply the method of HallRacineLi to reduce the impact of included irrelevant variables. This is similar to the approach proposed by FanLv2008, who suggested keeping 5% of the most correlated covariates and then using Lasso to refine the variable selection further. More importantly, we will show in Section (ref) that the method of HallRacineLi can eliminate the impact of certain covariates that may depend on $V$ but only through some other relevant covariates.
Choosing $\varsigma_{n}^{\ast}$ based on theoretical results
Based on the rates derived in Theorems (ref) and (ref), one can choose
for some positive constant $C$, which is equivalent to setting $\kappa_{n}=\sqrt{\log n}$. By doing so, we can identify relevant $Z_{l}$ satisfying $\mathcal{V}^{2}\left(V,Z_{l}\right)/S_{2}\left(V,Z_{l}\right)\gtrsim n^{-1/2}\kappa_{n}\sqrt{\log n}=n^{-1/2}\log n$ with high probability. Note that this choice results in a loss of power up to the order of $\log n$.
A practical way of choosing $\varsigma_{n}^{\ast}$
We propose choosing the threshold following the idea of FanLv2008. Specifically, we calculate $\mathcal{\hat{T}}_{nl}$ for all $l=1,2,...,p_{n}$ and then rank them from highest to lowest. We then retain the $Z$s with the highest $\mathcal{\hat{T}}_{nl}$ values while ensuring that we do not include too many variables to the extent that the method proposed in HallRacineLi becomes impractical. Denote the number of retained $Z$s as $\tilde{p}$. Note that $\tilde{p}$ is a fixed integer that does not change with $n$. Without loss of generality, assume that $Z_{1},Z_{2},...,Z_{\tilde{p}}$ are the $\tilde{p}$ variables with the highest $\mathcal{\hat{T}}_{nl}$ values. We then proceed to estimate the conditional density using the method presented in HallRacineLi, which effectively reduces the dimension of the conditional set and automatically eliminates the impact of irrelevant variables, as demonstrated in their work.
Here, we briefly overview the procedure proposed in HallRacineLi with more details in Appendix (ref). We recommend referring to the original paper for a more comprehensive understanding of the details and theoretical properties. Throughout our presentation, we will assume that all $Z$s are continuously distributed for brevity. First, we define some technical terms:
where $K_{h}\left(x\right)=h^{-1}K\left(x/h\right)$ and $K(\cdot)$ is a standard kernel function. The conditional density of $V$ on $\left(Z_{1},Z_{2},...,Z_{\tilde{p}}\right)$ is estimated as
The CV criterion is defined as
where $\hat{I}_{n1}$ and $\hat{I}_{n2}$ are computed as \[ \hat{I}_{n1}\left(h_{1},...,h_{\tilde{p}}\right)=\frac{1}{n}\sum_{i=1}^{n}\frac{\hat{G}_{-i}\left(z_{1i},...,z_{\tilde{p}i}\right)}{\hat{f}_{-i}\left(z_{1i},...,z_{\tilde{p}i}\right)^{2}}\text{ and }\hat{I}_{n2}\left(h_{V},h_{1},...,h_{\tilde{p}}\right)=\frac{1}{n}\sum_{i=1}^{n}\frac{\hat{f}_{-i}\left(v_{i},z_{1i},...,z_{\tilde{p}i}\right)}{\hat{f}_{-i}\left(z_{1i},...,z_{\tilde{p}i}\right)} \] with $\hat{G}_{-i}$ and $\hat{f}_{-i}$ being the leave-one-out (leave the $i$-th observation out) estimator. Note that we omit the weighting function used in HallRacineLi for conciseness.
We can compute the optimal bandwidth
and then plug $\hat{\boldsymbol{h}}$ into equation ((ref)) to obtain $\hat{f}\left(v|z_{1},...,z_{\tilde{p}}\right)$.
As demonstrated in Theorem 2 of HallRacineLi, the CV procedure ((ref)) can select asymptotically optimal bandwidths for relevant variables and push the bandwidths for irrelevant variables to their upper limits, effectively making them appear uniform in the conditional set and thus less informative to $V$, with the probability approaching one as $n\rightarrow\infty$.
We provide the convergence rate of the post-selection conditional density estimation in Theorem (ref), which is an immediate result from Theorem 2 of HallRacineLi and our Theorem (ref). The following technical conditions are necessary, akin to those in HallRacineLi. One note is that we need to assume $\tilde{p}\geq p^{\ast}$ to ensure we do not miss any relevant $Z$s. Although we assume all covariates are continuous below, our procedure can also accommodate discrete variables. Demonstrating the validity of this extension is straightforward, given HallRacineLi.
Note that we omit some minor technical details required for the above theorem, specifically, conditions (22)--(24) in HallRacineLi. The main implication of Theorem (ref) is that irrelevant variables have no impact on the convergence rate of $\hat{f}\left(v|z_{1},...,z_{\tilde{p}}\right)$ since the fastest convergence rate achievable for $\hat{f}\left(v|z_{1},...,z_{\tilde{p}}\right)$ would be $n^{-\frac{r}{p^{\ast}+1+2r}}$ if we only included the $p^{\ast}$ truly relevant variables in the estimation. When $r>\left(p^{\ast}+1\right)/2$, the nonparametric estimation has an error on the order of $o_{P}\left(n^{-1/4}\right)$, which meets the minimum requirement for parametric estimation with nuisance nonparametric estimation as the first step. The $o_{P}\left(n^{-1/4}\right)$ holds uniformly (with some loss of power of $\log n$), which can be demonstrated easily as in LiRacine. Since $p^{\ast}$ is unknown, practitioners may want to choose a high-order kernel with order $r>\left(\tilde{p}+1\right)/2.$
By utilizing the method proposed in HallRacineLi, we do not need a perfect variable selection during the initial screening for conditional density estimation. We may be reasonably conservative and include some potentially irrelevant variables in the first-step screening, but they will have no impact on the asymptotic convergence rate. It is worth noting that HallRacineLi does not eliminate these irrelevant variables from the estimation. Instead, their method reduces the influence of irrelevant variables by assigning them large bandwidths selected through the CV procedure ((ref)).
Consider the following illustrative scenario:
where $g(\cdot,\cdot)$ is an unknown continuous function. Apparently, $V\perp\left(Z_{2},Z_{3},...,Z_{s}\right)|Z_{1}$, but $V\not\perp Z_{1}|(Z_{2},...,Z_{s})$ due to the presence of $U_{1}$. Thus, the conditional independence is unambiguous in this case. However, it is possible that $\left(Z_{2},Z_{3},...,Z_{s}\right)$ and $V$ are dependent through $U^{\ast}$. If such dependence exists, then ((ref)) violates Assumption (ref)(2) (and hence condition (19) in HallRacineLi). Proposition (ref) presented below demonstrates that the procedure proposed in Section (ref) can be applied to scenario ((ref)), yielding the same convergence rate as obtained in Theorem (ref). This broadens the applicability of the theory of HallRacineLi.
Note that in our exposition, we use ((ref)) as an illustrating example, where only a scalar $Z_{1}$ enters $g(\cdot,\cdot)$. We intentionally adopt this simplified setup to facilitate understanding and avoid tedious discussions on the potential ambiguity of conditional (in)dependence that can arise with different sub-vectors of $(Z_{1},...,Z_{s})$. In practice, when the conditional independence is ambiguous, the CV criterion may have multiple local minima, and we may be unable to identify the global minimum. Extending our approach to cases with multiple relevant $Z$s in $g(\cdot,\cdot)$ is straightforward, albeit more technically involved. Due to space constraints, we do not pursue this direction formally in this paper.
To understand the significance of the above corollary, note that it obtains the same rate as Theorem (ref) when $p^{\ast}=1$.
However, our method is unsuitable in scenarios where all noise variables are highly dependent on $V$ (e.g., through $Z_{1}$), as our screening procedure fails to exclude these noise variables.
{0pt} We conclude this section by summarizing our proposed procedure for conditional density estimation with dimension reduction:
{1em}
In this section, we explore high-dimensional binary choice models. We start by using the feasible estimator $\hat{f}(v|\boldsymbol{z}_{1}^{\tilde{p}})$, as derived in Section (ref), to approximate $f\left(v|\boldsymbol{z}_{1}^{p_{n}}\right)$. This leads to the feasible $\tilde{Y}$ in equation ((ref)). We then apply the special regressor approach in high-dimensional settings. The analysis in this section follows existing methods but incorporates an estimated $\tilde{Y}$.
To facilitate understanding in this section, consider viewing $\hat{f}(v|\boldsymbol{z}_{1}^{\tilde{p}})$ from the previous section as a “regular” nonparametric estimator and ignore that it is a post-selection estimator. Then, the estimators introduced in this section are treated as building upon existing methods but incorporating a nonparametric first-step plug-in. We emphasize that we have accounted for the selection issue related to the estimation of $f(v|\boldsymbol{z}_{1}^{p_{n}})$ when deriving all the theoretical results. This treatment is demonstrated in Appendices (ref) and (ref).
We begin by reviewing the special regressor estimator and presenting the conditions required to validate the moment conditions defined in equation ((ref)). With a feasible estimator for the conditional density, these moment conditions enable the estimation of the model parameters. We can then apply either an LS or a GMM estimator, incorporating appropriate regularization techniques commonly used for linear models. To illustrate this, we explore the classic variable selection problem in Section (ref) and discuss how to choose valid moment conditions from a large set of candidate IVs in Section (ref). Other scenarios, such as selecting optimal instruments (see, e.g., Belloni_et_al2012), can be similarly studied but may involve more complex technical details.
It is important to note that the problems discussed in Sections (ref) and (ref) are fundamentally distinct: The notations used in each section are specific to their context and may have different meanings, even though they might appear identical.
The following technical conditions placed in Lewbel2000 are required for the validity of the special regressor estimator.
Note that in the case where all candidate covariates $\boldsymbol{X}_{1}^{p_{n}}$ are exogenous, $\boldsymbol{Z}_{1}^{p_{n}}$ is simply $\boldsymbol{X}_{1}^{p_{n}}$, and hence Assumption (ref)(2) can be expressed as $V\perp\varepsilon|\boldsymbol{X}_{1}^{p_{n}}$. Assumption (ref)(3) is the “large support” condition for $V$. This condition is similar to the “overlap” condition for the average treatment effects estimator, which assumes that the propensity score is bounded away from zero and one. The “large support” condition can be restrictive. For instance, the special regressor “age” in DongLewbel2015 and “precipitation” in XueYangZhou have bounded support, which excludes the inclusion of $\varepsilon$ with larger support. Lastly, Assumption (ref)(4) is a technical condition that is necessary to prevent the “irregular slower-than-$\sqrt{n}$ convergence” property, as explained in KhanTamer.
The special regressor approach comes at the cost of imposing strong assumptions on $V$. Specifically, these assumptions require $V$ to have strong exogeneity and large support. However, despite these constraints, these assumptions yield the following useful identity: \[ \mathbb{E}[\tilde{Y}-X_{1}\beta_{1}^{\ast}-X_{2}\beta_{2}^{\ast}-...-X_{s^{\ast}}\beta_{s^{\ast}}^{\ast}\left\vert \boldsymbol{Z}_{1}^{p_{n}}\right.]=\mathbb{E}\left[\varepsilon\left\vert \boldsymbol{Z}_{1}^{p_{n}}\right.\right] \] with
which implies the moment conditions that resemble those for the linear IV models: For $j\in\{1,...,p_{n}\}$, \[ \mathbb{E}[Z_{j}(\tilde{Y}-X_{1}\beta_{1}^{\ast}-X_{2}\beta_{2}^{\ast}-...-X_{s^{\ast}}\beta_{s^{\ast}}^{\ast})]=0, \] if $\mathbb{E}\left[\varepsilon|Z_{j}\right]=0$ is satisfied.\footnote{To see this, consider
} We refer interested readers to Lewbel2000 for a detailed discussion on the technical conditions, the proofs, and the estimators. In subsequent sections, we explore the applications of the special regressor estimators in high-dimensional settings.
In this section, we attempt to replicate some findings regarding the adaptive Lasso in a high-dimensional setting. The adaptive Lasso, introduced by Zou2006, is a widely used penalized regression method in econometrics and statistics. This approach was extended to high-dimensional settings in subsequent studies (HuangEtal2008 and LinEtal2009). HuangEtal2008 examined scenarios with an extremely high dimension, where the number of covariates increases exponentially. In these cases, it was assumed that both the distributions of covariates and error terms decay at an exponential rate to manage the extremely high-dimensional covariates. On the other hand, LinEtal2009 eased the tail (or moment) conditions set by HuangEtal2008 and concentrated on moderately high-dimensional situations, where the number of covariates increases at polynomial rates.
In typical economic applications, while economists might handle numerous covariates, their number is generally far less than that of observations. Additionally, many economic variables, such as wages, exhibit heavy-tailed distributions, making the assumption of light tails rather restrictive. Given these considerations, we focus on replicating results for the scenario described in LinEtal2009, as it presents a more practical case for economic applications.
Suppose the true model is given by
where $(X_{1},\dots,X_{s^{*}})$ are the true signals among all candidate covariates $\boldsymbol{X}_{1}^{p_{n}}$. Researchers do not know which components of $\boldsymbol{X}_{1}^{p_{n}}$ correspond to $(X_{1},\dots,X_{s^{*}})$ until the data reveal this. We assume $\boldsymbol{X}_{1}^{p_{n}}$ to be exogenous, and thus it plays the role of $\boldsymbol{Z}_{1}^{p_{n}}$ in Assumption (ref). We use the notation $X$ instead of $Z$ in this section to avoid confusion. The true parameter values are collected into a $p_{n}\times1$ vector $\boldsymbol{\beta}^{*}=(\beta_{1}^{*},\beta_{2}^{*},\dots,\beta_{s^{*}}^{*},0,0,\dots,0)^{\prime}$, corresponding to the true signal variables, followed by zeros for the remaining candidate covariates. It is worth noting that the true signals, $(X_{1},\dots,X_{s^{*}})$, may or may not overlap with the $X$s relevant to $V$, because the two issues (namely, “true signals in the binary choice model” and “relevant $X$s to $V$”) are generally independent. However, we do not reorder $X$ to distinguish these two issues in this section for notational convenience.
We denote the feasible estimator of $\tilde{y}_{i}$ by
where $\hat{f}(v_{i}|\boldsymbol{x}_{1i}^{\tilde{p}})$ is obtained using the procedure outlined in Section (ref), but with a $2r$-th order kernel function (refer to Theorem (ref) and the subsequent discussion for a detailed explanation).
The adaptive Lasso estimator for our case is then obtained by:
where $\boldsymbol{x}_{i}\equiv(x_{1i},\dots,x_{s^{*}i},\dots,x_{p_{n}i})^{\prime}$, and \[ \varpi_{j}=|\tilde{\beta}_{j}|^{-\gamma},\quad j=1,\dots,p_{n}, \] are adaptive weights used for penalizing different coefficients in the $L_1$ penalty, with $\gamma>0$, and $\tilde{\beta}_{j}$ being an initial estimator, as discussed in Section (ref).
We make the following assumptions, which are similar to those placed in LinEtal2009.
While LinEtal2009 assumed $\mathbf{X}_{n} \equiv (\boldsymbol{x}_{1},\ldots,\boldsymbol{x}_{n})'$ to be non-random, we extend their work to accommodate random $\mathbf{X}_{n}$. Assumption (ref)(1) is the full rank condition, the same as in LinEtal2009. Assumption (ref)(2) assumes the existence of an initial estimator converging to the true parameters at a polynomial rate and imposes the so-called “beta-min" condition, which is crucial for deriving the “oracle" property and conducting inference. Otherwise, well-known impossibility results for the uniform validity of post-selection estimators, such as those obtained in LeebPotscher2006,LeebPotscher2008, apply. Assumption (ref)(4) is standard considering $\Omega_{n}$ is a variance matrix. We now discuss Assumption (ref)(3) in detail. We impose a slightly stricter condition than LinEtal2009 due to the relaxation of their setting to one with random $\mathbf{X}_{n}$. The trade-off is that the rate of $s^{*}$ imposes some restrictions on the rates of tuning parameters as described in ((ref)). For intuition, by setting $C_{s}=0$, we need $1/2-\varrho\gamma<C_{\lambda}<1/2-\alpha\gamma$, and this set is not empty because $\alpha<\varrho$. For $\gamma$, we only require $\gamma>0$. It is clear that the larger the $\gamma$, the wider the valid range of $C_{\lambda}$. When $C_{s}$ is positive, the range of valid $C_{\lambda}$ is narrower, and we need $\gamma$ to be sufficiently large to ensure the range of valid $C_{\lambda}$ is not empty.
With these technical conditions in place, we establish theoretical results for the adaptive Lasso estimator, akin to Theorem 1 of LinEtal2009. It is important to note that our goal is not to improve existing results on Lasso estimation, which is beyond the scope of this paper. We define \[ \Sigma_{n}\equiv\left(\left[\mathbb{E}\left(\boldsymbol{X}_{1}^{s^{*}}\boldsymbol{X}_{1}^{s^{*}\prime}\right)\right]^{-1}\Omega_{n}\left[\mathbb{E}\left(\boldsymbol{X}_{1}^{s^{*}}\boldsymbol{X}_{1}^{s^{*}\prime}\right)\right]^{-1}\right), \] which is useful for the asymptotics of $\hat{\boldsymbol{\beta}}$. $\Omega_n$ is defined in ((ref)).
In this theorem, we require that $r>\frac{p^{*}+1}{2}$ so that \[ \hat{f}(v|\boldsymbol{x}_{1}^{\tilde{p}})-f(v|\boldsymbol{x}_{1}^{p^{*}})=o_{P}\left(n^{-1/4}\right) \] as shown in Theorem (ref). Furthermore, we adopt a $2r$-th order kernel for $\hat{f}(v|\boldsymbol{x}_{1}^{\tilde{p}})$ to ensure that the bias term is of order $o\left(n^{-1/2}\right)$, which is required for the $\sqrt{n}$-convergence of the parameters, as explained in Lewbel2000.
The above theorem presents the “oracle” property, implying the “super efficiency” that the procedure can effectively eliminate irrelevant regressors with probability approaching one, and the resulting estimator has the same asymptotic properties as one obtained using only the true signals. We can obtain this result thanks to the “beta-min” condition. A relaxation of this condition and the related error bounds for the Lasso estimator can be found in Appendix (ref).
We recommend that practitioners select $\lambda_n$ by 10-fold CV and experiment with a range of $\gamma$ values, choosing the value at which the results stabilize. For further guidance and illustrating examples, we refer readers to Sections (ref) and (ref).
A final note is that we use the U-statistics technique to derive the asymptotic property of our estimator $\hat{\boldsymbol{\beta}}$. As an anonymous referee pointed out, an interesting alternative way is to find the Neyman-orthogonal score with respect to $\hat{f}(v|x)$ and apply empirical processes techniques (e.g., Andrews1994) afterwards. We note that finding the Neyman-orthogonal score in our case is challenging and requires considerable effort. Unlike in ChernozhukovEtal2022 and many related papers, where the first-step estimation involves $X$, our first step is to obtain the dependent variable $\widehat{\tilde{Y}}$ for the final estimation. This distinction inhibits us from directly applying the results in existing literature. To apply empirical processes techniques, we need to verify that $\hat{f}(v|x)$ belongs to a Donsker class. Note that $\hat{f}(v|x)$ involves bandwidths as tuning parameters. Typically, different technical conditions are required to handle $\hat{f}(v|x)$ in various scenarios. For example, see MammenEtal2011, Escanciano2014, and SeoOtsu2018. We leave this direction for future work.
We recommend two options for the initial estimator: Ordinary Least Squares (OLS) and the Lasso estimator, with the latter detailed in Appendix (ref). If the number of covariates is considerably smaller than the number of observations and the Gram matrix of covariates ($n^{-1}\mathbf{X}_{n}'\mathbf{X}_{n}$) is full rank, we recommend using the OLS estimator for its simplicity. However, if the number of covariates is large or the Gram matrix is nearly singular, we suggest the Lasso estimator, as outlined in Appendix (ref).
In this section, we cite a result from LinEtal2009, which demonstrates that the OLS estimator can yield consistent estimates that may serve as initial estimators. For details on the Lasso estimator, we refer readers to Appendix (ref). Let $\widehat{\mathbf{\tilde{Y}}} \equiv (\widehat{\tilde{y}}_{1}, \ldots, \widehat{\tilde{y}}_{n})'$. By definition, the OLS estimator is given by \[ \tilde{\boldsymbol{\beta}} = \left(\mathbf{X}_{n}'\mathbf{X}_{n}\right)^{-1}\mathbf{X}_{n}'\widehat{\mathbf{\tilde{Y}}}. \] The following result is Proposition 2.1 of LinEtal2009, if we adopt $\mathbf{\tilde{Y}}$ in $\tilde{\boldsymbol{\beta}}$. The generalization of using $\widehat{\mathbf{\tilde{Y}}}$ instead of $\mathbf{\tilde{Y}}$ is straightforward given the proof of Theorem (ref). We omit the proof for brevity.
Clearly, if we can find a $\varrho>0$ such that \[ \rho_{n}^{-1}n^{-\frac{1-C_{p}}{2}}\apprle n^{-\varrho}, \] then $\tilde{\boldsymbol{\beta}}$ is a valid candidate for the initial estimator.
Assuming certain moment conditions hold, Liao2013 examined the selection of valid moment conditions from a fixed set of candidates for GMM estimation. ChengLiao2015 expanded this approach by selecting moment conditions that are both valid and relevant while allowing for an increasing number of potential IVs. For simplicity, we focus solely on validity in this section, assuming all IVs are relevant. While it is possible to consider both validity and relevance, doing so introduces additional complexity. In what follows, we adopt the approach of ChengLiao2015 (ChengLiao2015), allowing for an increasing number of potential IVs, and aim to replicate their results under similar technical conditions.
The model remains unchanged: \[ Y=\boldsymbol{1}\left(V+X_{1}\beta_{1}^{\ast}+X_{2}\beta_{2}^{\ast}+...+X_{s^{\ast}}\beta_{s^{\ast}}^{\ast}+\varepsilon>0\right). \] Following Liao2013 and ChengLiao2015, we assume that we have prior knowledge of the first $k^{\ast}$ ($k^{\ast}\geq s^{\ast}$) IVs, denoted as $\boldsymbol{Z}^{\ast}\equiv\left(Z_{1},...,Z_{k^{\ast}}\right)^{\prime}$, all of which are valid for identifying the parameters $\boldsymbol{\beta}^{\ast}\equiv\left(\beta_{1}^{\ast},\beta_{2}^{\ast},...,\beta_{s^{\ast}}^{\ast}\right)^{\prime}$.\footnote{Note that $\boldsymbol{Z}^{\ast}$ (known valid IVs) may or may not overlap with $\boldsymbol{Z}_{1}^{p^{\ast}}$ ($Z$s relevant to $V$) defined in Section (ref).} We assume $k^{\ast}$ is fixed to simplify analysis. There are other valid IVs, denoted by $\underline{\boldsymbol{Z}}_{A}=(Z_{k^{*}+1},Z_{k^{*}+2},...,Z_{k^{*}+d_{A}})^{\prime}$, which are mixed with invalid IVs, denoted by \[ \underline{\boldsymbol{Z}}_{B}=(Z_{k^{*}+d_{A}+1},Z_{k^{*}+d_{A}+1},...,Z_{k^{*}+d_{A}+d_{B}})^{\prime}. \] We do not know which IVs are valid or invalid. We use $A$ to denote the indices of the valid instruments ($\underline{\boldsymbol{Z}}_{A}$), $B$ to denote the indices of the invalid instruments ($\underline{\boldsymbol{Z}}_{B}$), and $D=A\cup B$ to denote the indices for all candidate IVs to be selected. The objective is to choose the valid IVs from $\underline{\boldsymbol{Z}}_{D}=(\underline{\boldsymbol{Z}}_{A}^{\prime},\underline{\boldsymbol{Z}}_{B}^{\prime})'$. We refer to Liao2013 and ChengLiao2015 for more information and applications related to this setting. We denote \[ p_{n}\equiv k^{\ast}+d_{A}+d_{B}. \] With this notation, $\boldsymbol{Z}_{1}^{p_{n}}=\left(\boldsymbol{Z}^{\ast\prime},\underline{\boldsymbol{Z}}_{A}^{\prime},\underline{\boldsymbol{Z}}_{B}^{\prime}\right)^{\prime}$. To further simplify notations, we denote $\boldsymbol{X}\equiv\left(X_{1},...,X_{s^{\ast}}\right)^{\prime}$ and $\boldsymbol{\beta}\equiv\left(\beta_{1},...,\beta_{s^{\ast}}\right)^{\prime}$. For the same reason as in the last section, the two issues (namely, “valid IV” and “relevant $Z$s to $V$”) are generally independent. We do not reorder $Z$ to distinguish these two issues in this section for ease of notation.
To implement the moment selection, we introduce the auxiliary parameter $\boldsymbol{\eta}$ and its true value $\boldsymbol{\eta}$$^{\ast}$ as \[ \mathbb{E}[\underline{\boldsymbol{Z}}_{D}(\tilde{Y}-\boldsymbol{X}^{\prime}\boldsymbol{\beta})]=\underline{\boldsymbol{\eta}}\text{ and }\mathbb{E}[\underline{\boldsymbol{Z}}_{D}(\tilde{Y}-\boldsymbol{X}^{\prime}\boldsymbol{\beta}^{\ast})]=\underline{\boldsymbol{\eta}}^{\ast} \] with $\tilde{Y}$ defined in ((ref)). By definition, $\eta$$_{j}^{\ast}=0$ if $j\in A,$\ and $\eta$$_{j}^{\ast}\neq0$ if $j\in B$. We denote the parameters to be estimated as $\boldsymbol{\theta}\equiv\left(\boldsymbol{\beta}^{\prime},\underline{\boldsymbol{\eta}}^{\prime}\right)^{\prime}$. The moment conditions can now be expressed as \[ \boldsymbol{m}\left(\boldsymbol{\theta}\right)=\left[
\right]. \] The feasible sample analog is given by \[ \overline{\boldsymbol{\hat{m}}}_{n}\left(\boldsymbol{\theta}\right)=\left[
\right], \] where $\underline{\boldsymbol{z}}_{Di}$ is a $\left(d_{A}+d_{B}\right)\times1$ vector representing the realization of $\underline{\boldsymbol{Z}}_{D}$ for the observation $i$, $\boldsymbol{z}_{i}^{\ast}$ is a $k^{\ast}\times1$ vector, defined similarly to $\underline{\boldsymbol{z}}_{Di}$, and \[ \widehat{\tilde{y}}_{i}=\frac{y_{i}-\boldsymbol{1}\left(v_{i}>0\right)}{\hat{f}(v_{i}|\boldsymbol{z}_{1i}^{\tilde{p}})} \] with $\hat{f}(v_{i}|\boldsymbol{z}_{1i}^{\tilde{p}})$ being obtained using the procedure outlined in Section (ref), but with a $2r$-th order kernel function (for the same reason as listed in the last section).
We perform the selection and estimation using weighted $L_{1}$ penalty proposed as in ChengLiao2015. Alternative methods employing different penalty functions can be handled similarly. Specifically, we define a penalized GMM estimator for $\boldsymbol{\theta}$ as
where $\mathbf{W}_{n}$ is a $p_{n}\times p_{n}$ positive definite weighting matrix and \[ \varpi_{j}= |\underline{\tilde{\eta}}_{j} |{}^{-\gamma},\quad j=1,...,d_{A}+d_{B}, \] where $\gamma>0,$ $\underline{\tilde{\eta}}_{j}$ is some initial estimator, and one candidate is from ((ref)). We impose the following technical conditions, which are essentially the same as those placed in ChengLiao2015, including restrictions on tuning parameters.
Parts (1), (4), and (5) of Assumption (ref) are standard. We place the classic rank condition in Assumption (ref)(2), which implicitly assumes $k^{\ast}\geq s^{\ast}$ since $\mathbb{E}\left(\boldsymbol{Z}^{\ast}\boldsymbol{X}^{\prime}\right)$ is a $k^{\ast}\times s^{\ast}$ matrix. The rest of Assumption (ref)(2) defines valid and invalid IVs. In Assumption (ref)(3), $\Omega_{n}$ is the variance of the moment conditions. Note that the dimension of $\boldsymbol{\theta}_{B}^{\ast}$ (defined below in ((ref))) might be diverging, so we consider a linear combination of $\boldsymbol{\hat{\theta}}_{B}\boldsymbol{-\theta}_{B}^{\ast}$ for the asymptotics, as did in ChengLiao2015. We require $\Omega_{n}$ to be of full rank because we need $\Omega_{n}^{-1}$ to construct the weight of this linear combination. In Assumption (ref)(6), we assume $k^{\ast}$ is fixed and we need $\min_{j\in B}\{|\underline{\eta}_{j}^{\ast}|\}$ to be big enough so that our procedure can effectively distinguish valid IVs from invalid ones.
Assumption (ref)(7) places restrictions on the rate of initial estimator $\underline{\tilde{\eta}}_{j}$ that can be obtained as
which is the GMM estimator without the penalty term. ChengLiao2015 has shown that $\max_{j=1,...,d_{A}+d_{B}} |\underline{\tilde{\eta}}_{j}-\underline{\eta}_{j}^{*} |=O_{P} (\sqrt{\left.p_{n}\right/n} )$ (see their discussion after Assumption 5.1). $p_{n}\ll n^{1/3}/\log n$ is needed to verify the Lindeberg condition for asymptotic normality. Note that ChengLiao2015 essentially assumes that $a_{n}=C>0$. When $a_{n}=C>0$, we can, for example, take $\lambda_{n}=1/\sqrt{np_{n}\log n}$ and any $\gamma\geq1,$ due to $p_{n}\ll n^{1/3}/\log n.$
Since $\underline{\boldsymbol{\hat{\eta}}}_{A}=\boldsymbol{0}$ will be shown to occur with high probability, regarding the asymptotic distribution, the parameters of our interest are
We denote the corresponding estimator as $\boldsymbol{\hat{\theta}}_{B}\equiv(\boldsymbol{\hat{\beta}}^{\prime},\underline{\boldsymbol{\hat{\eta}}}_{B}^{\prime})^{\prime}$ and the true value as $\boldsymbol{\theta}_{B}^{\ast}\equiv(\boldsymbol{\beta}^{\ast\prime},\underline{\boldsymbol{\eta}}_{B}^{\ast\prime})^{\prime}$. The partial derivative of the moment conditions with respect to $\boldsymbol{\theta}_{B}\equiv(\boldsymbol{\beta}^{\prime},\underline{\boldsymbol{\eta}}_{B}^{\prime})^{\prime}$ is denoted as
$\Gamma_{\boldsymbol{\theta}_{B}}^{\prime}\Gamma_{\boldsymbol{\theta}_{B}}$ is of rank $s^{\ast}+d_{B},$ due to the assumption that $\mathbb{E}\left(\boldsymbol{Z}^{\ast}\boldsymbol{X}^{\prime}\right)$ is of rank $s^{\ast}.$ We define the following matrix: \[ \Sigma_{n}\equiv\left(\Gamma_{\boldsymbol{\theta}_{B}}^{\prime}\mathbf{W}_{n}\Gamma_{\boldsymbol{\theta}_{B}}\right)^{-1}\Gamma_{\boldsymbol{\theta}_{B}}^{\prime}\mathbf{W}_{n}\Omega_{n}\mathbf{W}_{n}\Gamma_{\boldsymbol{\theta}_{B}}\left(\Gamma_{\boldsymbol{\theta}_{B}}^{\prime}\mathbf{W}_{n}\Gamma_{\boldsymbol{\theta}_{B}}\right)^{-1}, \] which is useful for the asymptotics.
We reproduce the results in Liao2013 and ChengLiao2015 in the following theorem.
For the same reason as for Theorem (ref), we impose restrictions on the adopted kernel function, as specified in the theorem. As in the last section, we recommend to choose $\lambda_n$ by 10-fold CV and try a range of $\gamma$. For practical examples and further guidance, see Sections (ref) and (ref).
Finally, an interesting direction would be to permit $s^{*}$ to diverge while simultaneously considering variable selection and moment selection. We leave this topic for future research.
This section assesses the finite sample performance of the procedures proposed in Section (ref) through Monte Carlo experiments. For all designs studied in this section, we consider sample sizes $n=500,1000, 2000$ and base our findings on 500 independent replications conducted using the R programming language and MATLAB. All simulation results are reported in tables collected in Appendix (ref).
We first investigate the adaptive Lasso method as discussed in Section (ref). We explore three distinct simulation designs based on the following data-generating process: \[ Y=\mathbf{1}(V+\beta_{0}+\beta_{1}X_{1}+\cdots+\beta_{p}X_{p_n}+\varepsilon>0). \] Here, we consider $p_n\in\{14,29,49\}$, corresponding to binary choice models with $p_{n}=15,30$, and 50 covariates, respectively. The true parameters are set as follows: $\beta_{0}=\beta_{1}=0$, $\beta_{2}=\beta_{3}=1$, and $\beta_{l}=0$ for all $l=4,\ldots,p_n$. $V$ is defined as $V=X_{1}+X_{1}X_{2}+\mathbf{1}(X_{2}>0)+e_{v}$, where $e_{v}\sim\text{Logistic}(0,2)$. The three designs differ in the distributions of $\varepsilon$ and $\boldsymbol{X}_{1}^{p_n}=(X_{1},\ldots,X_{p_n})'$:
Design 1 is a benchmark design, where the Probit model is correctly specified up to a scale of $\pi/\sqrt{3}$. We anticipate that both the standard Probit and our adaptive Lasso estimators are consistent, with the former being more efficient. Design 2 is a heteroskedastic Probit model, where $\varepsilon$ has the same standard deviation as Design 1. We expect our adaptive Lasso regression to give consistent estimates in this design, but the Probit estimator will be biased. Design 3 examines the effectiveness of our procedure in handling the more general dependence structure as in ((ref)), where $V$ and $X_{3},...,X_{p_n}$ are dependent but only through their dependence on $(X_{1},X_{2})$. It is worth noting that by construction, $V$ has a stronger dependence on $X_{1}$ than on $X_{2}$. For instance, in Designs 1 and 2, $\mathrm{Corr}(V,X_{1})=0.25$, while $\mathrm{Corr}(V,X_{2})=0.1$.
For each simulation design, we estimate the parameters $\boldsymbol{\beta}\equiv(\beta_{0},\beta_{1},\ldots,\beta_{p_n})'$ using two methods: adaptive Lasso regression and Probit regression. The former is implemented with the following procedure:
Probit results are obtained using R's built-in function $\mathtt{glm}()$. Note that $\mathtt{glm}()$ assumes a standard deviation of 1 for $\varepsilon$ (default scale normalization) and estimates the coefficient on $V$, denoted as $\rho$.\footnote{As a result, $\mathtt{glm}()$ fits the following Probit model: \[ Y=\mathbf{1}(\rho V+\beta_{0}^{Probit}+\beta_{1}^{Probit}X_{1}+\cdots+\beta_{p}^{Probit}X_{p_n}+e>0), \] where $e=\sqrt{3}\varepsilon/\pi\sim N(0,1)$, $\rho=\sqrt{3}/\pi$, and $\beta^{Probit}=(\beta_{0}^{Probit},\beta_{1}^{Probit},\cdots,\beta_{p_n}^{Probit})'=\sqrt{3}\beta/\pi$.} To facilitate comparison, we compute the estimates of ratios $\beta_{l}/\rho$ for the Probit estimation.
The simulation results are displayed in tables presented in Appendix G.1, beginning with “Table 1" for Design 1, “Table 2" for Design 2, and so on. For each design, tables labeled “A", “B", and “C" correspond to the variable selection results in Steps 1 and 4 of our proposed procedure, the outcomes of the adaptive Lasso regression, and the Probit estimation, respectively. For example, in Table (ref), we report $\Pr(\{1\} \in \tilde{\mathcal{A}})$ and $\Pr(\{2\} \in \tilde{\mathcal{A}})$, representing the probabilities that $X_{1}$ and $X_{2}$ are selected as variables relevant to $V$. The column labeled “Correct (%)" in Table (ref) displays the average number (percentage) of the true zero coefficients correctly identified as zero, and the column labeled “Incorrect" shows the average count of the two true nonzero coefficients incorrectly set to zero. Table (ref) summarizes the adaptive Lasso estimator's performance for Design 1, presenting mean bias (MEANB), root mean squared errors (RMSE), median bias (MEDB), and median absolute deviation (MAD) for $(\hat{\beta}_{2}, \hat{\beta}_{3})$. Table (ref) contains the same set of statistics for the Probit estimator of $(\beta_{2}/\rho,\beta_{3}/\rho)$ as in Table (ref).
The variable selection outcomes are encouraging, as presented in Tables (ref)--(ref). The screening selects $(X_{1},X_{2})$ as relevant variables for $V$ with probability approaching one. Notably, for $X_{1}$, which has moderate dependence with $V$, the $\Pr(\{1\}\in\tilde{\mathcal{A}})$ attains 1 even when $n\leq 500$ across all three designs. For $X_{2}$, which has weak dependence on $V$, the $\Pr(\{2\}\in\tilde{\mathcal{A}})$ approaches 1 when $n\geq1000$ and is more sensitive to the number of candidate variables ($p$) in smaller sample sizes.\footnote{It is worth highlighting that the screening performs slightly better in Design 3 compared to Designs 1 and 2. This improvement can be attributed to the enhanced correlation between $X_{1}$ and $X_{2}$ in Design 3, strengthening their association with $V$. In fact, in Design 3, $\text{Corr}(V,X_{1})\approx 0.30$ and $\text{Corr}(V,X_{2})\approx 0.21$.} Additionally, the adaptive Lasso regression can automatically exclude irrelevant regressors (e.g., $X_{1}$) from the regression with increasing probability as the sample size increases, which is evidenced in the “Correct (%)" columns of Tables (ref)--(ref).
As expected, when the Probit model is correctly specified (up to scale), the Probit estimator is consistent and performs better than our adaptive Lasso estimator, as shown in Tables (ref) and (ref). However, in the presence of heteroskedasticity, observed in Tables (ref) and (ref), the Probit estimator maintains a bias of approximately 10%–20% for $\beta_{3}/\rho$ that persists as the sample size increases. Our adaptive Lasso estimator, which uses an inverse density as a plug-in term, is prone to generating extreme estimates, resulting in comparatively larger RMSEs, especially in smaller sample sizes. However, focusing on the MAD, which is robust to outliers, we observe that our adaptive Lasso estimator retains its $\sqrt{n}$-consistency across all three designs, aligning with our asymptotic theory, as illustrated in Tables (ref)–(ref). Finally, although the adaptive Lasso variable selection is influenced by the choice of the tuning parameter $\gamma$, the resulting estimates appear to be pretty robust to this choice.
In this section, we examine the penalized GMM approach proposed in Section (ref) through three simulation designs (Designs 4--6) based on the following binary choice model: \[ Y=\boldsymbol{1}(V+\beta_{0}+\beta_{1}X+\varepsilon>0), \] where the true parameters $\boldsymbol{\beta}=(\beta_{0},\beta_{1})'=(0,1)$ and $X$ is an endogenous regressor. To estimate $\boldsymbol{\beta}$, we use instrumental variables $\boldsymbol{Z}^{*}=(1,Z^{*})'$, $\underline{\boldsymbol{Z}}_{A}=(Z_{1},...,Z_{d_{A}})'$, and $\underline{\boldsymbol{Z}}_{B}=(Z_{d_{A}+1},...,Z_{d_{A}+d_{B}})'$ with $(d_{A},d_{B})\in\{(6,7),(14,14),(24,24)\}$, corresponding to cases with $p_{n}=15,30,$ and 50, respectively. Designs 4--6 differ in how $V$ and $\boldsymbol{Z}_{1}^{p_{n}}=(\boldsymbol{Z}^{*\prime},\underline{\boldsymbol{Z}}_{A}^{\prime},\underline{\boldsymbol{Z}}_{B}^{\prime})'$ are dependent and how variables in $\underline{\boldsymbol{Z}}_{B}$ are correlated with $\varepsilon$. Specifically, letting $\boldsymbol{e}=(e_{1},...,e_{4})'\sim\text{MVN}(\boldsymbol{0},\boldsymbol{I}_{4})$, $\boldsymbol{u}=(u_{1},...,u_{d_{A}+d_{B}+1})'\sim\text{MVN}(\boldsymbol{0},\boldsymbol{I}_{d_{A}+d_{B}+1})$, and $\boldsymbol{e}\perp\boldsymbol{u}$, we set
Design 4 is a benchmark design where $V$ is solely dependent on $(Z^{*},Z_{1})$ with different levels among all candidate IVs. In this design, it is evident from the data-generating process that all candidate IVs are relevant to the endogenous regressor $X$. Note that $Z^{*}$ and $\underline{\boldsymbol{Z}}_{A}$ are valid IVs for $X$, while $\underline{\boldsymbol{Z}}_{B}$ are invalid IVs. In Design 5, we modify Design 4 by replacing $e_{3}$ with $e_{2}$ to examine the effectiveness of the DC screening procedure under a more general dependence structure, as in ((ref)). Consequently, $V$ becomes dependent on $Z_{2},...,Z_{d_{A}}$ as well, but only through their correlation with $(Z^{*},Z_{1})$, so that $V\perp(Z_{2},...,Z_{d_{A}})|(Z^{*},Z_{1})$ and this conditional independence is not ambiguous.\footnote{This conditional independence doesn't alter the identification strength of $\boldsymbol{Z}_{1}^{p_{n}}$.} Design 6 aims to strengthen the correlation between $\underline{\boldsymbol{Z}}_{B}$ and $\varepsilon$ to assess the sensitivity of the penalized GMM method in identifying invalid IVs.
To implement the penalized GMM method, we first obtain a feasible estimate $\hat{f}(v|z_{1}^{\tilde{p}})$ for $f(v|z_{1}^{p_{n}})$. This follows the same Steps 1–3 as in Designs 1–3, with $\tilde{p}=4$, but substituting $X$s with $Z$s in this context. After computing the feasible $\hat{f}(v|z_{1}^{\tilde{p}})$ estimates, we proceed to obtain $\hat{\boldsymbol{\theta}}$ from ((ref)). We solve the minimization problem in the penalized GMM estimation using the optimal projected gradient algorithm proposed by schmidt2010graphical, setting the tuning parameters as $\lambda_{n} = 1/\sqrt{n p_n}$ and $\gamma\in\{1,2,3\}$. Throughout all designs, we use the identity matrix as the weighting matrix, i.e., \(\mathbf{W}_{n} = \boldsymbol{I}_{p_{n}}\). Our simulation studies reveal that the simple identity matrix often outperforms the theoretically more efficient optimal GMM weighting matrix. An intuitive explanation for this observation is that, in finite samples, the efficiency gains from using the optimal weighting matrix are often insufficient to outweigh the errors introduced by its estimation, particularly given that our estimator incorporates a nonparametrically estimated plug-in term. Thus, we recommend practitioners use the identity matrix as the default weighting matrix.
We present the simulation results for Designs 4–6 in Tables 4–6, respectively. The outcomes for variable screening and moment selection are shown in tables labeled “A", and the performance of the penalized GMM estimation is illustrated in tables labeled “B". Specifically, let $\tilde{\mathcal{Z}}$ denote the set of relevant $Z$s selected by DC screening. Table (ref) reports $\Pr(Z^{*} \in \tilde{\mathcal{Z}})$ and $\Pr(Z_{1} \in \tilde{\mathcal{Z}})$. Additionally, the “Correct (%)" column in Table (ref) presents the average count (percentage) of truly valid IVs (other than $Z^{*}$) that are successfully selected, and the “Incorrect (%)" column reports the average number (percentage) of invalid IVs mistakenly categorized as valid. Table (ref)--(ref) provides the same set of performance metrics as Tables (ref)–(ref) for $(\hat{\beta}_{0}, \hat{\beta}_{1})$ in Designs 4--6.
As shown in Tables (ref)--(ref), the DC screening procedure consistently performs well. When $n\geq500$, it demonstrates an almost certain ability to select all variables relevant to $V$. For IV (moment) selection, the penalized GMM procedure effectively distinguishes valid IVs from invalid ones with increasing accuracy as the sample size grows, as indicated by the “Correct (%)" and “Incorrect (%)" columns. The choice of tuning parameter $\gamma$ presents a trade-off: for a given $\lambda_n$, a larger (smaller) $\gamma$ improves (reduces) the likelihood of correctly identifying valid IVs but also raises (lowers) the risk of mistakenly classifying some invalid IVs as valid. We recommend that practitioners try a range of $\gamma$ values and carefully examine their IV selection outcomes to avoid inconsistent estimates.
Tables (ref)--(ref) present the performance of the penalized GMM estimator. The estimator shows noticeable bias in small samples, but as the sample size increases, this bias quickly diminishes, consistent with our asymptotic theory. Interestingly, the RMSE of the penalized GMM estimator does not consistently decrease with sample size. This may be due to the use of the estimated inverse density function as a plug-in term; when the estimated density approaches zero, “outliers” can appear in $\widehat{\tilde{Y}}$, leading to estimates that deviate significantly from the true values. However, focusing on the MAD, which is less sensitive to extreme values, we observe that the estimator converges at approximately a parametric rate, aligning with our Theorem (ref). Finally, as the result of stronger IVs, the penalized GMM estimator performs better in Design 6 compared to Designs 4 and 5.
China's economic reforms since the late 1970s have significantly transformed the country’s economy. In the early 1980s, the government began to relax restrictions on population mobility. Over time, rural residents were granted more freedom to leave their villages and seek employment in larger cities for higher wages. In the following analysis, we explore factors influencing the migration intentions of rural residents by applying our method to a dataset drawn from the Rural-Urban Migration in China (RUMiC) project.\footnote{RUMiC consists of three components: the Urban Household Survey, the Rural Household Survey, and the Migrant Household Survey. It was initiated by researchers from the Australian National University, the University of Queensland, and Beijing Normal University, with support from the Institute for the Study of Labor (IZA). RUMiC received funding from the Australian Research Council, the Australian Agency for International Development, the Ford Foundation, IZA, and the Chinese Foundation of Social Sciences. Further details on the survey are available at \url{https://datasets.iza.org/dataset/58/longitudinal-survey-on-rural-urban-migration-in-china}.} The survey focuses on individuals moving from rural areas to major Chinese cities. Participants (both migrants and workers) answered a wide range of questions. For more detailed information on the survey design and variable construction, refer to the survey website and Gongetal2008. Currently, data from the 2008 wave are publicly available.
Previous studies have investigated the factors influencing rural residents' migration intentions using various methodologies, yielding mixed results. Zhao1999, Zhao2003, and Mullan2011 applied logistic regression models to cross-sectional data, while XueYangZhou utilized the special regressor approach within a panel data framework that included interactive fixed effects.
Thanks to the distribution-free nature of our method, we relax the assumption that the error term follows a logistic distribution, which is imposed in Zhao1999, Zhao2003, and Mullan2011. Additionally, we explicitly examine the dependence of the special regressor on other included covariates using our screening procedure, whereas XueYangZhou relied on intuition for such assumption. However, unlike their work, we focus on cross-sectional data, as our methods are not designed for panel data settings.
Departing from aforementioned works that considered only a limited set of explanatory variables--thereby risking omitted variable bias--we significantly broaden the scope of variables included in our analysis and avoid specifying the regression model based solely on the researcher's subjective judgment. Instead, we employ our proposed data-driven methods to identify the key factors influencing migration intentions. However, a limitation of this approach is that the final estimation results may lack clear causal interpretation. We emphasize that the primary purpose of this application is to illustrate the use of our methodology and provide insights for future research that can adopt more focused and rigorously designed studies.
The migration intention is modeled as follows:
where \( Y \) represents the binary migration intention. We use data from the 2008 wave of the RUMiC survey. Following XueYangZhou, we employ the negative logarithm of the average daily precipitation between April and August from two and three years prior as the special regressor, $V$. The model includes 51 explanatory variables (i.e., \( p = 51 \)). By applying the procedure proposed in Section (ref), we aim to identify predictors relevant to migration decisions. After excluding observations with missing values, the final sample consists of 3,787 observations. Tables (ref) and (ref) in Appendix G.2 provide the definitions and summary statistics for these variables, respectively.
We standardize all variables, including $V$, to have a mean of 0 and variance of 1 before estimation. We apply the adaptive Lasso estimator proposed in Section (ref), specifically the one defined in equation ((ref)). The initial estimator is obtained using OLS, as outlined in Section (ref). Data-driven bandwidths are determined according to the procedure in Section (ref).
For the screening procedure, we set $\tilde{p} = 4$, meaning we retain the four covariates with the highest values of $\mathcal{\hat{T}}_{n}$. After applying the CV method from HallRacineLi on these covariates, we found that the bandwidth of one covariate reached its upper limit, indicating that it has no impact on the conditional density. This suggests that $\tilde{p} = 4$ is sufficiently large for this application. The remained covariates are “Height”, “Income”, and “Age”.
Next, we discuss the tuning parameters for the adaptive Lasso. For a given \(\gamma\) in the weight \(\varpi_j\), the parameter \(\lambda_n\) in equation ((ref)) is selected via 10-fold CV. After experimenting with \(\gamma\) values ranging from 1 to 15, we found that setting \(\gamma\) between 2 and 12 resulted in selecting the same five covariates. Based on this, we conduct the variable selection by setting \(\gamma\) within this range.
As the last step, we perform post-selection OLS and make inferences using the asymptotics from Theorem (ref), constructing confidence intervals by estimating the sample counterpart of \(\Sigma_n\). In \(\Sigma_n\), \(\Omega_n\) represents the variance of the influence term, the same as in Lewbel2000 for a fixed dimension setting. A plug-in estimate of \(\Omega_n\) can be constructed as follows. First, estimate \(\hat{u}_i = \widehat{\tilde{y}}_i - \boldsymbol{x}_i' \boldsymbol{\hat{\beta}}\) from the post-selection OLS. Then, compute \[ \hat{\boldsymbol{q}}_i = \hat{u}_i \boldsymbol{x}_i + \widehat{\mathbb{E} \left( \hat{u}_i \boldsymbol{x}_i \,\middle|\, \boldsymbol{x}_{1i}^{p^{*}} \right)} - \widehat{\mathbb{E} \left( \hat{u}_i \boldsymbol{x}_i \,\middle|\, \boldsymbol{x}_{1i}^{p^{*}}, v_i \right)}, \] where \(\widehat{\mathbb{E}(\cdot|\cdot)}\) denotes standard nonparametric kernel regression. In this case, we include “Height", “Income", and “Age" as \(\boldsymbol{x}_{1i}^{p^{*}}\). Finally, compute \[ \hat{\Omega}_n = \frac{1}{n}\sum_{i=1}^{n} \hat{\boldsymbol{q}}_i \hat{\boldsymbol{q}}'_i. \] With \(\hat{\Omega}_n\), \(\hat{\Sigma}_n\) can be straightforwardly obtained. Note \(\hat{\Omega}_n\) can be similarly calculated for the scenario in Section (ref).
To ensure a fair and meaningful comparison, we apply the adaptive Logistic Lasso, using the initial estimator from logistic regression. This estimation is implemented using the standard R package glmnet. The tuning parameter $\lambda_n$ is also selected via 10-fold CV, with a fixed value of $\gamma$ for the weight $\varpi_j$. However, we find that the selection results are sensitive to the choice of $\gamma$. Specifically, when $\gamma \leq 3$, more than 19 covariates are selected, for $\gamma = 4$, six covariates are selected, whereas with $\gamma = 5$, only one covariate is selected. Based on insights from our method, we set $\gamma = 4$, as both methods (our proposed approach and adaptive Logistic Lasso) select a comparable number of covariates. The final estimates are then obtained by applying post-selection logistic regression.
We summarize the key findings as follows:
First, our results indicate that $V$ (the average precipitation) is not independent of all the covariates considered. Using the CV method by HallRacineLi, we find that $V$ is dependent on “Height", “Income", and “Age". This dependence is reasonable given that northern China is typically drier but less economically developed than southern China, and individuals from the north tend to be taller on average, explaining the association between $V$ and both “Height" and “Income". However, we do not have a clear explanation for the relation between $V$ and “Age". These findings emphasize the importance of the screening procedure, as a model-free data analysis can reveal significant (possibly nonlinear) relationships that may otherwise be overlooked. Second, $V$ is consistently selected by the adaptive Logistic Lasso procedure and is highly significant with the expected sign in the post-selection logistic regression. This supports the decision to include $V$ in the model and justifies the normalization of its coefficient to 1 in the special regressor approach.
We report the post-selection estimation results for both methods in Table (ref). Standard errors are computed based solely on the post-selection estimates. Both methods select several common regressors, including a10 (“Height”), “Age”, and a23 (“Medical Expenses Reimbursed”). Post-selection results from both approaches suggest that a23 is insignificant. However, there are some discrepancies. For example, the signs of the coefficients for “Height” differ between the two methods. Our approach uniquely selects nold (“No. of Elderly”), whereas the adaptive Logistic Lasso selects deduc5, e35, and \texttt{e36}. These discrepancies may arise from biased estimates in the Logistic estimation, possibly due to model mis-specifications or heteroskedasticity in the error term.
We conclude this section by discussing the implications of our empirical results. The positive effect of nold on migration intentions is likely driven by financial pressures, as elderly individuals in rural China generally lack access to pensions, pushing younger family members to migrate for better financial support. The negative effect of “Age" can be explained by younger individuals being less settled and more willing to take risks, making them more inclined to migrate. The negative coefficient for c18_2 (“Monthly Bonus and Allowance") suggests that better working conditions discourage migration, as individuals receiving higher compensation may be less motivated to leave their current jobs. However, we do not have a clear causal interpretation of the effect of “Height". As previously noted, one challenge with including many covariates in the model is that the results may lack clear causal interpretations.
In this paper, we introduce a new estimation procedure for semiparametric binary choice models in high-dimensional settings. This is achieved by combining an innovative dimension reduction method for conditional density estimation with the special regressor approach. We study the classic variable and moment selection problems in the binary choice model context. Monte Carlo simulations illustrate the finite sample properties of our methods. We illustrate our proposed approaches by studying migration intentions of rural residents in China. Our empirical findings indicate that the special regressor used in XueYangZhou (“Precipitation”) is related to several covariates and therefore may not satisfy the strong exogeneity assumption. This underscores the value of our data-driven approach in guiding practitioners through empirical model specification.
For future research, one promising avenue is to integrate our novel conditional density estimator with other estimators that engage with high-dimensional conditional density estimates.
We thank the editor, Xiaohong Chen, an associate editor, and three anonymous referees for their helpful comments, which have substantially improved the paper. We are grateful for the valuable feedback and discussions provided by seminar participants at the University of Melbourne, the University of Queensland, and the University of Sydney, as well as conference attendees at the 17th International Symposium on Econometric Theory and Applications (SETA 2023) and the 31st Australia New Zealand Econometric Study Group Meeting (ANZESG 2023). Any remaining errors are our responsibility.