EconBase
← Back to paper

Post-Selection Inference in Three-Dimensional Panel Data

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.

60,098 characters · 12 sections · 56 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Post-Selection Inference in Three-Dimensional Panel Data

abstract{5.8mm} Three-dimensional panel models are widely used in empirical analysis. Researchers use various combinations of fixed effects for three-dimensional panels. When one imposes a parsimonious model and the true model is rich, then it incurs mis-specification biases. When one employs a rich model and the true model is parsimonious, then it incurs larger standard errors than necessary. It is therefore useful for researchers to know correct models. In this light, LuMiaoSu2018 propose methods of model selection. We advance this literature by proposing a method of post-selection inference for regression parameters. Despite our use of the lasso technique as means of model selection, our assumptions allow for many and even all fixed effects to be nonzero. Simulation studies demonstrate that the proposed method is less biased than under-fitting fixed effect estimators, is more efficient than over-fitting fixed effect estimators, and allows for as accurate inference as the oracle estimator. \\ {\bf Keywords:} post-selection inference, three-dimensional panel data. \\ {\bf JEL Code:} C23

Introduction

Matyas1997 suggests the three-dimensional panel model

align[align omitted — 146 chars of source]

for $(i,j,t) \in \{1,...,N\} \times \{1,...,M\} \times \{0,...,T\}$, where $y_{ijt}$ denotes an outcome variable of unit $(i,j)$ at time $t$, ${x}_{ijt}$ denotes $k$-dimensional explanatory variables of unit $(i,j)$ at time $t$, and $\alpha_i$, $\gamma_j$, and $\lambda_t$ are fixed effects associated with indices $i$, $j$, and $t$, respectively. To fix our ideas, consider the gravity model Tinbergen1962 from the empirical trade literature where $y_{ijt}$ denotes the logarithm of the volume of exports from country $i$ to country $j$ in year $t$, and the $k$-dimensional covariates ${x}_{ijt}$ contain observed characteristics of the trade pair $(i,j)$ in year $t$, including the log GDP of country $i$ in year $t$ ($GDP_{it}$), the log GDP of country $j$ in year $t$ ($GDP_{jt}$), the log distance between countries $i$ and $j$ ($DIST_{ij}$), and the dummy variable of a bilateral trade agreement between countries $i$ and $j$ ($TA_{ij}$), among others. The fixed effects $\alpha_i$, $\gamma_j$, and $\lambda_t$ represent the unobserved exporting country effects, destination country effects, and year effects, respectively. Researchers are often interested in the coefficient of $DIST_{ij}$ interpreted as the trade elasticity or the trade cost. Another important parameter of empirical interest is the coefficient of $TA_{ij}$ interpreted as the effect of bilateral trade agreements on trade volumes. See HeadMayer2014 for a comprehensive review of gravity models.

To date, variants of the three-dimensional panel model ((ref)) have been extensively used in empirical analysis of international trade (see BaltagiEggerErhardt2017 for a survey), housing (see BaltagiBresson2017 for a survey), migration (see Ramos2017 for a survey), and consumer price. In these analyses, researchers employ various combinations of fixed effects, including (I) $\alpha_i + \gamma_j$, (II) $\alpha_i + \gamma_j + \lambda_t$, and (III) $\alpha_{it} + \gamma_{jt}$, among others.\footnote{ Parameters $\beta$ of certain types of controls are not identified under more general combinations of fixed effects. For example, the coefficients of $GDP_{it}$ and $GDP_{jt}$ are not identified under the fixed effect model (III) due to the collinearity. However, the coefficients of $DIST_{ij}$ and $TA_{ij}$ would be identifiable under any of the three models. In empirical analysis of bilateral trade flows, the latter two coefficients are of more common interest. In fact substituting fixed effects (such as $\alpha_{it}$ and $\gamma_{jt}$) for observed proxies (such as $GDP_{it}$ and $GDP_{jt}$) is “now common practice and recommended by major empirical trade economists” HeadMayer2014. } See BalazsiMatyasWansbeek2017 for a comprehensive list of empirical papers and their specifications of the combinations of fixed effects. Researchers in general do not know which combination of fixed effects correctly specifies the model of their interest. If the true model is parsimonious and a researcher erroneously assumes a rich specification, then na\"ive fixed effect estimators generally entail exacerbated variances. On the other hand, if the true model is rich and a researcher erroneously assumes a parsimonious specification, then na\"ive fixed effect estimators generally entail mis-specification biases. The lack of knowledge of the true model specification therefore leads to undesired econometric results in any event.

A recent paper by LuMiaoSu2018 develops a method of model selection. Their method serves as a useful guideline for empirical researchers to choose a correct combination of fixed effects in three-dimensional panel models. When a researcher uses a selected model to compute estimates of $\beta$ and their standard errors, it is also important that she takes into account the statistical effects of the model selection. To our knowledge, the existing literature does not provide a method of post-selection inference for three-way panel models. In this light, we extend the frontier of this existing econometric literature LuMiaoSu2018 by providing a method of inference for $\beta$ accounting for the effect of the model selection. We make use of the lasso technique along with de-biasing to this end, but our method does not require exactly sparse fixed effects. In other words, our assumptions do allow for many and even all of the fixed effects to be nonzero in a general combination of fixed effects.

{\bf Related Literature} A three-dimensional panel model was suggested by Matyas1997. The literature on multi-dimensional panels is extensive today, and is surveyed in the book of article collections edited by Matyas2017. Its chapter written by BalazsiMatyasWansbeek2017 provides a comprehensive list of empirical research papers employing multi-dimensional panel data.

Methods of model selection in three-dimensional panels are developed by LuMiaoSu2018, and this paper was motivated by LuMiaoSu2018. As stated earlier, we aim to extend this frontier of the literature by developing a post-selection inference for the regression parameters.

We use the lasso technique for model selection and post-selection inference, but our assumptions do allow for all fixed effects to be nonzero. This is because we rely on the approximate sparsity condition as opposed to the conventional sparsity. Post-selection inference via lasso is studied by an extensive body of the literature in various contexts. This literature includes, but are not limited to, BelloniChenChernozhukovHansen2012 for IV models, and BelloniChernozhukovHansen2014, JavanmardMontanari2014, vandeGeeretal2014, and ZhangZhang2014 for linear regression models.

Lasso estimation for panel models are suggested by Koenker2004, Lamarche2010, Kock2013, CanerHan2014, LuSu2016, LiQianSu2016, QianSu2016, CanerHanLee2018, HardingLamarche2019, among others. Classification and estimation by lasso for panel models are proposed by SuShiPhillips2016 -- also see LuSu2017, SuJu2018, and SuWangJin2017. For post-selection inference with panel data using lasso, BelloniChernozhukovHansenKozbur2016 work with de-meaned fixed effect models with high-dimensional controls using post-double-selection estimator. Kock2016 and KockTang2018 work with correlated random effect panel models and dynamic panel models with sparse fixed effects via de-biased lasso, respectively. We extend this frontier of the literature to three-dimensional panels. Besides the different framework of three-dimensional panels as opposed to two-dimensional ones, this paper is different from Kock2016 and KockTang2018 in the following four technical points. First, we extend the theory of nodewise lasso by allowing for different convergence rates to incorporate a larger class of fixed effect models. Second, we use a different proof strategy with the sparsity requirement of $ss_l(\log (p\vee (NM)))^2/(N \wedge M)=o(1)$ inspired by BelloniChenChernozhukovHansen2012, whereas an adaptation of the proof strategies of Kock2016\footnote{See Assumption A3 (b) of Kock2016.} and KockTang2018\footnote{See Assumption 5 (c) of KockTang2018.} to our framework would require $ss^2_l (\log(p\vee (NM)))^2/(N \wedge M)=o(1)$. This feature further extends the class of models that can be handled under our framework. Third, the sub-gaussianity assumption of covariates, which is assumed by the majority of papers in the de-biased lasso literature, is not required. Fourth, we allow for non-sparse coefficients based on the notion of approximate sparsity following that of BelloniChenChernozhukovHansen2012 instead of the $L^v$ sparsity for $0<v<1$ as in KockTang2018.

With all these technical relations to the existing literature, we once again emphasize that our main contribution is the robust inference method for three-dimensional panels. Unlike two-dimensional panels, there are a number of alternative combinations of fixed effect specifications in three-dimensional panels, and hence model selection is more important in these models LuMiaoSu2018. We apply and extend state-of-the-art technology BelloniChenChernozhukovHansen2012,Kock2016,KockTang2018 to this three-dimensional panel framework which concerns many empirical researchers.

{\bf Organization:} The rest of this paper is organized as follows. We introduce the model framework in Section (ref). An overview of our proposed method is presented in Section (ref). The main theoretical result is presented in Section (ref), followed by sufficient conditions discussed in Section (ref). We discuss the key assumption in the context of gravity analysis of international trade in Section (ref). We conduct simulation studies in Section (ref). Section (ref) concludes the paper.

The Model Framework

Consider the following representation of a general class of three-dimensional panel models with large $N$ and large $M$.

align[align omitted — 428 chars of source]

This representation consists of a $k$-dimensional parameter vector $\beta$, $N$-dimensional parameter vector $\alpha_{[N]}=(\alpha_1,...,\alpha_{N})'$, $M$-dimensional parameter vector $\gamma_{[M]}=(\gamma_1,...,\gamma_{M})'$, $T$-dimensional parameter vector $\lambda_{[T]}=(\lambda_1,...,\lambda_{T})'$, $NT$-dimensional parameter vector $\alpha_{[NT]}=(\alpha_{11},...,\alpha_{NT})'$, and $MT$-dimensional parameter vector $\gamma_{[MT]}=(\gamma_{11},...,\gamma_{MT})'$. In total, there are $k+N+M+T+NT+MT$ parameters involved in this representation ((ref)).

Recall that conventional fixed effect models include

enumerate[(I)] • $\alpha_i + \gamma_j$, • $\alpha_i + \gamma_j + \lambda_t$, and • $\alpha_{it} + \gamma_{jt}$,

among others. Model (I) entails $k+N+M$ of possibly nonzero parameters $(\beta',\alpha_{[N]}',\gamma_{[M]}')'$, while the rest of the $T+NT+MT$ parameters $(\lambda_{[T]}',\alpha_{[NT]}',\gamma_{[MT]}')'$ are all zero. Similarly, Model (II) entails $k+N+M+T$ of possibly nonzero parameters $(\beta',\alpha_{[N]}',\gamma_{[M]}',\lambda_{[T]})'$, while the rest of the $NT+MT$ parameters $(\alpha_{[NT]}',\gamma_{[MT]}')'$ are all zero. Likewise, Model (III) entails $k+NT+MT$ of possibly nonzero parameters $(\beta',\alpha_{[NT]}',\gamma_{[MT]}')'$, while the rest of the $N+M+T$ parameters $(\alpha_{[N]}',\gamma_{[M]}',\lambda_{[T]}')'$ are all zero. Furthermore, the representation ((ref)) includes many other combinations than these three models.

When Model (I) is true for example, then the representation ((ref)) has $T+NT+MT$ redundant parameters and hence estimating the model ((ref)) generally yields much larger standard errors for the parameters $\beta$ of interest than necessary. This motivates the need of model selection. We propose to use the lasso to select such redundant fixed effect parameters out of the representation ((ref)), and then conduct inference robustly accounting for the statistical effects of the model selection.

For ease of conducting econometric analysis, we further rewrite the representation ((ref)) as

align[align omitted — 201 chars of source]

where $\mathbf{x}_{ijt} = ({x}_{ijt}', \mathbbm{1}_{t=1}, ..., \mathbbm{1}_{t=T})'$ and $\boldsymbol{\beta} = (\beta', \lambda_1 , ..., \lambda_T)'$ are of dimension $k_0=k+T$, $\mathbf{d}_{1,it} = (\mathbbm{1}_{i=1}, ..., \mathbbm{1}_{i=N}, \mathbbm{1}_{i=1}\mathbbm{1}_{t=1}, ..., \mathbbm{1}_{i=N}\mathbbm{1}_{t=T})'$ and $\boldsymbol{\alpha} = (\alpha_{[N]},\alpha_{[NT]})'$ are of dimension $N_0=N+NT$, and $\mathbf{d}_{2,jt} = (\mathbbm{1}_{j=1}, ..., \mathbbm{1}_{j=N}, \mathbbm{1}_{j=1}\mathbbm{1}_{t=1}, ..., \mathbbm{1}_{j=M}\mathbbm{1}_{t=T})'$ and $\boldsymbol{\gamma} = (\gamma_{[M]},\gamma_{[MT]})'$ are of dimension $M_0=M+MT$.

Suppose that we can decompose the fixed effects $\boldsymbol{\alpha}$ into $\overline{\boldsymbol{\alpha}}$ and $\boldsymbol{\alpha} - \overline{\boldsymbol{\alpha}}$ and decompose the fixed effects $\boldsymbol{\gamma}$ into $\overline{\boldsymbol{\gamma}}$ and $\boldsymbol{\gamma} - \overline{\boldsymbol{\gamma}}$ such that

align[align omitted — 393 chars of source]

and

align[align omitted — 393 chars of source]

hold, where $\left \|\cdot\right \|_0$ denotes the support cardinality (the $L^0$ norm).\footnote{With this said, we emphasize that this decomposition is merely theoretical, and a researcher need not implement such a decomposition in practice. Precise requirements for the decomposition are stated in Assumptions (ref) and (ref) (4) ahead, followed by discussions in the context of our motivating application ((ref)) in Remark (ref). In Section (ref), we use world trade data to argue that these assumptions are plausible in the application ((ref)). } Such a decomposition is constructed for example by setting $\overline{\boldsymbol{\alpha}}_\ell$ equal to $\boldsymbol{\alpha}_\ell$ for those coordinates $\ell$ for which $\left\vert \boldsymbol{\alpha}_\ell \right\vert$ is large and setting $\overline{\boldsymbol{\alpha}}_\ell$ equal to zero for those coordinates $\ell$ for which $\left\vert \boldsymbol{\alpha}_\ell \right\vert$ is small, and similarly for $\boldsymbol{\gamma}$. Consequently, we can further rewrite the representation ((ref)) as

align*[align* omitted — 195 chars of source]

where $r_{ijt}$ is the approximation error defined by

align*[align* omitted — 195 chars of source]

and it satisfies

align*[align* omitted — 214 chars of source]

Stacking the three-dimensional panel data across the $NMT$ observations, we in turn construct the matrix representation

align[align omitted — 205 chars of source]

where $Y=(y_{111},...,y_{NMT})'$, $R=(r_{111},...,r_{NMT})'$, and $\varepsilon=(\varepsilon_{111},...,\varepsilon_{NMT})'$, are vectors of dimension $NMT$, $X=(\mathbf{x}_{111},...,\mathbf{x}_{NMT})'$ is a matrix of size $NMT \times k_0$, $D_1=(\mathbf{d}_{1,11},...,\mathbf{d}_{1,NT})'$ is a matrix of size $NMT \times N_0$, $D_2=(\mathbf{d}_{2,11},...,\mathbf{d}_{2,MT})'$ is a matrix of size $NMT \times M_0$, $Z=[X \ D_1 \ D_2]$, and $\overline{\boldsymbol{\eta}} = [\boldsymbol{\beta}' \ \overline{\boldsymbol{\alpha}}' \ \overline{\boldsymbol{\gamma}}']'$ is a vector of dimension $p=k_0+N_0+M_0$.

If the true model is parsimonious, like Model (I), then a large number of the elements of the high-dimensional parameters, $\boldsymbol{\alpha}$ and $\boldsymbol{\gamma}$, will be zero. Thus, a large number of the elements of $\overline{\boldsymbol{\alpha}}$ and $\overline{\boldsymbol{\gamma}}$ will be zero. Furthermore, for those coordinates of $\boldsymbol{\alpha}$ and $\boldsymbol{\gamma}$ that are small in absolute value, the corresponding coordinates of $\overline{\boldsymbol{\alpha}}$ and $\overline{\boldsymbol{\gamma}}$ are set to zero in the decomposition in light of the relatively smaller approximation errors caused by setting them to zero. We propose to use the lasso technique to select such redundant parameters in $\overline{\boldsymbol{\alpha}}$ and $\overline{\boldsymbol{\gamma}}$ out of this high-dimensional model as means of model selection for the purpose of obtaining smaller standard errors. Furthermore, accounting for the statistical effects of this model selection, we then conduct robust inference for the main parameters $\boldsymbol{\beta}$ in the panel model. Section (ref) illustrates an overview of our proposed method. A formal theoretical analysis will then follow in Sections (ref) and (ref).

Overview of the Method

Our proposed method consists of four steps. The first step is a lasso estimation of the parameter vector $\boldsymbol{\eta}$ entailing a model selection. The second step is an auxiliary step to calculate an approximate inverse of the Gram matrix to be used in the subsequent two steps. The third step de-biases the regularized lasso estimate from the first step. The fourth step is a calculation of the asymptotic variance of each coordinate of the de-biased lasso estimator of $\boldsymbol{\beta}$. \\ {\bf Step 1:} For the representing equation ((ref)), define the lasso estimator

align[align omitted — 173 chars of source]

where $\mu \in [0,\infty)$ is a regularization tuning parameter and the penalty function $P$ is defined by

align*[align* omitted — 252 chars of source]

for some diagonal normalization matrix $\widehat\Upsilon_\ell$ for each $\ell \in \{1,2,3\}$.\footnote{See Remark (ref) in Appendix (ref).} In practice, the regularization tuning parameter $\mu$ can be chosen using a cross validation via software packages. \\ {\bf Step 2:} The next step is an auxiliary process to obtain a $p \times p$ matrix $\widehat\Theta$ of approximate inverse of the Gram matrix to be used in Step 3. We define the nodewise lasso estimator

align[align omitted — 255 chars of source]

of the $\ell$-th column $Z^{\ell}$ on all the other $(p-1)$ columns $Z^{-\ell}$ for each $\ell \in \{1,...,p\}$, where $\mu_{\text{node}}^\ell \in [0,\infty)$ is a regularization tuning parameter, $\widehat\Upsilon_{\text{node}}^\ell$ is some diagonal normalization matrix for each $\ell \in \{1,...,p\}$, and $S_{-\ell}$ is the $(p-1) \times (p-1)$ matrix obtained by removing the $\ell$-th row and the $\ell$-th column of

align*[align* omitted — 140 chars of source]

In practice, the regularization tuning parameter $ \mu_{\text{node}}^\ell$ can be chosen using a cross validation via software packages.

Once the nodewise lasso estimates $\widehat\phi^\ell$ are obtained, a $p \times p$ matrix $\widehat\Theta$ approximating the inverse Gram matrix can be constructed by

align[align omitted — 700 chars of source]

with $\widehat\tau_\ell$ given by

align*[align* omitted — 232 chars of source]

for each $\ell \in \{1,...,p\}$ and $\widehat\phi^\ell_l$ denoting the $l$-th coordinate of the nodewise lasso estimate $\widehat\phi^\ell$ for each $\ell \in \{1,...,p\}$ and $l \in \{1,...,p-1\}$. \\ {\bf Step 3:} The shrinkage by the regularization $\mu P(\boldsymbol{\eta})$ forces a sub-vector of the lasso estimates $\widehat{\boldsymbol{\eta}}$ to be zero, and this mechanism serves as means of model selection. Since this regularization biases the second-stage lasso estimator $\widehat{\boldsymbol{\eta}}$, we further `de-bias' it according to

align[align omitted — 177 chars of source]

for each $\ell \in [p]$, where $\widehat\Theta_\ell$ is the $\ell$-th column of $\widehat\Theta$ and $\widehat\Theta$ is the $p \times p$ approximate inverse Gram matrix constructed in Step 2. The sub-vectors of $\widetilde{\boldsymbol{\eta}}$ will be denoted by $\widetilde{\boldsymbol{\eta}} = \left( \widetilde{\boldsymbol{\beta}}', \widetilde{\boldsymbol{\alpha}}', \widetilde{\boldsymbol{\gamma}} \right)'$. \\ {\bf Step 4:} The asymptotic variance of $\sqrt{NM}\left(\widetilde{\boldsymbol{\beta}}_\ell - \boldsymbol{\beta}_\ell\right)$ for $\ell \in \{1,...,k_0\}$ is approximated by

align*[align* omitted — 94 chars of source]

where $\widehat\Theta_\ell$ is defined in Step 3,

align*[align* omitted — 186 chars of source]

and $\widehat\varepsilon_{ijt}$ is the residual from the lasso in Step 1.

The Main Theory

Define the de-biased lasso estimator by

align[align omitted — 158 chars of source]

where $P'$ denotes the sub-gradient of $P$. Recall that the sub-vectors of $\widetilde{\boldsymbol{\eta}}$ are denoted by $\widetilde{\boldsymbol{\eta}} = [ \widetilde{\boldsymbol{\beta}}', \widetilde{\boldsymbol{\alpha}}', \widetilde{\boldsymbol{\gamma}} ]'$, corresponding to $\overline{\boldsymbol{\eta}} = [\boldsymbol{\beta}' \ \overline{\boldsymbol{\alpha}}' \ \overline{\boldsymbol{\gamma}}']'$. This section presents a general limit distribution result for each coordinate of the de-biased lasso estimator $\widetilde{\boldsymbol{\beta}}$ for the coefficients of $x_{ijt}$. We focus on short panels with fixed $T$ and large $(N,M)$, although an extension to large $T$ cases may be feasible with alternative assumptions. While we maintain high-level assumptions in the current section for the sake of generality, we will follow up with lower-level sufficient conditions in Section (ref). Define the $p \times p$ rate-adjusted Gram matrix

align[align omitted — 303 chars of source]

Let $[n] = \{1,...,n\}$ for any $n \in \mathbb{N}$. With these notations, consider the following assumption.

assumption[Asymptotic Normality] For all $(N,M)$, there exists a column random vector $\widehat \Theta_l$ such that the following conditions hold for an $(N,M)$-dependent choice of $\mu$ as $N,M\rightarrow\infty$. \begin{enumerate}[(i)] • $\max_{l\in[k_0]} \left\vert \sqrt{NM}(\widehat\Theta_l' Q\bar\Psi Q -e_l')(\widehat {\boldsymbol{\eta}} - \overline{\boldsymbol{\eta}} ) \right\vert = o_p(1)$. • $\max_{l\in[k_0]} \left\vert \widehat \Theta_l' Z'R/\sqrt{NM} \right\vert=o_p(1)$. • For each $l\in[k_0]$, there exists $V_{ll} \in (0,\infty)$ that can depend on $(N,M)$ such that \begin{align*} V_{ll}^{-1/2}\widehat\Theta_l' Z'\varepsilon /\sqrt{NM} \leadsto N(0,1). \end{align*} \end{enumerate}

In the current general theoretical discussions, Assumption (ref) merely requires an existence of some $\widehat \Theta_l$ satisfying the three conditions, and does not say how it should be constructed. Recall that the overview of the method in Section (ref) suggests a concrete way to construct such $\widehat \Theta_l$. Section (ref) ahead will discuss lower-level sufficient conditions to guarantee that such a concrete construction of $\widehat \Theta_l$ satisfies the three high-level conditions in Assumption (ref).

theorem[Asymptotic Normality] Suppose that Assumption (ref) (i)--(ii) are satisfied. Then, \begin{align*} \widetilde {\boldsymbol{\eta}}_l - \overline{\boldsymbol{\eta}}_l = \frac{1}{NM} \widehat \Theta'_l Z'\varepsilon+o_p \left( 1/\sqrt{NM} \right) \end{align*} for each $l \in [p]$. Furthermore, if Assumption (ref) (iii) is satisfied in addition, then we have \begin{align*} \sqrt{NM}(\widetilde {\boldsymbol{\beta}}_l - {\boldsymbol{\beta}}_l)\leadsto N(0,V_{ll}) \end{align*} for each $l\in [k_0]$.

A proof is found in Appendix (ref).

remarkThe de-biased lasso estimator $\widetilde {\boldsymbol{\eta}}_l = \widehat {\boldsymbol{\eta}}_l -\frac{\mu }{NM}\widehat \Theta'_l P'(\widehat {\boldsymbol{\eta}})$ can be also rewritten by replacing $\mu P'(\widehat {\boldsymbol{\eta}})$ by $-Z'(Y-Z\widehat {\boldsymbol{\eta}})$ following the K.K.T. condition, i.e., $\widetilde {\boldsymbol{\eta}}_l = \widehat {\boldsymbol{\eta}}_l +\frac{1}{NM}\widehat \Theta_l' Z'(Y-Z\widehat {\boldsymbol{\eta}})$. This representation yields the de-biased lasso formula proposed in ((ref)).

Sufficient Conditions and Variance Estimation

In this section, we propose lower-level sufficient conditions for the high-level general statements in Assumption (ref). These conditions provide a theoretical guarantee for the concrete practical procedure of Section (ref) to work. While the general limit distribution result in Theorem (ref) did not specify a concrete form of the asymptotic variance $V_{ll}$, the current section also provides a formula for it under these sufficient conditions. Furthermore, we propose an analog variance estimator $\widehat V_{ll}$, and show its consistency under these sufficient conditions.

Throughout this section, we will assume $\widehat \Upsilon=I_p$ and $\widehat \Upsilon_{\text{node},l}=I_{p-1}$ for all $l\in[k_0]$ for simplicity, although these restrictions are not essential at all. We use the following notations for the parameter supports: $J_1 = \text{supp}(\boldsymbol{\beta}),$ $J_2 = \text{supp}(\overline{\boldsymbol{\alpha}}),$ $J_3 = \text{supp}(\overline{\boldsymbol{\gamma}}),$ and $J = \text{supp}(\overline{\boldsymbol{\eta}})$, Their cardinalities are denoted by $s_1 = |J_1|,$ $s_2 = |J_2|,$ $s_3 = |J_3|,$ and $s = |J|.$ We note that $s$ is non-decreasing in $N$ and/or $M$. Similarly to the decomposition ((ref)) for the main regression model, we also consider the decomposition

align[align omitted — 115 chars of source]

for each coordinate $l\in [k_0]$ of the regressors.

Sufficient Conditions

We present sufficient conditions as five modules, Assumptions (ref), (ref), (ref), (ref), and (ref), listed below.

assumption[Approximate Sparsity] (1) $\| \overline{\boldsymbol{\eta}} \|\le K$. (2) $\|R\|\le c_s \lesssim \sqrt{s }$ with probability $1-o(1)$. (3) $ \| Z'R \| = o_p\Big( \sqrt{NM} \Big). $

Recall that the fixed effects $\boldsymbol{\alpha}$ are decomposed into $\overline{\boldsymbol{\alpha}}$ and $\boldsymbol{\alpha} - \overline{\boldsymbol{\alpha}}$ such that ((ref)) is satisfied, and the fixed effects $\boldsymbol{\gamma}$ are decomposed into $\overline{\boldsymbol{\gamma}}$ and $\boldsymbol{\gamma} - \overline{\boldsymbol{\gamma}}$ such that ((ref)) is satisfied. These conditions ((ref)) and ((ref)) are imposed to satisfy Assumption (ref) (1) and (2). Assumption (ref) (3) can be relaxed to a weaker condition,\footnote{For example, $\sup_{\substack{\|\xi\|=1 \\ \|\xi\|_0 = Cs}} \| \xi' Z'R \| = o_p(\sqrt{NM})$ for some finite positive $C$.} but we present the current condition for its better interpretation.

remark[Discussion of the Approximate Sparsity Condition] We emphasize that the approximate sparsity condition of Assumption (ref) (together with Assumption (ref) (4) to be stated below) allows for many and even all the fixed effects (i.e., $\boldsymbol{\eta}$ as opposed to $\overline{\boldsymbol{\eta}}$) to be nonzero. The assumption should be interpreted as a requirement for how the fixed effects can be decomposed into the sparse components ($\overline{\boldsymbol{\alpha}}$ and $\overline{\boldsymbol{\gamma}}$) and the remaining components ($\boldsymbol{\alpha} - \overline{\boldsymbol{\alpha}}$ and $\boldsymbol{\gamma} - \overline{\boldsymbol{\gamma}}$) generating $ R = D_1 \left(\boldsymbol{\alpha} - \overline{\boldsymbol{\alpha}}\right) + D_2 \left(\boldsymbol{\gamma} - \overline{\boldsymbol{\gamma}}\right). $ Indeed, the assumption implicitly imposes a non-trivial restriction on sampling procedures. For example, an i.i.d. sampling of fixed effects is not accommodated, although this feature does not contradict with our sampling assumption to be stated below as Assumption (ref). With this said, the same limitations apply to all the preceding papers (cf. Section (ref)) that employ (approximate) sparsity conditions on fixed effects in panel data. In fact, the approximate sparsity is a rather plausible assumption for the sampling process in the context of our motivating application ((ref)). In gravity analysis of trade, researchers initially used only the G7 countries, later added the OECD countries, and smaller economies have been added more recently. Nearly half of all import and export flows are determined by the top ten largest economies. Newly added countries to the sample tend to have very small trade volumes. This sampling process entails fixed effects taking smaller values as sample size increases, and it does not contradict with the approximate sparsity requirement. Section (ref) elaborates on the approximate sparsity in trade volumes based on the actual world trade data. $\triangle$
assumption[Moments] For each $(N,M)$, the random vectors $(Y'_{ij1},Z'_{ij1},...,Y'_{ijT},Z'_{ijT})'$, $(i,j)\in [N]\times [M]$, are independently distributed. Furthermore, there exist $q \in (4, \infty)$ and $K \in (0,\infty)$ not depending on $(N,M)$ such that the following conditions hold for all $l\in[k_0]$. \begin{enumerate}[(1)] • $ \Big( \frac{1}{NM}\sum_{i=1}^N \sum_{j=1}^M E\Big[\max_{ t\le T}\|X_{ijt}\|^{2q}_\infty\Big]\Big)^{1/2q}\le B_{NM} $ and $ \Big( E|X_{ijt,l}|^{2q} \Big)^{1/2q} \le K $ hold for all $i,j,t,l$, where $B_{MN}$ satisfies $B_{NM}\sqrt{\log (p \vee (NM))}\lesssim (NM)^{1/2-1/q}$; • $\|(D_1,D_2)\|_\infty=1$; and • $\frac{1}{NM}\sum_{i=1}^N \sum_{j=1}^M \sum_{t=1}^T E\varepsilon_{ijt}^{2q}\vee \frac{1}{NM}\sum_{i=1}^N \sum_{j=1}^M \sum_{t=1}^T E(\zeta^l_{ijt})^{2q}\le K^{2q} <\infty$. \end{enumerate}

For any squared matrix $A$, define the sparse eigenvalues by

align*[align* omitted — 212 chars of source]

With these notations, we state the following assumption of sparse eigenvalues for the rate-adjusted Gram matrix $\bar\Psi$ defined in ((ref)).

assumption[Sparse Eigenvalues] For any $C > 0$, there exist constants $0<\underline k < \overline k<\infty$, not depending on $(N,M)$, such that \begin{align*} \underline k \le \varphi_{\min} (\bar\Psi,Cs)\le\varphi_{\max} (\bar\Psi,Cs) \le \overline k \end{align*} with probability approaching one.

For each $(N,M)$, we write $\Psi = E\bar\Psi$ depending on $(N,M)$, With this notation, the auxiliary decomposition ((ref)) is made according to the following conditions.

assumption[Nuisance Parameters] The following conditions are satisfied. \begin{enumerate}[(1)] • $\max_{l \in [k_0]}\|\phi^l\|_0\le s_l$ and $\max_{l \in [k_0]}\|\phi^l\|+(s_l)^{-1/2}\|\phi^l\|_1\le K$; • For all $l\in [k_0]$, $\|r_l\|\le \sqrt{s_l}$; • For all $(N,M)$, $0<L<\Lambda_{\min}(\Psi)<\Lambda_{\max}(\Psi)<U<\infty$ for $L$, $U$ independent of $(N,M)$; • $\max_{l\in[k_0]}(s_l\vee s)\sqrt{\frac{(\log (p \vee (NM)))^2}{N\wedge M}}=o(1)$. \end{enumerate}

Accounting for the possible dependence, we define the cluster-robust variance matrix

align*[align* omitted — 184 chars of source]

For each $(N,M)$, we write $\Theta=(E[\frac{Z'Z}{NM}])^{-1}$ depending on $(N,M)$. Let $\Theta_l$ denote the $l$-th column of $\Theta$. We state the following assumption of finite and non-zero variance.

assumption[Variance] For any $(N,M)$ and for all $l\in[k_0]$, $ \|\Omega\| <\infty $ and $ \Theta_l'\Omega\Theta_l \ge\underline k >0 $ for a constant $\underline k $ which is independent of the sample size.
remarkNotice that the conditions above are imposed on the Gram matrices, $\bar \Psi$ and $\Psi$, re-weighted by effective sample size, rather than the original Gram matrices, $Z'Z/NM$ and $EZ'Z/NM$. Assumption (ref) is weaker than the common assumptions required in the literature, such as sub-gaussianity or uniform boundedness. Assumption (ref) is also assumed by BelloniChenChernozhukovHansen2012 and BelloniChernozhukovHansenKozbur2016. It requires some small sub-matrices of the big $p\times p$ re-weighted Gram matrix to be well-behaved. Lower level sufficient conditions are also possible by using Lemma P1 in BelloniChernozhukovChetverikovWei2018, but are not pursued here. Assumption (ref) (1) and (2) impose sparsity on the nodewise regression parameters and the approximation errors. Assumption (ref) (3) requires $\Psi$, the expectation of the re-weighted Gram matrix, to be positive definite uniform over $(N,M)$. These are rather standard in the literature. Assumption (ref) limits the models that can be handled in terms of their dimensionality and sparsity. Note that we need only $ss_l(\log (p\vee (NM)))^2/(N \wedge M)=o(1)$, whereas an adaptation of the proof strategies of Kock2016 and KockTang2018 to our framework would entail $ss^2_l(\log(p\vee (NM)))^2/(N \wedge M)=o(1)$. Finally, Assumption (ref) requires $\Omega$ in the sandwich form to be well-behaved.

The following proposition states that Assumptions (ref), (ref), (ref), (ref), and (ref) are sufficient for the high-level conditions in Assumption (ref), with a concrete variance formula motivating the practical guideline of Section (ref).

propositionAssumptions (ref), (ref), (ref), (ref), and (ref) imply Assumption (ref) with $V_{ll} = \Theta_l' \Omega \Theta_l$.

A proof is found in Appendix (ref). Combining Theorem (ref) and Proposition (ref) together, we state the following corollary.

corollary[Asymptotic Normality] If Assumptions (ref), (ref), (ref), (ref), and (ref) are satisfied, then \begin{align*} \sqrt{NM}(\widetilde {\boldsymbol{\beta}}_l - {\boldsymbol{\beta}}_l)\leadsto N(0,V_{ll}) \end{align*} for each $l\in [k_0]$, where $V_{ll} = \Theta_l' \Omega \Theta_l$.
remarkWe conjecture that one can further enhance the results of Corollary (ref) by showing the honesty property (uniform validity over a large set of parameters) of confidence intervals using the proposed procedure with no extra assumption by adapting the proof strategy of Theorem 3 of CanerKock2018 or Theorem 3 of KockTang2018 to our framework.

Asymptotic Variance Estimation

Based on the asymptotic variance formula presented in Proposition (ref), we suggest to compute the cluster-robust asymptotic variance of $\sqrt{NM}\left(\widetilde{\boldsymbol{\beta}}_\ell - \boldsymbol{\beta}_\ell\right)$ by

align*[align* omitted — 90 chars of source]

as suggested in Section (ref). This estimator is consistent in the current assumptions as formally stated in the following theorem.

theorem[Variance Estimator] If Assumptions (ref), (ref), (ref), (ref) and (ref) are satisfied, then $$ \max_{l\in [k_0]}|\widehat V_{ll} - V_{ll}|=o_p(1). $$

A proof is found in Appendix (ref).

Approximate Sparsity in Gravity Analysis of Trade

In this section, we discuss our key assumption, namely the assumption of approximate sparsity (Assumptions (ref) and (ref) (4) -- also see Remark (ref)), in the gravity model ((ref)) of international trade. The idea behind the approximate sparsity assumption is that only a small number of observations have large fixed effect values, and the remaining majority of observations have relatively modest fixed effect values that can be summarized into the approximation error term $r_{ijt}$. The assumption is likely satisfied in sampling processes where, after collecting observations with relatively large values of fixed effects (e.g., G7 and OECD countries), the remaining additions tend to have smaller values of fixed effects. We argue that this is plausible in common settings such gravity analysis in international trade.

To illustratea this point, we retrieved data from the World Integrated Trade Solution (WITS) Database, a common source of trade flows and trade costs used in gravity analysis.\footnote{This database was developed by the World Bank in conjunction with the United Nations Conference on Trade and Development (UNCTAD), the International Trade Center, United Nations Statistical Division (UNSD) and the World Trade Organization (WTO). The database combines information on trade flows from the UN Comtrade database, tariff and non-tariff barriers from the UN TRAINS database, and the both preferential and MFN tariffs from the WTO's Integrated Data Base.} We focus on country-specific import and export flows and aim to make two specific points. First, in any given year, trade is largely dominated by a few large countries. For instance, in 2015, the WITS database contains positive import flows for 237 countries and positive export flows for 232 countries. Nonetheless, nearly half of all import (respectively, export) flows are determined by the top 10 largest importers (respectively, exporters) alone. Not surprisingly, the largest importers are also the largest exporters. Second, the importance of these countries has remained stable over time, despite the fact that (a) world trade has grown exponentially over time and (b) WITS records exports and imports for a substantially larger number of countries in recent years than it did even a few years ago. In this sense, the `new' additions to trade databases tend to have very small trade flows.

Using a country's share of world imports (Table (ref)) as a measure of `importer' importance or a country's share of world exports (Table (ref)) as a measure of exporter importance, we document the 10 largest trading nations every 5 years starting in 1990. We note the following three attributes of standard trade data: (1) a small number of countries account for the large majority of world trade; (2) whether a country represents a large or small fraction of trade flows changes slowly over time; and (3) even though many developing countries have grown substantially since 1990, the average share of small countries has not changed very much. This last feature is largely due to the fact that the `new' countries which are added to world trade databases are nearly always very small. To make points (1) and (2) particularly clear, we would expect that a typical country would have an import/export share of roughly 0.5% for a sample of about 200 countries. However, in any given year, fewer than 40 countries have import or export shares of 0.5%. Of the countries which have been added to the import database since 1990, their average (median) import share was 0.07% (0.01%). Similarly, among the countries added to the export database since 1990, their average (median) export share was 0.09% (0.03%). Regardless of how to measure the size of these peripheral countries, their overall contribution to world trade is extremely small.

In summary, only a small number of observations have high trade volumes. The large majority of remaining observations have very modest and almost negligible trade shares. This pattern remains stable over time. Since researchers first collect observations with large volumes (e.g., G7 and OECD countries), new additions to the data thereafter entail relatively small volumes. This common sampling process in gravity analysis of international trade is compatible with our key assumption, namely the the assumption of approximate sparsity (Assumptions (ref) and (ref) (4) -- also see Remark (ref)).

Simulation Studies

Simulation Setting

Consider the following three fixed effect models of three-dimensional panel data.

align*[align* omitted — 325 chars of source]

Model (I) is nested by Model (II), and Model (II) is in turn nested by Model (III). Therefore, Model (I) is the most parsimonious and subject to under-fitting, whereas Model (III) is the richest and subject to over-fitting. If a researcher runs a fixed effect estimator under Model (I) when Model (II) or (III) is true, then the estimates generally suffers from mis-specification biases. If a researcher runs a fixed effect estimator under Model (III) when Model (I) or (II) is true, then the estimates generally suffers from larger standard errors than necessary.

We run simulations for varying sizes of $N$ and $M=N-1$, while the length of time is set to $T=5$ throughout. This setting follows from our asymptotic theory where $N$ and $M$ increases but $T$ does not. The $i$ and $j$ fixed effects are generated by $\alpha_i \sim N\left(m_\alpha, s_\alpha^2 \left/\left(\sqrt{i} \cdot (\log (i+1))^3\right)\right.\right)$ and $\gamma_j \sim N\left(m_\gamma, s_\gamma^2 \left/\left(\sqrt{j} \cdot (\log (j+1))^3\right)\right.\right)$ independently, where $m_\alpha=m_\gamma=0$ and $s_\alpha=s_\gamma=1$. The $t$ fixed effects are generated by $\lambda_t=0$ for all $t$ but for one year $t$ when a universal shock of $\lambda_t = 2$ is applied. The $it$ and $jt$ fixed effects are generated by $\alpha_{it} \sim N\left(m_{\alpha}, s_{\alpha}^2 \left/\left(\sqrt{i} \cdot (\log (i+1))^3\right)\right.\right)$, $\gamma_{jt} \sim N\left(m_{\gamma}, s_{\gamma}^2 \left/\left(\sqrt{j} \cdot (\log (j+1))^3\right)\right.\right)$, $m_{\alpha}=m_{\gamma}=0$, and $s_{\alpha}=s_{\gamma}=1$. We generate $X$ dependently on the fixed effects according to the mixture

align*[align* omitted — 100 chars of source]

where $m_x=0$, $s_x=2$, $\rho=0.5$, $\tilde x_{ijt} \sim N(0,1)$, and $F_{ijt}$ is the standardized sum of fixed effects for the unit $(i,j,t)$, i.e.,

align*[align* omitted — 587 chars of source]

for each $(i,j,t) \in \{1,...,N\} \times \{1,...,M\} \times \{1,...,T\}$. The error term is generated by $\varepsilon_{ijt} \sim N(m_\varepsilon,s_\varepsilon^2)$ independently where $m_\varepsilon=0$ and $s_\varepsilon=10$. The main coefficient of interest is set to $\beta = 1$. Each set of simulations consists of 10,000 Monte Carlo iterations of data generation, estimation, and inference.

We compare five methods of estimation and inference. These are the OLS without any individual fixed effects, the fixed effect estimator based on Model (I), the fixed effect estimator based on Model (II), the fixed effect estimator based on Model (III), and our proposed de-biased lasso estimator and post-selection inference. Note that the OLS is always under-fitting the true data generating model, and hence is expected to produce mis-specification biases. The fixed effect estimator based on Model (I) is correctly specified when the true data generating model is Model (I), but is under-fitting Model (II) and Model (III). The fixed effect estimator based on Model (II) is over-fitting Model (I), correctly specified when the true data generating model is Model (II), and under-fitting Model (III). The fixed effect estimator based on Model (III) is over-fitting Model (I) and Model (II), but is correctly specified when the true data generating model is Model (III).

Simulation Results

Table (ref) displays Monte Carlo simulation results under Model (I) (top panel), Model (II) (middle panel), and Model (III) (bottom panel) with the sample size $N=10$ ($NMT=450$). Similarly, Tables (ref) and (ref) display Monte Carlo simulation results with the sample sizes $N=15$ ($NMT=1050$) and $N=200$ ($NMT=1900$), respectively. The displayed statistics are the averages, biases, standard deviations, and root mean squared errors of estimates. Also displayed are the coverage frequencies of the true value of $\beta$ by the 95% confidence intervals. The first column of each table shows the OLS results without any individual fixed effects. The next three columns of each table show results of fixed effect estimators based on estimating equations of Model (I), Model (II), and Model (III). We shall call them FE-I, FE-II, and FE-III for succinctness. The last column of each table shows results of our proposed de-biased lasso estimator with valid post-selection inference. We shall call it POST for succinctness.

In the top panel of each table, where the true data generating model is Model (I), OLS is biased while FE-I, FE-II, and FE-III yield little biases. These results are consistent with the current simulation setting as OLS mis-specifies the true model while FE-I, FE-II, and FE-III correctly specify the true model. The bias of POST is in the middle between that of OLS and those of FE-I, FE-II, and FE-III. In other words, POST is de-biased to some extent but not to the full extent so that desired balances between the bias and variance are maintained. OLS yields a smaller standard deviation than FE-I or FE-II, and FE-III yields by far the largest standard deviation. These results are also consistent with the fact that OLS is the most parsimonious while FE-III is the most redundant in specification. POST yields an even smaller standard deviation than OLS. FE-I, as the oracle estimator, yields a smaller root mean square error than OLS, FE-II, or FE-III. Furthermore, POST yields an even smaller root mean square error than the oracle estimator, FE-I. The coverage frequency by FE-I, as the oracle estimator, is closer to the nominal level 95% than those of OLS, FE-II, or FE-III. Furthermore, POST yields the coverage frequency as close to the nominal level as the oracle estimator, FE-I. In summary, we observe that, when the true model is parsimonious, POST is more efficient than redundantly rich models and allows for as accurate inference as the oracle estimator.

In the middle panel of each table, where the true data generating model is Model (II), OLS and FE-I are biased while FE-II and FE-III yield little biases. These results are consistent with the current simulation setting as OLS and FE-I mis-specify the true model while FE-II and FE-III correctly specify the true model. The bias of POST is slightly larger than those of FE-II and FE-III, but much smaller than those of OLS and FE-I. In other words, POST is de-biased to a large extent but not to the full extent so that desired balances between the bias and variance are maintained. FE-II, as the oracle estimator, yields a smaller root mean square error than OLS, FE-I, or FE-III. Furthermore, POST yields an even smaller root mean square error than the oracle estimator, FE-II. The coverage frequency by FE-II, as the oracle estimator, is closer to the nominal level 95% than those of OLS, FE-I, or FE-III. POST yields the coverage frequency as close to the nominal level as the oracle estimator, FE-II. In summary, we observe that POST is more precise than biased parsimonious estimators, is more efficient than redundant estimators, and allows for as accurate inference as the oracle estimator.

In the bottom panel of each table, where the true data generating model is Model (III), OLS, FE-I, and FE-II are biased while FE-III yields a little bias. These results are consistent with the current simulation setting as OLS, FE-I, and FE-II mis-specify the true model while FE-III correctly specifies the true model. The bias of POST is in the middle between those of OLS, FE-I, and FE-II and that of FE-III. In other words, POST is de-biased to some extent but not to the full extent so that desired balances between the bias and variance are maintained. POST yields a smaller root mean square error than any other estimator, including the oracle estimator, FE-III. POST also yields the coverage frequency closer to the nominal level than any estimator, including the oracle estimator, FE-III. In summary, we observe that, when the true model is rich, POST is more precise than parsimonious estimators and allows for as accurate inference as the oracle estimator.

The simulation results reported above demonstrate that the proposed method (POST) can be used as a robustly applicable method of inference when a researcher does not know the correct fixed effect specification in practice. We also implemented many additional sets of simulations under alternative data generating parameters, and confirm that the qualitative pattern of these additional results remain the same as those of our baseline setting presented above. Specifically, we consistently observe that POST is more precise than biased parsimonious estimators, is more efficient than redundant estimators, and allows for as accurate inference as the oracle estimator.

Discussions

Three-dimensional panel models are widely used in empirical analysis of international trade, housing, migration, and consumer price, among others. Empirical researchers use various combinations of fixed effects for three-dimensional panels. When a researcher imposes a parsimonious model and the true model is rich, then estimation based on the assumed parsimonious model generally incurs mis-specification biases. When a researcher employs a rich model and the true model is parsimonious, then estimation based on the redundantly rich model generally incurs larger standard errors than necessary. It is therefore useful for researchers to know correct models for an application of interest. In this light, LuMiaoSu2018 propose methods of model selection in three-dimensional panel data. In this paper, we advance this literature by proposing a method of post-selection inference for regression parameters. We propose to use the lasso technique as means of model selection and to de-bias the lasso estimate, but our assumptions allow for many and even all fixed effects to be nonzero. Simulation studies demonstrate that the proposed method is more precise than biased estimators by parsimonious models, is more efficient than noisy estimators by redundant models, and allows for as accurate inference as the oracle estimator.

We suggest a couple of directions for future research. First, our model framework does not allow for $ij$ fixed effects, while $i$, $j$, $t$, $it$ and $jt$ fixed effects are allowed. Although allowing for $ij$ fixed effects is not of interest in our motivating example,\footnote{In gravity models for international trade, the main parameters of interest are the coefficient of $DIST_{ij}$, interpreted as the trade elasticity or trade cost, and the coefficient of $TA_{ij}$, interpreted as the effects of bilateral trade agreements on trade volume. These parameters will not be identified once $ij$ fixed effects enter the model.} it may be possible to allow for such fixed effects provided that the asymptotic setting allows for large $T$ as well as large $N$ and/or large $M$. Formal theoretical development for this case is left for future research. Second, we conjecture that our limit distribution result can be extended to establish honest (uniformly valid) confidence intervals, and formal theoretical investigation of the honesty property is left for future research.