Exact Post Model Selection Inference for Marginal Screening
Jason D Lee, Jonathan E Taylor
Introduction
where is an arbitrary function, and . Our goal is to perform inference on , which is the best linear predictor of . In the classical setting of , the least squares estimator
is a commonly used estimator for . Under the linear model assumption , the exact distribution of is
Using the normal distribution, we can test the hypothesis and form confidence intervals for using the z-test.
However in the high-dimensional setting, the least squares estimator is an underdetermined problem, and the predominant approach is to perform variable selection or model selection . There are many approaches to variable selection including AIC/BIC, greedy algorithms such as forward stepwise regression, orthogonal matching pursuit, and regularization methods such as the Lasso. The focus of this paper will be on the model selection procedure known as marginal screening, which selects the most correlated features with the response .
Marginal screening is the simplest and most commonly used of the variable selection procedures . Marginal screening requires only computation and is several orders of magnitude faster than regularization methods such as the Lasso; it is extremely suitable for extremely large datasets where the Lasso may be computationally intractable to apply. Furthermore, the selection properties are comparable to the Lasso . In the ultrahigh dimensional setting , marginal screening is shown to have the SURE screening property, ), that is marginal screening selects a superset of the truly relevant variables . Marginal screening can also be combined with a second variable selection procedure such as the Lasso to further reduce the dimensionality; our statistical inference methods extend to the Marginal Screening+Lasso method.
Since marginal screening utilizes the response variable , the confidence intervals and statistical tests based on the distribution in (3) are not valid; confidence intervals with nominal coverage may no longer cover at the advertised level:
Several authors have previously noted this problem including recent work in . A major line of work has described the difficulty of inference post model selection: the distribution of post model selection estimates is complicated and cannot be approximated in a uniform sense by their asymptotic counterparts.
In this paper, we describe how to form exact confidence intervals for linear regression coefficients post model selection. We assume the model (1), and operate under the fixed design matrix setting. The linear regression coefficients constrained to a subset of variables is linear in , for some . We derive the conditional distribution of for any vector , so we are able to form confidence intervals and test regression coefficients.
In Section 2 we discuss related work on high-dimensional statistical inference, and Section 3 introduces the marginal screening algorithm and shows how z intervals may fail to have the correct coverage properties. Section 4 and 5 show how to represent the marginal screening selection event as constraints on , and construct pivotal quantities for the truncated Gaussian. Section 6 uses these tools to develop valid hypothesis tests and confidence intervals.
Although the focus of this paper is on marginal screening, the “condition on selection” framework, first proposed for the Lasso in , is much more general; we use marginal screening as a simple and clean illustration of the applicability of this framework. In Section 7, we discuss several extensions including how to apply the framework to other variable/model selection procedures and to nonlinear regression problems. Section 7 covers
marginal screening+Lasso, a screen and clean procedure that first uses marginal screening and cleans with the Lasso,
Related Work
Most of the theoretical work on high-dimensional linear models focuses on consistency. Such results establish, under restrictive assumptions on , the Lasso is close to the unknown and selects the correct model . We refer to the reader to for a comprehensive discussion about the theoretical properties of the Lasso.
There is also recent work on obtaining confidence intervals and significance testing for penalized M-estimators such as the Lasso. One class of methods uses sample splitting or subsampling to obtain confidence intervals and p-values . In the post model selection literature, the recent work of proposed the POSI approach, a correction to the usual t-test confidence intervals by controlling the familywise error rate for all parameters in any possible submodel. The POSI approach will produce valid confidence intervals for any possible model selection procedure; however for a given model selection procedure such as marginal regression, it will be conservative. In addition, the POSI methodology is extremely computationally intensive and currently only applicable for .
A separate line of work establishes the asymptotic normality of a corrected estimator obtained by “inverting” the KKT conditions . The corrected estimator has the form where is a subgradient of the penalty at and is an approximate inverse to the Gram matrix . The two main drawbacks to this approach are 1) the confidence intervals are valid only when the M-estimator is consistent, and thus require restricted eigenvalue conditions on , 2) obtaining is usually much more expensive than obtaining , and 3) the method is specific to regularized estimators, and does not extend to marginal screening, forward stepwise, and other variable selection methods.
Most closely related to our work is the “condition on selection” framework laid out in for the Lasso. Our work extends this methodology to other variable selection methods such as marginal screening, marginal screening followed by the Lasso (marginal screening+Lasso), orthogonal matching pursuit, and non-negative least squares. The primary contribution of this work is the observation that many model selection methods, including marginal screening and Lasso, lead to “selection events” that can be represented as a set of constraints on the response variable . By conditioning on the selection event, we can characterize the exact distribution of . This paper focuses on marginal screening, since it is the simplest of variable selection methods, and thus the applicability of the “conditioning on selection event” framework is most transparent. However, this framework is not limited to marginal screening and can be applied to a wide a class of model selection procedures including greedy algorithms such as matching pursuit and orthogonal matching pursuit. We discuss some of these possible extensions in Section 7, but leave a thorough investigation to future work.
A remarkable aspect of our work is that we only assume is in general position, and the test is exact, meaning the distributional results are true even under finite samples. By extension, we do not make any assumptions on and , which is unusual in high-dimensional statistics . Furthermore, the computational requirements of our test are negligible compared to computing the linear regression coefficients.
Our test assumes that the noise variance is known. However, there are many methods for estimating in high dimensions. A data splitting technique is used in , while proposes a method that computes the regression estimate and an estimate of the variance simultaneously. We refer the reader to for a survey and comparison of the various methods, and assume is known for the remainder of the paper.
Marginal Screening
We will assume that is in general position and has unit norm columns. The algorithm estimates via Algorithm 1.
The marginal screening algorithm chooses the variables with highest absolute dot product with , and then fits a linear model over those variables. We will assume . For any fixed subset of variables , the distribution of is
We will use the notation , where is indexing a variable in the set . The z-test intervals for a regression coefficient are
and each interval has coverage, meaning . However if is chosen using a model selection procedure that depends on , the distributional result (5) no longer holds and the z-test intervals will not cover at the level. It is possible that
Similarly, the test of the hypothesis will not control type I error at level , meaning \Pr\left(\text{rejectH_{0}}|H_{0}\right)>\alpha.
We will illustrate empirically that the z-test intervals do not cover at when is chosen by marginal screening in Algorithm 1.
For this experiment we generated from a standard normal with and . The signal vector is sparse with , , and . The confidence intervals were constructed for the variables selected by the marginal screening algorithm. The z-test intervals were constructed via (6) with , and the adjusted intervals were constructed using Algorithm 3. The results are described in Figure 1. The y-axis plots the coverage proportion or the fraction of times the true parameter value fell in the confidence interval. Each point represents independent trials. The x-axis varies the SNR parameter over the values . From the figure, we see that the z intervals can have coverage proportion drastically less than the nominal level of , and only for SNR= does the coverage tend to . This motivates the need for intervals that have the correct coverage proportion after model selection.
Representing the selection event
Since Equation (5) does not hold for a selected when the selection procedure depends on , the z-test intervals are not valid. Our strategy will be to understand the conditional distribution of and contrasts (linear functions of ) , then construct inference conditional on the selection event . We will use to represent a random variable, and to represent an element of the range of . In the case of marginal screening, the selection event corresponds to the set of selected variables and signs :
for some matrix and vector can be taken to be for marginal screening, but this extra generality is needed for other model selection methods. We will use the selection event and the selected variables/signs pair interchangeably since they are in bijection.
The space is partitioned by the selection events,
The vector can be decomposed with respect to the partition as follows
The previous equation establishes that is a different constrained Gaussian for each element of the partition, where the partition is specified by a possible subset of variables and signs . The above discussion can be summarized in the following theorem.
The distribution of conditional on the selection event is a constrained Gaussian,
The event is in bijection with a pair , and is unconditionally Gaussian. Thus the conditional y\big{|}\{A(S,s)y\leq b(S,s)\} is a Gaussian constrained to the set . ∎
Truncated Gaussian test
This section summarizes the recent tools developed in for testing contrastsA contrast of is a linear function of the form . of a constrained Gaussian . The results are stated without proof and the proofs can be found in .
The conditioning set can be rewritten in terms of as follows:
Moreover, are independent of .
The geometric picture gives more intuition as to why and are independent of . Without loss of generality, we assume and (otherwise we could replace by ). Now we can decompose into two independent components, a 1-dimensional component and an -dimensional component orthogonal to :
The case of is illustrated in Figure 2. Since the two components are independent, the distribution of is the same as . If we condition on , it is clear from Figure 2 that in order for to lie in the set, it is necessary for , where and are functions of .
The distribution of conditioned on is a (univariate) Gaussian truncated to fall between and , i.e.
where . is the normal distribution truncated to lie between and .
In Figure 3, we plot the density of the truncated Gaussian, noting that its shape depends on the location of relative to as well as the width relative to .
The following pivotal quantityThe distribution of a pivotal quantity does not depend on unobserved parameters. follows from Corollary 5.2 via the probability integral transform.
Let denote the CDF of a random variable, and let denote the CDF of , i.e.:
Then is a pivotal quantity, conditional on :
where and are defined in (10) and (11).
Inference for marginal screening
In this section, we apply the theory summarized in Sections 4 and 5 to marginal screening. In particular, we will construct confidence intervals for the selected variables.
To summarize the developments so far, recall that our model (1) says that . The distribution of interest is , and by Theorem 4.1, this is equivalent to , where . By applying Theorem 5.3, we obtain the pivotal quantity
for any , where and are defined in (10) and (11).
In this section, we describe how to form confidence intervals for the components of . The best linear predictor of that uses only the selected variables is , and is an unbiased estimate of . In this section, we propose hypothesis tests and confidence intervals for . If we choose
then , so the above framework provides a method for inference about the variable in the model . This choice of is not fixed before marginal screening selects , but it is measurable with respect to the -algebra generated by the partition. Since it is measurable, is constant on each partition, so the pivot is uniformly distributed on each element of the partition, and thus uniformly distributed for all .
If we assume the linear model for some , , and is full rank, then by the following computation :
In , the screening property for the marginal screening algorithm is established under mild conditions. Thus under the screening property, our method provides hypothesis tests and confidence intervals for .
By applying Theorem 5.3, we obtain the following (conditional) pivot for :
The quantities and are both random through , a quantity which is fixed after conditioning, therefore Theorem 5.3 holds even for this choice of .
Consider testing the hypothesis . A valid test statistic is given by , which is uniformly distributed under the null hypothesis and . Thus, this test would reject when or .
The test of that accepts when
Under , we have , so by (17) F_{\beta_{j},\ \sigma^{2}||\eta_{j}||^{2}}^{[{\cal V}^{-},{\cal V}^{+}]}(\eta_{j}^{T}y)\big{|}\{\hat{E}(y)=E\} is uniformly distributed. Thus
and the type 1 error is exactly . Under , but not conditional on selection event , we have
For each element of the partition , the conditional (on selection) hypothesis test is level , so by summing over the partition the unconditional test is level . ∎
Our hypothesis test is not conservative, in the sense that the type 1 error is exactly ; also, it is non-asymptotic, since the statement holds for fixed and . We summarize the hypothesis test in this section in the following algorithm.
2 Confidence intervals for selected variables
Next, we discuss how to obtain confidence intervals for . The standard way to obtain an interval is to invert a pivotal quantity . In other words, since
one can define a (conditional) confidence interval for as
In fact, is monotone decreasing in , so to find its endpoints, one need only solve for the root of a smooth one-dimensional function. The monotonicity is a consequence of the fact that the truncated Gaussian distribution is a natural exponential family and hence has monotone likelihood ratio in .
We now formalize the above observations in the following result, an immediate consequence of Theorem 5.3.
Let be defined as in (16), and let and be the (unique) values satisfying
Then is a confidence interval for , conditional on :
The confidence region of is the set of such that the test of accepts at the level. The function is monotone in , so solving for and identify the most extreme values where is still accepted. This gives a confidence interval. ∎
In relation to the literature on False Coverage Rate (FCR) , our procedure also controls the FCR.
Furthermore, the FCR of the intervals is .
By (20), the conditional coverage of the confidence intervals are . The coverage holds for every element of the partition , so
We summarize the algorithm for selecting and constructing confidence intervals below.
3 Experiments on Diabetes dataset
In Figure 1, we have already seen that the confidence intervals constructed using Algorithm 3 have exactly coverage proportion. In this section, we perform an experiment on real data where the linear model does not hold, the noise is not Gaussian, and the noise variance is unknown. The diabetes dataset contains diabetes patients measured on baseline variables . The baseline variables are age, sex, body mass index, average blood pressure, and six blood serum measurements, and the response is a quantitative measure of disease progression measured one year after the baseline. The goal is to use the baseline variables to predict , the measure of disease progression after one year, and determine which baseline variables are statistically significant for predicting .
Extensions
The purpose of this section is to illustrate the broad applicability of the condition on selection framework. This framework was first proposed in to form valid hypothesis tests and confidence intervals after model selection via the Lasso. However, the framework is not restricted to the Lasso, and we have shown how to apply it to marginal screening. For expository purposes, we focused the paper on marginal screening where the framework is particularly easy to understand. In the rest of this section, we show how to apply the framework to marginal screening+Lasso, orthogonal matching pursuit, and non-negative least squares. This is a non-exhaustive list of selection procedures where the condition on selection framework is applicable, but we hope this incomplete list emphasizes the ease of constructing tests and confidence intervals post-model selection via conditioning.
The marginal screening+Lasso procedure was introduced in as a variable selection method for the ultra-high dimensional setting of . Fan et al. recommend applying the marginal screening algorithm with , followed by the Lasso on the selected variables. This is a two-stage procedure, so to properly account for the selection we must encode the selection event of marginal screening followed by Lasso. This can be done by representing the two stage selection as a single event. Let be the variables and signs selected by marginal screening, and the be the variables and signs selected by Lasso . In Proposition 2.2 of , it is shown how to encode the Lasso selection event as a set of constraints The Lasso selection event is with respect to the Lasso optimization problem after marginal screening., and in Section 4 we showed how to encode the marginal screening selection event as a set of constraints . Thus the selection event of marginal screening+Lasso can be encoded as . Using these constraints, the hypothesis test and confidence intervals described in Algorithms 2 and 3 are valid for marginal screening+Lasso.
2 Orthogonal Matching Pursuit
Orthogonal matching pursuit (OMP) is a commonly used variable selection method. At each iteration, OMP selects the variable most correlated with the residual , and then recomputes the residual using the residual of least squares using the selected variables. The description of the OMP algorithm is given in Algorithm 4.
Similar to Section 4, we can represent the OMP selection event as a set of linear constraints on .
The selection event encodes that OMP selected a certain variable and the sign of the correlation of that variable with the residual, at steps to . The primary difference between the OMP selection event and the marginal screening selection event is that the OMP event also describes the order at which the variables were chosen. The marginal screening event only describes that the variable was among the top most correlated, and not whether a variable was the most correlated or most correlated.
Since the selection event can be represented as constraints on , the hypothesis test and confidence intervals described in Algorithms 2 and 3 are valid for OMP selected .
3 Nonnegative Least Squares
Non-negative least squares (NNLS) is a simple modification of the linear regression estimator with non-negative constraints on :
Under a positive eigenvalue conditions on , several authors have shown that NNLS is comprable to the Lasso in terms of prediction and estimation errors. The NNLS estimator also does not have any tuning parameters, since the sign constraint provides a natural form of regularization. NNLS has found applications when modeling non-negative data such as prices, incomes, count data. Non-negativity constraints arise naturally in non-negative matrix factorization, signal deconvolution, spectral analysis, and network tomography; we refer to for a comprehensive survey of the applications of NNLS.
We show how our framework can be used to form exact hypothesis tests and confidence intervals for NNLS estimated coefficients. The primal dual solution pair is a solution iff the KKT conditions are satisfied,
Let . By complementary slackness , where is the complement to the “active” variables chosen by NNLS. Given the active set we can solve the KKT equation for the value of ,
which is a linear contrast of . The NNLS selection event is
The selection event encodes that for a given the NNLS optimization program will select a subset of variables . Similar to the case in OMP and marginal screening, we can use Algorithms 2 and 3, since the selection event is represented by a set of linear constraints .
Conclusion
Due to the increasing size of datasets, marginal screening has become an important method for fast variable selection. However, the standard hypothesis tests and confidence intervals used in linear regression are invalid after using marginal screening to select important variables. We have described a method to perform hypothesis and form confidence intervals after marginal screening. The conditional on selection framework is not restricted to marginal screening, and also applies to OMP, marginal screening + Lasso, and NNLS.
Acknowledgements
Jonathan Taylor was supported in part by NSF grant DMS 1208857 and AFOSR grant 113039. Jason Lee was supported by a NSF graduate fellowship, and a Stanford Graduate Fellowship.