Confidence Intervals and Hypothesis Testing for High-Dimensional Regression

Introduction

It is widely recognized that modern statistical problems are increasingly high-dimensional, i.e. require estimation of more parameters than the number of observations/samples. Examples abound from signal processing [LDSP08], to genomics [PZB+10], collaborative filtering [KBV09] and so on. A number of successful estimation techniques have been developed over the last ten years to tackle these problems. A widely applicable approach consists in optimizing a suitably regularized likelihood function. Such estimators are, by necessity, non-linear and non-explicit (they are solution of certain optimization problems).

The use of non-linear parameter estimators comes at a price. In general, it is impossible to characterize the distribution of the estimator. This situation is very different from the one of classical statistics in which either exact characterizations are available, or asymptotically exact ones can be derived from large sample theory [VdV00]. This has an important and very concrete consequence. In classical statistics, generic and well accepted procedures are available for characterizing the uncertainty associated to a certain parameter estimate in terms of confidence intervals or $p$-values [Was04, LR05]. However, no analogous procedures exist in high-dimensional statistics.

In this paper we develop a computationally efficient procedure for constructing confidence intervals and $p$-values for a broad class of high-dimensional regression problems. The salient features of our procedure are:

Our approach guarantees nearly optimal confidence interval sizes and testing power.

It is the first one to achieve this goal under essentially no assumptions beyond the standard conditions for high-dimensional consistency.

It allows for a streamlined analysis with respect to earlier work in the same area.

For the sake of clarity, we will focus our presentation on the case of linear regression, under Gaussian noise. Section 4 provides a detailed study of the case of non-Gaussian noise. A preliminary report on our results was presented in NIPS 2013 [JM13a], which also discusses generalizations of the same approach to generalized linear models, and regularized maximum likelihood estimation.

In the classic setting, $n\gg p$ and the estimation method of choice is ordinary least squares yielding $\widehat{\theta}^{\rm OLS}=({\mathbf{X}}^{{\sf T}}{\mathbf{X}})^{-1}{\mathbf{X}}^{{\sf T}}Y$. In particular $\widehat{\theta}^{\rm OLS}$ is Gaussian with mean $\theta_{0}$ and covariance $\sigma^{2}({\mathbf{X}}^{{\sf T}}{\mathbf{X}})^{-1}$. This directly allows to construct confidence intervalsFor instance, letting $Q\equiv({\mathbf{X}}^{{\sf T}}{\mathbf{X}}/n)^{-1}$, $\widehat{\theta}^{\rm OLS}_{i}-1.96\sigma\sqrt{Q_{ii}/n},\widehat{\theta}^{\rm OLS}_{i}+1.96\sigma\sqrt{Q_{ii}/n}]$ is a $95\%$ confidence interval [Was04]..

In case the right hand side has more than one minimizer, one of them can be selected arbitrarily for our purposes. We will often omit the arguments $Y$, ${\mathbf{X}}$, as they are clear from the context.

where we use the notation $[p]=\{1,\dotsc,p\}$. We further let $s_{0}\equiv|S|$. A copious theoretical literature [CT05, BRT09, BvdG11] shows that, under suitable assumptions on ${\mathbf{X}}$, the LASSO is nearly as accurate as if the support $S$ was known a priori. Namely, for $n=\Omega(s_{0}\log p)$, we have $\|\widehat{\theta}^{n}-\theta_{0}\|_{2}^{2}=O(s_{0}\sigma^{2}(\log p)/n)$.

We will prove in Section 2.1 that $\widehat{\theta}^{u}$ is approximately Gaussian, with mean $\theta_{0}$ and covariance $\sigma^{2}(M\widehat{\Sigma}M)/n$, where $\widehat{\Sigma}=({\mathbf{X}}^{{\sf T}}{\mathbf{X}}/n)$ is the empirical covariance of the feature vectors. This result allows to construct confidence intervals and $p$-values in complete analogy with classical statistics procedures. For instance, letting $Q\equiv M\widehat{\Sigma}M$, $[\widehat{\theta}^{u}_{i}-1.96\sigma\sqrt{Q_{ii}/n},\widehat{\theta}^{u}_{i}+1.96\sigma\sqrt{Q_{ii}/n}]$ is a $95\%$ confidence interval. The size of this interval is of order $\sigma/\sqrt{n}$, which is the optimal (minimum) one, i.e. the same that would have been obtained by knowing a priori the support of $\theta_{0}$. In practice the noise standard deviation is not known, but $\sigma$ can be replaced by any consistent estimator $\widehat{\sigma}$ (see Section 3 for more details on this).

From a technical point of view, our proof starts from a simple decomposition of the de-biased estimator $\widehat{\theta}^{u}$ into a Gaussian part and an error term, already used in [vdGBRD13]. However –departing radically from earlier work– we realize that $M$ need not be a good estimator of $\Sigma^{-1}$ in order for the de-biasing procedure to work. We instead set $M$ as to minimize the error term and the variance of the Gaussian term. As a consequence of this choice, our approach applies to general covariance structures $\Sigma$. By contrast, earlier approaches applied only to sparse $\Sigma$, as in [JM13b], or sparse $\Sigma^{-1}$ as in [vdGBRD13]. The only assumptions we make on $\Sigma$ are the standard compatibility conditions required for high-dimensional consistency [BvdG11]. A detailed comparison of our results with the ones of [vdGBRD13] can be found in Section 2.3.

Our presentation is organized as follows.

considers a general debiased estimator of the form $\widehat{\theta}^{u}=\widehat{\theta}^{n}+(1/n)\,M{\mathbf{X}}^{\sf T}(Y-{\mathbf{X}}\widehat{\theta}^{n})$. We introduce a figure of merit of the pair $M,{\mathbf{X}}$, termed the generalized coherence parameter ${\mu}_{*}({\mathbf{X}};M)$. We show that, if the generalized coherence is small, then the debiasing procedure is effective (for a given deterministic design), see Theorem 2.3.

We then turn to random designs, and show that the generalized coherence parameter can be made as small as $\sqrt{(\log p)/n}$, though a convex optimization procedure for computing $M$. This results in a bound on the bias of $\widehat{\theta}^{u}$, cf. Theorem 2.5: the largest entry of the bias is of order $(s_{0}\log p)/n$. This must be compared with the standard deviation of $\widehat{\theta}^{u}_{i}$, which is of order $\sigma/\sqrt{n}$. The conclusion is that, for $s_{0}=o(\sqrt{n}/\log p)$, the bias of $\widehat{\theta}^{u}$ is negligible.

applies these distributional results to deriving confidence intervals and hypothesis testing procedures for low-dimensional marginals of $\widehat{\theta}_{0}$. The basic intuition is that $\widehat{\theta}^{u}$ is approximately Gaussian with mean $\theta_{0}$, and known covariance structure. Hence standard optimal tests can be applied.

We prove a general lower bound on the power of our testing procedure, in Theorem 3.5. In the special case of Gaussian random designs with i.i.d. rows, we can compare this with the upper bound proved in [JM13b], cf. Theorem 3.6. As a consequence, the asymptotic efficiency of our approach is constant-optimal. Namely, it is lower bounded by a constant $1/\eta_{\Sigma,s_{0}}$ which is bounded away from , cf. Theorem 3.7. (For instance $\eta_{{\rm I},s_{0}}=1$, and $\eta_{\Sigma,s_{0}}$ is always upper bounded by the condition number of $\Sigma$.)

uses the a central limit theorem for triangular arrays to generalize the above results to non-Gaussian noise.

illustrates the above results through numerical simulations both on synthetic and on real data.

Note that our proofs require stricter sparsity $s_{0}$ (or larger sample size $n$) than required for consistent estimation. We assume $s_{0}=o(\sqrt{n}/\log p)$ instead of $s_{0}=o(n/\log p)$ [CT07, BRT09, BvdG11]. The same assumption is made in [vdGBRD13], on top of additional assumptions on the sparsity of $\Sigma^{-1}$.

It is currently an open question whether successful hypothesis testing can be performed under the weaker assumption $s_{0}=o(n/\log p)$. We refer to [JM13c] for preliminary work in that direction. The barrier at $s_{0}=o(\sqrt{n}/\log p)$ is possibly related to an analogous assumption that arises in Gaussian graphical models selection [RSZZ13].

The problem of quantifying statistical significance in high-dimensional parameter estimation is, by comparison, far less understood. Zhang and Zhang [ZZ14], and Bühlmann [Büh13] proposed hypothesis testing procedures under restricted eigenvalue or compatibility conditions [BvdG11]. These papers provide deterministic guarantees but –in order to achieve a certain target significance level $\alpha$ and power $1-\beta$– they require $|\theta_{0,i}|\geq c\,\max\{\sigma s_{0}\log p/\,n,\sigma/\sqrt{n}\}$. The best lower bound [JM13b] shows that any such test requires instead $|\theta_{0,i}|\geq c(\alpha,\beta)\sigma/\sqrt{n}$. (The lower bound of [JM13b] is reproduced as Theorem 3.6 here, for the reader’s convenience.)

Lockhart et al. [LTTT13] develop a test for the hypothesis that a newly added coefficient along the LASSO regularization path is irrelevant. This however does not allow to test arbitrary coefficients at a given value of $\lambda$, which is instead the problem addressed in this paper. These authors further assume that the current LASSO support contains the actual support ${\rm supp}(\theta_{0})$ and that the latter has bounded size.

Belloni, Chernozhukov and collaborators [BCH11, BCW13] consider inference in a regression model with high-dimensional data. In this model the response variable relates to a scalar main regressor and a $p$-dimensional control vector. The main regressor is of primary interest and the control vector is treated as nuisance component. Assuming that the control vector is $s_{0}$-sparse, the authors propose a method to construct confidence regions for the parameter of interest under the sample size requirement $(s_{0}^{2}\log p)/n\to 0$. The proposed method is shown to attain the semi-parametric efficiency bounds for this class of models. The key modeling assumption in this paper is that the scalar regressor of interest is random, and depends linearly on the $p$-dimensional control vector, with a sparse coefficient vector (with sparsity again of order $o(\sqrt{n/\log p})$. This assumption is closely related to the sparse inverse covariance assumption of [vdGBRD13] (with the difference that only one regressor is tested).

After the present paper was submitted for publication, we became aware that Bühlmann and Dezeure [DB13] had independently worked on similar ideas.

2 Preliminaries and notations

In this section we introduce some basic definitions used throughout the paper, starting with simple notations.

We let $\widehat{\Sigma}\equiv{\mathbf{X}}^{\sf T}{\mathbf{X}}/n$ be the sample covariance matrix. For $p>n$, $\widehat{\Sigma}$ is always singular. However, we may require $\widehat{\Sigma}$ to be nonsingular for a restricted set of directions.

In the following, we shall drop the argument $\widehat{\Sigma}$ if clear from the context. Note that a slightly more general definition is used normally [BvdG11, Section 6.13], whereby the condition $\|\theta_{S^{c}}\|_{1}\leq 3\|\theta_{S}\|_{1}$, is replaced by $\|\theta_{S^{c}}\|_{1}\leq L\|\theta_{S}\|_{1}$. The resulting constant $\phi(\widehat{\Sigma},S,L)$ depends on $L$. For the sake of simplicity, we restrict ourselves to the case $L=3$.

The sub-gaussian norm of a random variable $X$, denoted by $\|X\|_{\psi_{2}}$, is defined as

The sub-exponential norm of a random variable $X$, denoted by $\|X\|_{\psi_{1}}$, is defined as

Compensating the bias of the LASSO

In this section we present our characterization of the de-biased estimator $\widehat{\theta}^{u}$ (subsection 2.1). This characterization also clarifies in what sense the LASSO estimator is biased. We discuss this point in subsection 2.2.

For notational simplicity, we shall omit the arguments $Y,{\mathbf{X}},M,\lambda$ unless they are required for clarity. The quality of this debiasing procedure depends of course on the choice of $M$, as well as on the design ${\mathbf{X}}$. We characterize the pair $({\mathbf{X}},M)$ by the following figure of merit.

Note that the minimum coherence can be computed efficiently since $M\mapsto{\mu}_{*}({\mathbf{X}};M)$ is a convex function (even more, the optimization problem is a linear program).

The motivation for our terminology can be grasped by considering the following special case.

The quantity (9) is known as the coherence parameter of the matrix ${\mathbf{X}}/\sqrt{n}$ and was first defined in the context of approximation theory by Mallat and Zhang [MZ93], and by Donoho and Huo [DH01].

Assuming, for the sake of simplicity, that the columns of ${\mathbf{X}}$ are normalized so that $\|{\mathbf{X}}e_{i}\|_{2}=\sqrt{n}$, a small value of the coherence parameter ${\mu}_{*}({\mathbf{X}};{\rm I})$ means that the columns of ${\mathbf{X}}$ are roughly orthogonal. We emphasize however that ${\mu}_{*}({\mathbf{X}};M)$ can be much smaller than its classical coherence parameter ${\mu}_{*}({\mathbf{X}};{\rm I})$. For instance, ${\mu}_{*}({\mathbf{X}};{\rm I})=0$ if and only if ${\mathbf{X}}/\sqrt{n}$ is an orthogonal matrix. On the other hand, ${\mu}_{\rm min}({\mathbf{X}})=0$ if and only if ${\mathbf{X}}$ has rankOf course this example requires $n\geq p$. It is the simplest example that illustrates the difference between coherence and generalized coherence, and it is not hard to find related examples with $n<p$. $p$.

The following theorem is a slight generalization of a result of [vdGBRD13]. Let us emphasize that it applies to deterministic design matrices ${\mathbf{X}}$.

Further, assume that ${\mathbf{X}}$ satisfies the compatibility condition for the set $S={\rm supp}(\theta_{0})$, $|S|\leq s_{0}$, with constant $\phi_{0}$, and has generalized coherence parameter $\mu_{*}=\mu_{*}({\mathbf{X}};M)$, and let $K\equiv\max_{i\in[p]}({\mathbf{X}}^{{\sf T}}{\mathbf{X}}/n)_{ii}$. Then, letting $\lambda=\sigma\sqrt{(c^{2}\log p)/n}$, we have

Further, if $M=M_{\rm min}({\mathbf{X}})$ minimizes the convex cost function $|M\widehat{\Sigma}-{\rm I}|_{\infty}$, then $\mu_{*}$ can be replaced by ${\mu}_{\rm min}({\mathbf{X}})$ in Eq. (11).

The above theorem decomposes the estimation error $(\widehat{\theta}^{*}-\theta_{0})$ into a zero mean Gaussian term $Z/\sqrt{n}$ and a bias term $\Delta/\sqrt{n}$ whose maximum entry is bounded as per Eq. (11). This estimate on $\|\Delta\|_{\infty}$ depends on the design matrix through two constants: the compatibility constant $\phi_{0}$ and the generalized coherence parameter $\mu_{*}({\mathbf{X}};M)$. The former is a well studied property of the design matrix [BvdG11, vdGB09], and assuming $\phi_{0}$ of order one is nearly necessary for the LASSO to achieve optimal estimation rate in high dimension. On the contrary, the definition of $\mu_{*}({\mathbf{X}};M)$ is a new contribution of the present paper.

The next theorem establishes that, for a natural probabilistic model of the design matrix ${\mathbf{X}}$, both $\phi_{0}$ and $\mu_{*}({\mathbf{X}};M)$ can be bounded with probability converging rapidly to one as $n,p\to\infty$. Further, the bound on $\mu_{*}({\mathbf{X}},M)$ hold for the special choice of $M$ that is constructed by Algorithm 1.

Then there exists $c_{*}\leq 2000$ such that the following happens. If $n\geq\nu_{0}\,s_{0}\log(p/s_{0})$, $\nu_{0}\equiv 4c_{*}(C_{\rm max}\kappa^{4}/C_{\rm min})$, $\phi_{0}=C_{\rm min}^{1/2}/2$, and $K\geq 1+20\kappa^{2}\sqrt{(\log p)/n}$, then

For $a>0$, ${\cal G}_{n}={\cal G}_{n}(a)$ be the event that the problem (LABEL:eq:optimization) is feasible for ${\mu}=a\sqrt{(\log p)/n}$, or equivalently

Then, for $n\geq a^{2}C_{\min}\log p/(4e^{2}C_{\max}\kappa^{4})$

The proof of this theorem is given in Section 6.2 (for part $(a)$) and Section 6.3 (part $(b)$).

The proof that event $\mathcal{E}_{n}$ holds with high probability relies crucially on a theorem by Rudelson and Zhou [RZ13, Theorem 6]. Simplifying somewhat, the latter states that, if the restricted eigenvalue condition of [BRT09] holds for the population covariance $\Sigma$, then it holds with high probability for the sample covariance $\widehat{\Sigma}$. (Recall that the restricted eigenvalue condition is implied by a lower bound on the minimum singular valueNote, in particular, at the cost of further complicating the last statement, the condition $\sigma_{\rm min}(\Sigma)=\Omega(1)$ can be further weakened., and that it implies the compatibility condition [vdGB09].)

Finally, by putting together Theorem 2.3 and Theorem 2.4, we obtain the following conclusion.

Consider the linear model (1) and let $\widehat{\theta}^{u}$ be defined as per Eq. (5) in Algorithm 1, with $\mu=a\sqrt{(\log p)/n}$. Then, setting $Z=M{\mathbf{X}}^{{\sf T}}W/\sqrt{n}$, we have

Further, under the assumptions of Theorem 2.4, and for $n\geq\max(\nu_{0}s_{0}\log(p/s_{0}),\nu_{1}\log p)$, $\nu_{1}=\max(1600\kappa^{4},a/4)$, and $\lambda=\sigma\sqrt{(c^{2}\log p)/n}$, we have

Finally, the tail bound (17) holds for any choice of $M$ that is only function of the design matrix ${\mathbf{X}}$, and satisfies the feasibility condition in Eq. (LABEL:eq:optimization), i.e. $|M\widehat{\Sigma}-{\rm I}|_{\infty}\leq{\mu}$.

Assuming $\sigma,C_{\rm min}$ of order one, the last theorem establishes that, for random designs, the maximum size of the ‘bias term’ $\Delta_{i}$ over $i\in[p]$ is:

On the other hand, the ‘noise term’ $Z_{i}$ is roughly of order $\sqrt{[M\widehat{\Sigma}M^{\sf T}]_{ii}}$. Bounds on the variances $[M\widehat{\Sigma}M^{\sf T}]_{ii}$ will be given in Section 3.3 showing that, if $M$ is computed through Algorithm 1, $[M\widehat{\Sigma}M^{\sf T}]_{ii}$ is of order one for a broad family of random designs. As a consequence $|\Delta_{i}|$ is much smaller than $|Z_{i}|$ whenever $s_{0}=o(\sqrt{n}/\log p)$. We summarize these remarks below.

Theorem 2.5 only requires that the support size satisfies $s_{0}=O(n/\log p)$. If we further assume $s_{0}=o(\sqrt{n}/\log p)$, then we have $\|\Delta\|_{\infty}=o(1)$ with high probability. Hence, $\widehat{\theta}^{u}$ is an asymptotically unbiased estimator for $\theta_{0}$.

A more formal comparison of the bias of $\widehat{\theta}^{u}$, and of the one of the LASSO estimator $\widehat{\theta}^{n}$ can be found in Section 2.2 below. Section 2.3 compares our approach with the related one in [vdGBRD13].

As it can be seen from the statement of Theorem 2.3 and Theorem 2.4, the claim of Theorem 2.5 does not rely on the specific choice of the objective function in optimization problem (LABEL:eq:optimization) and only uses the constraint on $\|\widehat{\Sigma}m-e_{i}\|_{\infty}$. In particular it holds for any matrix $M$ that is feasible. On the other hand, the specific objective function problem (LABEL:eq:optimization) minimizes the variance of the noise term ${\rm Var}(Z_{i})$.

2 Discussion: The bias of the LASSO

Theorems 2.3 and 2.4 provide a quantitative framework to discuss in what sense the LASSO estimator $\widehat{\theta}^{n}$ is asymptotically biased, while the de-biased estimator $\widehat{\theta}^{u}$ is asymptotically unbiased.

Given an estimator $\widehat{\theta}^{n}$ of the parameter vector $\theta_{0}$, we define its bias to be the vector

Note that, if the design is random, ${\sf Bias}(\widehat{\theta}^{n})$ is a measurable function of ${\mathbf{X}}$. If the design is deterministic, ${\sf Bias}(\widehat{\theta}^{n})$ is a deterministic quantity as well, and the conditioning is redundant.

Theorem 2.5 with high probability, $\|\Delta\|_{\infty}=O(s_{0}\log p/\sqrt{n})$. The next corollary establishes that this translates into a bound on ${\sf Bias}(\widehat{\theta}^{u})$ for all ${\mathbf{X}}$ in a set that has probability rapidly converging to one as $n$, $p$ get large.

Under the assumptions of Theorem 2.5, let $c_{1}$, $c_{2}$ be defined as per Eqs. (13), (15). Then we have

The proof of this corollary can be found in Appendix B.1.

This result can be contrasted with a converse result for the LASSO estimator. Namely, as stated below, there are choices of the vector $\theta_{0}$, and of the design covariance $\Sigma$, such that ${\sf Bias}(\widehat{\theta}^{n})$ is the sum of two terms. One is of order order $\lambda=c\sigma\sqrt{(\log p)/n}$ and the second is of order $\|{\sf Bias}(\widehat{\theta}^{u})\|_{\infty}$. If $s_{0}$ is significantly smaller than $\sqrt{n/\log p}$ (which is the main regime studied in the rest of the paper), the first term dominates and $\|{\sf Bias}(\widehat{\theta}^{n})\|_{\infty}$ is much larger than $\|{\sf Bias}(\widehat{\theta}^{u})\|_{\infty}$. If on the other hand $s_{0}$ is significantly larger than $\sqrt{n/\log p}$ then $\|{\sf Bias}(\widehat{\theta}^{n})\|_{\infty}$ is of the same order as $\|{\sf Bias}(\widehat{\theta}^{u})\|_{\infty}$. This justify referring to $\widehat{\theta}^{u}$ as to an unbiased estimator.

Notice that, since we want to establish a negative result about the LASSO, it is sufficient to exhibit a specific covariance structure $\Sigma$ satisfying the assumptions of the previous corollary. Remarkably it is sufficient to consider standard designs, i.e. $\Sigma={\rm I}_{p\times p}$.

In particular $\|{\sf Bias}(\widehat{\theta}^{u})\|_{\infty}\leq\lambda/3$ (which follows from $(s_{0}^{2}\log p)/n\leq(c/(3c_{**}))^{2}$) then we have

On the other hand, if $\|{\sf Bias}(\widehat{\theta}^{u})\|_{\infty}\geq\lambda$, then

A formal proof of this statement is deferred to Appendix B.2, but the underlying mathematical mechanism is quite simple and instructive. Recall that the KKT conditions for the LASSO estimator (3) read

Where $\widehat{\theta}^{*}$ a debiased estimator of the general form Eq. (7), for $M={\rm I}$. This suggest that ${\sf Bias}(\widehat{\theta}^{n})$ can be decomposed in two contributions as described above, and as shown formally in Appendix B.2,

3 Comparison with earlier results

In this Section we briefly compare the above debiasing procedure and in particular Theorems 2.3, 2.4 and 2.5 to the results of [vdGBRD13]. In the case of linear statistical models considered here, the authors of [vdGBRD13] construct a debiased estimator of the form (7). However, instead of solving the optimization problem (LABEL:eq:optimization), they follow [ZZ14] and use the regression coefficients of the $i$-th column of ${\mathbf{X}}$ on the other columns to construct the $i$-th row of $M$. These regression coefficients are computed –once again– using the LASSO (node-wise LASSO).

It useful to spell out the most important differences between our contribution and the ones of [vdGBRD13]:

The case of fixed non-random designs is covered by [vdGBRD13, Theorem 2.1], which should be compared to our Theorem 2.3. While in our case the bias is controlled by the generalized coherence parameter, a similar role is played in [vdGBRD13] by the regularization parameters of the nodewise LASSO.

The case of random designs is covered by [vdGBRD13, Theorem 2.2, Theorem 2.4], which should be compared with our Theorem 2.5. In this case, the assumptions underlying our result are significantly less restrictive. More precisely:

[vdGBRD13, Theorem 2.2, Theorem 2.4] assume ${\mathbf{X}}$ to have i.i.d. rows, while we only assume the rows to be independent.

[vdGBRD13, Theorem 2.2, Theorem 2.4] assume the rows inverse covariance matrix $\Sigma^{-1}$ be sparse. More precisely, letting $s_{j}$ be the number of non-zero entries of the $j$-th row of $\Sigma^{-1}$, [vdGBRD13] assumes $\max_{j\in[p]}s_{j}=o(n/\log p)$, that is much smaller than $p$. We do not make any sparsity assumption for $\Sigma^{-1}$, and $s_{j}$ can be as large as $p$.

(In fact [vdGBRD13, Theorem 2.4] also consider the assumption of ${\mathbf{X}}$ with bounded entries, but even stricter sparsity assumptions are made in that case.)

In addition our Theorem 2.5 provides the specific dependence on the maximum and minimum singular value of $\widehat{\Sigma}$.

Statistical inference

A direct application of Theorem 2.5 is to derive confidence intervals and statistical hypothesis tests for high-dimensional models. Throughout, we make the sparsity assumption $s_{0}=o(\sqrt{n}/\log p)$ and omit explicit constants that can be readily derived from Theorem 2.5.

As discussed above, the bias term $\Delta$ is negligible with respect to the random term $Z$ in the decomposition (16), provided the latter has variance of order one. Our first lemma establishes that this is indeed the case.

Let $M=(m_{1},\dotsc,m_{p})^{\sf T}$ be the matrix with rows $m_{i}^{\sf T}$ obtained by solving convex program (LABEL:eq:optimization) in Algorithm 1. Then for all $i\in[p]$,

Using this fact, we can then characterize the asymptotic distribution of the residuals $(\widehat{\theta}^{u}-\theta_{0,i})$. Theorem 2.5 naturally suggests to consider the scaled residual $\sqrt{n}(\widehat{\theta}^{u}_{i}-\theta_{0,i})/(\sigma[M\widehat{\Sigma}M^{\sf T}]_{i,i}^{1/2})$. In the next lemma we consider a slightly more general scaling, replacing $\sigma$ by a consistent estimator $\widehat{\sigma}$.

Consider the linear model (1) and let $\widehat{\theta}^{u}$ be defined as per Eq. (5) in Algorithm 1, with ${\mu}=a\sqrt{(\log p)/n}$ and $\lambda=\sigma\sqrt{(c^{2}\log p)/n}$, with $a,c$ large enough constants. Finally, let $\widehat{\sigma}=\widehat{\sigma}(y,{\mathbf{X}})$ an estimator of the noise level satisfying, for any ${\varepsilon}>0$,

The proof of this lemma can be found in Section 6.5. We also note that the dependence of $a,c$ on $C_{\rm min},C_{\rm max},\kappa$ can be easily reconstructed from Theorem 2.4.

The last lemma requires a consistent estimator of $\sigma$, in the sense of Eq. (30). Several proposal have been made to estimate the noise level in high-dimensional linear regression. A short list of references includes [FL01, FL08, SBvdG10, Zha10, SZ12, BC13, FGH12, RTF13, Dic12, FSW09, BEM13]. Consistency results have been proved or can be proved for several of these estimators.

In order to demonstrate that the consistency criterion (30) can be achieved, we use the scaled LASSO [SZ12] given by

This is a joint convex optimization problem which provides an estimate of the noise level in addition to an estimate of $\theta_{0}$.

The following lemma uses the analysis of [SZ12] to show that $\widehat{\sigma}$ thus defined satisfies the consistency criterion (30).

Under the assumptions of Lemma 3.2, let $\widehat{\sigma}=\widehat{\sigma}(\widetilde{\lambda})$ be the scaled LASSO estimator of the noise level, see Eq. (32), with $\widetilde{\lambda}=10\sqrt{(2\log p)/n}$. Then $\widehat{\sigma}$ thus satisfies Eq. (30).

The proof of this lemma is fairly straightforward and can be found in Appendix C.

2 Confidence intervals

In view of Lemma 3.2, it is quite straightforward to construct asymptotically valid confidence intervals. Namely, for $i\in[p]$ and significance level $\alpha\in(0,1)$, we let

Consider the linear model (1) and let $\widehat{\theta}^{u}$ be defined as per Eq. (5) in Algorithm 1, with ${\mu}=a\sqrt{(\log p)/n}$ and $\lambda=\sigma\sqrt{(c^{2}\log p)/n}$, with $a,c$ large enough constants. Finally, let $\widehat{\sigma}=\widehat{\sigma}(y,{\mathbf{X}})$ a consistent estimator of the noise level in the sense of Eq. (30). Then the confidence interval $J_{i}(\alpha)$ is asymptotically valid, namely

The proof is an immediate consequence of Lemma 3.2 since

3 Hypothesis testing

An important advantage of sparse linear regression models is that they provide parsimonious explanations of the data in terms of a small number of covariates. The easiest way to select the ‘active’ covariates is to choose the indexes $i$ for which $\widehat{\theta}_{i}^{n}\neq 0$. This approach however does not provide a measure of statistical significance for the finding that the coefficient is non-zero.

More precisely, we are interested in testing an individual null hypothesis $H_{0,i}:\theta_{0,i}=0$ versus the alternative $H_{A,i}:\theta_{0,i}\neq 0$, and assigning $p$-values for these tests. We construct a $p$-value $P_{i}$ for the test $H_{0,i}$ as follows:

The decision rule is then based on the $p$-value $P_{i}$:

where $\alpha$ is the fixed target Type I error probability. We measure the quality of the test $\widehat{T}_{i,{\mathbf{X}}}(y)$ in terms of its significance level $\alpha_{i}$ and statistical power $1-\beta_{i}$. Here $\alpha_{i}$ is the probability of type I error (i.e. of a false positive at $i$) and $\beta_{i}$ is the probability of type II error (i.e. of a false negative at $i$).

Consider the linear model (1) and let $\widehat{\theta}^{u}$ be defined as per Eq. (5) in Algorithm 1, with ${\mu}=a\sqrt{(\log p)/n}$ and $\lambda=\sigma\sqrt{(c^{2}\log p)/n}$, with $a,c$ large enough constants. Finally, let $\widehat{\sigma}=\widehat{\sigma}(y,{\mathbf{X}})$ a consistent estimator of the noise level in the sense of Eq. (30), and $\widehat{T}$ be the test defined in Eq. (39).

Then the following holds true for any fixed sequence of integers $i=i(n)$:

Theorem 3.5 is proved in Appendix 6.6. It is easy to see that, for any $\alpha>0$, $u\mapsto G(\alpha,u)$ is continuous and monotone increasing. Moreover, $G(\alpha,0)=\alpha$ which is the trivial power obtained by randomly rejecting $H_{0,i}$ with probability $\alpha$. As ${\gamma}$ deviates from zero, we obtain nontrivial power. Notice that in order to achieve a specific power $\beta>\alpha$, our scheme requires ${\gamma}\geq c_{\beta}(\sigma/\sqrt{n})$, for some constant $c_{\beta}$ that depends on $\beta$. This is because $\Sigma^{-1}_{i,i}\leq\sigma_{\max}(\Sigma^{-1})\leq(\sigma_{\min}(\Sigma))^{-1}=O(1)$.

The authors of [JM13b] prove an upper bound for the minimax power of tests with a given significance level $\alpha$, under random designs. For the readers’ convenience, we recall here this result. (The following is a restatement of [JM13b, Theorem 2.3], together with a standard estimate on the tail of chi-squared random variables.)

for any $\xi\in[0,(3/2)\sqrt{n-s_{0}+1}]$.

The intuition behind this bound is straightforward: the power of any test for $H_{0,i}:\,\theta_{0,i}=0$ is upper bounded by the power of an oracle test that is given access to the support of $\theta_{0}$, with the eventual exclusion of $i$. Namely, the oracle has access to ${\rm supp}(\theta_{0})\setminus\{i\}$ and outputs a test for $H_{0,i}$. Computing the minimax power of such oracle reduces to a classical hypothesis testing problem.

Let us emphasize that the last theorem applies to Gaussian random designs. Since this theorem establishes a negative result (an upper bound on power) it makes sense to consider this somewhat more specialized setting.

Using this upper bound, we can restate Theorem 3.5 as follows.

Consider a Gaussian random design model that satisfies the conditions of Theorem 3.5, and let $\widehat{T}$ be the testing procedure defined in Eq. (39), with $\widehat{\theta}^{u}$ as in Algorithm 1. Further, let

Under the sparsity assumption $s_{0}=o(\sqrt{n}/\log p)$, the following holds true. If $\{T_{i,{\mathbf{X}}}\}$ is any sequence of tests with $\lim\sup_{n\to\infty}\alpha_{i,n}(T)\leq\alpha$, then

In other words, the asymptotic efficiency of the test $\widehat{T}$ is at least $1/\eta_{\Sigma,s_{0}}$.

Hence, our test $\widehat{T}$ has nearly optimal power in the following sense. It has power at least as large as the power of any oter test $T$, provided the latter is applied to a sample size increased by a factor $\eta_{\Sigma,s_{0}}$.

Further, under the assumptions of Theorem 2.5, the factor $\eta_{\Sigma,s_{0}}$ is a bounded constant. Indeed

since $\Sigma^{-1}_{ii}\leq(\sigma_{\min}(\Sigma))^{-1}$, and $\Sigma_{i|S}\leq\Sigma_{i,i}\leq\sigma_{\max}(\Sigma)$ due to $\Sigma_{S,S}\succ 0$.

Note that $n$, $\gamma$ and $\sigma$ appears in our upper bound (44) in the combination $\gamma\sqrt{n}/\sigma$, which is the natural measure of the signal-to-noise ratio (where, for simplicity, we neglected $s_{0}=o(\sqrt{n}/\log p)$ with respect to $n$). Hence, the above result can be restated as follows. The test $\widehat{T}$ has power at least as large as the power of any oter test $T$, provided the latter is applied at a noise level augmented by a factor $\sqrt{\eta_{\Sigma,s_{0}}}$.

4 Generalization to simultaneous confidence intervals

In many situations, it is necessary to perform statistical inference on more than one of the parameters simultaneously. For instance, we might be interested in performing inference about $\theta_{0,R}\equiv(\theta_{0,i})_{i\in R}$ for some set $R\subseteq[p]$.

The simplest generalization of our method is to the case in which $|R|$ stays finite as $n,p\to\infty$. In this case we have the following generalization of Lemma 3.2. (The proof is the same as for Lemma 3.2, and hence we omit it.)

Under the assumptions of Lemma 3.2, define

where $(a_{1},\dots,a_{k})\leq(b_{1},\dots,b_{k})$ indicates that $a_{1}\leq b_{1}$,…$a_{k}\leq b_{k}$, and $\Phi_{k}(x)=\Phi(x_{1})\cdots\Phi(x_{k})$.

This lemma allows to construct confidence regions for low-dimensional projections of $\theta_{0}$, much in the same way as we used Lemma 3.2 to compute confidence intervals for one-dimensional projections in Section 3.2.

Then Lemma 3.8 implies (under the assumptions stated there) that $J_{R}(\alpha)$ is a valid confidence region

A more challenging regime is the one of large-scale inference, that corresponds to $|R(n)|\to\infty$ with $n$. Even in the seemingly simple case in which a correct $p$-value is given for each individual coordinate, the problem of aggregating them has attracted considerable amount of work, see e.g. [Efr10] for an overview.

In order to achieve familywise error control, we adopt a standard trick based on Bonferroni inequality. Given $p$-values defined as per Eq. (38), we let

Then we have the following error control guarantee.

Consider the linear model (1) and let $\widehat{\theta}^{u}$ be defined as per Eq. (5) in Algorithm 1, with ${\mu}=a\sqrt{(\log p)/n}$ and $\lambda=\sigma\sqrt{(c^{2}\log p)/n}$, with $a,c$ large enough constants. Finally, let $\widehat{\sigma}=\widehat{\sigma}(y,{\mathbf{X}})$ be a consistent estimator of the noise level in the sense of Eq. (30), and $\widehat{T}$ be the test defined in Eq. (54). Then:

The proof of this theorem is similar to the one of Lemma 3.2 and Theorem 3.5, and is deferred to Appendix D.

Non-Gaussian noise

As can be seen from the proof of Theorem 2.5, $Z=M{\mathbf{X}}^{\sf T}W/\sqrt{n}$, and since the noise is Gaussian, i.e., $W\sim{\sf N}(0,\sigma^{2}{\rm I})$, we have $Z|{\mathbf{X}}\sim{\sf N}(0,\sigma^{2}M\widehat{\Sigma}M^{\sf T})$. We claim that the distribution of the coordinates of $Z$ is asymptotically Gaussian, even if $W$ is non-Gaussian, provided the definition of $M$ is modified slightly. As a consequence, the definition of confidence intervals and $p$-values in Corollary 3.4 and (38) remain valid in this broader setting.

then $\sum_{j=1}^{n}\xi_{j}/\sqrt{n}|{\mathbf{X}}\overset{{\rm d}}{\longrightarrow}{\sf N}(0,1)$, from which we can build the valid $p$-values as in (38).

In order to ensure that the Lindeberg condition holds, we modify the optimization problem (LABEL:eq:optimization_mod) as follows:

Next theorem shows the validity of the proposed $p$-values in the non-Gaussian noise setting.

Let $M=(m_{1},\dotsc,m_{p})^{\sf T}$ be the matrix with rows $m_{i}^{\sf T}$ obtained by solving optimization problem (LABEL:eq:optimization_mod). Then under the assumptions of Theorem 2.5, and for sparsity level $s_{0}=o(\sqrt{n}/\log p)$, an asymptotic two-sided confidence interval for $\theta_{0,i}$ with significance $\alpha$ is given by $I_{i}=[\widehat{\theta}^{u}_{i}-\delta(\alpha,n),\widehat{\theta}^{u}_{i}+\delta(\alpha,n)]$ where

Further, an asymptotically valid $p$-value $P_{i}$ for testing null hypothesis $H_{0,i}$ is constructed as:

Numerical experiments

Regarding the regression coefficient, we consider a uniformly random support $S\subseteq[p]$, with $|S|=s_{0}$ and let $\theta_{0,i}=b$ for $i\in S$ and $\theta_{0,i}=0$ otherwise. The measurement errors are $W_{i}\sim{\sf N}(0,1)$, for $i\in[n]$. We consider several configurations of $(n,p,s_{0},b)$ and for each configuration report our results based on $20$ independent realizations of the model with fixed design and fixed regression coefficients. In other words, we repeat experiments over $20$ independent realization of the measurement errors.

We use the regularization parameter $\lambda=4\widehat{\sigma}\sqrt{(2\log p)/n}$, where $\widehat{\sigma}$ is given by the scaled LASSO as per equation (32) with $\widetilde{\lambda}=10\sqrt{(2\log p)/n}$. Furthermore, parameter $\mu$ (cf. Eq. (LABEL:eq:optimization)) is set to

This choice of $\mu$ is guided by Theorem 2.4 $(b)$.

Throughout, we set the significance level $\alpha=0.05$.

Confidence intervals. For each configuration, we consider $20$ independent realizations of measurement noise and for each parameter $\theta_{0,i}$, we compute the average length of the corresponding confidence interval, denoted by ${\rm Avglength}(J_{i}(\alpha))$ where $J_{i}(\alpha)$ is given by equation (33) and the average is taken over the realizations. We then define

We also consider the average length of intervals for the active and inactive parameters, as follows:

Similarly, we consider average coverage for individual parameters. We define the following three metrics:

False positive rates and statistical powers. Table 2 summarizes the false positive rates and the statistical powers achieved by our proposed method, the multisample-splitting method [MMB09], and the ridge-type projection estimator [Büh13] for several configurations. The results are obtained by taking average over $20$ independent realizations of measurement errors for each configuration. As we see the multisample-splitting achieves false positive rate 0 on all of the configurations considered here, making no type I error. However, the true positive rate is always smaller than that of our proposed method. By contrast, our method achieves false positive rate close to the pre-assigned significance level $\alpha=0.05$ and obtains much higher true positive rate. Similar to the multisample-splitting, the ridge-type projection estimator is conservative and achieves false positive rate smaller than $\alpha$. This, however, comes at the cost of a smaller true positive rate than our method. It is worth noting that an ideal testing procedure should allow to control the level of statistical significance $\alpha$, and obtain the maximum true positive rate at that level.

Here, we used the R-package hdi to test multisample-splitting and the ridge-type projection estimator.

Let $Z=(z_{i})_{i=1}^{p}$ denote the vector with $z_{i}\equiv\sqrt{n}(\widehat{\theta}^{u}_{i}-\theta_{0,i})/\widehat{\sigma}\sqrt{[M\widehat{\Sigma}M^{\sf T}]_{i,i}}$. Fig. 2 shows the sample quantiles of $Z$ versus the quantiles of the standard normal distribution for one realization of the configuration $(n,p,s_{0},b)=(1000,600,10,1)$. The scattered points are close to the line with unit slope and zero intercept. This confirms the result of Theorem 3.2 regarding the gaussianity of the entries $z_{i}$.

For the same problem, in Fig. 3 we plot the empirical CDF of the computed $p$-values restricted to the variables outside the support. Clearly, the $p$-values for these entries are uniformly distributed as expected.

2 Real data

As a real data example, we consider a high-throughput genomic data set concerning riboflavin (vitamin $B_{2}$) production rate. This data set is made publicly available by [BKM14] and contains $n=71$ samples and $p=4,088$ covariates corresponding to $p=4,088$ genes. For each sample, there is a real-valued response variable indicating the logarithm of the riboflavin production rate along with the logarithm of the expression level of the $p=4,088$ genes as the covariates.

Following [BKM14], we model the riboflavin production rate as a linear model with $p=4,088$ covariates and $n=71$ samples, as in Eq. (1). We use the ${\sf R}$ package ${\sf glmnet}$ [FHT10] to fit the LASSO estimator. Similar to the previous section, we use the regularization parameter $\lambda=4\widehat{\sigma}\sqrt{(2\log p)/n}$, where $\widehat{\sigma}$ is given by the scaled LASSO as per equation (32) with $\widetilde{\lambda}=10\sqrt{(2\log p)/n}$. This leads to the choice $\lambda=0.036$. The resulting model contains 30 genes (plus an intercept term) corresponding to the nonzero parameters of the lasso estimator.

We use Eq. (38) to construct $p$-values for different genes. Adjusting FWER to $5\%$ significance level, we find two significant genes, namely genes YXLD-at and YXLE-at. By contrast, the multisample-splitting method proposed in [MMB09] finds only the gene YXLD-at at the FWER-adjusted $5\%$ significance level. Also the Ridge-type projection estimator, proposed in [Büh13], returns no significance gene. (See [BKM14] for further discussion on these methods.) This indicates that these methods are more conservative and produce typically larger $p$-values.

In Fig. 4 we plot the empirical CDF of the computed $p$-values for riboflavin example. Clearly the plot confirms that the $p$-values are distributed according to uniform distribution.

Proofs

Substituting $Y={\mathbf{X}}\theta_{0}+W$ in the definition (7), we get

with $Z,\Delta$ defined as per the theorem statement. Further $Z$ is Gaussian with the stated covariance because it is a linear function of the Gaussian vector $W\sim{\sf N}(0,\sigma^{2}\,{\rm I}_{p\times p})$.

We are left with the task of proving the bound (11) on $\Delta$. Note that by definition (2.1), we have

By [BvdG11, Theorem 6.1, Lemma 6.2], we have, for any $\lambda\geq 4\sigma\sqrt{2K\log(pe^{t^{2}/2})/n}$

(More precisely, we consider the trivial generalization of [BvdG11, Lemma 6.2] to the case $({\mathbf{X}}^{T}{\mathbf{X}}/n)_{ii}\leq K$, instead of $({\mathbf{X}}^{T}{\mathbf{X}}/n)_{ii}=1$ for all $i\in[p]$.)

Substituting Eq. (67) in the last bound, we get

Finally, the claim follows by selecting $t$ so that $e^{t^{2}/2}=p^{c_{0}}$.

2 Proof of Theorem 2.4.(a)𝑎(a)

Note that the event ${\cal E}_{n}$ requires two conditions. Hence, its complement

We will bound separately the probability of ${\cal B}_{1,n}$ and the probability of ${\cal B}_{2,n}$. The claim of Theorem 2.4.$(a)$ follows by union bound.

It is also useful to recall the notion of restricted eigenvalue, introduced by Bickel, Ritov and Tsybakov [BRT09].

Rudelson and Zhou [RZ13] prove that, if the population covariance satisfies the restricted eigenvalue condition, then the sample covariance satisfies it as well, with high probability. More precisely [RZ13, Theorem 6], the following happens for some $c_{*}\leq 2000$, $m\equiv c_{*}s_{0}C_{\rm max}^{2}/\phi_{\rm RE}^{2}(\Sigma,s_{0},9)$, and every $n\geq 4c_{*}m\kappa^{4}\log(60ep/(m\kappa))$ we have

Note that $\phi_{\rm RE}(\Sigma,s_{0},9)\geq\sigma_{\rm min}(\Sigma)^{1/2}\geq C_{\rm min}^{1/2}$ and, by Cauchy-Schwartz $\min_{S:|S|\leq s_{0}}\phi(\widehat{\Sigma},S)\geq\phi_{\rm RE}(\widehat{\Sigma},s_{0},3)$. With the definitions in the statement (cf. Eq. (13)), we therefore have

By Bernstein-type inequality for centered subexponential random variables [Ver12], we get

Hence, for all ${\varepsilon}$ such that ${\varepsilon}/(e\kappa^{2})\in[\sqrt{(48\log p)/n},4]$,

3 Proof of Theorem 2.4.(b)𝑏(b)

and hence the statement follows immediately from the following estimate.

We have $\sigma_{\min}(\Sigma)\geq C_{\min}>0$, and $\sigma_{\max}(\Sigma)\leq C_{\max}<\infty$.

The rows of $X\Sigma^{-1/2}$ are sub-gaussian with $\kappa=\|\Sigma^{-1/2}X_{1}\|_{\psi_{2}}$.

Let $\widehat{\Sigma}=({\mathbf{X}}^{\sf T}{\mathbf{X}})/n$ be the empirical covariance. Then, for any constant $C>0$, the following holds true.

with $c_{2}=(a^{2}C_{\min})/(24e^{2}\kappa^{4}C_{\max})-2$.

Moreover, for any two random variables $X$ and $Y$, we have

Let $\kappa^{\prime}=2\sqrt{C_{\max}/C_{\min}}\kappa^{2}$. Applying Bernstein-type inequality for centered sub-exponential random variables [Ver12], we get

Choosing ${\varepsilon}=a\sqrt{(\log p)/n}$, and assuming $n\geq[a/(e\kappa^{\prime})]^{2}\log p$, we arrive at

The result follows by union bounding over all possible pairs $i,j\in[p]$. ∎

4 Proof of Theorem 2.5

be a shorthand for the bound on $\|\Delta\|_{\infty}$ appearing in Eq. (17). Then we have

where, in the firsr equation ${\cal A}^{c}$ denotes the complement of event ${\cal A}$ and the second inequality follows from Theorem 2.4. Notice, in particular, that the bound (13) can be applied for $K=3/2$ since, under the present assumptions $20\kappa^{2}\sqrt{(\log p)/n}\leq 1/2$.

Here the last inequality follows from Theorem 2.3 applied per given ${\mathbf{X}}\in\mathcal{E}_{n}(C_{\min}^{1/2}/2,s_{0},3/2)\cap{\cal G}_{n}(a)$ and hence using the bound (11) with $\phi_{0}=C_{\min}^{1/2}/2$, $K=3/2$, $\mu_{*}=a\sqrt{(\log p)/n}$.

5 Proof of Lemma 3.2

We will prove that, under the stated assumptions

A matching lower bound follows by a completely analogous argument.

which proves our claim. In order to prove Eq. (88), fix ${\varepsilon}>0$ and write

By taking the limit and using the assumption (30), we obtain

Since ${\varepsilon}>0$ is arbitrary, it is therefore sufficient to show that the limit on the right hand side vanishes for any ${\varepsilon}>0$.

Note that $[M\widehat{\Sigma}M^{\sf T}]_{i,i}\geq 1/(4\widehat{\Sigma}_{ii})$ for all $n$ large enough, by Lemma 3.1, and since ${\mu}=a\sqrt{(\log p)/n}\to 0$ as $n,p\to\infty$. We have therefore

where the last inequality follows from Eq. (17) since $s_{0}=o(\sqrt{n}/\log p)$ and hence $(16acs_{0}\log p)/(C_{\rm min}\sqrt{n})\leq{\varepsilon}/8$ for all $n$ large enough.

This completes the proof of Eq. (88). The matching lower bound follows by the same argument.

6 Proof of Theorem 3.5

We begin with proving Eq. (42). Defining $Z_{i}\equiv\sqrt{n}(\widehat{\theta}^{u}_{i}-\theta_{0,i})/(\widehat{\sigma}[M\widehat{\Sigma}M^{\sf T}]_{i,i}^{1/2})$, we have

where the last inequality follows from Lemma 3.2.

We next prove Eq. (43). Recall that $\Sigma^{-1}_{\cdot,i}$ is a feasible solution of (LABEL:eq:optimization), for $1\leq i\leq p$ with probability at least $1-2p^{-c_{2}}$, as per Lemma 6.2). On this event, letting $m_{i}$ be the solution of the optimization problem (LABEL:eq:optimization), we have

Therefore, by Borel-Cantelli (since we can make $c_{2}\geq 2$ by a suitable choice of $a$), we have, almost surely

This bound leads to a lower bound for the power. First of all, a straightforward manipulation yields as follows, letting $z_{*}\equiv\Phi^{-1}(1-\alpha/2)$:

Here $(a)$ follows from Eq. (101) and the fact $|\theta_{0,i}|\geq{\gamma}$.

7 Proof of Theorem 4.1

Under the assumptions of Theorem 2.5 and assuming $s_{0}=o(\sqrt{n}/\log p)$, we have

with $\|\Delta\|_{\infty}=o(1)$. Using Lemma 3.1, we have

The following lemma characterizes the limiting distribution of $Z_{i}|{\mathbf{X}}$ which implies the validity of the proposed $p$-value $P_{i}$ and confidence intervals.

A.J. is supported by a Caroline and Fabian Pease Stanford Graduate Fellowship. This work was partially supported by the NSF CAREER award CCF-0743978, the NSF grant DMS-0806211, and the grants AFOSR/DARPA FA9550-12-1-0411 and FA9550-13-1-0036.

Appendix A Proof of technical lemmas

Let $C_{i}({\mu})$ be the optimal value of the optimization problem (LABEL:eq:optimization). We claim that

To prove this claim notice that the constraint implies (by considering its $i$-th component):

The minimum over $m$ is achieved at $m=ce_{i}/2$. Plugging in for $m$, we get

Optimizing this bound over $c$, we obtain the claim (102), with the optimal choice being $c=2(1-{\mu})/\widehat{\Sigma}_{ii}$.

A.2 Proof of Lemma 6.3

where $\widetilde{W}_{j}=W_{j}/\sigma$ and the last limit follows by taking $a>(1/2-\beta)^{-1}$ as per the assumptions.

Using Lindenberg central limit theorem, we obtain $Z_{i}|{\mathbf{X}}$ converges weakly to standard normal distribution, and hence, ${\mathbf{X}}$-almost surely

What remains is to show that with high probability all the $p$ optimization problems in (LABEL:eq:optimization_mod) are feasible. In particular, we show that $\Sigma^{-1}_{i,\cdot}$ is a feasible solution to the $i$-th optimization problem, for $i\in[p]$. By Lemma 6.2, $|\Sigma^{-1}\widehat{\Sigma}-{\rm I}|_{\infty}\leq{\mu}$, with high probability. Moreover,

Using tail bound for sub-gaussian variables $\Sigma^{-1}_{i,\cdot}X_{j}$ and union bounding over $j\in[n]$, we get

for some constant $c>0$. Note that $s_{0}=o(\sqrt{n}/\log p)$ implies $p=e^{o(n^{2\beta})}$. Hence, eventually almost surely, $\Sigma^{-1}_{i,\cdot}$ is a feasible solution to optimization problem (LABEL:eq:optimization_mod), for all $i\in[p]$.

Appendix B Corollaries of Theorem 2.5

By Theorem 2.3, for any ${\mathbf{X}}\in\mathcal{E}_{n}(\sqrt{C_{\min}}/2,s_{0},3/2)\cap{\cal G}_{n}(a)$, we have

(This is obtained by setting $\phi_{0}=C_{\min}^{1/2}/2$, $K=3/2$, $\mu_{*}=a\sqrt{(\log p)/n}$ in Eq. (11). Hence

which coincides with Eq. (21). The probability estimate (22) simply follows from Theorem 2.4 using union bound.

B.2 Proof of Corollary 2.8

By Theorem 2.4.$(a)$, we have (setting $C_{\min}=C_{\max}=\kappa=1$):

Further, by Lemma 6.2, with $\widehat{\Sigma}\equiv{\mathbf{X}}^{{\sf T}}{\mathbf{X}}/n$, we have

Finally, by an obvious consequence of the proof of Theorem 2.4.$(a)$

we have the desired probability bound (24).

where $\widehat{\theta}^{n}=\widehat{\theta}^{n}(Y,{\mathbf{X}};\lambda)$ is the LASSO solution with $\lambda=\sigma\sqrt{(c^{2}\log p)/n}$. By Theorem 2.3, we have, for any ${\mathbf{X}}\in{\cal B}_{n}$

whence, proceeding as in the proof in the last section, we get, for some universal numerical constant $c_{**}$,

Note that $v(\widehat{\theta}^{n})_{i}=1$ whenever $\widehat{\theta}^{n}_{i}>0$ and, and $|v(\widehat{\theta}^{n})_{i}|\leq 1$, and therefore (letting $b_{0}\equiv c_{**}\sigma(s_{0}\log p)/n$)

with $\Phi(x)$ the standard normal distribution function, and in the last inequality we used the fact that $\max_{i\in[p]}\widehat{\Sigma}_{ii}\leq 3/2$ on ${\cal B}_{n}$. We then choose $\theta_{0}$ so that $\theta_{0,i}\geq b_{0}+\lambda+\sqrt{30\sigma^{2}/n}$, for $i\in[p]$ in the support of $\theta_{0}$. We therefore obtain

This finishes the proof of Eq. (23). Equations (26) and (27) are obtained by substituting $\lambda=c\sigma\sqrt{(\log p)/n}$ and using Eq. (23).

Appendix C Proof of Lemma 3.3

Let $\mathcal{E}_{n}=\mathcal{E}_{n}(\phi_{0},s_{0},K)$ be the event defined as per Theorem 2.4.$(a)$. In particular, we take $\phi_{0}=C_{\rm min}^{1/2}/2$, and $K\geq 1+20\kappa^{2}\sqrt{(\log p)/n}$ (for, instance $K=1.1$ will work for all $n$ large enough since $(s_{0}\log p)^{2}/n\to 0$, with $s_{0}\geq 1$, by assumption). Further note that we can assume without loss of generality $n\geq\nu_{0}\,s_{0}\log(p/s_{0})$, since $s_{0}=o(\sqrt{n}/\log p)$. Fixing ${\varepsilon}>0$, we have therefore

where $c_{1}>0$ is a constant defined as per Theorem 2.4.$(a)$.

where the last inequality follows for all $n$ large enough since $s_{0}=o(\sqrt{n}/\log p)$.

where we note that the right hand side is independent of $\theta_{0}$. The first term vanishes as $n\to\infty$ by a standard tail bound on the supremum of $p$ Gaussian random variables. The second term also vanishes because it is controlled by the tail of a chi-squared random variable [SZ12].

Appendix D Proof of Theorem 3.9

Since the second term vanishes as $n\to\infty$ by assumption Eq. (30), it is sufficient to consider the first term. Using Bonferroni inequality, letting z_{\alpha}({\varepsilon})\equiv(1-{\varepsilon})\Phi^{-1}\big{(}1-\frac{\alpha}{2p}\big{)}, we have

where, by Theorem 2.5, $\widetilde{Z}_{i}\sim{\sf N}(0,1)$ and $\Delta_{i}$ is given by Eq. (16). We then have

where in the first inequality, we used $[M\widehat{\Sigma}M^{\sf T}]_{i,i}\geq 1/(4\widehat{\Sigma}_{ii})$ for all $n$ large enough, by Lemma 3.1, and since ${\mu}=a\sqrt{(\log p)/n}\to 0$ as $n,p\to\infty$. Now the second term in the right hand side of Eq. (132) vanishes by Theorem 2.4.$(a)$, and the last term is zero by Theorem 2.5 since $\sqrt{n}/\log(p)\geq s_{0}\geq 1$. Therefore