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.
86,745 characters · 13 sections · 105 citation commands
Linear Regression for Panel with Unknown Number of Factors as Interactive Fixed Effects
\abstract{
}
{\bf Keywords:} Panel data, interactive fixed effects, factor models, perturbation theory of linear operators, random matrix theory.\\[4pt] {\bf JEL-Classification:} C23, C33
\linespread{1.3}
Panel data models typically incorporate individual and time effects to control for heterogeneity in cross-section and over time. While often these individual and time effects enter the model additively, they can also be interacted multiplicatively, thus giving rise to so called interactive effects, which we also refer to as a factor structure. The multiplicative form captures the heterogeneity in the data more flexibly, since it allows for common time-varying shocks (factors) to affect the cross-sectional units with individual specific sensitivities (factor loadings).\footnote{ The conventional additive model can be interpreted as a two factor interactive fixed effects model.} It is this flexibility that motivated the discussion of interactive effects in the econometrics literature, e.g. Holtz-Eakin, Newey and Rosen HoltzEakin-Newey-Rosen1988, Ahn, Lee and Schmidt AhnLeeSchmidt2001,AhnLeeSchmidt2013, Pesaran Pesaran2006, Bai Bai2009,Bai2013likelihood, Zaffaroni Zaffaroni2009, Moon and Weidner MoonWeidner2013, and Lu and Su LuSu2013.
Let $N$ be the number of cross-sectional units, $T$ be the number of time periods, $K$ be the number of regressors, and $R^0$ be the true number of interactive fixed effects. We consider a linear regression model with observed outcomes $Y$, regressors $X_k$, and unobserved error structure $\varepsilon$, namely
where $Y$, $X_k$, $\varepsilon$ and $e$ are $N\times T$ matrices, $\lambda^0$ is an $N \times R^0$ matrix, $f^0$ is a $T \times R^0$ matrix, and the regression parameters $\beta^0_k$ are scalars --- the superscript zero indicates the true value of the parameters. We write $\beta$ for the $K$-vector of regression parameters, and we denote the components of the different matrices by $Y_{it}$, $X_{k,it}$, $e_{it}$, $\lambda^0_{ir}$ and $f^0_{tr}$, where $i = 1, \ldots, N$, $t=1, \ldots, T$, and $r=1,\ldots,R^0$. It is convenient to introduce the notation $\beta \cdot X := \sum_{k=1}^{K} \, \beta_{k} \, X_k$. All matrices, vectors and scalars in this paper are real valued.
We consider the interactive fixed effect specification, i.e. we treat $\lambda^0$ and $f^0$ as nuisance parameters, which are estimated jointly with the parameters of interest $\beta$.\footnote{ When we refer to interactive fixed effects we mean that both factors and factor loadings are treated as non-random parameters. Ahn, Lee and Schmidt AhnLeeSchmidt2001 take a hybrid approach in that they treat the factors as non-random, but the factor loadings as random. The common correlated effects estimator of Pesaran Pesaran2006 was introduced in a context, where both the factor loadings and the factors follow certain probability laws, but it exhibits many properties of a fixed effects estimator.} The advantages of the fixed effects approach are for instance that it is semi-parametric, since no assumption on the distribution of the interactive effects needs to be made, and that the regressors can be arbitrarily correlated with the interactive effect parameters.
We study the least squares (LS) estimator of model (ref), which minimizes the sum of squared residuals to estimate the unknown parameters $\beta$, $\lambda$ and $f$.\footnote{ The LS estimator is sometimes called “concentrated” least squares estimator in the literature, and in an earlier version of the paper we referred to it as the “Gaussian Quasi Maximum Likelihood Estimator”, since LS estimation is equivalent to maximizing a conditional Gaussian likelihood function. Note also that for fixed $\beta$ the LS estimator for $\lambda$ and $f$ is simply the principal components estimator. } To our knowledge, this estimator was first discussed in Kiefer Kiefer1980. Under an asymptotic where $N$ and $T$ grow to infinity, the asymptotic properties of the LS estimator were derived in Bai Bai2009 for strictly exogeneous regressors, and extended in Moon and Weidner MoonWeidner2013 to the case of pre-determined regressors.
An important restriction of these papers is that the number of factors $R^0$ is assumed to be known. However, in many empirical applications there is no consensus about the exact number of factors in the data or in the relevant economic model. If $R^0$ is not known beforehand, then it may be estimated consistently,\footnote{See the discussion in Bai Bai2009supp, supplemental material, regarding estimation of $R^0$.} but difficulties in obtaining reliable estimates for the number of factors are well-documented in the literature (see, e.g., the simulation results in Onatski Onatski2010, and also our empirical illustration in Section (ref)). Furthermore, in order to use the existing inference results on $R^0$ one still needs a good preliminary estimator for $\beta$, so that working out the asymptotic properties of the LS-estimator for $R \geq R^0$ is still useful when taking that route.
We investigate the asymptotic properties of the LS estimator when the true number of factors $R^0$ is unknown and $R$ ($\geq R^0$) number of factors are used in the estimation.\footnote{ For $R<R^0$ the LS estimator can be inconsistent, since then there are interactive fixed effects in the model which can be correlated with the regressors but are not controlled for in the estimation. We therefore restrict attention to the case $R \geq R^0$. } We denote this estimator by $\widehat{\beta}_{R}$.
The main result of the paper, presented in Section (ref), is that under certain assumptions the LS estimator $\widehat{\beta}_{R}$ has the same limiting distribution as $\widehat{\beta}_{R^0}$ for any $R \geq R^0$ under an asymptotic where both $N$ and $T$ become large, while $R^0$ and $R$ are constant. This implies that the LS estimator $\widehat{\beta}_{R}$ is asymptotically robust towards inclusion of extra interactive effects in the model, and within the LS estimation framework there is no asymptotic efficiency loss from choosing $R$ larger than $R^0$. The important empirical implication of our result is that the number of factors $R^0$ need not be known or estimated accurately to apply the LS estimator.
To derive this robustness result, we impose more restrictive conditions than those typically assumed with known $R^0$. These include that the errors $e_{it}$ are independent and identically (iid) normally distributed and that the regressors are composed of a “low-rank” strictly stationary component, a “high-rank” strictly stationary component, and a “high-rank” pre-determined component.\footnote{ The pre-determined component of the regressors allows for linear feedback of $e_{it}$ into future realizations of $X_{k,it}$. } Notice that while some of these restrictions are necessary for our robustness result, some of them (e.g. iid normality of $e_{it}$) are imposed for technical reasons, because in the proof we use certain results from the theory of random matrices that are currently only available in that case (see the discussion in Section (ref)). In the Monte Carlo simulations in Section (ref), we consider DGPs that violate some technical conditions to demonstrate robustness of the result.
Under less restrictive assumptions we provide intermediate results that sequentially lead to the main result in Section (ref) and Appendix (ref) and (ref). In Section (ref) we show $\sqrt{\min(N,T)}$-consistency of the LS estimator $\widehat{\beta}_{R}$ as $N,T \rightarrow \infty$ under very mild regularity condition on $X_{it}$ and $e_{it}$, and without imposing any assumptions on $\lambda^0$ and $f^0$ apart from $R \geq R^0$. We thus obtain consistency of the LS estimator not only for unknown number of factors, but also for weak factors,\footnote{ See Onatski Onatski2010,Onatski2012 and Chudik, Pesaran and Tosetti ChudikPesaranTosetti2011 for a discussion of “strong” vs. “weak” factors in factor models.} which is an important robustness result.
In Section (ref) we derive an asymptotic expansion of the LS profile objective function that concentrates out $f$ and $\lambda$, for the case $R=R^0$. Given that the profile objective function is a sum of eigenvalues of a covariance matrix, its quadratic approximation is challenging because the derivatives of the eigenvalues with respect to $\beta$ are not generally known. We thus cannot use a conventional Taylor expansion, but instead apply the perturbation theory of linear operators to derive the approximation.
In Section (ref) we provide an example that satisfies the typical assumptions imposed with known $R^0$, so that $\widehat{\beta}_{R^0}$ is $\sqrt{NT}$ consistent, but we show that $\widehat{\beta}_{R}$ with $R > R^0$ is only $\sqrt{\min(N,T)}$ consistent in that example. This shows that stronger conditions are required to derive our main result.
In Appendix (ref) we show faster than $\sqrt{\min(N,T)}$-convergence of $\widehat{\beta}_{R}$ under assumptions that are less restrictive than those employed for the main result, in particular allowing for either cross-sectional or time-serial correlation of the errors $e_{it}$. In Appendix (ref) we provide an alternative version of our main result of asymptotic equivalence of $\widehat{\beta}_{R^0}$ and $\widehat{\beta}_{R}$, $R \geq R^0$, which is derived under high-level assumptions.
In Section (ref) we follow Kim and Oka KimOka2014 in employing the interactive fixed effects specification to study the effect of US divorce law reforms on divorce rates. This empirical example illustrates that the estimates for the coefficient $\beta$ indeed become insensitive to the choice of $R$, once $R$ is chosen sufficiently large, as expected from our theoretical results.
Section (ref) contains Monte Carlo simulation results for a static panel model. For the simulations we consider a DGP that violates the iid normality restriction of the error term. The simulation results confirm our main result of the paper even with a relatively small sample size (e.g. $N=100$, $T=10$) and non-iid-normal errors. In the supplementary appendix, we report the Monte Carlo simulation results of an AR(1) panel model. It also confirms the robustness result in large samples, but in finite samples it shows more inefficiency than the static case. In general, one should expect some finite sample inefficiency from overestimating the number of factors when the sample size is small or the number of overfitted factors is large.
A few words on notation. The transpose of a matrix $A$ is denoted by $A'$. For a column vectors $v$ its Euclidean norm is defined by $\| v \| = \sqrt{v^{\prime}v}$ . For an $m\times n$ matrix $A$ the Frobenius or Hilbert Schmidt norm is $\| A \|_{HS} = \sqrt{{\rm Tr} (AA^{\prime})}$, and the operator or spectral norm is $\| A \| = \max_{0 \neq v \in \mathbb{R}^n} \, \frac{ \| A v \|} {\| v\|}$. Furthermore, we use $P_A = A (A^{\prime}A)^{\dagger} A'$ and $M_A = \mathbbm{1} - A (A^{\prime}A)^{\dagger} A'$, where $\mathbbm{1}$ is the $m\times m$ identity matrix, and $(A^{\prime}A)^{\dagger}$ denotes some generalized inverse, in case $A$ is not of full column rank. For square matrices $B$, $C$, we use $B>C$ (or $B\geq C$) to indicate that $B-C$ is positive (semi) definite. We use “wpa1” for “with probability approaching one”.
In this section we provide a set of conditions under which the regression coefficient $\beta^0$, the interactive fixed effects $\lambda^0 f^{0 \prime}$, and the number of factors $R^0$ are determined uniquely by the data. Here, and throughout the whole paper, we treat $\lambda$ and $f$ as non-random parameters, i.e. all stochastics in the following are implicitly conditional on $\lambda$ and $f$. Let $x_k={\rm vec}(X_k)$, the $NT$-vectorization of $X_k$, and let $x=(x_1,\ldots,x_K)$, which is an $NT \times K$ matrix.
Assumption (ref)$(i)$ imposes existence of second moments. Assumption (ref)$(ii)$ is an exogeneity condition, which demands that $x_{it}$ and $e_{it}$ are not correlated contemporaneously, but allows for pre-determined regressors like lagged dependent variables. Assumption (ref)$(iv)$ imposes that the true number of factors $R^0 := {\rm rank}(\lambda^0 f^{0 \prime})$ is bounded by a non-negative integer $R$, which cannot be too large (e.g. the trivial bound $R=\min(N,T)$ is not possible), since otherwise Assumption (ref)$(iii)$ cannot be satisfied.
Assumption (ref)$(iii)$ is a non-collinearity condition, which demands that the regressors have significant variation across $i$ and over $t$ after projecting out all variation that can be explained by the factor loadings $\lambda^0$ and by arbitrary factors $F \in \mathbbm{R}^{T \times R}$. This generalizes the within variation assumption in the conventional panel regression with time invariant individual fixed effects, which in our notation reads $\mathbbm{E}[x' (M_{1_T} \otimes \mathbbm{1}_N) x] >0$.\footnote{ The conventional panel regression with additive individual fixed effects and time effects requires a non-collinearity condition of the form $\mathbbm{E}[x' (M_{1_T} \otimes M_{1_N}) x] >0$.} This conventional fixed effect assumption rules out time-invariant regressors. Similarly, Assumption (ref)$(iii)$ rules out more general “low-rank regressors”,\footnote{ We do not consider such “low-rank regressors” in this paper. Note also that Assumption A in Bai Bai2009 is the sample version of our Assumption (ref)$(iii)$.} see our discussion of Assumption (ref) below.
The estimator we investigate in this paper is the least squares (LS) estimator, which for a given choice of $R$ reads\footnote{ The optimal $\widehat \Lambda_R$ and $\widehat F_R$ in (ref) are not unique, since the objective function is invariant under right-multiplication of $\Lambda$ with a non-degenerate $R \times R$ matrix $S$, and simultaneous right-multiplication of $F$ with $(S^{-1})'$. However, the column spaces of $\widehat \Lambda_R$ and $\widehat F_R$ are uniquely determined. }
where $\|.\|_{HS}$ refers to the Hilbert Schmidt norm, also called Frobenius norm. The objective function $ \left\| Y \, - \, \beta \cdot X \, - \, \Lambda \, F' \right\|^2_{HS}$ is simply the sum of squared residuals. The estimator for $\beta^0$ can equivalently be defined by minimizing the profile objective function that concentrates out the $R$ factors and the $R$ factor loadings, namely
with\footnote{ The profile objective function ${\cal L}_{NT}^R(\beta)$ need not be convex in $\beta$ and can have multiple local minima. Depending on the dimension of $\beta$ one should either perform an initial grid search or try multiple starting values for the optimization when calculating the global minimum $\widehat \beta_R$ numerically. See also Section (ref) of the supplementary material. }
where, $\mu_r(.)$ is the $r$'th largest eigenvalue of the matrix argument. Here, we first concentrated out $\Lambda$ by use of its own first order condition. The resulting optimization problem for $F$ is a principal components problem, so that the the optimal $F$ is given by the $R$ largest principal components of the $T \times T$ matrix $\left(Y - \beta \cdot X\right)'\left(Y - \beta \cdot X\right)$. At the optimum the projector $M_F$ therefore exactly projects out the $R$ largest eigenvalues of this matrix, which gives rise to the final formulation of the profile objective function as the sum over its $T-R$ smallest eigenvalues.\footnote{ This last formulation of ${\cal L}_{NT}^R(\beta)$ is very convenient since it does not involve any explicit optimization over nuisance parameters. Numerical calculation of eigenvalues is very fast, so that the numerical evaluation of ${\cal L}_{NT}^R(\beta)$ is unproblematic for moderately large values of $T$. Since the model is symmetric under $N \leftrightarrow T$, $\Lambda \leftrightarrow F$, $Y \leftrightarrow Y'$, $X_k \leftrightarrow X_k'$ there also exists a dual formulation of ${\cal L}_{NT}^R(\beta)$ that involves solving an eigenvalue problem for an $N \times N$ matrix.} We write ${\cal L}_{NT}^0(\beta)$ for ${\cal L}_{NT}^{R^0}(\beta)$, the profile objective function obtained for the true number of factors. Notice that in (ref) the parameter set for $\beta$ is the whole Euclidean space $\mathbbm{R}^K$ and we do not restrict the parameter set to be compact.
\paragraph{Remarks}
Theorem (ref) follows from Theorem (ref) and Lemma (ref) in the appendix, whose proof is given in the supplementary material. The theorem guarantees that the asymptotic distribution of $\widehat \beta_{R}$, $R \geq R^0$, is identical to that of $\widehat \beta_{R^0}$ in (ref) below.
The limiting distribution of $\sqrt{NT}\big(\widehat \beta_{R^0} - \beta^0\big)$ with known $R^0$ is available in the existing literature. According to Bai Bai2009 and Moon and Weidner MoonWeidner2013,
where $W$ is the $K \times K$ matrix with elements $W_{k_1 k_2} = \frac 1 {NT} {\rm Tr}( M_{\lambda^0} X_{k_1} M_{f^0} X_{k_2}' )$ and $B$ is the $K$-vector with elements $B_k = \frac 1 N {\rm Tr}[ P_{f^0} \mathbbm{E}(e' X_k)]$.\footnote{ The asymptotic distribution in (ref) can also be derived from Corollary (ref) below under more general conditions than in Assumption (ref) (see Moon and Weidner MoonWeidner2013 for details). Here we have used the homoscedasticity of $e_{it}$ to simplify the structure of the asymptotic variance and bias. Bai Bai2009 finds further asymptotic bias in $\widehat \beta_{R^0}$ due to heteroscedasticity and correlation in $e_{it}$, which in our asymptotic result is ruled out by Assumption (ref)$(ii)$, but is studied in our Monte Carlo simulations below. Moon and Weidner MoonWeidner2013 work out the additional asymptotic bias in $\widehat \beta_{R^0}$ due to pre-determined regressors, which is allowed for in Theorem (ref). }
The result (ref) holds under the assumptions of Theorem (ref) and also assuming that $ \operatorname*{plim} W^{-1} B$ and $\operatorname*{plim} W^{-1}$ exist, where $\operatorname*{plim}$ refers to the probability limit as $N,T \rightarrow \infty$. Note that Assumption (ref) guarantees that $W$ is invertible asymptotically. The asymptotic bias in (ref) is an incidental parameter bias due to pre-determined regressors and is equal to zero for strictly exogenous regressors (for which $\mathbbm{E}(e' X_k)=0$); it generalizes the well-known Nickell Nickell1981 bias of the within-group estimator for dynamic panel models.
Estimators for $\sigma^2$, $W$ and $B$ are given by\footnote{The first factor in $ \widehat \sigma^2$ reflects the degree of freedom correction from estimating $\Lambda$, $F$ and $\beta$, but could simply be chosen as $1/NT$ for the purpose of consistency. Note also that $P_{\widehat F_R,t \tau} = {\cal O}_P(1/T)$, which explains why no $1/T$ factor is required in the definition of $\widehat B_{R,k}$.}
where $\widehat e_{R,it}$ denotes the $(i,t)^{th}$ element of $\widehat e_R = Y - \widehat \beta_R \cdot X- \widehat \Lambda_R \widehat F_R'$, and $P_{\widehat F_R,t \tau}$ denotes the $(t,\tau)^{th}$ element of $P_{\widehat F_R} = \mathbbm{1}_T - M_{\widehat F_R} = \widehat F_R (\widehat F_R' \widehat F_R)^\dagger \widehat F_R' $, and $M \in \{1,2,3,\ldots\}$ is a bandwidth parameter that also depends on the sample size $N,T$. Let $\widehat W_R$ and $\widehat B_R$ be the matrix and vector with elements $ \widehat W_{R,k_1 k_2}$ and $\widehat B_{R,k}$, respectively.
The next theorem establishes the consistency of these estimators. Let $\lambda^{\rm red} \in \mathbbm{R}^{N \times (R-R^0)}$ and $f^{\rm red} \in \mathbbm{R}^{T \times (R-R^0)}$ be the leading $R-R^0$ principal components obtained from the $N \times T$ matrix $M_{\lambda^0} e M_{f^0}$, i.e. $\lambda^{\rm red}$ and $f^{\rm red}$ minimize the objective function $\left\| M_{\lambda^0} e M_{f^0} - \lambda^{\rm red} \, f^{\rm red \prime} \right\|^2_{HS}$, analogous to $\widehat \Lambda_R$ and $\widehat F_R$ defined in (ref).\footnote{The superscript “red” stands for redundant, because it turns out that $\lambda^{\rm red}$ and $f^{\rm red}$ are asymptotically close to the $R-R^0$ redundant principal components that are estimated in (ref). }
Combining Theorems (ref) and (ref) and the asymptotic distribution in (ref) allows inference on $\beta$, for $R \geq R^0$. In particular, the bias corrected estimator $\widehat \beta^{\rm BC}_{R} = \widehat \beta_{R} + \frac 1 T \widehat W_R^{-1} \widehat B_R$ satisfies\footnote {Instead of estimating the bias analytically one can use the result that the bias is of order $T^{-1}$ and perform split panel bias correction as in Dhaene and Jochmans DhaeneJochmans2010, which instead of the conditions of Theorem (ref)(ii) only requires some stationary condition over time.} \[ \sqrt{NT}\big(\widehat \beta^{\rm BC}_{R} - \beta^0\big) \Rightarrow {\cal N}( 0 , \sigma^2 W^{-1} ). \]
\paragraph{Heuristic Discussion of the Main Result} \\ Intuitively, the inclusion of unnecessary factors in the LS estimation is similar to the inclusion of irrelevant regressors in an OLS regression. In the OLS case it is well known that if those irrelevant extra regressors are uncorrelated with the regressors of interest, then they have no effect on the asymptotic distribution of the regression coefficients of interest. It is therefore natural to expect that if the extra estimated factors in $\widehat F_R$ are asymptotically uncorrelated with the regressors, then the result of Theorem (ref) should hold. To explore this, remember that $\widehat F_R$ is given by the first $R$ principal components of the matrix $(Y-\widehat \beta_R \cdot X)' (Y-\widehat \beta_R \cdot X)$, and write
The strong factor assumption and the consistency of $\widehat \beta_R$ guarantee that the first $R^0$ principal components of $(Y-\widehat \beta_R \cdot X)' (Y-\widehat \beta_R \cdot X)$ are close to $f^0$ asymptotically, i.e. the true factors are correctly picked up by the principal component estimator. The additional $R-R^0$ principal components that are estimated for $R>R^0$ cannot pick up anymore true factors and are thus mostly determined by the remaining term $e - (\widehat \beta_R - \beta^0) \cdot X$. The key question for the properties of the extra estimated factors, and thus of $\widehat \beta_R$, is therefore whether the principal components obtained from $e - (\widehat \beta_R - \beta^0) \cdot X$ are dominated by $e$ or by $(\widehat \beta_R - \beta^0) \cdot X$. Only if they are dominated by $e$ can we expect the extra factors in $\widehat F_R$ to be uncorrelated with $X$ and thus the result in Theorem (ref) to hold. The result on $P_{\widehat F_R}$ in Theorem (ref) shows that the additional estimated factors are indeed close to $f^{\rm red}$, i.e. are mostly determined by $e$, but this result is far from obvious a priori, as the following discussion shows.
Under our assumptions we have $\| e\| = {\cal O}_P(\sqrt{N})$ and $\| X_k \| = {\cal O}_P( \sqrt{NT} )$ as $N$ and $T$ grow at the same rate. Thus, if the convergence rate of $\widehat \beta_R$ is faster than $\sqrt{N}$, i.e. $\| \widehat \beta_R - \beta^0 \| = o_P(\sqrt{N})$, then we have $\| e \| \gg \left\| (\widehat \beta_R - \beta^0) \cdot X \right\|$ asymptotically, and we expect the extra $\widehat F_R$ to be dominated by $e$. A crucial step in the derivation of Theorem (ref) is therefore to show faster than $\sqrt{N}$ convergence of $\widehat \beta_R$. Conversely, we expect counter examples to the main result to be such that the convergence rate of the estimator $\widehat \beta_R$ is not faster than $\sqrt{N}$, and we provide such a counter example -- which, however, violates Assumptions (ref) -- in Section (ref) below. Whether the intuition about “inclusion of irrelevant regressors” carries over to the “inclusion of irrelevant factors” thus crucially depends on the convergence rate of $\widehat \beta_R$.
Here we introduce key intermediate results for the proof of the main Theorem (ref) stated above. These intermediate results may be useful independently of the main result, e.g. Moon and Weidner MoonWeidner2013 and Moon, Shum, and Weidner MoonShumWeidner2014 crucially use the results established in Section (ref) for the case of known $R=R^0$. The assumptions introduced below are all implied by the low-level Assumptions (ref) above, see to Lemma (ref) in the appendix.
Here we present a consistency result for $\widehat \beta_R$ under an arbitrary asymptotic $N,T \rightarrow \infty$, i.e. without the assumption that $N$ and $T$ grow at the same rate, which is imposed everywhere else in the paper. In addition to Assumption (ref) we require the following high level assumptions to obtain the result.
\paragraph{Remarks}
To derive the limiting distribution of $\widehat \beta_R$, we study the asymptotic properties of the profile objective function ${\cal L}_{NT}^R(\beta)$ around $\beta^0$. The expression in (ref) cannot easily be discussed by analytic means, since no explicit formula for the eigenvalues of a matrix is available. In particular, a standard Taylor expansion of ${\cal L}^R_{NT}(\beta)$ around $\beta^0$ cannot easily be derived. Here, we consider the case of known $R=R^0$ and we perform a joint expansion of the corresponding profile objective function ${\cal L}^0_{NT}(\beta)$ in the regression parameters $\beta$ and in the idiosyncratic error terms $e$. To perform this joint expansion we apply the perturbation theory of linear operators (e.g., Kato Kato). We thereby obtain an approximate quadratic expansion of ${\cal L}_{NT}^0(\beta)$ in $\beta$, which can be used to derive the first order asymptotic theory of the LS estimator $\widehat \beta_{R^0}$, see Appendix (ref) for details. In addition to the $K \times K$ matrix $W$ already defined in Section (ref) we now also define
Let $C^{(1)}$ and $C^{(2)}$ be the $K$-vectors with elements $C^{(1)}_k$ and $C^{(2)}_k$, respectively.
The bound on remainder\footnote{ The expansion in Theorem (ref) contains a term that is linear in $\beta$ and linear in $e$ ($C^{(1)}$ term), a term that is linear in $\beta$ and quadratic in $e$ ($C^{(2)}$ term), and a term that is quadratic in $\beta$ ($W$ term). All higher order terms of the expansion are contained in the remainder term ${\cal L}_{NT}^{0,{\rm rem}}(\beta)$.} in Theorem (ref) is such that it has no effect on the first order asymptotic theory of $\widehat \beta_{R^0}$, as stated in the following corollary (see also Andrews Andrews1999).
Note that our assumptions already guarantee $C^{(2)}={\cal O}_P(1)$ and that $W$ is invertible with $W^{-1}={\cal O}_P(1)$, so this need not be explicitly assumed in Corollary (ref).
\paragraph{Remarks}
The results in Bai Bai2009 and Corollary (ref) above show that under appropriate assumptions the estimator $\widehat \beta_R$ is $\sqrt{NT}$-consistent for $R=R^0$. For $R>R^0$ we know from Theorem (ref) that $\widehat \beta_R$ is $\sqrt{N}$ consistent as $N$ and $T$ grow at the same rate, but we have not shown faster than $\sqrt{N}$ converge of $\widehat \beta_R$ for $R>R^0$, yet, which according to the heuristic discussion at the end of Section (ref) is a very important intermediate step to obtain our main result.\footnote{ One reason why $\widehat \beta_R$ might only converge at $\sqrt{N}$ rate, but not faster, are weak factors (both for $R>R^0$ and for $R=R^0$). A weak factor (see e.g. Onatski Onatski2010,Onatski2012 and Chudik, Pesaran and Tosetti ChudikPesaranTosetti2011) might not be picked up at all or might only be estimated very inaccurately by the principal components estimator $\widehat F_R$, in which case that factor is not properly accounted for in the LS estimation procedure. If this happens and the weak factor is correlated with the regressors, then there is some uncorrected weak endogeneity problem, and $\widehat \beta_R$ will only converge at $\sqrt{N}$ rate. We do not consider the issue of weak factors any further in this paper. } However, one might not obtain a faster than $\sqrt{N}$ convergence rate of $\widehat \beta_R$ for $R>R^0$ without imposing further restrictions, as the following example shows.
The proof of the last statement is provided in the supplementary material. The DGP in this example satisfies all the assumptions imposed in Corollary (ref) to derive the limiting distribution of the LS-estimator for $R=R^0$, including $\sqrt{NT}$-consistency of $\widehat \beta_R$ for $R=R^0$ (=0 in this example). It also satisfies all the regularity conditions imposed in Bai Bai2009.\footnote{See Section (ref) in the supplementary material for details.} The aspect that is special about this DGP is that $\lambda_x$ and $f_x$ feature both in $X_{it}$ and in the second moment structure of $e_{it}$. The heuristic discussion at the end of Section (ref) provides some intuition why this can be problematic, because the leading principal components obtained from only the error matrix $e$ will have a strong sample correlation with $X_{it}$ for this DGP.
In Appendix (ref), we summarize our results on faster than $\sqrt{N}$ convergence of $\widehat \beta_R$ for $R \geq R^0$. The above example shows that this requires more restrictive assumptions than those imposed for the analysis of the case $R=R^0$ above, but the assumptions that we impose for this intermediate results are still significantly weaker than the Assumption (ref) required for our main result above, in particular either cross-sectional correlation or time-serial correlation of $e_{it}$ are still allowed.
In that appendix we also provide one set of assumptions (Assumption (ref)) for faster than $\sqrt{N}$ convergence such that no additional conditions on $e$ are required, but where the regressors are restricted to essentially be lagged dependent variables in an AR(p) model with factors.
We establish the asymptotic equivalence of $\widehat \beta_R$ and $\widehat \beta_{R^0}$ in Theorem (ref) by showing that the LS objective function $\mathcal{L}_{NT}^{R}(\beta)$ can, up to a constant, be uniformly well approximated by $\mathcal{L}_{NT}^{0}(\beta)$ in shrinking neighborhoods around the true parameter. For this, we need not only the faster than $\sqrt{N}$ convergence rate of $\widehat \beta_R$, but also require the Assumption (ref) in Appendix (ref). This is a high-level assumption on the eigenvalues and eigenvectors of the random covariance matrices $E E'$ and $E' E$, where $E=M_{\lambda^0} e M_{f^0}$. The assumption essentially requires the eigenvalues of those matrices to be sufficiently separated from each other, as well as the eigenvectors of those matrices to be sufficiently uncorrelated with the regressors $X_k$, and with $e P_{f^0}$ and $P_{\lambda^0} e$.
We use the iid normality of $e_{it}$ to verify those high-level conditions in Section (ref) of the supplementary appendix. There are three reasons why we can currently only verify those conditions for iid normal errors:
In spite of these severe mathematical challenges, we believe that in principle our high-level Assumption (ref) could be verified for more general error distributions, implying that our main result of asymptotic equivalence of $\widehat \beta_R$ and $\widehat \beta_{R^0}$ holds more generally. This is also supported by our Monte Carlo simulations, where we explore non-independent and non-normal error distributions.
As an illustrative empirical example, we estimate the dynamic effects of unilateral divorce law reforms on the state-wise divorce rates in the US. The impact of the divorce law reform has been studied by many researches (e.g., Allen Allen1992, Peters Peters1986,Peters1992, Gray Gray1998, Friedberg Friedberg1998, Wolfers Wolfers2006, and Kim and Oka KimOka2014). In this section we revisit this topic, extending Wolfers Wolfers2006 and Kim and Oka KimOka2014 by controlling for interactive fixed effects and also a lagged dependent variable.
Let $Y_{it}$ denote the number of divorces per 1000 people in state $i$ at time $t$, and let $D_{i}$ denote the year in which state $i$ introduced the unilateral divorce law, i.e. before year $D_i$ state $i$ had a consent divorce law, while from $D_i$ onwards state $i$ had a unilateral “no-fault” divorce law, which loweres the barrier for divorce. The goal is to estimate the dynamic effects of this law change on the divorce rate. The empirical model we estimate is
where we follow Wolfers Wolfers2006 in defining the regressors as bi-annual dummies:
The dummy variable and quadratic trend specification $\alpha_{i}+ \gamma_i \, t + \delta_i \, t^{2} + \mu_t$ is also used in Friedberg Friedberg1998 and Wolfers Wolfers2006. The additional interactive fixed effects $\lambda _{i}^{\prime }f_{t}$ were added in Kim and Oka KimOka2014 to control for additional unobserved heterogeneity in the divorce rate, e.g. due to social, cultural or demographic factors. We extend the specification further by adding a lagged dependent variable $Y_{i,t-1}$ to control for state dependence of the divorce rate, but we also report results without $Y_{i,t-1}$ below. We use the dataset of Kim and Oka KimOka2014,\footnote{ The data is available from {\tt http://qed.econ.queensu.ca/jae/2014-v29.2/kim-oka/}} which is a balanced panel of $N=48$ states over $T=33$ years, leaving $T=32$ time periods if the lagged dependent variable is included.
For estimation we first eliminate $\alpha_{i}$, $\gamma_i$, $\delta_i$ and $\mu_t$ from the model by projecting the outcome variable and all regressors accordingly, e.g. $\widetilde Y = M_{1_N} Y M_{(1_T, {\bf t}, {\bf t}^2)}$, where $1_N$ and $1_T$ are $N$- and $T$-vectors, respectively, with all entries equal to one, and ${\bf t}$ and ${\bf t}^2$ are $T$-vectors with entries $t$ and $t^2$, respectively. The model after projection reads $\widetilde Y_{it}=\beta_{0} \, \widetilde Y_{i,t-1} + \sum_{k=1}^{8} \beta _{k} \widetilde X_{k,it} +\widetilde \lambda _{i}^{\prime } \widetilde f_{t} + \widetilde e_{it}$, which is exactly the model we have studied so far in this paper.\footnote{To construct $\widetilde Y_{i,t-1}$ we first apply the lag-operator and then apply the projections $M_{1_N}$ and $M_{(1_T, {\bf t}, {\bf t}^2)}$.} We use the LS estimator described above to estimate this model. The projection reduces the effective sample size to $N=48-1=47$ and $T=32-3=29$, which should be accounted for when calculating standard errors, e.g. in the formula for $\widehat \sigma_R^2$ above (degree of freedom correction). Our theoretical results are still applicable.\footnote{ If $e_{it}$ is iid normal, then $\widetilde e_{it}$ is not, but one can apply appropriate orthogonal rotations in $N$- and $T$-space such that $\widetilde e_{it}$ becomes iid normal again, although with sample size reduced to $N=47$ and $T=29$. The rotation has no effect on the LS estimator, i.e. it does not matter whether we work in the original or the rotated frame. }
We need to decide on a number of factors $R$ when implementing the LS estimator. As already mentioned in the last remark in Section (ref) above, we can can apply known techniques from the literature on factor models without regressors to obtain a consistent estimator of $R^0$. To do so we choose a maximum number of factors of $R_{\max}=9$ to obtain the preliminary estimate $\widehat \beta_{R_{\max}}$ and then calculate the residuals $\widehat u_{it} = \widetilde Y_{it} - \widehat \beta_{R_{\max},0} \, \widetilde Y_{i,t-1} - \sum_{k=1}^{8} \widehat \beta _{R_{\max},k} \widetilde X_{k,it}$. We then apply the IC, PC and BIC3 criteria of Bai and Ng BaiNg2002,\footnote{ Following Onatski Onatski2010 and Ahn and Horenstein AhnHorenstein2013 we report only BIC3 among the AIC and BIC criteria of Bai and Ng BaiNg2002. } the criterion described in Onatski Onatski2010, and the ER and GR criteria of Ahn and Horenstein AhnHorenstein2013 to $\widehat u$.\footnote{ To include $R=0$ as a possible outcome for the Ahn and Horenstein (2013) criterion, we use the mock eigenvalue used in their simulations.} Most of these criteria also require specification of $R_{\max}$, and we continue to use $R_{\max} = 9$. The corresponding estimation results for $R$ are presented in Table (ref). In addition, we also report the log scree plot, i.e. the sorted eigenvalues of $\widehat u' \widehat u$ in Figure 1.
The log scree plot already shows that it is not obvious how to decompose the eigenvalue spectrum into a few larger eigenvalues stemming from factors and the remaining smaller eigenvalues stemming from the idiosyncratic error term.\footnote{The first largest eigenvalue is 2.2 times larger than the second eigenvalue, the second is 1.6 times larger than the third, the third is 1.9 times larger than fourth. So the largest view eigenvalues are larger than the remaining ones, and the strong factor assumption might not be completely inappropriate here. However, deciding on a cutoff between factor and non-factor eigenvalues is difficult.} This problem is also reflected in the very different estimates for $R$ that one obtains from the various criteria. It might appear that IC1, IC3, PC1, PC2 and PC3 all agree on $\widehat R=9$, but this is simply $\widehat R=R_{\max}$, and if we choose $R_{\max}=10$, then all these criteria deliver $\widehat R=10$, so this should not be considered a reliable estimate.
On the other hand, our asymptotic theory suggests, that the exact choice of $R$ in the estimation of $\widehat \beta_R$ should not matter too much, as long as $R$ is chosen large enough to cover all relevant factors. Table (ref) contains the estimation results for the bias corrected $\widehat \beta_R$ for $R \in \{0,1,\ldots,9\}$. Table (ref) contains estimates if the lagged dependent variable is not included into the model.\footnote{ The result for $R=7$ in Table (ref) should be equal to column (6) in Table III of Kim and Oka KimOka2014. The discrepancy is explained by a coding error in their bias computation. Note also that the result for $R=0$ in Table (ref) does not match the one in Wolfers Wolfers2006, because he uses WLS with state population weights, while we use OLS for simplicity. Kim and Oka KimOka2014 estimate both WLS and OLS and find that the difference between the resulting estimates becomes insignificant, once a sufficient number of interactive fixed effects is controlled for.} For all reported estimates we perform bias correction and standard error estimation as described in Bai Bai2009 and Moon and Weidner MoonWeidner2013.\footnote{We correct for the biases due to heterscedasticity in both panel dimensions worked out in Bai Bai2009, as well as for the dynamic bias worked out in Moon and Weidner MoonWeidner2013. For the latter we use the formula for $ \widehat B_{R,k}$ above, with bandwidth $M=2$. For the standard error estimation we allow for heterscedasticity in both panel dimensions, also following Bai Bai2009 and Moon and Weidner MoonWeidner2013. The bias and standard error formulas in those paper assume $R=R^0$ known, but we strongly expect that those formulas are robust towards $R>R^0$, as partly justified by Theorem (ref) above. For the model without lagged dependent variable we also allow for serial correlation in $e_{it}$ when estimating the bias and standard deviation of $\widehat \beta_R$. }
When ignoring the lagged dependent variable coefficient, one finds that in both Table (ref) and Table (ref) the estimation results for $\widehat \beta_R$ and the corresponding t-values are quite sensitive to changes in $R$ for very small values of $R$, but become much more stable as $R$ increases, and actually do not change too much anymore from roughly $R=2$ onwards. These findings are very well in line with our asymptotic theory, and the dynamic effect of divorce law reform that we find are also similar to the findings in Wolfers Wolfers2006 and Kim and Oka KimOka2014. The effect of the law reform on the divorce rates initially increases over time, is certainly significant in year 3-4 after the reform, and declines and becomes insignificant afterwards.\footnote{ The magnitude of the estimates is smaller than those in Wolfers Wolfers2006, i.e. controlling for unobserved factors reduced the effect size, as already pointed out by Kim and Oka KimOka2014.}
In contrast, the estimated coefficient on the lagged dependent variable in Table (ref) is quite large and highly significant for small values of $R$, but decreases steadily with $R$, until it gets close to zero and insignificant for $R \geq 8$. A plausible interpretation of this finding is that the model that includes the lagged dependent variable is misspecified, and that the estimated value of $\beta_0$ for small values of $R$ does not correspond to a true state dependence of $Y_{it}$, but simply reflects the time-serial correlation of the error process being picked up by the autoregressive model. According to this interpretation, once we include more and more factors into the model we control for more and more serial dependence of the unobserved error term, thus uncovering the true insignificance of $\beta_0$ in the estimates for $R \geq 8$.
This empirical example shows that instead of relying on a single estimate $\widehat R$ for the number of factors and reporting the corresponding $\widehat \beta_{\widehat R}$ it can be very informative to calculate $\widehat \beta_{R}$ for multiple values of $R$. Whether the estimated coefficients become stable for sufficiently large $R$ values, as our asymptotic theory suggests, is a useful robustness check for the model. When reporting the final results, then, it is better, within a reasonable range, to choose an $R$ that is too large than one that is too small.
We also perform a Monte Carlo simulation that is tailored towards the empirical application. For this we use the static model without lagged dependent variable. To generate $Y_{it}$ in equation (ref) with $\beta_0=0$ we use the observed regressors $X_{k,it}$, as described above, and as true parameters we use the $\beta$ (bias corrected), $\alpha_i$, $\gamma_i$, $\delta_i$, $\mu_t$, $\lambda_i$ and $f_t$ obtained from the estimation with $R=4$ (i.e. $\beta_k$ as reported in the $R=4$ column of Table (ref)). We generate $e_{it}$ from an MA(1) model with $t(5)$ distributed innovations. Note that this error distribution violates the assumption (ref)$(ii)$.
In this “empirical Monte Carlo” we have $N=48$, $T=33$ and true number of factors $R^0=4$. We find that the bias corrected estimates for $\beta_k$, $k=1,\ldots,8$, are essentially unbiased when $R \geq R^0$ factors are used in the estimation, but for $R<R^0$ the coefficient estimates are often biased. For $\beta_k$, $k \geq 3$, there are only small changes in the standard deviation of the estimator between $R=4$ and $R=9$, but for $\beta_k$, $k=1,2$, we observe standard deviation inflation of up to $25 \%$ between $R=4$ and $R=9$. Given the relatively small sample size the difference between $R=9$ and $R^0=4$ is relatively large, and some finite sample inefficiency is not too surprising. The detailed results are available in the supplementary appendix.
In addition to the “empirical Monte Carlo” discussed above we now investigate the finite sample properties of $\widehat \beta_{R}$ and $\widehat \beta_{R}^{\rm BC}$ further. In the simulations in this section we use a generated regressor $X_{it}$ that is correlated with the interactive fixed effects. The serial correlation of the error term $e_{it}$ together with the data generating process (DGP) for $X_{it}$, $\lambda_i$ and $f_t$ are such that the naive LS estimator has an asymptotic bias. This allows to verify whether the bias is essentially unchanged for $R>R^0$ and whether bias correction works well for $R>R^0$ in finite sample. We also study various combinations of $N$ and $T$.
The model is a static panel model with one regressor ($K=1$), two factors ($R^0=2$), and the following DGP:
The random variables $\widetilde X_{it}$, $\lambda_{ir}$, $f_{tr}$, $\chi_{ir}$ and $v_{it}$ are mutually independent; with $\widetilde X_{it}$ and $f_{tr} \, \sim \, iid \, {\cal N}(0,1)$; $\lambda_{ir}$ and $\chi_{ir} \, \sim iid \, {\cal N}(1,1)$; and $v_{it} \, \sim \, iid \, t(5)$, i.e. $v_{it}$ has a Student's t-distribution with 5 degrees of freedom.
Note that this model satisfies Assumptions (ref), (ref), and (ref)(i), but not (ref)(ii). The error term $e_{it}$ is not distributed as $iid$ normal. The time series of $e_{it}$ follows an MA(1) process with innovations distributed as $t(5)$.
We choose $\beta^0=1$, and use $10,000$ repetitions in our simulation. The true number of factors is chosen to be $R^0=2$. For each draw of $Y$ and $X$ we compute the LS estimator $\widehat \beta_R$ according to equation (ref) for different values of $R$, namely $R \in \{0,1,2,3,4,5\}$.
Table (ref) reports bias and standard deviation of the estimator $\widehat \beta_R$ for different combinations of $R$, $N$ and $T$. For $R<R^0=2$ the model is misspecified and $\widehat \beta_R$ turns out to be severely biased. There is also bias in $\widehat \beta_R$ for $R \geq R^0$, due to time-serial correlation of $e_{it}$. This bias was worked out in Bai Bai2009, and bias correction is also discussed there.
Table (ref) reports various quantiles of the distribution of $\sqrt{NT}( \widehat \beta_R - \beta^0 )$ for $N=T=100$ and $N=T=300$, and different values of $R \geq R^0$. From these tables, we see that as $N,T$ increases the distribution of $\widehat \beta_R$ gets closer to that of $\widehat\beta_{R^0}$.
Table (ref) reports the size of a t-test with nominal size equal to $5 \%$ for $R \geq R^0$. We use the results in Bai Bai2009 to correct for the leading $1/N$ (not actually present in our DGP) and $1/T$ (present in our DGP) biases in $\widehat \beta_R$ before calculating the t-test statistics, allowing for heteroscedsticity in both panel dimensions and for time-serial correlation when estimating the bias and standard deviation of $\widehat \beta_R$. The finite sample size distortions are mostly due to residual bias after bias correction, but also partly due to some finite sample downward bias in the standard error estimates. The size distortions increase with $R$, but for all values of $R \geq R^0$ in Table (ref) the size distortions decrease rapidly as $T$ increases.
Monte Carlo Simulation results for an AR(1) model with factors can be found in Section (ref) of the supplementary material. Those additional simulations show that the finite sample properties (e.g. for $T=30$) of $\widehat \beta_{R^0}$ and $\widehat \beta_{R}$, $R>R^0$, can be quite different, but those differences vanish as $T$ becomes large, as predicted by our asymptotic theory. In general, we always expect some finite sample inefficiency from overestimating the number of factors.
We show that under certain assumptions the limiting distribution of the LS estimator of a linear panel regression with interactive fixed effects does not change when we include redundant factors in the estimation. The implication of this is that one can use an upper bound of the number of factors $R$ in the estimation without asymptotic efficiency loss. However, some finite sample efficiency loss from overestimating $R$ is likely, so that $R$ should not be chosen too large in actual applications. We impose $iid$ normality of the regression errors to derive the asymptotic result, because we require certain results on the eigenvalues and eigenvectors of random covariance matrices that are only known in that case. We expect that progress in the literature on large dimensional random covariance matrices will allow verification of our high-level assumptions under more general error distributions, and our Monte Carlo simulations suggest that the result also holds for non-normal and correlated errors. We also provide multiple intermediate asymptotic results under more general conditions.