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.
113,357 characters · 14 sections · 55 citation commands
Inference for Large Panel Data with Many Covariates
\onehalfspacing
\onehalfspacing
Our goal is the selection of a parsimonious sparse model from a large set of candidate covariates that explains a large dimensional panel. This problem is common in many social science applications, where a large number of potential covariates are available to explain the time-series of a large cross-section of units or individuals. An example is empirical asset pricing, where the literature has produced a “factor zoo” of potential risk factors to explain the large cross-section of stock returns. This problem requires a large panel, as a successful asset pricing model should explain the many available investment strategies, resulting in a large panel of test assets. At the same time, there is no consensus about what are the appropriate risk factors, which leads to a statistical selection problem from a large set of candidate covariates. So far, the literature has only provided solutions to one of the two subproblems, while keeping the dimensionality of the other problem small. Our paper closes this gap.
The inferential theory on a large panel with many covariates is a challenging problem. As a first step, we have to select a sparse set of covariates from a large pool of candidates with a regularized estimator. The challenge is to provide valid $p$-values from this estimation that account for the post-selection inference. Furthermore, researchers might want to impose economic priors on which variables should be more likely to be selected. The second challenge is that the panel cross-section results in a large number of $p$-values. Hence, some of them are inadvertently very small, which if left unaddressed leads to “$p$-hacking”. The multiple testing adjustment conditional on the selected subset of covariates from the first step is a novel problem, and requires to redesign what hypotheses should be tested jointly. A naive counting of all tests is overly conservative, and the test design and simultaneity counts should to be conditional on the covariate selection.
This paper proposes a new method for covariate selection in large dimensional panels, tackling all of the above challenges. We develop the inferential theory for large dimensional panel data with many covariates by combining post-selection inference with a new multiple testing method specifically designed for panel data. Our novel data-driven hypotheses are conditional on sparse covariate selections and valid for any regularized estimator. Based on our panel localization procedure, we control for family-wise error rates for the covariate discovery and can test unordered and nested families of hypotheses for large cross-sections. As an easy-to-use and practically relevant procedure, we propose Panel-PoSI, which combines the data-driven adjustment for panel multiple testing with valid post-selection $p$-values of a generalized LASSO, that allows to incorporate priors. In simulations and an empirical study we show that our selection method provides correct false detection control but has substantially higher power than existing approaches.
Our paper proposes the novel conceptual idea of data-driven hypotheses families for panels. This allows us to put forward a unifying framework of valid post-selection inference and multiple testing. Leveraging our data-driven hypotheses family, we adjust for multiple testing with a localized simultaneity count, which increases the power, while maintaining false discovery rate control. An essential step for a formal statistical test is to formulate the hypothesis. This turns out to be non-trivial for a large panel with a first stage selection step for the covariates. It is a fundamental insight of our paper, that the hypothesis of our test has to be conditional on the selected set of active covariates of the first stage. Once we have defined the appropriate hypothesis, we can deal with the multiple testing adjustment, which by construction is also conditional on the selection step.
Our method is a disciplined approach based on formal statistical theory to construct and interpret a parsimonious model. It goes beyond the selection of a sparse set of covariates as it also provides the inferential theory. This is important as it allows us to rank the covariates based on their statistical significance and can also be applied for relatively short time horizons, where cross-validation for tuning a regularization parameter might not be reliable. We answer the question which covariates are needed to explain the full panel jointly, and can also accommodate “weak” covariates or factors that only affect a small subset of the cross-sectional units.
Our data-driven hypothesis perspective exploits the geometric structure implied by the first stage selection step. Given valid post-selection $p$-values of a regularized sparse estimator from time-series regressions, we collect them across the large cross-section into a “matrix” of $p$-values. Only active coefficients, that are selected in the first stage, contribute $p$-value entries, whereas covariates that were non-active lead to “holes” in this matrix. We leverage the non-trivial shape of this matrix to form our adaptive hypotheses. This allows us to make valid multiple testing adjusted inference statements, for which we design a panel modified Bonferroni-type procedure that can control for the family-wiser error rate (FWER) in the discovery of the covariates. As one loosens the FWER requirements, the inferential thresholds admits more and more explanatory variables. Hence, the number of admitted covariates and the FWER control level form an “false-discovery control frontier”. We provide a method that allows us to traverse the inferential results and determine the least number of covariates that have to be included given a user-specified FWER level. In other words, we can make a statement on the number of covariates or factors needed to explain a panel based on a statistical significance requirement.
We propose the novel procedure Panel-PoSI, which combines the data-driven adjustment for panel multiple testing with valid post-selection $p$-values of a generalized LASSO. While our multiple testing procedure is valid for any sparsity constrained model, Panel-PoSI is an easy-to-use and practically relevant special case. We propose Weighted-LASSO for the first stage selection regression and provide valid $p$-values through post-selection inference (PoSI), which yields a truncated-Gaussian distribution for an adjusted LASSO estimator. This geometric perspective is less common in the LASSO literature, but has the advantage that it avoids the use of infeasible quantities, in particular the second moment of the large set of potential covariates. The Weighted-LASSO generalizes LASSO by allowing to put weights onto prior belief sets. For example, a researcher might have economic knowledge that she wants to include in her statistical selection method, and impose an infinite prior weight to include specific covariates in the sparse selection model. Our Weighted-LASSO makes several contributions. First, the expression for the truncated conditional distribution with weights become much more complex than for the special case of the conventional LASSO. Second, we provide a simple, easy-to-use and asymptotically valid conditional distribution in the case of an estimated noise variance.
We demonstrate in simulations and empirically that our inferential theory allows us to select better models. We compare different estimation approaches to select covariates and show that our approach better trades off false discovery and correct selections and hence results in a better out-of-sample performance. Our empirical analysis studies the fundamental problem in asset pricing of selecting a parsimonious factor model from a large set of candidate factors that can jointly explain the asset prices of a large cross-section of investment strategies. We consider a standard data set of 114 candidate asset pricing factors to explain 243 double sorted anomaly portfolios. We show that Panel PoSI selects 3 factors which form the best model to explain out-of-sample the expected returns and the variations of the test assets. The selected factors are economically meaningful and include the size and value factors of the Fama-French model. Hence, our statistical selection procedure confirms two of the most widely used asset pricing factors. Our findings contribute to the discussion about the number of asset pricing factors. We confirm that independent of the rotation of the covariates we select 3 factors for a 5% FWER control.
The rest of the paper is organized as follows. Section (ref) relates our work to the literature. Section (ref) introduces the model and the Weighted-LASSO. Section (ref) discusses the appropriate hypotheses to be considered for inference on the entire panel. Section (ref) proposes a joint unordered test for the panel using multiple testing adjustment so that we can maintain FWER control, and shows how to traverse this procedure to acquire the least covariate count associated with each FWER target. In section (ref), we consider the case of nested hypotheses, where the covariates follow a fixed ordering, which is of independent interest, and we propose a step-down procedure for this setting that maintains false discovery control. Section (ref) provides the results of our simulation and Section (ref) discusses our empirical analysis on a large asset pricing panel data set. Section (ref) concludes. The Appendix collects a detailed discussion about post-selection LASSO. The proofs and more technical details are available in the Online Appendix.
The problem of multiple testing is an active area of research with a long history. The statistical inference community has studied the problem of controlling the classical FWER since Bonf, and controlling for false-discover rate (FDR) going back to BH95 and BY01. Bonf allows for arbitrary correlations in the test statistics because its validity comes from a simple union bound argument, and is in fact the optimal test when statistics are “close to independent” under true sparse non-nulls. FDR control on the other hand requires a discussion about the estimated covariance in the test statistics. Recent developments include a stream of papers led by 15-AOS1337 and rssb.12265, which constructs a generative model to produce fake covariates and control for FDR. fithian2022conditional is a more recent work that iteratively adjusts the threshold for each hypothesis in the family to seek finite sample exact FDR control and dominates BH95 and BY01 in terms of power. Another notion on temporal false discovery control has been revived more recently by doi:10.1287/opre.2021.2135, who consider the industry practice of constantly checking $p$-values and provide an early stopping in line with SiegmundSeq that adjusts for a bias from sequentially picking favorable evidence, whereas we consider a static panel that is not an on-going experiment.
There are cases where the covariates warrant a natural order such that the hypothesis family possesses a special testing logic. A hierarchical structure in covariates arises when the inclusion of the next covariate only make sense if the previous covariates is included. An example is the use of principal component (PC) factors, where PCs are included sequentially from the dominating one to the least dominating one. We distinguish this from putting weights and assigning importance on features because this variant of family of hypotheses warrants a new definition of FWER. We propose a step-down procedure that can be considered as a panel extension of rssb.12122, relying on an approximation of the R\'enyi representation of $p$-values. The step-down control for nested FWER is based on 2336545, which along with Bonf can be seen as comparing sorted $p$-values against linear growth. Our framework contributes to estimating the number of principal component factors in a panel. There are have been many studies that provide consistent estimators for the number of PCs based on the divergence in eigenvalues of the covariance matrix, which include 1468-0262.00273, 40985808, ECTA8968 and PELGER201923. Another direction uses sequential testing procedures that presume correct nested family of hypotheses, which include jbes.2009.07239 and 16-AOS1536. In contrast, we characterize the least number of covariates (which can also be based on principal components), which should be expected when a FWER rate is provided. The nested version of our procedure is close in nature to a panel version of “when-to-stop” problem of a multiple testing procedure.
The problem of post-LASSO statistical testing for small dimensional cross-sections is studied in a stream of papers including 009053606000000281, rssb.12026, 14-AOS1221 and 10.1214/17-AOS1630, which consider inference statements by debiasing the LASSO estimator. An alternative stream of post-selection or post-machine learning inference literature includes annurev-economics-012315-015826, kuchibhotla2018valid and zrnic2020postselection, who provide non-parametric post-selection or post-regularization valid confidence intervals and $p$-values. These papers do not make conditional statements and presume that the researcher sets the hypotheses before seeing the data, which we will refer to as data agnostic hypothesis family. We follow a different train of thought that treats LASSO, among a family of conic maximum likelihood estimator, as a polyhedral constraint on the support of the response variable. This geometric perspective that provides inferential theory post-LASSO is pioneered by the work of lee2016exact and followed up by fithian2017optimal and tian2018selective, assuming Gaussian linear models. markovic2018unifying extend the results to LASSO with cross-validation, tian2017selective discuss a square-root LASSO variant that takes an unknown covariance into consideration. taylortibshirani2016inference and tian2017asymptotics study asymptotic results that allow to relax the assumptions of Gaussian errors or generalized linear models. This body of literature is often referred to as PoSI, and traverses the Karush-Kuhn-Tucker (KKT) condition of a LASSO optimization problem to show that the LASSO fit can be expressed as a polyhedral constraint on the support of the response variable. We extend this work by allowing to put weights onto prior belief sets, and by bringing it to the panel setting with multiple testing adjustment.
We consider a large dimensional panel data set $\bm{Y} \in \mathbb{R}^{T\times N}$ which we want to explain with a large number of potential covariates $\bm{X} \in \mathbb{R}^{T\times J}$. The panel data and explanatory variables are both observed over $T$ time periods.\footnote{Our setting and multiple testing results can be readily extended to the case of unbalanced panel, although we focus on the balanced panel case for now to highlight the core multiple testing insight of our method. We will further discuss on this once we introduce our main procedure in Section (ref)} The size of the cross-section $N$ and the dimension of the covariate candidate set $J$ are both large in our problem.
We assume a linear relationship between $\bm{Y}$ and $\bm{X}$:
which reads in matrix notation as
We refer to the coefficients $\bm{\beta}$ as loading matrix, where the $n$th column $\beta^{(n)}\in\mathbb{R}^{J}$ corresponds to the $n$th unit and $\beta^{(n)}_{j}$ denotes the loading of the $n$th unit on the $j$-th covariate. The remainder term $\bm{\epsilon}$ is unexplained noise. Throughout this paper, we use the superscript $\cdot^{(n)}$ to denote cross-sectional variables corresponding to the $n$th unit, subscript $\cdot_j$ for variables corresponding to the $j$th covariate, and subscript $\cdot_t$ for time-series corresponding to the $t$th time period.
We assume that a sparse linear model can explain jointly the full panel. Formally, a sparse linear model with $s$ active covariates is
where $s=|S|$ is the cardinality of the set of active covariates $S=\{j:\exists \beta^{(n)}_{j}\neq 0, n \in \{1,...,N\}$, that is, the set of covariates with non-zero loadings. $\bm{X}_S$ is the subset of covariates that belong to $S$. Our goal is to estimate this low dimensional model, that can explain the full panel, from a large number of candidate covariates, and provide a valid inferential theory.
Note that our sparse model formulation allows for two important properties. First, different units can be explained by different covariates with different loadings. This means that $\beta^{(n)}\neq \beta^{(m)}$ for $n \neq m$ is allowed. For example, a subset of the cross-sectional units might be modeled by different covariates than the remaining part of the panel. Second, we can accommodate “weak” covariates. A covariate is included in $S$ if it is required by at least one cross-sectional unit as explanatory variable. In other words, a sparse model can include covariates in $\bm{X}_S$ that explain only a very small subset of the panel $\bm{Y}$.
The first step is to estimate the sparse models over the time-series for each unit separately due to the heterogeneity in the loadings. In a second step, we provide the valid inferential theory for the loadings on the full panel. The time-series estimation requires an appropriate regularization to select a small subset of covariates that contains all the relevant covariates for each unit. We allow for prior belief weights $\omega_j\in (0,+\infty]$ on the $J$ candidate covariates, so that different $\bm{X}$ can have different relative penalizations, and a global $\lambda \in\mathbb{R}_+$ scalar penalty parameter. For the $n$th unit, we denote its $\beta^{(n)}$ regularized estimate as $\hat{\beta}^{(n)}$ and the active set $M^{(n)}=\{j:\hat{\beta}^{(n)}_j\neq 0\}$ as the set of $j$'s with non-zero loadings $\hat{\beta}^{(n)}_j$. A general regularized linear estimator solves the following optimization problem
for a penalty function $f$ and appropriate weights, where $Y^{(n)}$ is the vector of response variables of $n$th unit. In this paper, we consider the weighted-LASSO estimator with the regularization function
and weights $\omega_j>0$ for all $j\in \{1,...,J\}$ and $\sum_{j=1}^J\omega_j^{-1}=J$. We consider the penalty $\lambda$ as exogenously provided such that the set $\|\hat{\beta}^{(n)}\|_0=|M^{(n)}|$ is low dimensional.\footnote{In Appendix (ref) we discuss the case where variances of $\bm{\epsilon}$ are unknown and need to be estimated. We provide a specific discussion on $\lambda$'s rate conditions in terms of $J$ and $T$, as is typically required to ensure consistency in the LASSO literature.} Importantly, we do not need to assume that the selected set contains all “true” active covariates. Our goal is to provide a valid inferential theory conditional on the selected set. Our estimator generalizes the conventional LASSO with the $l_1$ regularization function of 2346178 by allowing for different relative weighting in the penalty. Importantly, we also allow for an infinite weight, which can be interpreted as a prior on a set of covariates. This allows researchers to take advantage of prior information and for example ensure that a specific set of covariates will always be included. The weighted-LASSO will be particularly relevant in our empirical study, where we can answer the question which risk factors should be added to a given set of economically motivated risk factors. Our weighted-LASSO formulation can also be interpreted as a Bayesian estimator with the canonical Laplacian prior.
Conventional regression theory will not provide correct inferential statements on the weighted-LASSO estimates. We face two challenges. First, regularized estimation results in a bias, which needs to be corrected. Second and more challenging, post-selection inference changes the distribution of the estimators. When we observe an active $\hat{\beta}^{(n)}_j$ from ((ref)), it would be incorrect to simply calculate its $p$-value from a conventional Student $t$-distribution. This invalidity stems from the fact that conditional on observing a LASSO output, ${\beta}^{(n)}_j$ must be large enough in magnitude for its $\hat{\beta}^{(n)}_j$ to be active. In other words, the probability distribution of the estimators is truncated.
The correct inference has to be conditional on the covariates being selected by the LASSO estimator. Hence, valid $p$-values have to be the tail probability conditional on being in the selection set. The key to quantify such styles of inference is to recognize that a sparsity constrained estimator is typically the result of solving Karush-Kuhn-Tucker (KKT) conditions, which can in turn be geometrically characterized as polyhedral constraints on the support of response variables. This is first established in lee2016exact, who provide the stylized results that Post-Selection Inference (PoSI) of debiased non-weighted LASSO estimators can be calculated as polyhedral truncation on $\bm{Y}$. This line of research is also referred to as Selective Inference, for example in taylor2015statistical. We extend this line of literature to allow for the Weighted-LASSO. We derive these results with assumptions common in the PoSI LASSO literature, detailed in Appendix (ref), and referred to as conventional regularity conditions for the ease of exhibition.
Theorem (ref) shows how we calculate $p$-values from the post-selection distribution of the debiased estimate $\bar{\beta}_j^{(n)}$. The regularized estimate $\hat{\beta}_j^{(n)}$ has a well-known bias. We debiase the LASSO estimate by a shifting argument. While we use a geometric argument to remove the bias, the bias adjustment takes the usual form in the LASSO literature as for example in 10.3150/11-BEJ410. The debiased LASSO estimator simply equals a standard OLS estimation on the subset $M^{(n)}$ selected by the Weighted-LASSO.
Theorem (ref) has two key elements. First, the distribution of the linear coefficients is not a usual Gaussian distribution, but it is truncated due to studying post-selection coefficients. This geometric perspective is less common in the LASSO literature, but provides several advantages. One advantage of the geometric approach is that it avoids the use of infeasible quantities, in particular the second moment of the large set of potential covariates. Second, conditional on the selection and under the null hypothesis that the coefficients of the active covariates for unit $n$ are zero, the post-selection $p$-values of the active coefficients follow a uniform distribution. This is important as it implies that we obtain valid post-selection $p$-values, which provide the correct Type-I error control for individual regressions. These individually valid $p$-values are the key for deriving the multiple testing adjustment in large panels.
Appendix (ref) provides the detailed information on constructing $\bar{\beta}$ and the definitions of $\eta, A, b(\omega)$ along with lemmas that lead up to this result. It also shows that, under the assumption of a known noise variance, the distribution result is not asymptotic in $T$, but also valid in finite samples. We can obtain these results because we make the stronger assumption that the noise is normally distributed. When using a consistent sample estimate of the noise variance, we need to require additionally that $T \rightarrow \infty$. Appendix (ref) clarifies the implications of different assumptions in the three Theorems (ref), (ref) and (ref). Theorem (ref) with estimated noise variance represents the explicit form of Theorem (ref), which we use for our empirical analysis. It is possible to relax the normality assumption of the noise and instead use a pivot convergence similarly to tian2017asymptotics to obtain asymptotically a truncated Gaussian distribution. However, this would not change the nature of our statement.
Our Weighted-LASSO results make several contributions. First, the expression for the truncated conditional distribution with weights become more complex than for the special case of the conventional LASSO. Second, we provide a simple, easy-to-use and asymptotically valid conditional distribution in the case of an estimated noise variance. Last but not least, we show the formal connection with alternative debiased LASSO estimators by showing that debiasing can be interpreted as one step in a Newton-Ralphson method of solving a constrained optimization.
Theorem (ref) allows us to obtain valid post-selection $p$-values for Weighted-LASSO coefficients. We obtain these $p$ values from the simulated cumulative distribution function of the truncated Gaussian distribution. Crucially, all results for multiple testing adjustment in panels that we study in the following sections neither require us to use a weighted Lasso estimator nor to use the $p$-values implied by Theorem (ref). We only require to have a set of valid post-selection $p$-values for sparsity constrained models. These can be obtained with any suitable regularized estimator and post-selection inference. The key element is the selection of a low dimensional subset with $p$-values conditional on this selection. We propose the weighted LASSO conditional inference results as an example of the type of sparsity constraint models we are interested in, and demonstrate a machinery with which we can obtain valid $p$-values for sparsity constrained models. In our empirical studies, we use Weighted-LASSO as our sparsity constrained model since we want to specify strong prior beliefs on a few covariates and it is common practice to use LASSO in the context of our empirical studies. Nonetheless, the testing methods in the next sections accommodate any sparse estimator, and can be detached from inference for Weighted-LASSO.
Our goal is to provide formal statistical tests that allow us to establish a joint model across a large cross-section with potentially weak covariates. This requires us to provide a form of statistical significance test with multiple testing adjustment that properly accounts for covariates that only explain a small subset of the cross-sectional units. This is important as in many problems in economic and finance, there is substantial cross-sectional variation in the explanatory power of covariates, and a model that simply minimizes an average error metric might neglect weaker covariates.
An essential step for a formal statistical test is to formulate the hypothesis. This turns out to be non-trivial for a large panel with a first stage selection step for the covariates. It is a fundamental insight of our paper, that the hypothesis of our test has to be conditional on the selected set of active covariates of the first stage. Once we have defined the appropriate hypothesis, we can deal with the multiple testing adjustment, which by construction is also conditional on the selection step.
Our hypothesis formulation and test construction only requires valid post-selection $p$-values from a first stage selection estimator as formalized in Section (ref). The results of the next two sections do not depend on a specific model for obtaining these $p$-values and the active set. The results are valid for any model including non-linear ones. The input to the analysis is a $N \times J$ matrix, which specifies which covariates are active for each unit and the corresponding post-selection $p$-values. The Weighted-LASSO is only one possible model, but it can be replaced by any regularized model. We have introduced the sparse linear model as it is the horse race model for many problems in economics and finance, and therefore of practical relevance.
We illustrate the concept of a data-driven hypothesis with a simple example, which we will use throughout this section. For simplicity we assume that we have $J=4$ covariates and want to explain $N=6$ cross-sectional units. In the first stage, we have estimated a sparse model and have obtained the post-selection valid $p$-values for each of the $N$ units. We collect the fitted sparse estimator $\bar{\beta}^{(n)}$ for the $n$th unit in the matrix $\bar{\bm\beta}$. Note, that this matrix has “holes” due to the sparsity for each $\bar{\beta}^{(n)}$. Figure (ref)(a) illustrates $\bar{\bm\beta}$ for this example.
Similarly, we collect the corresponding $p$-values in the matrix $\bm{P}$. For the $n$th unit, we only have $p$-values for those covariates that are active in the $n$th linear sparse model. Thus, Figure (ref)(b) also has white boxes showing the same pattern of unavailable $p$-values due to the conditioning on the output of the linear sparse model. These holes can appear at different positions for each unit, which makes this problem non-trivial. This non-trivial shape of either subplot (a) or (b) is completely data-driven and a consequence of linear sparse model selection. We show that the hypothesis should be formed around these non trivial shapes as well, which is why we name it the data-driven hypothesis family.
We want to test which covariates are jointly insignificant in the full panel. A data-agnostic approach would simply test if all covariates are jointly insignificant, independent of the data-driven selection step in the first stage. A data-agnostic hypothesis is unconditional as it does not depend on any model output. However, as we will show, this perspective is problematic for the high-dimensional panel setting with many covariates as it ignores the dimension reduction from the selection step. Therefore, an unconditional multiple testing adjustment accounts for “too many” tests, which severely reduces the power.
We propose to form the hypothesis conditional on the first stage selection step. The data-driven hypothesis only tests the significance of the covariates that were included in the selection, and hence can drastically reduce the number of hypothesis. However, given the non-trivial shape of the active set, the multiple testing adjustment for the data-driven hypothesis is more challenging.
Before formally defining the families of hypothesis, we illustrate them in our running example. The data-agnostic hypothesis $H_A$ for explaining the full panel takes the following form:
The data-driven hypothesis $H_D$ only includes the active set and hence equals
Clearly, $H_A$ has a larger cardinality of $|H_A|=24>|H_D|=14$. This holds in general, unless the first stage selects all covariates for each unit, in which case the two hypotheses coincide.
Formally, the data-agnostic family of hypothesis is defined as follows:
It is evident that $H_A$ does not need any model output or exploratory analysis, so it is indeed data-agnostic.
As soon as we use a sparsity constrained model that has censoring capabilities, we no longer observe $(\bm{Y},\bm{X})$ from its data generating process. Consequently, unless our hypotheses depend on how we built the model, or equivalently on how the data was censored, the data-agnostic hypotheses forgo power without any benefit in false discovery control. Therefore, we formulate the hypothesis on the $j$th covariate $H_{0,j}^{(n)}$ only if $j\in M^{(n)}$, that is, it is in the active set of the $n$th unit. Conditional on observing the model output, there is no inference statement to be made about $H_{0,j}^{(n)}$ if $j\notin M^{(n)}$, because its estimator is censored by the model.
We denote as $\mathcal{K}_j$ the set of units for which the $j$th covariate is active. We define the cross-sectional hypothesis for the $j$th covariate as:
By combining all covariates $\{j:\mathcal{K}_j\neq \emptyset\}$ that show up at least once in one of the active sets of our sparse linear estimators, we arrive at a data-driven hypothesis associated with our panel. This is defined as follows:
This demonstrates the non-trivial nature of writing down a hypothesis in high-dimensional panel: we can only collect $\mathcal{K}_j$ - the set of units for which the $j$th covariate is active - after seeing the sparse selection estimation result.
We show how to adjust for multiple testing of data-driven hypotheses. Given the the first stage selection of active covariates in $\mathcal{P}$, we form the data-driven hypothesis $H_D$. The only assumption that we require for the multiple testing adjustment is that we have valid post-selection $p$-values $p_j^{(n)}$ for covariate $j\in M^{(n)}$ and unit $n \in \mathcal{K}_j$. This is formalized in the following assumption:
Valid post-selection $p$-values abstract away from model-specific conditions on how the $p$-values are obtained. Post-selection $p$ values are trivially valid, if they follow a uniform distribution under the null hypothesis conditional on the selection event. This special case can be interpreted as exact valid post-selection $p$-values, since $p^{(n)}_j \stackrel{H_D|\mathcal{M},\omega}{\sim} \textrm{Unif }[0,1]$ implies $\mathbb{P}_{H_D|\mathcal{M},\omega}(p^{(n)}_j \leq x)= x$. The result, that valid post-selection p-values conditional on the selection event and weights are uniform under the null hypothesis, is first introduced in tibshirani2016exact for sparse regressions. In Theorem (ref), we extend this result to generic sparse regressions with weights. Hence, the post-selection $p$-values of our Weighted-LASSO satisfy Assumption (ref).
Our definition of valid post-selection $p$-values is natural as it simply states that we have correct conditional Type-I error control for each cross-sectional unit and each active covariate. This definition is in line with tibshirani2016exact and heard2018choosing, but more general. It also allows for more conservative $p$-values, that would have smaller size and worse power. The valid post-selection $p$-values provide the correct Type-I error control for each unit individually, and do not take the multiple testing issue into account. Hence, the $p$-values of Assumption (ref) are essentially the results of post-selection inference that is applied separately to each cross-sectional unit $n$. Given those values we show how to correct them to adjust for multiple testing.
Our derivations for the multiple testing adjustment only take advantage of the conditional valid distribution of the individual post-selection $p$-values, and hence we impose this property as the fundamental underlying assumption. Note that that our results do not require a linear model, but hold for any set of valid post-selection $p$-values.\footnote{This notion of valid $p$-values deviates from what is commonly used in regression analysis as its entire statement is based on a conditional distribution, highlighting its “post-selection” nature. In other words, the Assumption (ref) holds under the null $H_D$ and conditional on the selection event $\mathcal{M}$ and weights $\omega$, as opposed to the classical regression analysis that does not depend on the selection.} The assumptions on the data generating process and asymptotic regime are implicitly included in Assumption (ref). Theorem (ref) is a specific example that imposes Gaussian errors and $T \rightarrow \infty$. Hence, as long as a researcher has a selection estimation approach that provides valid post-selection $p$-values, our multiple testing results are applicable.
Our goal is to reject members of $H_D$ while controlling the Type I error, and the common way to measure such an error is the family-wise error rate. This is the same underlying logic that is used to define confidence intervals and determine significance of covariates in a conventional setup. The crucial difference is that we need to account for multiple testing given the large number of cross-sectional units. The family-wise error rate (FWER) is defined as follows:
Similar to the conventional definition, we simply count the number of Type I false rejections $V$, and define FWER as the probability of making at least one false rejection. Importantly, the FWER accounts for the fact that we might repeatedly test a specific covariate for multiple cross-sectional units rather than just for one unit. Our contribution to FWER control in the panel setting is thus to take into consideration both the multiplicities in units and covariates when we deal with the “matrix” of $p$-values $\bm{P}$. To achieve this goal, we propose a new simultaneity account for the $j$th covariate, calculated as
Figure (ref) illustrates the simultaneity counting for our running example with $N=6$ units and $J=4$ covariates. The blue boxes represent the active set for a specific covariate. The yellow boxes indicate the “co-active” covariates, which have to be accounted for in a multiple testing adjustment. In the case of the first covariate $j=1$, only the second unit $n=2$ has selected this covariate. This second unit has also selected covariate $j=3$ and $j=4$, which are jointly tested with the first covariates. Hence, they are “co-active”, and the simultaneity count equals $N_1=3$. Intuitively, $N_j$ represents all relevant comparisons for the $j$th covariate because it counts how many covariates are active with the $j$th covariate in the regressions. Hence, $N_j$ quantifies the number of “multiple tests” for each covariate.
In subplot (ref)(a), we see that $\mathcal{K}_1=\{2\}$ for the 1st covariate, indicated by the blue box, because it is only active in the second unit's regression. The multiple testing adjustment needs to consider all yellow boxes, and $N_1=3$ is thus the total count of 1 blue and 2 yellow boxes. Similarly, for the second covariate, $\mathcal{K}_2=\{1,3,5,6\}$, so we shade boxes yellow for the 2nd, 3rd and 5th units and obtain $N_2=9$. We can already see that our design of simultaneity counts takes all relevant pairwise comparisons into considerations, but avoids counting the white boxes - which would cause overcounting and result in over-conservatism.
Our multiplicity counting is a generalization of the classical Bonferroni adjustment for multiple testing. A conventional Bonferroni method for the data-agnostic hypothesis $H_A$ has a simultaneity count of $|H_A|=N\cdot J=24$ for testing each covariate. A direct application of a vanilla Bonferroni method to the panel of all selected units and the data-driven hypothesis $H_D$, would use a simultaneity count of $|H_D|=14$ for testing each covariate. Our proposed multiplicity counting is a refinement that leverages the structure of the problem, and takes the heterogeneity of the active sets for each covariate into account. Our count has only $N_1=3$, $N_2=9$ and $N_4=8$ for the covariates $j=1,2$ and $4$. Only for covariate $j=3$ is the simultaneity count the same as a vanilla Bonferroni count applied to $H_D$, i.e. $N_3=14$.
In addition to the simultaneity count of each covariate, we need an additional “global” metric for our testing procedure. We define a panel cohesion coefficient $\rho$ as a scalar that measures how sparse or de-centralized the proposed hypotheses family is:
The panel cohesion coefficient $\rho$ is conditional on the data-driven selection of the overall panel. It is straightforward to compute once we observe the sparse selection of the panel. This coefficient takes values between $J^{-1}$ and 1,\footnote{We prove this bound in the Online Appendix, without leveraging sparsity of first-stage models but rather as an algebraic result with intuitive interpretation.} where larger values of $\rho$ imply that the active set is more dependent in the cross-section. This can be interpreted as that the panel $Y$ has a stronger dependency due to the covariates $X$. Intuitively, in the extreme case when $\rho=J^{-1}$, the panel can be separated into $J$ smaller problems, each containing a subset of response units explained by only one covariate. Thus the panel would be very incohesive, and could be studied with $J$ separate tests. In the other extreme, if $\rho $ approaches 1, the first-stage models include the same active covariates for all units. We consider this as a very cohesive panel. If $\rho$ is between theses bounds, the panel is cohesive in a non-trivial way such that some units can be explained by some covariates and there is no clear separation of the panel into independent subproblems.
Figure (ref) illustrates the panel cohesion coefficient with examples. The subplots show four active sets that are different from our running example. The left subplot (ref)(a) shows the extreme case of $\rho=J^{-1}$, where the panel is the least cohesive. The right subplot (ref)(d) illustrates the other extreme for $\rho=1$, where the panel is the most cohesive. The middle subplots (ref)(b) and (c) correspond to the complex cases of a medium cohesion coefficient.
Our novel simultaneity count and cohesiveness measure are the basis for modifying a Bonferroni test for FWER-controlled inference. Theorem (ref) formally states the FWER control. The proof is in Appendix (ref).
Theorem (ref) is based on an algebraic union bound argument, that leverages the structure of the panel hypotheses $H_D$. The assumptions on the data generating process and asymptotic distribution are implicitly included in the valid post-selection $p$-values.
This completes the joint testing procedure. First, we calculate $p$-values after running a sparse linear time-series regression. Second, we use the output of the sparse linear estimation to write down a hypothesis and, third, we provide a FWER control inference procedure by combining the $p$-values across the cross-section and test the hypothesis.
The difference between a naive Bonferroni and our FWER control is particularly pronounced for weak covariates that affect only a subset of the cross-sectional units. Given a FWER control level of $\gamma$, the rejection threshold for a naive Bonferroni test is $\frac{\gamma}{J N}$ for every covariate. The rejection threshold for our FWER control is always higher, and differs in particular when $N_j$ is small and $\rho$ is large. This is the case for weak covariates in a cohesive panel.
As it is common in statistical inference, we focus on Type I error control. Type II error rates require the specification of alternatives. While we do not provide formal theoretical results for the power of our inference approach, we show comprehensively in the simulation and empirical part, that our approach has substantially higher power than conventional approaches.
We point out that the validity of our procedure holds for unbalanced panels as well. This is because even when there are different number of observations for the $n$th and $m$th units, i.e. $T_n\neq T_{m}$ for $n \neq m$, they can still be estimated separately in the first stage of the regularized regression. The hypothesis testing and selection of a parsimonious model only requires the matrix $\bm{P}$ of valid $p$-values, which can be based on samples of different sizes.
The typical logic of statistical inference is to determine which covariates we should admit from $X_M$, given a significance level $\gamma$. We use $K$ to denote the number of selected covariates. When $\gamma$ is specified as a lower quantity, we expect $K$ to decrease as well, that is, the rejection becomes harsher.
As the number of admitted covariates of our procedure is monotonically increasing in $\gamma$, we want to ask the following converse question: How do we need to set $\gamma$ such that we reject $K$ covariates? Concretely, we want to find:
Let $p_j=\min_{n\in \mathcal{K}_j}\{p_j^{(n)}\}$ be the $1$st order statistic for $j=1,...,J$. Then ((ref)) is simply the $K$-th order statistics of $N_j p_j / \rho$:
Since this minimization scan is monotone, we can determine how many covariates at least should be admitted, given a control level, which is similar to the “SimpleStop” procedure described in 16-AOS1536. The following corollary formalizes this inversion method that finds the least number of covariates to admit:
The statement simply states that the simplest linear model should have at least $K^*(\gamma)$ covariates for a given $\gamma$. Note that it is possible that, for example, $\gamma^*(5)$ and $\gamma^*(6)$ are both equal to $0.05$, while $\gamma^*(7)>0.05$. In this case the minimum number of covariates is $K^*(0.05)=6$ because it does not hurt FWER-wise to include 6 covariates in the model. Hence, we are making a slightly different statement than that there would be exactly $K^*(\gamma)$ covariates in the true linear model. The number of covariates is obviously conditional on the set of candidate covariates $\bm{X}$, and we can only make statements for this given set.
In our empirical study we consider candidate asset pricing factors $\bm{X}$ to explain the investment strategies $\bm{Y}$. More generally, the linear model that we consider is often referred to as a factor model. Therefore, we will also refer to the selected covariates as factors, and use these two expressions as synonyms moving forward. This directly links our procedure to the literature on estimating the number of factors to explain a panel. A common approach in this literature is to use statistics based on the eigenvalues of either $\bm{Y}$ or $\bm{X}$ to make statements about the underlying factor structure. Our approach is different, as it provides significance levels for the selected factors and FWER control for the number of factors.
Table (ref) illustrates the estimation of the number of factors and their ranking with our running example introduced in Figure (ref). We calculate the simultaneity counts $N_j$'s as given in ((ref)) and demonstrated in Figure (ref), and $p_j$ as the smallest $p$-values associated with the $j$th covariate. Then, the rejection rule in Theorem (ref) is based on whether a pre-specified level $\gamma$ satisfies $p_j<\frac{\rho \gamma}{N_j}$, which is equivalent to $\rho ^{-1}\cdot N_j\cdot p_j <\gamma$.
Thus, the natural ranking of the covariates is to sort all covariates in descending order of the $\rho ^{-1}\cdot N_j\cdot p_j $ values as shown in Table (ref). It is then trivial to determine $K^*(\gamma)$ for any choice of $\gamma$. For example, for $\gamma=1\%$, we would select factors 3 and 4, but not 1 or 2. On the other hand, for $\gamma>2\%$, we would include all four factors. Hence, the ranking of $\rho ^{-1}\cdot N_j\cdot p_j $ directly maps into $K^*(\gamma)$. Moreover, the ordered list of $\rho ^{-1}\cdot N_j\cdot p_j $ provides an importance ranking of the factors. Furthermore, the number $N_j$ reveals if significant factors are “weak”. In our case, factor 1 has $N_1=3$, which indicates that it affects only a small number of hypothesis. Its $p$-value $p_1$ is sufficiently small to still imply significance in terms of FWER control.
For comparison, Table (ref) also includes the corresponding analysis for the data-agnostic hypothesis and a conventional Bonferroni correction. The Bonferroni analysis uses the same $p$-values but a different multiple testing adjustment. In our case, the $p$ values would be multiplied by $J \cdot N=24$ as this corresponds to the total number of hypothesis tests. This will obviously make the inference substantially more conservative. Indeed, even for a FWER control of $\gamma=4\%$, we would only select factors 3 and 4. We would need to raise the FWER control to $\gamma=12\%$ to include factor 1. Hence, weak factors, like factor 1, are more likely to be discarded by the data-agnostic hypothesis with conventional multiple testing adjustment.
We emphasize that a data-agnostic hypotheses with conventional Bonferroni correction does provide correct FWER control, but it is overly conservative, and does not sufficiently leverage information already observable in LASSO estimation results. By construction, the data-agnostic Bonferroni approach will test a larger number of hypothesis, which means that the corresponding “significance levels” will always be lower or equal to our data-driven simultaneity count. Second, the data-agnostic Bonferroni approach does not differentiate the “strength” of the factors, while our approach provides a selection-based heterogeneous adjustment of the $p$-values. This is essential for detecting weak factors.
Having introduced all building blocks of our novel method to detect covariates, we put the entire procedure together as “Panel-PoSI”:
As we have now all results in place, we can summarize the advantages of our procedure. First, we want to clarify that our goals and results are different from just some form of optimal shrinkage selection. Selecting a shrinkage parameter with some form of cross-validation in a regularized estimator like LASSO does not provide the same insights and model that we do. A shrinkage estimator can either be applied to each unit separately, as we do it in our first step, or to the full panel in a LASSO panel regression. The separate covariate selection for each cross-sectional unit does not answer the question which covariates are needed to explain the full panel jointly. A shrinkage selection on the full panel for some form of panel LASSO can neglect weaker factors, as those receive a low weight in the cross-validation objective function. Second, tuning parameter selection with cross-validation requires a sufficiently large amount of data. Our approach is attractive as we can do the complete analysis on the same data. That means, an initial LASSO is used to first reduce the number of covariates, but this set is then further trimmed down using inferential theory. Hence, we can construct a parsimonious model even for data with a relatively short time horizon, but large cross-sectional dimension. Third, the statements that we can make are much richer than a simple variable selection. We can formally assess the relative importance of factors in terms of their significance. The model selection is directly linked to a form of significance level, which allows us to assess the relevance of including more factors. Last but not least, we can also make statements about the strength of factors. In summary, Panel-PoSI is a disciplined approach based on formal statistical theory to construct and interpret a parsimonious model.
So far, our hypothesis family $H_D$ has no hierarchy and consequently, we have not imposed a sequential structures on the admission order of covariates of $\bm{X}$. However, there are cases where the covariates or factors warrant a natural order such that the family possesses a special testing logic. A hierarchical structure in covariates arises when the inclusion of the next covariate only make sense if the previous covariates is included. One example would be if the next covariates refines a property of the previous covariate. Another case is the use of principal component (PC) factors. The conventional logic is to include PCs sequentially from the dominating one to the least dominating one. This is similar to the motivation for 16-AOS1536, but different from them, we treat the PCs as exogenous without taking the estimation of PCs explicitly into account. In this section, we will use exogenous PCs as hierarchical covariates, as this is the main example in our empirical study. However, all the results hold for any set of exogenous hierarchical covariates.
Without loss of generality, we presume $\bm{X}$ has the $k$th column as the $k$th nested factor. A $k$-order nested model is of the following form
where $[k]=\{1,...,k\}$ is the set that includes indices up to $k$. For example, a hierarchical 3-order model corresponds to the case where variables $\bm{X}_{\{1,2,3\}}$ are included, but not for the rest of the covariates in $\bm{X}$. When formulating our hypothesis family, we must represent the sequential testing structure, as reflected in our definition of nested families of hypotheses:
$H_{N,0}$ completes the case when no rejection on any factor is made. Whenever $H_{N,k}$ is true, then $H_{N,k'}$ is also true for $k<k'\leq J$. Moreover, in the cases where $\mathcal{K}_k=\emptyset$ but $\mathcal{K}_{k'}\neq\emptyset$ with $k<k'$, the notation ensures that the hypothesis $H_{N,k}$ is included in $H_N$ simply because $\mathcal{K}_{k'}$ is present. In other words, if a less dominating hypothesis $H_{N,k'}$ is suggested by data (that is, its active set is non-empty $\mathcal{K}_{k'}\neq\emptyset$), $H_N$ would automatically include all $H_{N,k}$ for $k\leq k'$.
The FWER control property needs to be adapted to the nested nature of this family. 16-AOS1536 argue that the proper measurement is to control for ordered factor count over-estimation with level $\gamma$, as follows:
Given the hierarchical logic embedded in the model, we need the following assumption, which is more restrictive than Assumption (ref):
Assumption (ref) only needs to hold for the tail hierarchical covariates, but requires them to be independent, whereas Assumption (ref) for the unordered tests does not require independence. In the case of PCs, it only applies to the lower order tail PC factors that should not be included for a given null hypothesis. For example, if the true model is $H_{N,5}$, we only need $p^{(n)}_{k'}\stackrel{iid}{\sim} \text{Unif}[0,1]$ for $k'>5$, which is a usual type of assumption in this literature such as in rssb.12122. Moreover, because the nested nature guarantees that the higher-order PCs are more likely to be null, a step-down procedure is expected to increase the power relative to a step-up procedure.
As our focus is to control for false discoveries, we also need to propose new simultaneity counts in the nested family setting. Concretely, we consider first taking a union to obtain the active unit set ${\mathcal{K}}_k^{\text{order}}$ and then calculate conservative simultaneity counts $N^{\text{order}}_k$, with $|M_j|$ as the number of units that the $j$th variable is active in:
It is possible for some $|M_k|$ to be 0 (for instance, the $k$th PC could be inactive for all units), but its ${N}_k^{\text{order}}$ would be 0 if and only if higher-order PCs all have $|M_{k'}|=0$ for $k'>k$. Note, that we are not assuming that $|M_j|$ is decreasing in $j$. Hence, our data-driven approach imposes only a nested structure in the hypotheses, but not a nested structure for the number of active covariates.
Figure (ref) illustrates the process of our step-down simultaneity count. From the left, we start with factor $k=4$ and move step-wise down to factor $k=1$ on the right. The dark blue columns present the active factors, while the light blue columns capture factors of higher-order. In the left-most sub-figure, we only need to account for the 4th PC, implying ${N}_4^{\text{order}}=3$, whereas in the mid-left sub-figure, the 3rd PC has ${N}_3^{\text{order}}=2+3=5$. Eventually, in the right-most sub-figure, we have swept through the entire panel and the 1st PC has a simultaneity count of ${N}_1^{\text{order}}=12$.
Now we can introduce a step-down procedure adapted to the nested structure of $H_{N}$:
This procedure will have FWER control at level $\gamma$ as stated in the following theorem:
The proof is deferred to the Online Appendix. This design extends Procedure 2 from rssb.12122 and “Rank Estimation” from 16-AOS1536, both of which focus on a single sequence of $p$-values rather than the panel setting.
In Step 2, we use Assumption (ref) to transform $p$-values into $\ln(p^{(n)}_k)$, which are i.i.d. standard exponential random variables. Since the family $H_{N}$ has $J$ members, we need to modify our simultaneity count and in a sense condense the panel into a sequence of statistics associated with the ordered covariates. We built a staircase sequence of conservative simultaneity counts ${N}_k^{\text{order}}$ in Step 1 to accumulate the number of $p$-values we use up to the $k$th ordered covariate, starting from the end. By the R\'enyi representation of Rnyi1953OnTT, the ${Z}_k^{\text{order}}$ of Step 2 approximate exponential order statistics and the ${q}_k^{\text{order}}$ approximate uniform order statistics. The nature of these approximations is to create a more conservative rejection, the technical details of which are examined in the proof in our Online Appendix. Finally, we run the order statistics through a step-down procedure proposed by 2336545 so that we find the $\hat{k}$ largest number of ordered covariates rejected by the data with FWER control. Also note that even if the global null, i.e. $H_{N,0}$, is true, and every linear sparse model active set is empty, that is ${N}^{\text{order}}_1=0$, the procedure in Step 3 is still valid because we do not reject $H_{N,1}$.
We demonstrate in simulations that our inferential theory allows us to select better models. We compare different estimation approaches to select covariates and show that our approach better trades off false discovery and correct selections and hence results in a better out-of-sample performance.
Table (ref) summarizes the benchmark models. Our framework contributes among three dimensions: the selection step for the sparse model, the construction of the hypothesis and the multiple testing adjustment. We consider variations for these three dimensions which yields in total six estimation methods. By varying the different elements of the estimators, we can understand the benefit of each component.
Our baseline model is Panel PoSI, which uses post-selection inference LASSO, and a simultaneity count for a data driven hypothesis. The first component that we modify is the selection of the sparse model. A simple OLS regression without shrinkage does not produce a sparse model. This gives us the methods Naive OLS and Bonferroni OLS. A conventional LASSO results in a sparse selection, but the p-values are not adjusted for the post-selection inference and the bias adjustment. The corresponding models are the Naive LASSO and the Bonferroni Naive LASSO. The second component is the hypothesis, which is agnostic for methods besides Panel PoSI. For the comparison models, we either consider no multiple testing adjustment or the conventional Bonferroni adjustment. Under the multiple testing adjustment we obtain the Bonferroni OLS, the Bonferroni Naive LASSO and the Bonferroni PoSI. The outcome of all the estimations are adjusted p-values for the covariates, which we use to select our model for a given target threshold. For a given value of $\gamma$ we include a covariate if its adjusted p-value is below the critical values summarized in the last column of Table (ref).
We simulate a simple and transparent model. Our panel follows the linear model
The covariates and errors are sampled independently as normally distributed random variables:
The noise is either generated as independent noise with covariance matrix $\Sigma= \sigma^2 I$ or as cross-sectionally dependent noise with non-zero off-diagonal elements $\Sigma_{ij}=\kappa$ and diagonal elements $\Sigma_{ii}=\sigma^2$. Note that our theorems for PoSI assume homogeneous noise, while dependent noise violates our assumptions. Hence, the dependent noise allows us to test how robust our method is to misspecification. We set $\sigma^2=2$ and $\kappa=1$, but the results are robust to many choices. The covariance $\Sigma$ as treated as unknown and hence has to be estimated.
We construct the active set based on the staircase structure depicted in Figure (ref). Of the $J$ covariates in $X$, we have $K=10$ active independent factors. Figure (ref) demonstrates the setting for the $10$ factors, where loadings are shaded based on whether they are active. The first factor affects all units, the 2nd factor affects $90\%$, and so on, and lastly the 10th factor affects $10\%$ of all units. This setting is relevant, and also challenging from a multiple testing perspective. It results in a large cohesion coefficient $\rho$, which makes the correct FWER control even more important. The loadings are sampled from a uniform distribution, if they are in the active set:
We simulate a panel of dimension $N=120$, $J=100$ and $T=300$ with $K=10$ active factors. The first half of the time-series observations is used for the in-sample estimation and selection, while the second half serves for the out-of-sample analysis. All results are averages of 100 simulations. We use the covariates selected on the in-sample data for regressions out-of-sample. Our focus is on the inferential theory, and not on the bias correction for shrinkage. Hence, we first use the inferential theory on the in-sample data to select our set of covariates. Second, we use the selected subset of covariates in an OLS regression on the in-sample data to obtain the loadings. Last but not least, we apply the estimated loadings of the selected subset to the out-of-sample data to obtain the model fit. Note that this procedure helps a Naive LASSO, which in contrast to PoSI LASSO does not have a bias correction. The out-of-sample explained variation is measured by $R^2$, which is the sum of explained variation normalized by the total variation. The rejection FWER is set to $\gamma=5\%$ or $\gamma=1\%$. The LASSO shrinkage penalty $\lambda$ is selected by 5-fold cross-validation on the in-sample data.
Table (ref) compares the selection results for the different methods. For each method we report the number of selected covariates, the number of falsely selected covariates and the number of correctly selected covariates. We also report the out-of-sample $R^2$. The upper panel shows the results for independent noise, while the lower panel collects the results for cross-sectionally dependent noise.
PanelPoSI clearly dominates all models. It provides the best trade-off between correct and false selection, which results in the best out-of-sample performance. In the case of $\gamma=5\%$ and independent noise, Panel PoSI selects 10.8 factors in a model generated by 10 factors. 7.9 of these factors are correct. A simple Bonferroni correction is overly conservative. The Bonferroni PoSI selects only 4.7 correct factors. While this overly conservative selection protects against false discovery, it omits over half of the relevant factors which lowers the out-of-sample performance. Using post-selection inference is important, as a naive lasso provides wrong p-values which makes the overly conservative selection even worse. The other extreme is to have neither shrinkage nor multiple testing adjustment. As expected the naive OLS has an extreme number of false selections with a correspondingly terrible out-of-sample performance.
As expected, tightening the FWER control to 1% lowers the number of false rejections, but also the number of correct selections. It reveals again that Panel PoSI provides the best inferential theory among the benchmark models. Panel PoSI selects 7.5 correct covariates, while it controls the false rejections at 1.1. The overly conservative Bonferroni methods select even fewer correct covariates, which further deteriorates the out-of-sample performance. The gap in OOS $R^2$ between Panel PoSI and Bonferroni PoSI widens to 5.4%. All the other approaches cannot be used for a meaningful selection.
Panel PosI performs well, even when some of the underlying assumptions are not satisfied. The lower panel of Table (ref) shows the results for dependent noise. As the dependence in the noise is relatively strong, it can be interpreted as omitting a relevant factor in the set of candidate covariates $X$. Even thought the PoSI theory is developed for homogeneous noise, Panel PoSI continues to perform very well. In contrast, the comparison methods perform even worse, and the Bonferroni approaches select even fewer correct covariates.
Our empirical analysis studies a fundamental problem in asset pricing. We select a parsimonious factor model from a large set of candidate factors that can jointly explain the asset prices of a large cross-section of investment strategies. Our data is standard and obtained from the data libraries of Kenneth French and HouEtAl.
We consider monthly excess returns from January 1967 to December 2021, which results in a time dimension of $T=660$. Our test assets are the $N = 243$ double-sorted portfolios of Kenneth French's data library summarized in Table (ref) in the Appendix. The candidate factors are $J=114$ univariate long-short factors based on the data of HouEtAl. We include all univariate portfolio sorts from their data library that are available for our time period, and construct top minus bottom decile factor portfolios. In addition, we include the five Fama-French factors of FAMA20151 from Kenneth French's data library.
Our analysis projects out the excess return of the market factor. We are interested in the question which factors explain the component that is orthogonal to market movements. Hence, we regress out the market factor from the test assets and use the residuals as test assets. We also do not include a market factor in the set of long-short candidate factors. The original test assets have a market component as they are long only portfolios. Our results are essentially the same when we include the market component in the test assets, with the only difference that we would need to include the market factor as an additional factor in our parsimonious models. The market factor would always be selected by all models as significant, but this by itself is neither a novel nor interesting result.
We present in-sample and out-of-sample results. The in-sample analysis uses the first 330 observations (January, 1967 to June, 1994), while the out-of-sample results are based on the second 330 observations (July, 1994 to December, 2021). As in the simulation, we first use the inferential theory on the in-sample data to select our set of covariates. Second, we use the selected subset of covariates in an OLS regression on the in-sample data to obtain the loadings. Last but not least, we use the estimated loadings on the selected subset of factors for the out-of-sample model. The LASSO penalty $\lambda$ is selected via 5-fold cross-validation on the in-sample data.\footnote{Our cross-validation follows the one-standard-deviation rule for selecting parsimonious models, that is, the largest choice of $\lambda$ within 1 standard error of minimizing the squared errors. This is the default setting of popular implementations like glmnet and argued for in \S3.4 of hastie2009elements. We select $\lambda$ from the grid $\exp(a)\cdot \log J/\sqrt{T}$ with $a=-8,...,8$. This grid choice satisfies the Assumptions in chatterjee2014assumptionless and hence Assumption A.4.} Hence, LASSO represents a first-stage dimension reduction tool, and we need the inferential theory to select our final sparse model.
We allow our selection to impose a prior on two of the most widely used asset pricing models. More specifically, we estimate models without a prior, and two specific priors that impose an infinite weight on the Fama-French 3 factors (FF3) and the Fama-French 5 factors (FF5). This prior as part of PoSI LASSO enforces that the FF3 and FF5 factors are included in the active set. Note that because we work with data orthogonal to the market return, we do not include the market factor in the prior, but only the size and value factors for FF3 and in addition the investment and profitability factor for FF5. We denote these weights by $\omega_{\text{FF3}}$ and $\omega_{\text{FF5}}$. This is an example where the researcher has economic knowledge that she wants to include in her statistical selection method.
We evaluate the models with standard metrics. The root-mean-squared error (RMSE) is based on the squared residuals relative to the estimated factor models. Hence, in-sample the models are estimated to minimize the RMSE. The pricing error is the economic quantity of interest. It is the time-series mean of the residual component of the factor model, and corresponds to the mean return that is not explained by the risk premia and exposure to the factors. In summary, we obtain the residuals as $\hat \epsilon =Y_{t,n}- X_S \hat \beta_S$ for the selected factors, where the loadings are estimated on the in-sample data. The metrics are the RMSE and mean absolute pricing error (MAPE):
In addition to Panel PoSI without and with the FF3 and FF5 priors, we consider the benchmark methods of Table (ref). We compare Panel PoSI (P-PoSI), Panel PoSI with infinite priors on FF3 and FF5 (P-PoSI $\omega_{\text{FF3}}$ respectively $\omega_{\text{FF5}}$), Bonferroni Naive LASSO (B-LASSO), Naive LASSO (N-LASSO), Bonferroni OLS (B-OLS) and Naive OLS (N-OLS). Our main analysis sets the FWER control to the usual $\gamma=5\%$.
Panel PoSI selects parsimonious factor models with the best out-of-sample performance among the benchmarks. For the FWER rate of $\gamma=5\%$ the number of factors differs substantially among the different methods. Panel PoSI selects 3 factors. Imposing infinite priors on FF3 or FF5 results in 4 and 5 factors for P-PoSI $\omega_{\text{FF3}}$ respectively $\omega_{\text{FF5}}$. In contrast, the alternative approaches select too many factors. Bonferroni Naive LASSO includes 10, Naive Lasso 70, Bonferroni OLS 107 and Naive OLS 114. These over-parametrized models lead to overfitting of the in-sample data.
Figure (ref) shows in-sample and out-of-sample RMSE for each set of double-sorts. The composition of the double sorts is summarized in Table (ref) in the Appendix. The in-sample performance in the left subfigure has the expected result that more factors mechanically decrease the RMSE. The important findings are in the right subfigure with the out-of-sample RMSE. The uniformly best performing model is Panel PoSI without any priors. In fact, imposing a prior on the Fama-French factors increases the out-of-sample RMSE. The conventional LASSO and OLS estimates have substantially higher RMSE, which can be more than twice as large.
The Panel PoSI models also explain the average returns the best. In Figure (ref), we compare the mean absolute pricing errors among the benchmarks for each set of double sorts. Importantly, the pricing errors are not used as in objective function of the estimation, and hence the fact that the models with the smallest RMSE explain expected returns is an economic finding supporting arbitrage pricing theory. Our Panel PoSI has the smallest out-of-sample pricing errors, which can be up to six times smaller compared to the OLS estimates. Including the Fama-French factors as a prior does not improve the models, except for the profitability and investment double sort, which uses the same information as two of the Fama-French factors.
The Panel PoSI models select economically meaningful factors. Table (ref) reports the ranking of factors based on their FWER bound without prior and infinite prior weights on the Fama-French 3 and 5 factors. The rows are ordered based on sorted ascending $\rho^{-1} N_jp_j$, which corresponds to the FWER bound. It allows us to infer the number of factors for different levels of FWER control values. Setting $\gamma=5\%$ leads to 3, 4 and respectively 5 factors, while a $\gamma=1\%$ results in 2, 4 and 5 factors, respectively.
In addition to their significance, we can infer the relative importance of factors. The baseline PoSI with $\gamma=5\%$ selects a size, dollar trading volume and value factor. The size and value factors are among the most widely used asset pricing factors. Their selection is in line with their economic importance and confirms the Fama-French 3 factor model. The dollar trading volume factor is less conventional, but is correlated with many assets in our cross-sections. The size factor is the most important as measured by the FWER bound, that is, the product of the number of relevant assets and its minimum p-value are the smallest. The short term reversal factor is less important and would require a FWER control of 10% to be included.
Imposing a prior affects the p-values of PoSI and the simultaneity count. For example, the cohesiveness coefficient increases from $\rho=0.16$ for no priors to $\rho=0.18$ in the case of the two priors. Hence, the FWER bounds of all factors can change when we impose a prior. The FF3 prior increases the significance of the short-term reversal factor, which is widely used in asset pricing. Interestingly, even for a FF5 prior, the profitability and investment factors remain insignificant.
Our method contributes to the discussion about the number of asset pricing factors. Many popular asset pricing models suggest between three and six factors. Our approach allows a disciplined estimate for the number of factors based on inferential theory. The level of sparsity of a linear model also depends on the rotation of the covariates. Therefore, we also study the principal components (PCs) of the covariates $X$ as candidate factors. In this case, we use the step-down procedure, which we refer to as “Ordered PoSI” or O-POSI for short.
Figure (ref) shows the number of factors for different FWER rates $\gamma$. The factor count is obtained by traversing $K^*(\gamma)$ equal to 0.01, 0.02, 0.05 and 0.1. Panel PoSI without priors selects 2 factors for $\gamma=0.01$ and 3 for $\gamma=0.05$. Once, we impose an infinite weight on the Fama-French 3 factors, we select 4 factors for all FWER levels, while the prior on the Fama-French 5 factors results in a 5 factor model for all FWER levels. The Ordered PoSI with PCA rotated factors selects 3 factors for all FWER levels. In summary, our results confirm that depending on the desired significance, the number of asset pricing factors for a good model seems to be between 2 and 4. Note that our analysis is orthogonal to the market factor, which would also be added to the final model. Thus, the final model would have between 3 and 5 factors.
Table (ref) further confirms our findings about the number of asset pricing factors. We compare the number of factors for $\gamma=5\%$ selected either from the univariate high-minus-low factors (HL), their PCA rotation or the combination of the high-minus-low factors and their PCs. Panel PoSI selects consistently 3 factors from the long-short factors and their PCs. When combined, PoSI selects 4 factors, which is plausible as the optimal sparse model can be different for this larger set of candidate factors. The Bonferroni PoSI is overly conservative and selects only 2 HL factors. The models based on Naive LASSO or OLS select excessively many factors independent of the rotation. Overall, the findings support that parsimonious asset pricing models can be described by three to four factors. Of course, any discussion about the number of asset pricing factors is always subject to the choice of test assets and candidate factors.
This paper proposes a new method for covariate selection in large dimensional panels. We develop the conditional inferential theory for large dimensional panel data with many covariates by combining post-selection inference with a new multiple testing method specifically designed for panel data. Our novel data-driven hypotheses are conditional on sparse covariate selections and valid for any regularized estimator. Based on our panel localization procedure, we control for family-wise error rates for the covariate discovery and can test unordered and nested families of hypotheses for large cross-sections. We provide a method that allows us to traverse the inferential results and determine the least number of covariates that have to be included given a user-specified FWER level.
As an easy-to-use and practically relevant procedure, we propose Panel-PoSI, which combines the data-driven adjustment for panel multiple testing with valid post-selection p-values of a generalized LASSO, that allows to incorporate weights for priors. In an empirical study, we select a small number of asset pricing factors that explain a large cross-section of investment strategies. Our method dominates the benchmarks out-of-sample due to its better control of false rejections and detections.