When Does the First Spurious Variable Get Selected by Sequential
Regression Procedures?

arXiv:1708.03046v1 [math.ST] 10 Aug 2017

Weijie Su
Department of Statistics, University of Pennsylvania, Philadelphia, PA 19104, USA
August 9, 2017
Abstract
Applied statisticians use sequential regression procedures to produce a ranking of explanatory variables and, in settings of low correlations between variables and strong true effect sizes,
expect that variables at the very top of this ranking are true. In a regime of certain sparsity
levels, however, three examples of sequential procedures—forward stepwise, the lasso, and least
angle regression—are shown to include the first spurious variable unexpectedly early. We derive
a rigorous, sharp prediction of the rank of the first spurious variable for the three procedures,
demonstrating that the first spurious variable occurs earlier and earlier as the regression coefficients get denser. This counterintuitive phenomenon persists for independent Gaussian random
designs and an arbitrarily large magnitude of the true effects. We further gain a better understanding of the phenomenon by identifying the underlying cause and then leverage the insights
to introduce a simple visualization tool termed the “double-ranking diagram” to improve on
sequential methods.
As a byproduct of these findings, we obtain the first provable result certifying the exact
equivalence between the lasso and least angle regression in the early stages of solution paths
beyond orthogonal designs. This equivalence can seamlessly carry over many important model
selection results concerning the lasso to least angle regression.

Keywords. Linear model; Lasso; Least angle regression; Forward stepwise regression; False variable; Signal-to-noise ratio; Bias.

1

Introduction

Consider observing an n-dimensional response vector y that is generated by a linear model
y = Xβ + z,
where X ∈ Rn×p is a design matrix, β ∈ Rp is a vector of regression coefficients, and z ∈ Rn
is a noise term. To find explanatory variables that are associated with the response y, especially
in the setting where p > n, three sequential regression procedures are frequently used: forward
stepwise regression, the lasso (Tibshirani, 1996), and least angle regression (Efron et al., 2004).
These popular methods build a model by sequentially adding or removing variables based upon
some criterion. In a very natural way, a sequential method ranks explanatory variables according
to when the variables enter the solution path. With this ranking of variables in place, a routine

1

practice for forming the final model is to select all variables ranked earlier than a certain cutoff and
discard the rest.
When running a sequential procedure, a practitioner often wishes to understand where along the
solution path noise variables (regressors with zero regression coefficients) start to enter the model.
In particular, when is the first noise variable selected? A better understanding of this problem is
desirable from at least two perspectives. First, the rank of the first noise variable sheds light on the
difficulty of consistent model selection, offering guidelines for selecting important variables. More
precisely, if the rank is about the same size as the sparsity (the total number of nonzero regression
coefficients), we could obtain a model retaining most of the important variables without causing
many false selections using a sequential method, whereas a small rank implies that signal variables
(regressors with nonzero regression coefficients) and noise variables are interspersed early on in
the solution path and, as a result, a false selection must occur long before the power reaches one.
Second, the empirical performance of numerous tools for post-selection inference in linear regression
is, to a large extent, contingent upon whether the first noise variable occurs early or not (Lockhart
et al., 2014; G’Sell et al., 2016; Tibshirani et al., 2016). Insights into the first false variable would
be valuable for improving these tools and developing new ones.
However, despite an extensive body of work on these sequential methods, the literature remains
relatively silent on questions of the first false variable. Existing results address these questions
in a limited setting, mostly characterizing under what conditions all the signal variables precede
the first noise variable, that is, perfect support recovery or, put more simply, selecting the exactly
correct model. Specifically, this set of results guarantees perfect support recovery using a certain
sequential method provided sufficiently strong effect sizes compared to the noise level and a form
of local orthogonality of the design matrix. These results can be found for both the lasso (Zhao
and Yu, 2006; Bickel et al., 2009b; Wainwright, 2009) and forward stepwise (Tropp, 2004; Zhang,
2009; Cai and Wang, 2011).
Figure 1 illustrates a simulation study that examines when the first noise variable gets selected
by the lasso. The design matrix X is of size 2000 × 1800 consisting of independent N (0, 1/2000)
entries, the noise term z is comprised ofpindependent standard normals, and the regression coefficients are set to β1 = · · · = βk = 100 2 log(1800) = 387.2 and βj = 0 for all j > k, with the
sparsity k varying from 10 to 320. Note that the true effect sizes can be practically thought of as
infinitely strong and the sample correlations between the regressors are small due to the independence. Fig. 1 shows that, in the low sparsity regime, the pairs (sparsity, rank) lie close to the 45◦
line (precisely, it is the line y = x + 1). This behavior is equivalent to saying that all the signal
variables are selected prior to any false variables, which is in perfect agreement with a copious body
of theoretical results available in the literature.
Strikingly, once the sparsity exceeds a certain level (around 140 in the example), a phenomenon
that is not explained by existing theory occurs: the average rank of the first noise variable becomes
substantially smaller than the sparsity k and, more surprisingly, the rank keeps decreasing as the
sparsity increases. This phenomenon clearly demonstrates the impossibility of perfect support
recovery in this non-extreme sparsity regime using the lasso, even though it is under high signalto-noise ratios and low correlations. Presumably, as the signal β gets amplified by setting more
components to a large magnitude, one might instinctively anticipate that a sequential method such
as the lasso tends to include more signal variables at the beginning and, thus, would imagine that
the first noise variable would get selected later and later. Unfortunately, the counterintuitive results
as shown in Fig. 1 falsify this belief. We remark that a similar phenomenon is observed earlier in

2

300
250
200
Rank

150

●

● ● ● ●

100

●
●

●
●

● ●
● ●

●
●
●

50

●

● ●
●

● ●
● ● ●
●

●

●
●
●
●

0

●

0

50

100

150

200

250

300

Sparsity

Figure 1: Rank of the first spurious variable along the lasso path. Recall that the rank equals
one plus the number of signal variables preceding the first spurious variable. We plot averages
◦
from 500 independent
p replicates in dots. The 45 dashed line is shown for comparison. The
solid line plots exp( (2n log p)/k − n/(2k) + log(n/(2p log p))) as a function of k, starting from
n/(2 log p).

Su et al. (2017), although they did not provide any justification for the observation.
Thus, concrete predictions and explanations are urged to better understand and improve sequential methods in this sparsity regime. In response, we derive an analytical prediction that is
asymptotically exact for the first noise variable. The prediction applies to the three methods under
our consideration, namely forward stepwise, the lasso, and least angle regression, and potentially
to other sequential methods. Denote by T the rank of the first noise variable. Informally, the
prediction states that, in the setting of strong effect sizes and independent regressors as in Fig. 1,
the three sequential procedures in the non-extreme sparsity regime all satisfy
r
2n log p
n
n
log T ≈
−
+ log
.
k
2k
2p log p
The formal statement of this result is given in Theorem 2 in § 2.
The prediction of T is additionally presented in Fig. 1, showing excellent agreement between
the predicted and observed behaviors. To better appreciate this result, note that the quantity as
an approximation to log T can be equivalently written as
r
r 2

p
2n log p
n
n
n
n
−
+ log
=−
log p −
+ log
.
k
2k
2p log p
2k
2 log p
Once the sparsity k is greater than n/(2 log p), this predicted quantity decreases as the sparsity k
increases. Put differently, the first noise variable is bound to occur earlier as the signal vector β gets
denser, successfully predicting the phenomenon shown in Fig. 1. While problems in selecting the
true model by the lasso have been empirically documented in earlier work (Fan and Song, 2010),
such sharp and analytical predictions are not available in the literature, perhaps due to technical
difficulties.
3

This result has several implications. First, once the underlying signals go beyond the very sparse
regime, using sequential procedures would inevitably lead to a very low power with familywise error
rate control, which is the probability of selecting one or more noise variable, no matter how large
the effect sizes are. Taking a simple example in which both n and p are set to be equal and large
and k = p for some fixed 0 <  < 1, the prediction asserts that the first false variable comes no
later than
h
p
i
h
i
p
exp (1 + o(1))
2(log p)/ − 1/(2) − log(2 log p) = exp (1 + o(1)) 2(log p)/
p
2h log p/ is the leading
i component for a large p). For a fixed , however, the
p
predicted rank exp (1 + o(1)) 2 log p/ only accounts for a vanishing fraction of the k = p
p
signal variables, which can be gleaned from the fact that 2 log p/ = o(log p). In other words,
the three sequential methods being considered yield vanishing power if no noise variable is allowed
to be included, even in the noiseless case (z = 0). In particular, these negative results are derived
under Gaussian designs with independent columns, which have vanishing sample correlations and
satisfy some conditions believed to be favorable for model selection, including restricted isometry
properties (Candès and Tao, 2005) and restricted eigenvalue conditions (Bickel et al., 2009a). Thus,
the negative results are likely to carry over to a much broader class of design matrices. In fact,
extensive simulations carried out in § 3 demonstrate that problems of the first false variable are
only exacerbated in more general settings.
Another implication yielded by this prediction is that the three sequential regression methods
seem to behave similarly in ranking variables, at least in the independent random design setting.
Compared with forward stepwise, the lasso and least angle regression, along with their infinitesimal
version forward stagewise regression (see, for example, Efron et al. (2004)), are long-time considered less greedy because at each step they gradually blend in a new variable instead of adding it
discontinuously (Efron et al., 2004). To be more precise, the forward stepwise selects the predictor
with the largest absolute correlation with the residual vector and then aggressively takes a large
step in the direction of the selected predictor, whereas the latter three proceeds in a more democratic manner along a direction equiangular between the set of selected predictors (the lasso and
forward stagewise regression bear certain restrictions on this equiangular approach). This critical
distinction between the two strategies is anticipated—or, at least wished—to lead to contrasting
model selection performance. Interestingly, this is not the case; these two strategies yield the same
behavior of selecting the first noise variables in our setting. As a byproduct, we obtain Theorem
3 for the lasso and least angle regression, which, to the best of our knowledge, is the first mathematically provable result certifying the exact equivalence between early solution paths of these two
procedures beyond orthogonal designs.
In the non-extreme sparsity regime, what are underlying these different sequential methods
responsible for selecting the first false variable so early? Taking a closer look at the derivation of
the prediction, we can identify the cause, which is, loosely speaking, due to the greedy nature of
these sequential regression methods. Moreover, the equiangular strategy adopted by the lasso and
least angle regression fails to alleviate greediness from the perspective of when the first noise variable
gets selected. To shed light on this cause, recall that all the three methods at each time includes
a variable that roughly has the largest absolute inner product with the current residuals. As the
regression coefficients get denser, the solution at the beginning of the path is overwhelmingly biased
and the residual vector absorbs many of the true effects contributed by the nonzero components of

steps (note that

4

β. As a result, some irrelevant variable would exhibit high correlations with the residuals and are
hence selected incorrectly and early. That being said, it requires several novel ideas to precisely
characterize what we describe here.
With this underlying cause in mind and to improve on sequential methods, we introduce the
double-ranking diagram to identify early false variables along solution paths. In slightly more
detail, this diagram contrasts the rank of each variable given by a sequential procedure (horizontal
axis) with that given by a low-biased estimator (vertical axis) such as the least-squares estimator.
In spite of a significant horizontal rank, an early noise variable might be revealed by its possibly
less significant vertical rank. Related ideas have appeared in recent variable screening work, for
instance, Wang and Leng (2016). We demonstrate the usefulness of this diagram via a mix of
theoretical and empirical results.

2
2.1

Understanding the Phenomenon
Predicting the first spurious variable

We consider a sequence of problems indexed by (kl , nl , pl ), where kl , nl , and pl are all assumed to
grow to infinity as l → ∞ in asymptotic statements. The subscript l is often omitted when clear
from the context. Letters ci and Ci in various settings denote positive constants that do not depend
on the problem index l. Below we formalize our working hypothesis concerning the linear model
y = Xβ + z.
Assumption 1. The design X
of independent N (0, σ 2 ) errors.
vector β has k fixed components
c1 p/ logc2 p ≤ n ≤ c3 p and c4 n
c1 , c2 , c3 , c4 , and c5 .

∈ Rn×p has independent N (0, 1/n) entries and z ∈ Rn consists
We further assume X and z are independent. The coefficient
equal to some M 6= 0 and the rest are all zero. Last, we assume
≤ k ≤ min{0.99p, c5 n log0.99 p} for arbitrary positive constants

This assumption implies that each column of X is approximately normalized, having about unit
Euclidean norm. The nonrandom parameters σ ≥ 0 and M both implicitly depend on the index l
and thus are allowed to vary freely. In particular, the noiseless case σ = 0 is not excluded, in which
the signal-to-noise ratio is essentially infinite. The constraint on (k, n, p) includes some popular
models studied in the literature as special examples, for instance, the linear sparsity framework
where k/p and n/p converges to some constants (Bayati and Montanari, 2012). For completeness,
the number 0.99 can be replaced by any positive number smaller than 1.
Before presenting our main results Theorems 1 and 2, we give a brief overview of the three
methods for ease of reading. In broad outlines, least angle regression increases the coefficients of
included variables in their joint least squares direction until an unselected variable has as much
inner product with the residuals, which is included in the next step. Least angle regression stops
when the residuals are zero or all variables are included. If a nonzero coefficient is removed from
the active set whenever it hits zero, this adjustment leads to the lasso, which is better-known as
the minimizer of the convex program 21 ky − Xbk22 + λkbk1 over b ∈ Rp for λ ranging from infinity
to zero. Forward stepwise is described in detail in § 8.5 of Weisberg (1980). In our setting, the
intercept is not included and normalization is not applied to any columns of X. Last, recall that
T denotes the rank of the first noise variable, and oP (1) denotes a sequence of random variables
converging to zero in probability.

5

Theorem 1. Under Assumption 1, the first spurious variable selected by each of forward stepwise,
the lasso, and least angle regression satisfies
i
hp
2n(log p)/k − n/(2k) + log(n/(2p log p)) .
log T ≤ (1 + oP (1))
If the signal magnitude M is not sufficiently large compared with σ, the logarithm of T would
be much smaller than the upper bound appearing in the display above, meaning that the problems
of the first false variables could be worse. Interestingly, this bound is sharp when M is sufficiently
large compared with σ, as demonstrated in the theorem below.
Theorem 2. Under Assumption 1 and in addition provided that σ/M → 0, then each of the three
sequential methods obey
hp
i
log T = (1 + oP (1))
2n(log p)/k − n/(2k) + log(n/(2p log p)) .
Provided in the Appendix, the proofs of both theorems involve some techniques that are likely
to extend beyond the three sequential
p methods. The condition concerning the ratio between σ and
M can be relaxed to |M |/σ  n/k. Setting k ≈ n log0.99 p as in Assumption 1, for example,
Theorem 2 follows if M/σ is bounded away from 0. An immediate consequence of this theorem is
as follows.
Corollary 2.1. Under Assumption 1, each of the three methods in the noiseless case (σ = 0) obeys
hp
i
log T = (1 + oP (1))
2n(log p)/k − n/(2k) + log(n/(2p log p)) .
In addition to predicting the phenomenon observed in Fig. 1, Theorem 2 together with Corollary
2.1 demonstrates that having an even stronger signal magnitude does not affect T much as long as
it exceeds a certain level.
The theorems presented here differ from results that are found extensively in the literature
claiming a high probability of selecting the exactly correct model, mainly due to assuming different sparsity regimes of the regression coefficients β. Explicitly, the former assumes c4 n ≤ k ≤
min{0.99p, c5 n log0.99 p} whereas the latter often, if not always, assumes the extreme sparsity regime
k = O(n/ log p) or k  n/ log p. On the one hand, an interesting question for future research is
whether the lower bound c4 n can be improved to n/(2 log p). On the other hand, a number of cases
leading to k = o(p) satisfy Assumption 1, for instance, n = c1 p/ logc2 p and k = c4 n. Indeed, under
Assumption 1, it is even unrealistic to expect perfect model selection using sequential methods:
below a simple corollary of Theorem 1 shows the number of signal variables before the first false
variable only accounts for an insignificant fraction of the total number of signal variables.
Corollary 2.2. Under Assumption 1, each of the three methods satisfies
T
−→ 0 in probability.
k
To better appreciate Corollary 2.2, consider the scenario where k/p →  and n/p → δ for some
positive constants  < 1p
and δ. Theorem 1 shows that, up to a vanishing fraction,
the logarithm
p
of T is no larger than 2δ(log p)/ − δ/(2) + log(δ/(2 log p))
=
(1
+
o(1))
2δ(log
p)/. This
p
expression for approximating log T yields T ≤ exp((1 + o(1)) 2δ(log p)/)  p = k, confirming
Corollary 2.2 in this linear sparsity regime. We summarize the finding in the corollary below.
6

Corollary 2.3. Under Assumption 1 and additionally provided that k/p →  and n/p → δ for
arbitrary positive constants  < 1 and δ, each of the three methods satisfies
√
T ≤ e(1+oP (1)) 2δ(log p)/ .
This regime of linear sparsity is previously employed in Su et al. (2017), which studies limitations
of the lasso for the false discovery rate control. The techniques developed there are not applicable
to studying the first noise variable, which is a much finer problem.

2.2

Equivalence between lasso and least angle regression

In contrast to the other two methods, the lasso would drop a selected variable if its coefficient
hits zero. This irregularity of the lasso path might lead to ambiguity in interpreting the rank T
in Theorems 1 and 2. Fortunately, as a byproduct of the above, the theorem below rules out the
possibility of such ambiguity.
Theorem 3. Assume X has independent N (0, 1/n) entries. Then, with probability at least 1−p−2 ,
no drop-out occurs before the first
nl p
m o
min c n/ log p , p
variables along the lasso path are selected, where dxe denotes the least integer greater than or equal
to x and c > 0 is some universal constant.
Note that Theorem 3 only requires the normality of X, as opposed to additional conditions
imposed on z, β, and (k, n, p) in Assumption 1. As seen from its proof in the Appendix, the validity
of the theorem depends on the design matrix X basically only through its restricted isometry
property. Thus, this result can seamlessly carry over to other matrix ensembles with an appropriate
restricted isometry property constant,psuch as Bernoulli random
and Tao, 2005).
p matrices (Candès
p
Under Assumption 1, log(min{dc n/ log pe, p}) = logdc n/ log pe  2n(log p)/k − n/(2k).
Consequently, Theorem 3 together with Theorem 1 ensures that the first noise variable selected by
the lasso is not preceded by any drop-out with probability approaching one.
This byproduct provides new insights into the lasso path and is a contribution of independent
interest to high-dimensional statistics. The lasso is known to coincide exactly with least angle
regression until the first time the lasso drops a selected variable (Efron et al., 2004; Tibshirani and
Taylor, 2011). To our knowledge, however, the question of where along the path the lasso and least
angle regression differ has not been addressed in prior research, perhaps due to technical difficulties.
By confirming the equivalence between these two procedures, Theorem 3 allows us to carry over
well-known results on the lasso for model selection to least angle regression.

2.3

Heuristics and insights

In this section we give an informal derivation of Theorems 1 and 2. Although our discussion below
lacks rigor, nevertheless the goal is to gain insights into this counterintuitive phenomenon. For the
full proofs, see the Appendix.
We focus on the noiseless case z = 0 in Assumption 1, which is presumably the most ideal
scenario for model selection. Denote by S = {j : βj 6= 0} the support of the signals and by βb an
7

estimate given by any of the three sequential methods somewhere along the solution path. Write
j1 ∈
/ S for the index off the support having the largest inner product in magnitude with the residual
b and j2 ∈ S for the index on the support having the T − 1st largest inner
y − X βb = X(β − β),
product in magnitude with the residual. By using some technical arguments found in the Appendix,
we get
r
r
2k
log(p
−
k)
b ≈M
b ≈ M + M 2k log(k/T ) .
Xj>1 X(β − β)
, Xj>2 X(β − β)
(2.1)
n
n
Above and henceforth, X > denotes the transpose of X. Recognizing that the sequential methods rank variables essentially according to the correlations with the residual, where in our case
correlations are roughly equivalent to inner products since the columns of X are approximately
normalized, from (2.1) we must have
r
r
2k log(p − k)
2k log(k/T )
M
≈M +M
n
n
at the point where the first false variable is just about to enter the model. In the linear sparsity
regime k/p → , n/p → δ, this yields
r
2δ log p
log T ≈
.

The exposition above suggests that an early spurious variable is mainly due to a large inner
b which would not be the case if βb was a low-biased estimator of β. However,
product Xj>1 X(β − β),
until a significant proportion of the variables have been selected, a solution βb provided by a sequential method is overwhelmingly biased. Another way to formalize this point is that the residual
b still contains a significant amount of true effects, largely contributed by presently unseX(β − β)
lected variables. This bias acts as if it were noise and, as a consequence, some irrelevant variables
happen to correlate highly with the residual vector, leading to false variables selected early. This
is not a matter of the signal-to-noise ratio; an increasing signal magnitude would enlarge the bias
as well and, hence, noise variables always occur early. Other examples of pseudo noise caused
by bias have been observed in previous work (Bayati and Montanari, 2012). To be complete, we
remark that this phenomenon no longer appears in regimes of extreme sparsity (see, for example,
Wainwright (2009)).

3
3.1

Illustrations
Numerical examples

We present simulation experiments to illustrate the first false variable of the three sequential methods, along with the predictions given by Theorems 1 and 2. Specifically, we numerically examine
three studies concerning the effect of design matrix shapes, signal magnitudes, and correlations
between the columns of X, respectively, on the first spurious variable. Two scenarios are experimented for each study.
Study 1. In the first experiment (square design) the design X of size 1000×1000 has independent
N (0, 1/1000) entries, the signals βj = 100 for j ≤ k and βj = 0 for j ≥ k + 1, and each noise
8

component zi follows N (0, 1) independently. In the second experiment (fat design) the design
X is
√
changed to be of 800 × 1200 and has
√ independent Bernoulli entries, which take value 1/ 500 with
probability half and otherwise −1/ 500, while all the other assumptions remain the same. Results
of the two experiments are shown in Figure 2(a) and (b), respectively.
Study 2. In both experiments, the 500 × 1000 design matrix X consists of independent
N (0, 1/500) entries and each zi is independently distributed as N (0, 1). For the first experiment
(one mixture), we set βj = M for j = 1, . . . , 80 and βj = 0 for j √
= 81, . . . , 1000. For the second
one (two mixtures), we set βj = M for j = 1, . . . , 40, βj = M 2 /(10
√ 2 log p) for√j = 41, . . . , 80 and
βj = 0 for j = 81, . . . , 1000. The parameter M is varied
from 0.2 2 log p to 10 2 log p. Note that
√
the two mixtures take the same value when M = 10 2 log p = 37.17. Results are shown in Fig. 2(c)
and (d).
√
Study 3. This scenario uses β obeying βj = 100 2 log p for j ≤ 80 and βj = 0 otherwise. The
noise z consists of independent standard normals. The 500 × 1000 design matrix X has each row
independently drawn from N (0, Σ). For the 1000 × 1000 covariance matrix Σ, the first experiment
(equi correlation) assumes Σij = ρ/n if i 6= j and Σjj = 1/n. In the second one (decaying
correlation), Σij = ρ|i−j| /n. Results are shown in Fig. 2(e) and (f).
Both the lasso and least angle regression closely matches our predictions as indicated in red
lines. Notably, the two procedures yields exactly the same ranks of the first noise variables, hence
supporting Theorem 3. On the other hand, forward stepwise exhibits larger departures from the
theoretical predictions mainly due to the slow convergence to the asymptotics, while as well showing
a decreasing rank once the sparsity exceeds a cutoff.
As shown in Fig. 2(a) and (b), the first false variable occurs earlier as n decreases while p gets
larger. In particular, the behaviors of the methods under Bernoulli random designs as in Fig. 2(b)
closely resemble that under Gaussian random designs. In Fig. 2(c), the rank of the first false
variable increases as the signal magnitude M is amplified. While this increasing rank is expected,
Fig. 2(d) in contrast illustrates a rather surprising
√ phenomenon: the rank drops after M exceeds
a certain level. More precisely, given M ≥ 3.4 2 log p = 12.64, the lasso selects the first false
variable earlier and earlier even though the sparsity is fixed and each signal gets strengthened,
and the phenomenon is even more transparent for forward stepwise. Intuitively, this is because
the effective sparsity in the case of a moderately large M is smaller than the nominal sparsity
80. To see this, observe that the ratio of the signals of the first 40 components and the next 40
components is M/(M 2 /37.17) = 37.17/M , which is noticeably larger than 1. Put another way, the
first 40 components act like the main signals and, hence loosely speaking, the effective sparsity is
smaller than 80. In the presence of significant correlations between columns of X, Fig. 2(e) and
(f) clearly show that the problem of early false variables is further exacerbated.

3.2

HIV data

As a real data example, we consider the HIV-1 data introduced by Rhee et al. (2006) to study the
genetic basis of HIV-1 resistance to several drugs. Also used in a number of other works (Barber
and Candès, 2015; G’Sell et al., 2016; Janson and Su, 2016), this data set in particular contains
genotype information X ∈ R634×463 of 634 HIV-1-infected individuals across 463 locations after
removing duplicate and missing values. The columns of X are further standardized to have zero
mean and unit Euclidean norm. The response y is synthetically generated by assigning an effect of

9

100
●

●

●

Rank

●

●

●

●

●

●

50

●

●

●

●

●

●

●

●

●

50

Rank

100

150

(b)

150

(a)

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

0

●

0

●

0

50

100

150

0

50

100

Sparsity

150

Sparsity

60
50
●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

40

Rank
●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●
●

●

●●

0

0

●

●

●

●
●

●●

●
●

10

●

●

20

30
20
10

●●
●●

30

40

50

60

70

(d)

70

(c)

0

10

20

30

0

10

20

M

30

M

30
●

20

Rank

30

40

(f)

40

(e)

20

●

●

●

Rank

●

●

●

Rank

●

●
●

●

●

●

●

●

●
●

●
●

●

●

10

10

●
●

●
●
●

●

●
●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

0

●

0

●

0.00

0.02

0.04

0.06

0.08

0.10

0.0

Rho

0.2

0.4

0.6

0.8

Rho

Figure 2: Rank of the first spurious variable in three studies. Averaged over 500 replicates,
the ranks of forward stepwise, the lasso, and least angle regression are marked with triangles,
dots, and crosses, respectively (the dots and the crosses overlap exactly so they look like solid
squares). The solid lines indicate the predictions given by Theorem 2. Note that for (c,d,e,f)
the prediction is constant since k = 80 is fixed.

√
100 2 log p to each of k uniformly randomly chosen columns of X and setting a noise level σ set
to 1.
Table 1 reports the results averaged over 500 replicates. The three methods start to have a
decreasing rank around k = 25, which is much smaller than n/(2 log p) = 51.6. In addition, for
each level of sparsity, the first spurious variable is included much earlier than the predictions. This
10

gap is not surprising given that the predictions are tailored to independent Gaussian designs while
the design X from the HIV-1 data has strongly correlated columns. To be more precise, about
4600 column pairs of X have correlations greater than 10%.
Table 1: Rank of the first selected noise predictor averaged over 500 runs, with standard
errors given in parentheses. The predictions for sparsity no larger than n/(2 log p) = 51.6 are
just given as k + 1 and otherwise are given by Theorems 1 and 2. The last row presents the
predictions.
Lasso
Least angle regression
Forward stepwise
Gaussian designs

4

10
10.4 (2.1)
10.4 (2.1)
10.6 (1.9)
11.0

25
15.4 (9.6)
15.4 (9.6)
18.8 (10.0)
26.0

40
10.3 (9.2)
10.3 (9.2)
18.5 (15.9)
41.0

55
6.8 (6.0)
6.8 (6.0)
13.3 (16.0)
51.3

70
5.5 (4.6)
5.5 (4.6)
10.0 (11.7)
45.7

85
4.4 (3.8)
4.4 (3.8)
7.5 (8.2)
38.3

100
3.7 (3.5)
3.7 (3.5)
7.0 (7.6)
31.8

Visualizing Early Noise Predictors

As discussed in § 2.3, the three sequential procedures are marginal correlation-based at the beginning of their solution paths, picking variables essentially according to the correlations with the
residuals. In light of this viewpoint, an unbiased or low-biased estimator of the signals β might provide sequential methods with complementary information for variable selection. The least-squares
estimator, if available, is a natural candidate.
We introduce the double-ranking diagram to bring together the strengths of sequential methods
and low-biased estimators such as the least-squares estimator βbLS . Figure 3 presents two instances
of this diagram: one is in the same setting as Fig. 1 except for a different size 200 × 180 and a
fixed sparsity k = 50, and another is in the same setting as Table 1 with a fixed sparsity k = 60.
For each variable, the horizontal axis represents its rank by a sequential method, and the vertical
axis represents its rank by a low-biased estimator. For
p example, the horizontal rank of the jth
LS
b
variable is given according to the magnitude of |βj |/ [(X > X)−1 ]jj : the larger this statistic is,
the smaller the rank is. Equivalently, the variables can be ranked using the t-values.
The double-ranking diagram can serve as a simple data visualization tool to assist the identification of early false variables for sequential regression methods. Intuitively, an important variable
would presumably possess both a small horizontal rank and a small vertical rank, hence appearing
in the bottom-left corner of the diagram with a good chance. In light of this intuition, we screen
out variables that are selected early by a sequential method but have unusually large vertical ranks,
which in the case of least squares amount to small t-values or insignificant p-values. As seen from
Fig. 3, the first five false variables in each instance have much larger vertical ranks compared with
their horizontal ranks. In particular, these false variables are placed far away from the signal variables in the diagram. In view of this example, to use this diagram, one can set some threshold
for the vertical rank and only select variables that are below the threshold and in addition have
significant horizontal ranks. On the other hand, in the low signal-to-noise ratio regime the diagram may not give a clear-cut separation between false and true predictors, and its use requires
some caution. The following simple theorem states that the diagram can perfectly separate the
first spurious variable from all the true variables using the least-squares estimator under certain
conditions.

11

50

100

300
0

100

Vertical rank

150
100
0

50

Vertical rank

0

●

200

●
●●
●
● ● ●
● ●
●
●
●
●
●●●
●
●
●● ●
● ●●
● ●
●
●
● ●
●
●●●
●
●
●
● ●
●
● ●
●●
● ●
●
●
●●
●
●
●
●
●
●
●
●● ●
● ●
●
● ●
●
●
●
●
●
●●●
●
●
●
●
●● ●
●
●
●
●
●
●
●
●
●
●
● ●
●
●
●
● ●
●
●
● ● ●
●
●
● ● ●●
● ●● ● ●●
● ● ●
●
●●
●●
●●
●
●
● ● ●
●
●●
●
● ●●
● ● ●● ●●
● ● ●
●
●
●
● ● ●
●●●
●
● ●●●● ●
●●
●●● ● ●
●
●●●

400

●

●

150

●
● ●
●●
●●
●
●
●
● ●
● ●
●
●
●● ●●● ●
●
● ●
● ●
●●
●
● ●●
●●
● ● ●● ● ●
●
●
●
●
●
● ●
● ●
●●●
●● ●
●
●
● ●●
●
● ●
●
●●
● ●
● ● ● ●●
● ●●●●●
●
● ●
●
●
●
●
●●
● ●●
●●
● ● ●
●
●● ● ●
●
●
●
● ●
●
●
●
●
●
●
●
●
●
●
●●
●
● ● ● ● ●● ● ●
●
● ●
●
● ●● ● ●
●
●
●
●
●
●
●● ● ● ● ● ●
●
● ●● ●
●
● ●●
●
● ●
●●●
●
●
●
●
●
●
●●
●
●●
●
●
●
●
●
● ●
● ● ●
● ●
● ● ●●
●●
●
●
● ● ● ●● ● ● ● ● ●
●
●
●
●
●
● ● ●●●
●●
● ●● ● ● ● ●
●
●
●
●
●
●
●
● ●
●
● ● ● ●●● ●
●
●
●
●
●●
●● ●
● ●● ● ●
●
●
●
● ●● ●
●
● ●
●
●
●
●●
● ●
● ● ● ● ●● ●
●
●
● ●● ●●
● ●● ● ●● ●●
●
●
●
● ●●
●
●
●
●
●
●
● ● ●●
●
●●
●● ●
●● ● ● ● ● ●
●
● ●● ● ●
●● ●
●
●
●
● ●
●●
●
●●● ●
● ● ●●
●
● ●
●
●
●
●● ●● ● ● ● ● ● ●● ●●
●
● ● ●
●
●● ● ●
● ● ● ●● ● ● ● ● ●
●
●●
●
●
●●
●●●
●● ●
●●●
●●
●
●
●
● ●●●
●●
●
●
● ●● ●
●
●
●
●●
●
●
● ●●
●●
●
●●
●
●●
●●●●
●
● ● ●

0

●
●●
●● ●

100

Horizontal rank

200

300

400

Horizontal rank

Figure 3: Double-ranking diagrams. Left: under the same setting as Fig. 1. Right: under
the same setting as the HIV data example. Vertical rankings are given by the least-squares
estimators and horizontal rankings are given by the least angle regression. The first five noise
variables along the solution path of least angle regression are marked with crosses.

Theorem 4. Under Assumption 1 and provided that n > δp and
r
M
2δ log p
>3
σ
δ−1
for some constant δ > 1, then in the double-ranking diagram the first noise variable has a greater
vertical rank than all the true variables.
The main ingredient behind this tool is a blend of new ideas and old ones found in the literature.
On the one hand, our discussion in § 2.3 demonstrates that, while sequential methods work well
in very sparse settings, as the signals get denser, the pseudo noise can quickly accumulate and
thus may dwarf some true signals no matter how strong the corresponding coefficients are. On
the other hand, the method of least squares favors the case of dense signals since the estimator
variances basically stay the same as the sparsity of the signals increases. In particular, variables with
sufficiently strong effects can always stand out in the least-squares estimators even in the presence
of significantly correlated columns of designs matrices. This property of the least-squares estimator
and its variants plays a pivotal role for a number of variable screening procedures appearing in
recent work (Wasserman and Roeder, 2009; Pokarowski and Mielniczuk, 2015; Wang and Leng,
2016).

5

Discussion

In the regime of non-extreme sparsity, the common intuition that sequential regression procedures
find a significant portion of all important variables before the first false variable merits some skepticism. We have developed sharp predictions that disprove this intuition for forward stepwise, the
lasso, and least angle regression under independent Gaussian designs, which satisfy certain desirable properties for model selection. Additionally, the predictions hold irrespective of how strong
12

the effect sizes are. Thus, the first noise variable are likely to occur very early in more general
settings. Our numerical results are in agreement with this viewpoint.
In light of the above, more caution is required when using these sequential methods, unless the
true regression coefficients are very sparse. Useful information for identifying early noise variables
can be provided by low-biased methods such as the least-squares estimators. The double-ranking
diagram is a simple tool that unifies the strengths of the two groups of methods.
Avenues for further investigation are in order. First, it is of interest to improve the predictions
for forward stepwise and extend the predictions to more sequential methods such as backward
stepwise and forward-backward stepwise. The simulation studies imply that the lower bound c4 n
on the sparsity k in Assumption 1 could be possibly relaxed to n/(2 log p). Second, in the highdimensional setting where p > n, which low-biased estimator should we choose for the doubleranking diagram to yield the vertical ranking? Candidates worth considering include the lasso with
a small penalty, ridge regression with a small penalty, generalized least-squares estimators (see, for
example, Wang and Leng (2016)), and some recently proposed ranking procedures such as in Ke
and Yang (2017). Last, as seen from Table 1, the rank of the first noise variable has relatively large
standard errors. A question of practical relevance is to characterize this large variation.

References
Richard Baraniuk, Mark Davenport, Ronald DeVore, and Michael Wakin. A simple proof of the restricted
isometry property for random matrices. Constructive Approximation, 28(3):253–263, 2008.
Rina Foygel Barber and Emmanuel Candès. Controlling the false discovery rate via knockoffs. The Annals
of Statistics, 43(5):2055–2085, 2015.
Mohsen Bayati and Andrea Montanari. The LASSO risk for Gaussian matrices. IEEE Transactions on
Information Theory, 58(4):1997–2017, 2012.
Peter Bickel, Ya’acov Ritov, and Alexandre Tsybakov. Simultaneous analysis of Lasso and Dantzig selector.
The Annals of Statistics, pages 1705–1732, 2009a.
Peter J Bickel, Ya’acov Ritov, and Alexandre B Tsybakov. Simultaneous analysis of lasso and dantzig
selector. The Annals of Statistics, pages 1705–1732, 2009b.
T Tony Cai and Lie Wang. Orthogonal matching pursuit for sparse signal recovery with noise. IEEE
Transactions on Information Theory, 57(7):4680–4688, 2011.
T Tony Cai, Jianqing Fan, and Tiefeng Jiang. Distributions of angles in random packing on spheres. Journal
of Machine Learning Research, 14(1):1837–1864, 2013.
Emmanuel Candès and Terence Tao. Decoding by linear programming. IEEE Transactions on Information
Theory, 51(12):4203–4215, Dec. 2005.
Laurens de Haan and Ana Ferreira. Extreme value theory: An introduction. Springer Science & Business
Media, 2007.
Bradley Efron, Trevor Hastie, Iain Johnstone, and Robert Tibshirani. Least angle regression. The Annals
of Statistics, 32(2):407–499, 2004.
Jianqing Fan and Rui Song. Sure independence screening in generalized linear models with NPdimensionality. The Annals of Statistics, 38(6):3567–3604, 2010.

13

Max Grazier G’Sell, Stefan Wager, Alexandra Chouldechova, and Robert Tibshirani. Sequential selection
procedures and false discovery rate control. Journal of the Royal Statistical Society: Series B (Statistical
Methodology), 78(2):423–444, 2016.
Daniel Hsu, Sham Kakade, and Tong Zhang. A tail inequality for quadratic forms of subgaussian random
vectors. Electronic Communications in Probability, 17, 2012.
Lucas Janson and Weijie Su. Familywise error rate control via knockoffs. Electronic Journal of Statistics,
10(1):960–975, 2016.
Zheng Tracy Ke and Fan Yang. Covariate assisted variable ranking. arXiv preprint arXiv:1705.10370, 2017.
Richard Lockhart, Jonathan Taylor, Ryan J Tibshirani, and Robert Tibshirani. A significance test for the
lasso. The Annals of Statistics, 42(2):413, 2014.
Piotr Pokarowski and Jan Mielniczuk. Combined `1 and greedy `0 penalized least squares for linear model
selection. Journal of Machine Learning Research, 16(5), 2015.
Soo-Yon Rhee, Jonathan Taylor, Gauhar Wadhera, Asa Ben-Hur, Douglas L Brutlag, and Robert W Shafer.
Genotypic predictors of human immunodeficiency virus type 1 drug resistance. Proceedings of the National
Academy of Sciences, 103(46):17355–17360, 2006.
Weijie Su, Malgorzata Bogdan, and Emmanuel Candès. False discoveries occur early on the Lasso path. The
Annals of Statistics, to appear, 2017.
Robert Tibshirani. Regression shrinkage and selection via the lasso. Journal of the Royal Statistical Society.
Series B (Methodological), pages 267–288, 1996.
Ryan Tibshirani, Jonathan Taylor, Richard Lockhart, and Robert Tibshirani. Exact post-selection inference
for sequential regression procedures. Journal of the American Statistical Association, 111(514):600–620,
2016.
Ryan J Tibshirani and Jonathan Taylor. The solution path of the generalized lasso. The Annals of Statistics,
39(3):1335–1371, 2011.
Joel A Tropp. Greed is good: Algorithmic results for sparse approximation. IEEE Transactions on Information Theory, 50(10):2231–2242, 2004.
Roman Vershynin. Introduction to the non-asymptotic analysis of random matrices. In Compressed Sensing,
Theory and Applications, 2012.
Martin J Wainwright. Sharp thresholds for high-dimensional and noisy sparsity recovery using `1 -constrained
quadratic programming (Lasso). IEEE transactions on information theory, 55(5):2183–2202, 2009.
Xiangyu Wang and Chenlei Leng. High dimensional ordinary least squares projection for screening variables.
Journal of the Royal Statistical Society Series B, 78(3):589–611, 2016.
Larry Wasserman and Kathryn Roeder. High dimensional variable selection. The Annals of Statistics, 37
(5A):2178, 2009.
Sanford Weisberg. Applied linear regression. John Wiley & Sons, 1980.
Tong Zhang. On the consistency of feature selection using greedy least squares regression. Journal of Machine
Learning Research, 10(Mar):555–568, 2009.
Peng Zhao and Bin Yu. On model selection consistency of Lasso. Journal of Machine Learning Research, 7
(Nov):2541–2563, 2006.

14

A
A.1

Proofs
Theorem 1

We state some preparatory lemmas for the proof of Theorem 1. Throughout the Appendix, we
assume M > 0 and adopt the following notations:
β> X > y
kXβk2 + z > Xβ
√
Γ= √
≡
,
kM kyk
kM kyk

kyk
D≡ √ ,
k

where k · k denotes the usual `2 norm k · k2 .
Lemma A.1. Under Assumption 1, for an arbitrary constant c > 0, we have
−c ≤ Γ ≤ 1 + c
with probability tending to one.
Lemma A.2. Let ζ1 , . . . , ζm be independent standard normals and ζ(1) ≥ · · · ≥ ζ(m) be the order
statistics. For any (deterministic) sequence {im } such that im /m → 0 as m → ∞, we have
r
log log im
m
m
ζ(im ) = 2 log
.
− (1 + oP (1)) q
im
2 2 log m
im

The proof of Lemma A.2 is omitted. Interested readers can find its proof in Chapter 2 of
de Haan and Ferreira (2007).
Lemma A.3. Under Assumption 1, we have
r
max |Xj> y| ≤ 2D
j

2k log p
n

with probability converging to one.
Lemma A.4 (for the lasso case). Fit y on the true support XS using the lasso and denote by
βbS (λ) the lasso solution with penalty λ. Then, under Assumption 1, there exists a constant C such
that
√
C kkβbS (λ)k0 log p
>
S
b
·D
max Xi X β (λ) ≤
n
bS (λ)=0
i:β
i
holds uniformly for all λ > 0 with probability tending to one, where we make the (unusual) convention that max ∅ = 0.
Above, k · k0 equals the number of nonzero components of a vector, and βbS is a p-dimensional
vector which takes 0 on S := {1, . . . , p} \ S. Next, we proceed to a definition that is the subject of
Lemma A.6.
Definition A.5. Let Ic (γ) be the largest positive integer I such that
r
r
2k log(k/I)
2k log(p − k)
γ+
>
− c,
n
n
15

and let Jc (γ) denote the largest positive integer J such that
r
r
2k log(k/J)
2k log(p − k)
−γ +
>
− c.
n
n
Set Ic (γ) or Jc (γ) to 0 if it does not exist.
In applying this definition, while c is often fixed and small, γ can be a sequence indexed by
(n, p, k).
Lemma A.6. Let c be a constant. Under Assumption 1, each of the following cases holds.
(a) Assume c > 0. If −c/2 < γ ≤ 0, then Jc (γ) ≥ Ic (γ), and
#
"
r
2n log p (c − γ)2 n
n
.
−
+ log
log Jc (γ) = (1 + o(1)) (c − γ)
k
2k
2p log p
p
(b) Assume c > 0. If 0 < γ ≤ 1.01 and 2k log(p − k)/n > c + γ + c0 for some constant c0 > 0
(which is used to guarantee that Ic (γ) 6= k), then we have Ic (γ) ≥ Jc (γ) and
"
#
r
2n log p (c + γ)2 n
n
log Ic (γ) = (1 + o(1)) (c + γ)
−
+ log
.
k
2k
2p log p
(c) Assume −0.1 < c < 0. If 0.11 < γ ≤ 1.01, then Ic (γ) ≥ Jc (γ) and
#
"
r
2n log p (c + γ)2 n
n
−
+ log
.
log Ic (γ) = (1 + o(1)) (c + γ)
k
2k
2p log p
With these preparatory lemmas in place, we are ready to prove Theorem 1 in the lasso and
least angle regression cases. The first part of the proof offers a simple representation of the order
statistics of XS> y, and this representation serves as useful ingredients in proofs of other results.
Proof of Theorem 1 in the lasso and least angle regression cases. Without loss of generality, assume
that the true support set S = {1, . . . , k}. Let (1), (2), . . . , (k) be a permutation of 1, . . . , k such
that
>
>
>
X(1)
y ≥ X(2)
y ≥ · · · ≥ X(k)
y.
(A.1)
Conditional on z and XS βS ≡ µ, the k exchangeable random variables X1> y, . . . , Xk> y are
jointly normal with means all equal to µ> y/(kM ) and the equicorrelated covariance that has
(k − 1)kyk2 /(kn) on the diagonal and −kyk2 /(kn) off the diagonal. Let ξ be normally distributed
with mean 0 and variance kyk2 /(kn), and further assume ξ to be independent of the random design
X and the noise z. Then, X1> y + ξ, . . . , Xk> y + ξ are independent normals each with conditional
mean µ> y/(kM ) and conditional variance kyk2 /n. Thus, conditional on µ and z, we see that
µ> y kyk
+ √ ζi
kM
n
!
r
k
=D Γ+
ζi ,
n
d

Xi> y + ξ =

16

for independent standard normals ζ1 , . . . , ζk . Note that X1> y + ξ, . . . , Xk> y + ξ keeps the same
ordering as in (A.1).
Consider the first time that the full lasso (regressing on X) is just about to include some variable
among X(Ic +1) , X(Ic +2) , . . . , X(k−Jc ) , where Ic = Ic (Γ) and Jc = Jc (Γ) as in Definition A.5. Call
this variable X(L) and the penalty at that time λ∗ . Denote by A the event that all the selected
variables preceding X(L) are true variables. On A (the complement of A), there are at most Ic + Jc
variables selected before the first false variable. Hence,
T ≤ Ic + Jc + 1
on A, where c > 0 is an arbitrary constant. By Lemma A.6, we get
log(Ic + Jc + 1) = log max{Ic , Jc } + O(1)
"
r

#
2n log p (c + |Γ|)2 n
n
+ O(1)
= (1 + oP (1)) (c + |Γ|)
−
+ log
k
2k
2p log p
"
#
r
2n log p (c + |Γ|)2 n
n
= (1 + oP (1)) (c + |Γ|)
−
+ log
k
2k
2p log p
#
"
r
2n log p (1 + 3c/2)2 n
n
−
+ log
≤ (1 + oP (1)) (1 + 3c/2)
k
2k
2p log p

with probability approaching one. Above, we make use of Lemma A.1 in the last inequality. Since
c > 0 is arbitrary, we get
"r
#
2n log p
n
n
log T ≤ log(Ic + Jc + 1) ≤ (1 + oP (1))
−
+ log
.
k
2k
2p log p
on A.
The rest is devoted to proving P(A) → 0. To proceed, we point out an observation on the lasso:
the lasso running on the full design X is the same as regressing on XS until the T th variable gets
selected. In light of this observation, we perform the lasso on the true support XS . Since Ic + Jc is
much smaller than the bound given by Theorem 3, we can assume no drop-out has happened until
X(L) arrives and, hence, least angle regression is the same as the lasso. Because Ic +1 ≤ L ≤ k −Jc ,
we get
>
>
>
y.
(A.2)
X(I
y ≥ X(L)
y ≥ X(k−J
c)
c +1)
By Lemma A.2, the left-hand side of (A.2) obeys
"
#
r
√
2k
log(k/(I
+
1))
k
log
log(k/(I
+
1))
ξ
c
c
>
X(I
− (1 + oP (1)) p
−
y =D Γ+
c +1)
n
D
2 2n log(k/(Ic + 1))
A little analysis reveals that Assumption 1 implies
√
k log log(k/(Ic + 1))
p
= oP (1).
2 2n log(k/(Ic + 1))

17

Hence, we get
"

r

2k log(k/(Ic + 1))
ξ
>
y =D Γ+
X(I
− oP (1) −
c +1)
n
D
"
#
r
2k log(k/(Ic + 1))
=D Γ+
+ oP (1)
n
#
"r
2k log(p − k)
− c + oP (1) .
≤D
n
Similarly, the right-hand side of (A.2) obeys
"
r

2k log(k/(Jc + 1))
>
X(k−J
y =D Γ−
+ oP (1)
c)
n
#
" r
2k log(p − k)
+ c + oP (1)
≥D −
n
"r
#
2k log(p − k)
= −D
− c + oP (1) .
n

#

#

Thus, from (A.2) it follows that
"r
#
"r
#
2k log(p − k)
2k
log(p
−
k)
>
D
− c + oP (1) ≥ X(L)
y ≥ −D
− c + oP (1) ,
n
n
yielding
"r
>
X(L)
y ≤D

#
2k log(p − k)
− c + oP (1) .
n

(A.3)

The KKT conditions of the lasso at λ∗ give


>
X(L)
y − X βbS = λ∗ .
> Xβ
bS = oP (D), yields
The equality above together with (A.3) and Lemma A.4, which gives X(L)

∗

λ =

>
X(L)



bS

y − Xβ



"r

2k log(p − k)
>
≤D
− c + oP (1) + |X(L)
X βbS |/D
n
"r
#
2k log(p − k)
=D
− c + oP (1) .
n

#
(A.4)

Now, we seek a contradiction to get P(A) → 0. Since none of Xk+1 , . . . , Xp has been included
at λ∗ on the event A, we get
"r
#


2k
log(p
−
k)
max Xj> y − X βbS ≤ λ∗ ≤ D
− c + oP (1) .
(A.5)
k+1≤j≤p
n
18

Recognizing the independence between all Xj , j ≥ k + 1 and y − X βbS , we get


max Xj> y − X βbS ≥ max Xj> y − max Xj> X βbS
k+1≤j≤p

k+1≤j≤p

=

max

k+1≤j≤p

k+1≤j≤p

Xj> y − oP (D).

√
Now we turn to focus on maxk+1≤j≤p n|Xj> y|/kyk, which is distributed as the maximum absolute
values of p − k iid standard normals. Hence, Lemma A.2 gives
√
n|Xj> y| p
(0.5 + oP (1)) log log(p − k)
p
max
.
= 2 log(p − k) −
k+1≤j≤p
kyk
2 log(p − k)
Consequently,
max

k+1≤j≤p



Xj> y − X βbS ≥

max Xj> y − oP (D)
#
r "
(0.5 + oP (1)) log log(p − k)
k p
p
=D
2 log(p − k) −
− oP (D)
n
2 log(p − k)
r
r
2k log(p − k)
k (0.5 + oP (1)) log log(p − k)
p
=D
−D
− oP (D)
n
n
2 log(p − k)
k+1≤j≤p

By assumption,
r
q
q
log log(p − k)
log log p
k log log(p − k)
p
. log0.99 p · p
 log0.99 p · √
= o(1).
n 2 log(p − k)
2 log p
2 log(p − k)
Therefore, we get
max

k+1≤j≤p

Xj>




y − X βb ≥ D

"r

#
2k log(p − k)
+ oP (1) ,
n

(A.6)

contradicting (A.5). Thus, P(A) = o(1).
Next, we aim to prove Theorem 1 in the case of forward stepwise. Below are two preparatory
lemmas.
Lemma A.7. Given 1 ≤ m ≤ p, for any subset F ⊂ {1, . . . , p} of cardinality at most m, i ∈
/ F,
and j ∈ F , regress Xi on XF and denote by α
bji,F the least-squares coefficient of Xj in this fit.
Then, under Assumption 1, with probability approaching one, we have
r
5m log p
i,F
|b
αj | ≤
n
uniformly over all i, j, and F of cardinality at most m.

19

Lemma A.8. Under Assumption 1, we get
r
r
log p
log p
≤ kXi k ≤ 1 + 3
.
1−3
n
n
with probability approaching one uniformly for all 1 ≤ i ≤ p.
This lemma is a simple consequence of a well-known concentration inequality of the form
 2

χn
2
P
− 1 ≥ t ≤ 2e−nt /8
n
for all 0 ≤ t ≤ 1. The proof is thus omitted. Now we move to prove Theorem 1 for forward stepwise
regression.
Proof of Theorem 1 in the forward stepwise case. Here, we use the same proof strategy as in the
lasso and least angle regression cases. Define A and X(L) as earlier. To complete the proof, it
suffices to show that P(A) → 0.
We also assume that forward stepwise is performed on the true support XS . Denote by m∗
the number of variables selected before X(L) , that is, X(L) is selected in the m∗ + 1st step. Take
some m satisfying m ≤ pc for some c < 1/3. It is easy to show that m∗  m with probability
approaching one.
Denote by F ∗ the set of selected variables before X(L) . On the event A, the definition of forward
stepwise yields
min
ky − Xbk2 ≤
min
ky − Xbk2 .
supp(b)=F ∗ ∪(L)

supp(b)=F ∗ ∪l,l∈S

Recall that S = {k + 1, . . . , p} stands for the complement of S. In particular, we have
min

supp(b)=F ∗ ∪(L)

ky − Xbk2 ≤

min

supp(b)=F ∗ ∪j 0

ky − Xbk2 ,

(A.7)

f(L) be residual vector by regressing X(L) on XF ∗ . Then,
where j 0 = argmaxj∈S |Xj> y|. Let X
min

ky − Xbk2 −



f> y
X
(L)

2

f(L) k
kX



supp(b)=F ∗ ∪(L)

= −
"
=−

> −
(X(L)

P

i∈F ∗

αi Xi> )y

f(L) k
kX

20

min

supp(b)=F ∗

#2
,

ky − Xbk2

where αi ’s are the least-squares coefficients. By Lemma A.7, we get
X

αi Xi> y ≤

i∈F ∗

√
m∗ max∗ kαi Xi> yk
i∈F

√
≤ m max∗ kαi Xi> yk
i∈F
r
√
5m log p
≤ m
max∗ kXi> yk
i∈F
n
r
r
2
5m log p
2k log p
≤
· 2D
n
n
s
40m2 k log2 p
=D
n2
s
40m2 log2.99 p
.D
n
= o(D),

which makes use of Lemma A.3. Therefore, we get
!
X
>
>
>
X(L) −
αi Xi y = X(L)
y + oP (D).
i∈F ∗

Together with Lemma A.8, the equality above gives
#2
"
#2 "
P
>
> −
> y + o (D)
X(L)
(X(L)
P
i∈F ∗ αi Xi )y
P
≤
f(L) k
kX(L) k − i∈F ∗ |αi |kXi k
kX


≤


q
1 − 3 logn p

= (1 + O(
Similarly, we have
"
#2
P
> −
>
(X(L)
i∈F ∗ αi Xi )y
f(L) k
kX

2

> y + o (D)
X(L)

P


q
q

− m 5m nlog p 1 + 3 logn p

(A.8)


2
p
>
m3 n−1 log p)) X(L)
y + oP (D) .

≥ (1 − O(


2
p
>
m3 n−1 log p)) X(L)
y + oP (D) .

(A.9)

For the right-hand side of (A.7), we have
− (1 + O(


2
p
m3 n−1 log p)) Xj>0 y + oP (D)
≤

min

supp(b)=F ∗ ∪j 0

ky − Xbk2 −

≤ −(1 − O(
21

min

supp(b)=F ∗

ky − Xbk2


2
p
m3 n−1 log p)) Xj>0 y + oP (D)
(A.10)

with probability approaching one. Combining (A.7), (A.8), (A.9), and (A.10) gives
2


2
p
p
>
(1 + O( m3 n−1 log p)) X(L)
y + oP (D) ≥ (1 − O( m3 n−1 log p)) Xj>0 y + oP (D)

(A.11)

on the event A.
To complete the proof, we shall show that (A.11) holds with probability tending to zero. On
the one hand, from the earlier proof for the case of the lasso, we see that

2
p
p
>
(1 + O( m3 n−1 log p)) X(L)
y + oP (D) ≤ (1 + O( m3 n−1 log p))D2
"r
=D
which makes use of the fact that O(
tions. On the other hand,

2

"r

#2
2k log(p − k)
− c + oP (1)
n
#2

2k log(p − k)
− c + oP (1)
n

,

p
p
m3 n−1 log p)) · 2k log(p − k)/n = o(1) under our assump-


2
p
p
(1 − O( m3 n−1 log p)) Xj>0 y + oP (D) = (1 − O( m3 n−1 log p))D2
"r
=D

2

"r

#2
2k log(p − k)
+ oP (1)
n

#2
2k log(p − k)
+ oP (1)
n

Hence, (A.11) yields
"r
D2

"r
#2
#2
2k log(p − k)
2k
log(p
−
k)
− c + oP (1) ≥ D2
+ oP (1) ,
n
n

which is a clear contradiction. Hence, P(A) → 0.
To conclude this section, we present proofs of Lemma A.1 through Lemma A.7, respectively.
Proof of Lemma A.1. Recall the notation µ = Xβ. We first list some facts that are constantly
used in the proof. First,
√
kXβk = (1 + oP (1))M k.
This is because nkXβk2 /(kM 2 ) is just a χ2 random variable with
√ n degrees of freedom. Thus,
nkXβk2 /(kM 2 ) = (1 + oP (1))n, yielding kXβk = (1 + oP (1))M k. Hence, we get
kβ > X > kkyk
kXβk
β> X > y
Γ= √
≤ √
= √
≤ 1 + oP (1).
kM kyk
kM kyk
kM
Next we turn to prove that Γ ≥ −c. In fact, we have
Γ≡

kµk2 + µ> z
µ> z
√
≥√
.
kM kyk
kM kyk
22

(A.12)

For kyk, note that

p
√
kyk = (1 + oP (1)) kM 2 + nσ 2 ≥ (1 + oP (1))σ n.

And conditional on kµk, the numerator µ> z is normally distributed with mean 0 and variance
σ 2 kµk2 , implying
√
µ> z = OP (σkµk) = OP (M σ k).
Combining the results above, particularly (A.12), gives
√
√
OP (M σ k)
Γ ≥ (1 + oP (1)) √
√ = OP (1/ n) = oP (1).
kM σ n
Hence, Γ ≥ −c with probability approaching one for an arbitrary constant c > 0.
Proof of Lemma A.3. If βj = 0, then Xj is independent of y. This means conditional on y the
distribution of Xj> y is N (0, kyk2 /n). Thus,
r
kyk p
2k log p
>
2 log p = D
|Xj y| ≤ √
n
n
with probability at least 1 − o(1/p). Taking a union bound yields
r
r
2k log p
2k log p
>
< 2D
max |Xj y| ≤ D
n
n
j ∈S
/
with probability 1 − o(1). Now consider a j such that βj 6= 0. Note that
Xj> y = Xj> Xj βj + Xj> (X−j β−j + z) = M Xj> Xj + Xj> (X−j β−j + z).

(A.13)

The first term obeys



2
P |Xj> Xj − 1| ≥ t = P |χ2n /n − 1| ≥ t ≤ 2e−nt /8
1

for any 0 < t < 1. Setting t to n− 4 gives
max |Xj> Xj − 1| = oP (1).
j

In addition, we have
r
D

Recognizing that

r
r
2k log p
kM 2 + nσ 2 2k log p
= (1 + oP (1))
n
k
n
r
r
2
kM
2k log p
≥ (1 + oP (1))
k
n
r
2k log p
≥ (1 + oP (1))M
.
n

p
2k(log p)/n → ∞, we get
r
M Xj> Xj

= (1 + oP (1))M = oP (1) · D
23

2k log p
.
n

(A.14)

For the second term of (A.13), we get
max |Xj> (X−j β−j
j

maxj kX−j β−j + zk p
√
2 log p = (1 + oP (1))D
+ z)| ≤
n

r

2k log p
,
n

which makes use of the independence between Xj and (X−j β−j + z). Plugging the last display
and (A.14) into (A.13) gives
r
r
2k
log
p
2k log p
max |Xj> y| ≤ oP (1)D
+ (1 + oP (1))D
,
j∈S
n
n
p
which is smaller than 2D 2k(log p)/n with probability tending to one.
p
b 0 ≤ 0.25 n/ log p. Since
Proof of Lemma A.4. First, we consider the case where kβk
√

p
√
k × n/ log p/4 × log p
kkβbS k0 log p
·D ≥
·D
n
n
r
1 k log p
=
· D.
4
n

and
Xj> X βbS ≤ Xj> y + Xj> (y − X βbS )
≤ 2 max Xi> y
i
r
2k log p
≤ 4D
n
for all j, which makes use of Lemma A.3. In this case, the proof is simply as follows:
max
j:βbjS (λ)=0

Xj> X βbS (λ) ≤ max Xj> X βbS (λ)
1≤j≤p

r

2k log p
n
r
√ 1 k log p
= 16 2 ·
·D
4
n
√
√
kkβbS k0 log p
· D.
≤ 16 2 ·
n
≤ 4D

(A.15)

p
n/ log p/4. Theorem 6 of Cai et al. (2013) ensures
r
4 log p
>
max |Xj Xi | ≤ (1 + oP (1))
1≤i<j≤p
n

b 0≤
Now we move to the case where kβk

24

under Assumption 1. This result implies
X
Xj> X βbS ≤
|βbiS ||Xi> Xj |
i

=

X

|βbiS ||Xi> Xj |

i:βbiS 6=0

r
X

≤ (1 + oP (1))

|βbiS |

i:βbiS 6=0

4 log p
n

(A.16)

r
4 log p
S
S
b
b
≤ (1 + oP (1))kβ k0 max |βi |
i
n
for all j such that βbjS = 0. Let j ? be the index that |βbjS? | is the largest. Taking
r
|βbjS? |

≤ (8 + oP (1))D

2k log p
n

(A.17)

as given for the moment, from (A.16) we get
Xj> X βbS

r
4 log p
S
S
b
b
≤ (1 + oP (1))kβ k0 max |βi |
i
n
r
r
2k
log
p
4 log p
≤ (1 + oP (1))kβbS k0 × 8D
×
n
n
√
S
√
kkβb k0 log p
= (1 + oP (1))16 2 ·
·D
n

√
Summarizing from the two cases, the lemma holds if we take any constant C > 16 2.
The rest of the proof aims to verify (A.17). Note that, on the one hand,
λ = Xj>? (y − X βbS )
≤ λmax := max |Xj> y|
1≤j≤p
r
2k log p
≤ 2D
.
n

25

(A.18)

On the other hand, we have
Xj>? (y − X βbS ) ≥ Xj>? X βbS − Xj>? y
≥ Xj>? X βbS − max |Xj> y|
1≤j≤p
X
2
S
≥ kXj ? k |βbj ? | −
|βbiS ||Xi> Xj ? | − max |Xj> y|
1≤j≤p

i6=j ? ,βbiS 6=0

r
2

≥ kXj ? k |βbjS? | −

X

(1 + oP (1))

i6=j ? ,βbiS 6=0

4 log p bS
|βi | − max |Xj> y|
1≤j≤p
n

r

4 log p bS
|βj ? | − max |Xj> y|
1≤j≤p
n
i6=j ? ,βbiS 6=0
"
#
r
4
log
p
kβbS k0 |βbjS? | − max |Xj> y|
≥ kXj ? k2 − (1 + oP (1))
1≤j≤p
n


1
+ oP (1) |βbjS? | − max |Xj> y|
≥
1≤j≤p
2

≥ kXj ? k

2

|βbjS? |

−

X

(1 + oP (1))

Hence,


r
r

1
2k
log
p
2k log p
+ oP (1) |βbjS? | ≤ 2D
+ max |Xj> y| ≤ 4D
1≤j≤p
2
n
n

with probability tending to one. We see that (A.17) is an immediate consequence.
Proof of Lemma A.6. For the first case, by definition we have
r
r
2k log(k/J)
2k log(p − k)
−γ +
>
−c
n
n
and

r
2k log(k/(J + 1))
2k log(p − k)
≤
− c,
−γ +
n
n
where we write J for Jc (γ). These two inequalities are equivalent to
r
2n log(p − k) (c − γ)2 n
k
log Jc (γ) < (c − γ)
−
+ log
k
2k
p−k
r
2
2n log(p − k) (c − γ) n
k
log(Jc (γ) + 1) ≥ (c − γ)
−
+ log
.
k
2k
p−k
r

Under Assumption 1, it is easy to check that
r
2n log(p − k) (c − γ)2 n
k
(c − γ)
−
+ log
k
2k
p−k
#
"
r
2n log p (c − γ)2 n
n
= (1 + o(1)) (c − γ)
−
+ log
k
2k
2p log p
→ ∞.
26

(A.19)

Thus, from (A.19) we get
"

r

log Jc (γ) = (1 + o(1)) (c − γ)

#
2n log p (c − γ)2 n
n
.
−
+ log
k
2k
2p log p

The rest two cases follow from similar reasoning and, thus, their proofs are omitted.
Proof of Lemma A.7. Here, the least-squares estimate is
b i,F = (XF> XF )−1 XF> Xi ,
α
which, conditional on XF , is distributed as
N (0, (XF> XF )−1 /n).
Denote by Am the event that

1
max kXF> XF − Ik ≤ ,
2

|E|≤m

where k · k denotes the matrix spectral norm. Note that we have




p
p
αji,F | ≥ 5m(log p)/n, Am
P max |b
αji,F | ≥ 5m(log p)/n ≤ P(Am ) + P max |b


X
p
≤ P(Am ) +
P |b
αji,F | ≥ 5m(log p)/n, Am
i∈F,j∈F,|F
/
|≤m



p
1
i,F
>
P |b
αj | ≥ 5m(log p)/n, kXF XF − Ik ≤
≤ P(Am ) +
2
i∈F,j∈F,|F
/
|≤m


X
p
1
i,F
>
≤ P(Am ) +
,
P |b
αj | ≥ 5m(log p)/n kXF XF − Ik ≤
2
X

i∈F,j∈F,|F
/
|≤m

where the max operator is taken over all triples (i, j, F ) such that |F | ≤ m, i ∈
/ F , and j ∈ F .
>
>
−1
Conditional on kXF XF − Ik ≤ 1/2, all the diagonal elements of (XF XF ) are no greater than
2. As a result,
!
p


p
p
5m(log
p)/n
1
p
P |b
αji,F | ≥ 5m(log p)/n kXF> XF − Ik ≤
≤ P |N (0, 1)| ≥
= 2Φ(− 2.5m log p).
2
2/n
Hence,
X
i∈F,j∈F,|F
/
|≤m



p
1
≤
P |b
αji,F | ≥ 5m(log p)/n kXF> XF − Ik ≤
2

X

p
2Φ(− 2.5m log p)

i∈F,j∈F,|F
/
|≤m

 
p
p−m+1 p
≤ pm ·
· 2Φ(− 2.5m log p)
p − 2m + 1 m
p
p−m+1 m
≤ pm ·
p · 2Φ(− 2.5m log p)
p − 2m + 1
p
m+O(1)
≤p
Φ(− 2.5m log p)
1
≤ pm+O(1) √
e−1.25m log p
2.5m log p
→ 0.
27

Together with P(Am ) → 0, which follows from the proof of Theorem 3, the last display gives


p
P max |b
αji,F | ≥ 5m(log p)/n → 0.

A.2

Theorem 2

As earlier in Section A.1, we first state some lemmas before turning to the proof of Theorem 2.
Lemma A.9. Fit the response y on the true support XS using the lasso and denote by βbS (λ) the
lasso solution. Under Assumption 1, further assume that 0.5 < Γ < 1.1 with probability tending
to one. Then, for an arbitrary constant c > 0, with probability approaching one, all the variables
X(1) , X(2) , . . . , X(I−c ) enter the model before any of
X(Ic +1) , X(Ic +2) , . . . , X(k)
along the lasso path.
Lemma A.10. Under Assumption 1, further assume
r
σ n
→ 0.
M k
Then, Γ = 1 + oP (1).
With Lemmas A.9 and A.10 in place, now we give the proof of Theorem 2 where the procedure is
the lasso or least angle regression. Note that for large γ and sufficiently small c, we have Jc (γ) = 0.
Proof of Theorem 2 in the lasso and least angle regression cases. The proof starts by making use
of Lemma A.10, which concludes
1−c≤Γ≤1+c
with probability
tending to one for any constant c > 0. To see this, note that Assumption 1 along
p
√
ensures n/k ≤ 1/ c4 = O(1) and, as a result, we get
r
σ n
→0
M k
given the condition that σ/M → 0.
As in the proof of Theorem 1, here it also suffices to only focus on the lasso case. Consider
fitting the lasso on XS . Denote by X(I  ) the last variable among X(1) , X(2) , . . . , X(I−c ) that enters
the model and let λ be the lasso penalty when X(I  ) is just about to enter the lasso path. Note
S (λ ) = 0 and I
that βb(I
)
−c = I−c (Γ). The KKT conditions give


>
bS = λ .
X(I
) y − X β

28

(A.20)

Note that

"

>
X(I
)y

r

2k log(k/(I−c ))
ξ
=D Γ+
≥
+ oP (1) −
n
D
"
#
r
2k log(k/(I−c ))
=D Γ+
+ oP (1)
n
"r
#
2k log(p − k)
=D
+ c + oP (1) + oP (1)
n
"r
#
2k log(p − k)
=D
+ c + oP (1) .
n

#

>
y
X(I
−c )

(A.21)

From Lemma A.9 we see that, with probability approaching one, by the time of λ none of the
variables X(Ic +1) , X(Ic +2) , . . . , X(k) has been included in the lasso path. Hence, the lasso model
consists of no more than Ic variables at any time before λ , and it is easy to check that


n
S
kβb (λ)k0 ≤ Ic (Γ) ≤ oP √
k log p
for all λ > λ . Consequently, we get
√
C kkβbS (λ)k0 log p
→0
n
> Xβ
bS | = oP (D). Hence, taking (A.20) and
uniformly for all λ > λ . Thus, Lemma A.4 shows |X(I
)
(A.21) together shows that λ obeys
#
"r
2k
log(p
−
k)
(A.22)
+ c + oP (1) .
λ ≥ D
n

(Similarly, given a small Ic , Theorem 3 ensures that by that time no drop-out has ever happened
and, consequently, the lasso and least angle regression are equivalent in our discussion.) If we can
show that


Xj> y − X βbS (λ) < λ
(A.23)
for all λ > λ and all j ∈
/ S, then the full lasso (on X) selects variables only from XS before time
λ . Call this event B. Thus, on B the first I−c selected variables along the lasso path are all signal
variables. As for I−c , Lemma A.6 yields
log T ≥ log(I−c (Γ) + 1)
"
= (1 + oP (1)) (−c + Γ)

r

#
2n log(p − k) (−c + Γ)2 n
n
−
+ log
.
k
2k
2p log p

Provided 1 − c ≤ Γ ≤ 1 + c with probability tending to one, we get
r
r
2n log(p − k) (−c + Γ)2 n
2n log(p − k) (1 − 2c)2 n
(−c + Γ)
−
≥ (1 − 2c)
−
.
k
2k
k
2k
29

Hence, we have
"

r

log T ≥ (1 + oP (1)) (1 − 2c)

2n log(p − k) (1 − 2c)2 n
n
−
+ log
k
2k
2p log p

for all constant c > 0. Setting c → 0+ gives
"r

2n log(p − k)
n
n
log T ≥ (1 + oP (1))
−
+ log
k
2k
2p log p
"r
#
2n log p
n
n
= (1 + oP (1))
,
−
+ log
k
2k
2p log p

#

#

as desired.
Establishing (A.23) is the subject of the remaining proof. That is, prove P(B) → 1. For j ∈
/ S,
observe that


(A.24)
Xj> y − X βbS ≤ |Xj> y| + |Xj> X βbS |.
Making use of the independence between y and Xj , the first term Xj> y obeys
r
max |Xj> y|
j ∈S
/

≤ kyk

2 log(p − k)
n

with probability approaching one. The second term, |Xj> X βbS |, as shown earlier by Lemma A.4,
satisfies that
√
C kkβbS k0 log p
>
S
b
max |Xj X β (λ)| ≤
· D = oP (D)

n
j ∈S,λ>λ
/
with probability tending to one. Hence, (A.24) yields
r


2 log(p − k)
+ oP (D)
Xj> y − X βbS ≤ kyk
n
r
2k log(p − k)
=D
+ oP (D)
n
"r
#
2k log(p − k)
+ oP (1)
=D
n
"r
#
2k log(p − k)
<D
+ c + oP (1)
n
≤ λ
< λ,
where the second last inequality follows from (A.22). Therefore, P(B) → 1.
Next, we turn to prove Theorem 2 in the case of forward stepwise regression.

30

Lemma A.11. Fit the response y on the true support XS using forward stepwise. Under Assumption 1, further assume that 0.5 < Γ < 1.1 with probability tending to one. Then, for an arbitrary
constant c > 0, with probability approaching one, all the variables X(1) , X(2) , . . . , X(I−2c ) enter the
model before any of
X(I−c +1) , X(I−c +2) , . . . , X(k)
by forward stepwise.
Proof of Theorem 2 in the forward stepwise case. Consider running forward stepwise on the design
matrix XS . Denote by XI  the last variable among X(1) , . . . , X(I−2c ) that gets selected. Let m? − 1
be the number of variables selected prior to X(I  ) . From Lemma A.11 we see that, with probability
approaching to one, by then none of X(I−c +1) , . . . , X(k) has been selected, that is, m? ≤ I−c . The
proof would be completed once we show that by then no noise variables would be selected if we
perform forward stepwise on X instead of XS .
Denote by Fm the set of variables selected in the first m steps. We would like to show that
ky − Xbk2 ≤

min
supp(b)=Fm

ky − Xbk2

min

(A.25)

supp(b)=Fm−1 ∪l,l∈S

for all m ≤ m? . Denote by X(jm ) the variable selected in the mth step. To this end, note that from
the proof of Theorem 1 for the forward stepwise case, we see that for any fixed l ∈ S,
min

supp(b)=Fm−1 ∪l

ky − Xbk2 −

min

ky − Xbk2

supp(b)=Fm

>
= (X(j
y + oP (D))2 − (Xl> y + oP (D))2 .
m)

p
Note that since l ∈ S, we have |Xl> y| ≤ D( 2k log(p − k)/n + oP (1)), and since 1 ≤ jm ≤ I−c , we
get
"r
#
2k log(p − k)
>
>
X(jm ) y ≥ X(I−c ) y = D
+ c + oP (1) .
n
An immediate consequence is the following:
min

supp(b)=Fm−1 ∪l

"r
≥D

2

ky − Xbk2 −

min

ky − Xbk2

supp(b)=Fm

#2
"r
#2
2k log(p − k)
2k log(p − k)
2
+ c + oP (1) − D
+ oP (1)
n
n

> 0,
which certifies (A.25) provided the arbitrariness of l ∈ S. Therefore, with probability approaching
one, we get
log T ≥ log(I−2c (Γ) + 1).
Setting c → 0+ completes the proof. Details are the same as the proof for the case of the lasso and
least angle regression.
To conclude this section, we prove Lemmas A.9, A.10, and A.11.

31

Proof of Lemma A.9. Consider the first time along the lasso path that a variable among X(Ic +1) , X(Ic +2) , . . . , X(k)
is just about to enter the lasso model and denote by X(L) this variable, where Ic + 1 ≤ L ≤ k.
Specifically, writing λ∗ for the lasso penalty at this point, we know that
βb(L) (λ∗ ) = 0
and
βb(L) (λ) 6= 0
if λ∗ − c0 < λ < λ∗ for some c0 > 0. For a proof of this lemma by contradiction, assume that X(l)
has not yet entered the model at λ∗ for some l satisfying 1 ≤ l ≤ I−c . Under this assumption, our
discussion below considers the lasso solution βb at λ∗ . First of all, the KKT conditions of the lasso
give




>
>
y − X βb = λ∗ .
X(l)
y − X βb ≤ λ∗ , X(L)
Consequently,




>
>
X(l)
y − X βb ≤ X(L)
y − X βb .

(A.26)

Thus, the proof of the present lemma would be finished if we show that, with probability tending
to one, (A.26) cannot be satisfied.
> y satisfies
The rest part of the proof is devoted to disproving (A.26). Note that X(l)
"

r

2k log(k/I−c )
ξ
>
>
X(l)
y ≥ X(I
y =D Γ+
+ oP (1) −
−c )
n
D
#
"
r
2k log(k/I−c )
+ oP (1)
=D Γ+
n
"r
#
2k log(p − k)
=D
+ oP (1) ,
n

#

which, together with Lemma A.3, implies that
#
"r


2k
log(p
−
k)
>
>
b
+ oP (1) − |X(l)
X β|/D
X(l)
y − X βb ≥ D
n
"r
#
2k log(p − k)
=D
+ oP (1) .
n

(A.27)

Similarly, from Ic + 1 ≤ L ≤ k we get
>
>
>
X(I
y ≥ X(L)
y ≥ X(k)
y.
c +1)

The quantities appearing on the left- and right-hand sides obey, respectively,
#
"r
2k
log(p
−
k)
>
− c + oP (1)
X(I
yD
c +1)
n

32

(A.28)

and

"r
>
y
X(k)

= −D

#
2k log(p − k)
− Γ + oP (1) .
n

Hence, from (A.28) it follows that
"r
"r
#
#
2k log(p − k)
2k
log(p
−
k)
>
D
y ≥ −D
− c + oP (1) ≥ X(L)
− Γ + oP (1) ,
n
n
yielding
>
X(L)



"r

2k log(p − k)
>
b
X β|/D
− c + oP (1) + |X(L)
n
"r
#
2k log(p − k)
=D
− c + oP (1)
n


y − X βb ≤ D

#
(A.29)

Last, combining (A.27) and (A.29) shows that, with probability tending to one, (A.26) cannot
be satisfied, as desired.
Proof of Lemma A.10. Recall that the definition
(Xβ)> y
Γ= √
.
kM kyk
√
Observing that kyk = (1 + oP (1)) kM 2 + nσ 2 , the present lemma is equivalent to


Xβ
y
√
, √
= 1 + oP (1),
kM 2 + nσ 2
kM

(A.30)

where h, i denotes the usual inter product of vectors. Since
y
Xβ
√
and √
kM 2 + nσ 2
kM
are vectors of unit norm asymptotically, (A.30) boils down to claiming that the angle between the
two vectors are asymptotically zero. Recognizing that Xβ and z are asymptotically orthogonal, a
vanishing angle between Xβ and y = Xβ + z is equivalent to
kzk
→ 0.
kXβk
The display above is a direct consequence of the condition
√
nσ
√
→ 0,
kM
which is provided in the assumptions. Hence, (A.30) holds.

33

Proof of Lemma A.11. Consider the first time a variable among X(I−c +1) , . . . , X(k) enters the
model. Denote by m the rank of this variable. As earlier, call this variable X(jm ) . Note that, by
definition, I−c + 1 ≤ jm ≤ k. Suppose on the contrary that by the time X(jm ) is selected, at least
one variable among X(1) , . . . , X(I−2c ) , say X(l) , has not been included. Denote by C this event, on
which we must have
ky − Xbk2 ≤

min
supp(b)=Fm

min

supp(b)=Fm −1 ∪(l)

ky − Xbk2 .

(A.31)

From the proof of Theorem 1 for the forward stepwise case, we know that
min

supp(b)=Fm −1 ∪(l)

ky − Xbk2 −

min


supp(b)=Fm

ky − Xbk2
(A.32)

>
>
y + oP (D))2 .
y + oP (D))2 − (X(l)
= (X(j
m )

Since jm ≥ I−c + 1, we get
"r
>
X(j
y
m )

≤

>
X(I
y
−c )

=D

#
2k log(p − k)
+ c + oP (1) ,
n

(A.33)

and from jm ≤ k we get
>
>
X(j
y ≥ X(k)
y
m )
"
r

=D Γ−
"r
> −D

2k log k
+ oP (1)
n

#
(A.34)
#

2k log(p − k)
+ c + oP (1) ,
n

which makes use of the fact that 0.5 < Γ < 1.1. Similarly, we have
#
"r
2k
log(p
−
k)
>
>
+ 2c + oP (1) .
X(l)
y ≥ X(I
y≥D
−2c )
n

(A.35)

Plugging (A.33), (A.34), and (A.35) into (A.32) yields
min

supp(b)=Fm −1 ∪(l)

"r
≤D

2

ky − Xbk2 −

min


supp(b)=Fm

2k log(p − k)
+ c + oP (1)
n

ky − Xbk2

#2

"r
−D

2

#2
2k log(p − k)
+ 2c + oP (1)
n

< 0,
which contradicts (A.31). This implies the event C happens with probability vanishing to zero.

34

A.3

Theorem 3

Proof of Theorem 3. Write
m = min

m o
nl p
c n/ log p , p

for some constant c to be determined later. A variable might be dropped by the lasso only if its
fitted coefficient crosses zero. By using this fact, it suffices to show that the fitted coefficients of
the selected variables shall never cross zero in the first m selected variables, except for a rare event
with probability no more than 1/p2 .
Let λ0 be the first time along the lasso path a previously selected variable is just about to drop
out of the model. (The proof shall only focus on the event that such λ0 exists; otherwise no variable
drops out.) Denote by j the number of variables selected just before the first dropout and Sb the
set of these j variables. For the sake of contradiction, assume that j ≤ m − 1.
Pick an λ that is (slightly) smaller than λ0 and at which j variables have been included in the
model. Observe the following partial KKT condition for the lasso solution:
b b) + λ sgn(βb b) = 0.
− XS>
b βS
b (y − XS
S

(A.36)

In a componentwise manner, sgn(·) above returns 1 if the corresponding component is positive and
−1 if negative. From (A.36) it follows that
−1 >
>
−1
βbSb(λ) = (XS>
sgn(βbSb(λ)).
b) XS
b)
b XS
b y − λ(XS
b XS

(A.37)

−1 sgn(β
b b) has the
Now let us gradually decrease λ to λ0 , moving along the lasso path. If (X >
b)
b XS
S
S
b = j coordinates, then in this process of moving λ down to λ0 ,
same sign as sgn(βbSb) across all |S|
b
we see all positive coefficients of β(λ)
get larger while all negative coefficients become even smaller.

This implies that no coefficient will cross zero, a contradiction to the assumption. Therefore,
−1 sgn(β
b b) differs from sgn(βb b) in the sign of at least one coordinate, yielding
(X >
b)
b XS
S
S
S

−1
(XS>
sgn(βbSb) − sgn(βbSb) ≥ 1.
b)
b XS

(A.38)

On the other hand, we have
−1
−1
(XS>
sgn(βbSb) − sgn(βbSb) ≤ (XS>
−I
b)
b)
b XS
b XS

sgn(βbSb)
√
−1
−I
≤ (XS>
m−1
b)
b XS
≤ XS>
b−I
b XS

−1
(XS>
b)
b XS

(A.39)
√

m − 1.

Let θ be the restricted isometry constant for (m − 1)-sparse vectors, which is defined as the smallest
θ such that
(1 − θ)kbk2 ≤ kXbk2 ≤ (1 + θ)kbk2
for all (m − 1)-sparse b ∈ Rp (we use the notation θ instead of the conventional δ since the later
has been reserved for denoting the ratio n/p). It is not hard to see that this definition yields the
following more amenable expression
θ=

max kXU> XU − Ik,

|U |≤m−1

35

(A.40)

where the maximum is taken over all subsets of {1, 2, . . . , p} with cardinality smaller than m. In
the case where m ≥ 2 (if m = 1 then this theorem is trivially true), a well known result regarding
this constant states
r
(m − 1) log(p/(m − 1))
θ≤C
(A.41)
n
with probability at least 1 − 1/p2 , where C > 0 is a universal constant (see e.g. Baraniuk et al.
(2008); note that θ is random due to its dependence on X).
Now we prove that on the even (A.41), the assumption j ≤ m − 1 cannot hold with a suitable
of the constant c. To this end, we start by observing that (A.40) together with the assumption
b = j ≤ m − 1 ensures that all the eigenvalues of X > X b are between 1 − θ and 1 + θ. Thus, from
|S|
b
S
S
(A.39) we get
√
−1
−1
(XS>
sgn(βbSb) − sgn(βbSb) ≤ XS>
(XS>
m−1
b)
b−I
b)
b XS
b XS
b XS
≤

θ √
m − 1,
1−θ

Substituting (A.38) into the above display, under the assumption that j ≤ m − 1 we see that
q
C (m−1) log(p/(m−1))
√
θ √
n
q
· m−1≥
m−1≥1
(A.42)
1−θ
(m−1) log(p/(m−1))
1−C
n
holds on the event (A.41). If we can show that
q
C (m−1) log(p/(m−1))
√
n
q
· m − 1 < 1,
1 − C (m−1) log(p/(m−1))
n

(A.43)

meaning that (A.42) cannot be satisfied, then the assumption that j ≤ m − 1 must be violated
and, consequently, the first m variables are included along the lasso path without any drop-out
with probability at least 1 − 1/p2 (note that the event (A.41) happens with probability at least
1 − 1/p2 ). Indeed, (A.43) is true if we set (note that we have excluded the trivial case n < log p)


1
1
c = min 1,
,
.
4C 2 2C
This concludes the proof.

A.4

Theorem 4

Proof of Theorem 4. Note that the least-squares estimator takes the form:
βbLS = (X > X)−1 X > y
= β + (X > X)−1 X > z,
which is, conditional on X, distributed as
βbjLS
p
∼ N (0, 1)
σ [(X > X)−1 ]jj
36

for all j ∈
/ S and
βbjLS
p
∼ N (µj , 1)
σ [(X > X)−1 ]jj
for j ∈ S. Above,
M
µj = p
.
σ [(X > X)−1 ]jj
Take the following result as given for the moment:
h
i
max (X > X)−1
≤
1≤j≤p

jj

n
+c
n−p

(A.44)

with probability approaching one for an arbitrary constant c > 0. Then, we can show that
min
j∈S

βbjLS
βbjLS
p
> max p
.
j ∈S
/
σ [(X > X)−1 ]jj
σ [(X > X)−1 ]jj

To see this, note that
min
j∈S

p
βbjLS
M
p
≥ min p
− 2 log k
j∈S σ [(X > X)−1 ]jj
σ [(X > X)−1 ]jj

and
max
j ∈S
/

βbjLS
σ

p

[(X > X)−1 ]jj

≤

p

2 log(p − k)

with probability approaching one. Hence, (A.45) follows if one can show that
p
p
M
min p
− 2 log k − 2 log(p − k) > 0.
j∈S σ [(X > X)−1 ]jj
In fact, from (A.44) we have
p
p
M
min p
− 2 log k − 2 log(p − k)
j∈S σ [(X > X)−1 ]jj
p
M
> min p
− 2 2 log p
j∈S σ [(X > X)−1 ]jj
p
M
≥ q
− 2 2 log p
n
σ n−p
+c
p
M
≥ q
− 2 2 log p
δ
σ δ−1
+c
q
log p
3 2δδ−1
p
≥q
− 2 2 log p
δ
δ−1 + c
s
p
p
δ/(δ − 1)
= 3 2 log p ·
− 2 2 log p,
δ/(δ − 1) + c
37

(A.45)

which is positive by setting c sufficiently small. Hence, (A.45) holds, meaning that all the true
variables are vertically ranked higher than any of the false variables. In particular, the first false
variable is vertically ranked lower than all of the true variables.
In the remaining part of the proof, our aim is to prove (A.44). Write the singular value decomposition of X as
X = U DV > ,
where U ∈ Rn×n and V ∈ Rp×p are orthogonal matrices, and D ∈ Rn×p is diagonal. Then, the
left-hand side of (A.44) can be expressed as
h
i
>
−1 >
= e>
(X > X)−1
j V (D D) V ej ,
jj

where ej denotes the jth canonical basis vector in Euclidean space. Using this fact, we get


h
i
n
>
−1
P max (X X)
≥
+c
1≤j≤p
n−p
jj


n
>
>
−1 >
+c
= P max ej V (D D) V ej ≥
1≤j≤p
n−p


n
>
>
−1 >
= E P max ej V (D D) V ej ≥
+c D
1≤j≤p
n−p


 

p
X
n
>
−1 >
+ c D ,1
≤ E min
P e>
j V (D D) V ej ≥


n−p
j=1
 
 
n
>
−1 >
= E min p P e>
V
(D
D)
V
e
≥
+
c
D
,1 ,
1
1
n−p
where the last step makes use of the exchangeability of V > e1 , . . . , V > ep given D. Therefore, writing
η = V > e1 , it suffices to prove that
 
 
n
>
>
−1
E min p P η (D D) η ≥
+ c D , 1 → 0.
(A.46)
n−p
Recognizing that X has independent N (0, 1/n) entries, it is known that η is independent of
D and is uniformly distributed on the unit sphere in Rp . In particular, η assumes the following
representation
(ζ1 , . . . , ζp )>
ζ
η=q
=q
,
ζ12 + · · · + ζp2
ζ12 + · · · + ζp2
where ζ1 , . . . , ζp are independent N (0, 1) random variables. Using this representation, we get


n
P η > (D > D)−1 η ≥
+c D
n−p

 !
n
= P ζ > (D > D)−1 ζ ≥ kζk2
+c D
(A.47)
n−p
!


n
>
>
−1
0.75
≤ P ζ (D D) ζ ≥ (p − p )
+ c D + P(kζk2 ≤ p − p0.75 ).
n−p
38

√
p, Assumption 1 ensures that the following inequality must eventually hold:






√
δ
n
n
c
0.75
√
+ (2 pt + 2t) ×
(p − p )
+c ≥p
+
+c .
n−p
n−p 2
( δ − 1)2

Setting t =

Plugging the display above into (A.47) gives


n
>
>
−1
+c D
P η (D D) η ≥
n−p

 !


√
c
n
δ
>
>
−1
+ (2 pt + 2t) × √
+
+ c D + P(kζk2 ≤ p − p0.75 ).
≤ P ζ (D D) ζ ≥ p
2
n−p 2
( δ − 1)
Next, denote by

A=

>

tr((D D)

−1


)≥p

n
c
+
n−p 2



>

−1

or k(D D)


δ
k> √
+c
( δ − 1)2

and write 1A for the indicator function of A defined as
(
1 if the event A holds,
1A :=
0 otherwise.
Then, it is easy to see that
>

>

−1

P ζ (D D)


ζ≥p

c
n
+
n−p 2

√





+ (2 pt + 2t) ×

 !
δ
√
+c D
( δ − 1)2
!

√
≤ 1A + P ζ > (D > D)−1 ζ ≥ tr((D > D)−1 ) + (2 pt + 2t)k(D > D)−1 k D
≤ 1A + e−t ,

where the last step follows from a Gaussian concentration inequality (see, for example, Hsu et al.
(2012)). As a result, it suffices to show that

E min p1A + pe−t + p P(kζk2 ≤ p − p0.75 ), 1 → 0.
(A.48)
To this end, first note that from
δ
1
p
+c< √
+c
2
( δ − 1)2
(1 − p/n)
we get
!


δ
1
p
P k(D > D)−1 k > √
+ c ≤ P k(D > D)−1 k >
+ c → 0,
( δ − 1)2
(1 − p/n)2

(A.49)

where we
fact that the smallest singular value of the Wishart matrix X > X is concentrated
puse the
2
at (1 − p/n) with probability tending to one (see, for example, Vershynin (2012)). Second, it is
known that
np
,
tr((D > D)−1 ) = (1 + oP (1))
n−p
39

which implies



n
c
>
−1
→ 0.
P tr((D D) ) ≥ p
+
n−p 2

(A.50)

Taking (A.49) and (A.50) together, we get
P(A) → 0.
Recognizing that t =
(A.51) we have

√

(A.51)

p, we get pe−t → 0. In addition, p P(kζk2 ≤ p − p0.75 ) → 0. Hence, from

min p1A + pe−t + p P(kζk2 ≤ p − p0.75 ), 1 → 0

in probability. In addition, note that

min p1A + pe−t + p P(kζk2 ≤ p − p0.75 ), 1
is upper bounded by 1, which is an integrable function. Hence, by the dominated convergence
theorem, we get

E min p1A + pe−t + p P(kζk2 ≤ p − p0.75 ), 1 → 0,
as desired.

40

