Multiple Testing and Variable Selection along the path of the Least Angle Regression

TL;DR

This paper develops exact statistical testing procedures for variable selection along the LARS path in high-dimensional Gaussian noise settings, controlling error rates and false discoveries with theoretical guarantees and practical demonstrations.

Contribution

It introduces the Generalized t-Spacing tests and a false negatives test with exact non-asymptotic level, extending previous work to unknown noise levels and general designs.

Findings

Exact non-asymptotic level control for new tests

Perfect fit of theoretical results in numerical experiments

Control of false discovery rate under orthogonal design

Abstract

We investigate multiple testing and variable selection using the Least Angle Regression (LARS) algorithm in high dimensions under the assumption of Gaussian noise. LARS is known to produce a piecewise affine solution path with change points referred to as the knots of the LARS path. The key to our results is an expression in closed form of the exact joint law of a $K$-tuple of knots conditional on the variables selected by LARS, namely the so-called post-selection joint law of the LARS knots. Numerical experiments demonstrate the perfect fit of our findings. This paper makes three main contributions. First, we build testing procedures on variables entering the model along the LARS path in the general design case when the noise level can be unknown. These testing procedures are referred to as the Generalized $t$-Spacing tests (GtSt) and we prove that they have an exact non-asymptotic level (i.e., the Type I error is exactly controlled). This extends work of (Taylor et al., 2014) where the spacing test works for consecutive knots and known variance. Second, we introduce a new exact multiple false negatives test after model selection in the general design case when the noise level may be unknown. We prove that this testing procedure has exact non-asymptotic level for general design and unknown noise level. Third, we give an exact control of the false discovery rate under orthogonal design assumption. Monte Carlo simulations and a real data experiment are provided to illustrate our results in this case. Of independent interest, we introduce an equivalent formulation of the LARS algorithm based on a recursive function.

Equations

$$Y=X{\beta }^0+\eta \sim {\mathcal{N}}_n(X{\beta }^0,{\sigma }^2{\mathrm{Id}}_n),$$

$$\widehat{\imath}_{k}:=\overline{\imath}_{k}+p\,\Big{(}\frac{1-\varepsilon_{k}}{2}\Big{)}\in[2p],$$

$${\lambda }_k={\beta }_{(k)}^0,$$

$$\mathcal{E}:=\big{\{}\widehat{\imath}_{1}=\imath_{1},\ldots,\widehat{\imath}_{K}=\imath_{K},\lambda_{K+1}\big{\}}.$$

$${\mathrm{Q}}^{-1}_{(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K},\lambda_{K+1})}\,\Big{(}\prod_{k=1}^{K}\varphi_{m_{k},v_{k}^{2}}(\ell_{k})\Big{)}\,\mathds{1}_{\{\ell_{1}\geq\ell_{2}\geq\cdots\geq\ell_{K}\geq\lambda_{K+1}\}}\text{ at point }(\ell_{1},\ell_{2},\ldots,\ell_{K}),$$

$$\{0\}=:H_0\subset H_1\subset \cdots \subset H_k:=\mathrm{Span}(X_{{\overline{ı}}_1},\dots ,X_{{\overline{ı}}_k})\subset \cdots \subset H_K.$$

$$\forall k\in [K],m_k=c_k{\epsilon }_k⟨X_{{\overline{ı}}_k},P_{k-1}^{\perp }(X{\beta }^0)⟩\text{and}{\rho }_k^2=d_k\sin \measuredangle (X_{{\overline{ı}}_k},H_{k-1}),$$

$$\text{mathds}{H}_0:\text{‘}X{\beta }^0\in H_{a_0}\text{’},$$

$$\text{mathds}{H}_0:\text{‘‘}X{\beta }^0\in H_a\text{”},$$

$${\hat{\alpha }}_{abc}={\hat{\alpha }}_{abc}({\lambda }_a,{\lambda }_b,{\lambda }_c,{\hat{ı}}_1,\dots ,{\hat{ı}}_{c-1}):=1-\frac{\text{mathds}{F}_{abc}({\lambda }_b)}{\text{mathds}{F}_{abc}({\lambda }_a)}\sim \mathcal{U}(0,1),$$

$${\mathcal{S}}_{abc}:=\text{mathds}{1}_{\{{\hat{\alpha }}_{abc}\le \alpha \}},$$

$$\text{mathds}{H}_0:\text{‘}X{\beta }^0\in H_a\text{’},$$

$${\hat{\beta }}_{abc}={\hat{\beta }}_{abc}({\Lambda }_a,{\Lambda }_b,{\Lambda }_c,{\hat{ı}}_1,\dots ,{\hat{ı}}_K):=1-\frac{{\widetilde{\text{mathds}F}}_{abc}({\Lambda }_b)}{{\widetilde{\text{mathds}F}}_{abc}({\Lambda }_a)}\sim \mathcal{U}(0,1),$$

$$\text{mathds}{H}_0:\text{‘}X{\beta }^0\in H_{a_0}\text{’},$$

$${\hat{\alpha }}_{ab(c+1)}\preccurlyeq {\hat{\alpha }}_{abc}\text{ and }{\hat{\alpha }}_{a(b-1)c}\preccurlyeq {\hat{\alpha }}_{abc}\text{ and }{\hat{\alpha }}_{(a-1)bc}\preccurlyeq {\hat{\alpha }}_{abc},$$

$$\widehat{\imath}_{k}:=\overline{\imath}_{k}+p\,\Big{(}\frac{1-\varepsilon_{k}}{2}\Big{)},$$

$$Z:=(\overline{Z},-\overline{Z})=(X^{\top }Y,-X^{\top }Y),$$

$${\mu }^0:=(\overline{R}{\beta }^0,-\overline{R}{\beta }^0)=(X^{\top }X{\beta }^0,-X^{\top }X{\beta }^0)=({\overline{\mu }}^0,-{\overline{\mu }}^0),$$

$$R=\left[\begin{array}[]{cc}\overline{R}&-\overline{R}\\ -\overline{R}&\overline{R}\end{array}\right]=\left[\begin{array}[]{cc}X^{\top}X&-X^{\top}X\\ -X^{\top}X&X^{\top}X\end{array}\right].$$

$$i_{k}:=j_{k}+p\,\Big{(}\frac{1-s_{k}}{2}\Big{)}.$$

$$S^{0}:=\big{\{}k\in[p]\ :\ \beta^{0}_{k}\neq 0\big{\}}.$$

$$\underset{{\overline{S}}^1}{\underbrace{\{{\overline{ı}}_1\}}}\subset \underset{{\overline{S}}^2}{\underbrace{\{{\overline{ı}}_1,{\overline{ı}}_2\}}}\subset \cdots \subset \underset{{\overline{S}}^a}{\underbrace{\{{\overline{ı}}_1,{\overline{ı}}_2,\dots ,{\overline{ı}}_a\}}}\subset \cdots \subset \underset{{\overline{S}}^K}{\underbrace{\{{\overline{ı}}_1,{\overline{ı}}_2,\dots ,{\overline{ı}}_K\}}},$$

$$\hat{S}={\overline{S}}^{\hat{m}}.$$

$$\{0\}=:H_0\subset H_1\subset \cdots \subset H_a:=\mathrm{Span}(X_{{\overline{ı}}_1},\dots ,X_{{\overline{ı}}_a})\subset \cdots \subset H_K,$$

$$K\text{ is fixed and such that }1\le K<\min (n,r)\text{ where }r=\mathrm{rank}(X).$$

$$\mathcal{A}_{K}:=\big{\{}(i_{1},\ldots,i_{K})\in[2p]^{K}\ :\ j_{1},\ldots,j_{K}\text{ are pairwise different}\big{\}}.$$

$$\forall k\in[K],\ \forall j\notin{T^{k}}:=\{j_{1},\ldots,j_{k}\},\quad X_{j}^{\top}X_{T^{k}}\big{(}X_{T^{k}}^{\top}X_{T^{k}}\big{)}^{-1}s_{k}<1,$$

$$\forall S\subset[p]\ \text{s.t.}\ \#S\leq K,\quad\max_{j\in[p]\setminus S}\max_{||v||_{\infty}\leq 1}X_{j}^{\top}X_{S}\big{(}X_{S}^{\top}X_{S}\big{)}^{-1}v<1,$$

$$\text{For }K\text{ defined by }\text{eqref}e:K,({\hat{ı}}_1,\dots ,{\hat{ı}}_K)\text{ satifies }\text{eqref}hyp:IrrAlongThePath.$$

$$\forall(i_{1},\ldots,i_{k})\in[2p]^{k},\quad\theta_{j}(i_{1},\ldots,i_{k}):=\big{(}R_{j,i_{1}}\cdots R_{j,i_{k}}\big{)}M^{-1}_{i_{1},\ldots,i_{k}}(1,\ldots,1),$$

Full text

This paper makes three main contributions. First, we build testing procedures on variables entering the model along the LARS path in the general design case when the noise level can be unknown. These testing procedures are referred to as the Generalized $t$-Spacing tests (GtSt) and we prove that they have an exact non-asymptotic level (i.e., the Type I error is exactly controlled). This extends work of Tibshirani et al., (2016) where the spacing test works for consecutive knots and known variance. Second, we introduce a new exact multiple testing procedure after model selection in the general design case when the noise level may be unknown. We prove that this testing procedure has exact non-asymptotic level for general design and unknown noise level. Third, we prove exact control of the false discovery rate under orthogonal design assumption. Monte Carlo simulations and a real data experiment are provided to illustrate our results in this case. Of independent interest, we introduce an equivalent formulation of the LARS algorithm based on a recursive function. Multiple Testing; False Discovery Rate; High-Dimension; Selective Inference. \bUnif 2000 Math Subject Classification: Primary 62E15, 62F03, 60G15, 62H10, 62H15; secondary 60E05; 60G10; 62J05; 94A08

1 Introduction

In the past decades, statistical problems have become increasingly high-dimensional, i.e., they require estimation of more parameters than the number of available samples/observations. Some examples range from signal processing (Chen et al.,, 1998; Candès et al.,, 2006) to genomics (Rhee et al.,, 2006; Barber et al.,, 2015). Some successful techniques of estimation have been developed and a popular approach is based on optimizing a suitable regularized likelihood function. Most models of statistical parameters are well approximated by sparse vectors; and sparsity promoting regularizations, such as the $\ell_{1}$-norm, are now well recognized to tackle high-dimensional problems. Recent advances have focused on a deeper understanding of the law of the estimates of $\ell_{1}$-regularization procedures in high-dimension. One goal is to quantify the uncertainty of some linear statistic of the outcomes of sparse regression estimation. Such estimators are non-linear and non-explicit. They are defined as the minimum of some optimization program, or as the outcomes of some greedy method. Most of them estimate some set of relevant parameters, i.e., a small number of parameters that may explain the observation. This non-linear framework makes it impossible to characterize the distribution of the estimator. One possibility is to look at some conditional distribution of the estimator and this is the scope of the so-called selective inference, which produces an uncertainty quantification conditional on the set of indices of nonzero estimated parameters, referred to as the selection event. Selective inference aims at building some confidence intervals and some testing procedures on the estimates (see (van de Geer,, 2016, Chapter 6) and references therein), or controlling the false discovery rate, e.g., Barber et al., (2015) for instance.

One of the most popular regularized estimation procedure in high-dimensions is LASSO (Chen et al.,, 1998) and its asymptotic de-biased version referred to as the debiased LASSO. Controlling the FDR (resp., confidence intervals (CI)) built upon the debiased LASSO procedure has been studied in Javanmard et al., (2019) (resp., Javanmard and Montanari, (2014)) which provides an FDR with asymptotic control (resp., the CI with asymptotic control of the confidence level) for designs with some independent sub-Gaussian rows. The LASSO is based on $\ell_{1}$-norm regularization and one of its offsprings is the sorted-$\ell_{1}$ regularization, referred to as the SLOPE, which achieves minimax rate of prediction and estimation. Controlling the FDR for SLOPE with the Benjamini-Hochberg (BH) selection procedure has been achieved in Bogdan et al., (2015) for orthogonal designs.

Inference after model selection has been studied in several papers, such as Fithian et al., (2014); Taylor and Tibshirani, (2015) (resp., Tian et al., (2018)) for selective inference (resp., for a joint estimate of the noise level). These works give the non-asymptotic law of any linear statistics, i.e., any linear combination of the estimates of the parameters, conditional on the selection event. For the first time, this paper provides the non-asymptotic joint law of several linear statistics conditional on the selection event. These linear statistics are given by the knots of the LARS procedure. One may note that, conditional on the selection event, the law of three consecutive knots has been studied by Lockhart et al., (2014) who refer to it as the spacing test (ST) (Tibshirani et al.,, 2015). The article Azaïs et al., (2018) proved that the spacing test is unbiased and introduce a studentized version of this test. In the same direction, inference after model selection has been studied in several papers, such as Fithian et al., (2014); Taylor and Tibshirani, (2015) and respectively Tian et al., (2018) for selective inference and respectively a joint estimate of the noise level.

In the present paper, our test is based on the conditional joint law of three, not necessarily consecutive, knots. In this way, we extend the work from Tibshirani et al., (2016) where the spacing test works for consecutive knots. We refer to these new tests as the generalized spacing tests (GSt). Furthermore, the exact formulation of the spacing test of the pioneering work of Tibshirani et al., (2016) requires extra computations of the term denoted by $M^{+}$ in (Tibshirani et al.,, 2016, Lemma 5). They proved that the spacing test is asymptotically equivalent to the conservative spacing test. We remove this restriction and we prove that it suffices to check wether the so-called Irrepresentable Check Condition holds to get a non-asymptotic equivalence between the Spacing test and the conservative Spacing test. Finally, we theoretically prove that working with non-consecutive knots can render the testing procedure more powerful.

1.1 Joint law of LARS knots in Post-Selection Inference

In this paper, we consider linear models in high-dimensions where the number of observations $n$ may be less than the number of predictors $p$. We denote by $Y\in\mathbb{R}^{n}$ the response variable and we assume that

[TABLE]

where $\eta\sim\mathcal{N}_{n}(0,\sigma^{2}\mathrm{Id}_{n})$ is a Gaussian noise, the noise level $\sigma>0$ may be known or may have to be estimated (depending on the context), and $X\in\mathds{R}^{n\times p}$ has rank $r>0$. We consider the LARS and denote by $(\lambda_{k})_{k\geq 1}$ the sequence of knots and by $(\overline{\imath}_{k},\varepsilon_{k})_{k\geq 1}$ the sequence of variables $\overline{\imath}_{k}\in[p]$ and signs $\varepsilon_{k}\in\{\pm 1\}$ that enter the model along the LARS path. We encode by

[TABLE]

both the variables $\overline{\imath}_{k}\in[p]$ and the signs $\varepsilon_{k}\in\{\pm 1\}$, calling them the ‘signed variables’​. Section 5.1 recalls LARS (Algorithm 2) and present equivalent formulations in Algorithm 3 (using orthogonal projections) and Algorithm 4 (using a recursion). In particular, Algorithm 4 consists in three lines, applying the same function recursively, see Section 5.2. As far as we know, Algorithm 4 is new.

For a short moment, consider the simplest linear model, where one observes the target vector $\beta^{0}=(\beta^{0}_{1},\ldots,\beta^{0}_{p})\in\mathds{R}^{p}$, namely there is no noise and the design $X=\mathrm{Id}_{p}$ is the identity. In this case, LASSO and LARS give the same knots $\lambda_{1},\lambda_{2},\ldots$ and the estimate of the LASSO is the outcome of the proximal operator of the $\ell_{1}$-norm at point $\beta^{0}\in\mathds{R}^{p}$, see for instance (Tibshirani et al.,, 2015, Chapter 2). In this simple case, we deduce that the knots are

[TABLE]

where we have considered the reordering $\beta^{0}_{(1)}\geq\beta^{0}_{(2)}\geq\ldots$ of the entries of the target. Obviously, this is no longer true for general designs in high-dimensions with noise, but one may ask:

• [Q1]

What is the joint law of the LARS knots $\lambda_{1},\lambda_{2},\ldots$ and how do they relate to the target $\beta^{0}$?

We will answer $\bf{[Q1]}$ in high-dimensions under the assumption of Gaussian noise in Section 3.1.3 and Section 3.2. Working with the so-called ‘Irrepresentable Check’ Condition111See Section 2.3 for a definition and detailed comments on this assumption. (${\mathcal{A}}_{\mathrm{Irr.}}$), which can be efficiently checked in practice, we are able to provide the joint law of the LAR’s knots conditional on the so-called ‘selection event’ defined by

[TABLE]

This selection event states that the signed variable $\imath_{k}$ has been selected by the LARS algorithm at its $k^{\text{th}}$ step for $k=1,\ldots,K$. This is the cornerstone of the paper, showing that the conditional joint distribution of the LARS knots is a mixture of Gaussian order statistics, as presented in the next theorem.

Theorem 1.1** (Conditional Joint Law of the LARS knots).**

Let $(\lambda_{1},\ldots,\lambda_{K},\lambda_{K+1})$ be the first knots of the LARS and let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first variables entering along the LARS path. If $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$), then, conditional on the selection event $\{\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K},\lambda_{K+1}\}$, the vector $(\lambda_{1},\ldots,\lambda_{K})$ obeys a law with the following density (w.r.t. Lebesgue measure)

[TABLE]

where ${\mathrm{Q}}_{(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K},\lambda_{K+1})}$ is a normalizing constant, $\varphi_{m_{k},v_{k}^{2}}$ is the standard Gaussian density with mean $m_{k}$ and variance $v_{k}^{2}:=\sigma^{2}\rho_{k}^{2}$, are explicitly given by (3.12) and (3.13).

The proof of this theorem is given in Section 3.2.1. Now, let us describe the dependency between $(m_{k},\rho_{k}^{2})$ and $\overline{\mu}^{0}:=X^{\top}X\beta^{0}$. For a design matrix $X$ with columns $(X_{j})_{j=1}^{p}$, we denote by222Recall that the selected variables $\widehat{\imath}_{k}\in[2p]$ are decomposed into $\widehat{\imath}_{k}:=\overline{\imath}_{k}+p\,(\frac{1-\varepsilon_{k}}{2})$.

[TABLE]

By (3.12) and (3.13), one has

[TABLE]

where $c_{k},d_{k}>0$ are constants that depends only on $X_{\overline{\imath}_{1}},\ldots,X_{\overline{\imath}_{k-1}}$, $P_{k-1}^{\perp}$ denotes the orthogonal projection onto the orthogonal of $H_{k-1}$, $\varepsilon_{k}$ is the sign of the $k^{\text{th}}$ variable entering the LARS path, and $\measuredangle(X_{\overline{\imath}_{k}},H_{k-1})$ is the angle between $X_{\overline{\imath}_{k}}$ and $H_{k-1}$.

1.2 The Generalized $t$-Spacing test (GtSt)

This paper introduces a class of exact tests built from $\ell_{1}$-minimization regression in high-dimensions. More precisely, we design a testing procedure for a null hypothesis of the form

[TABLE]

where $H_{a_{0}}:=\mathrm{Span}(X_{\overline{\imath}_{1}},\ldots,X_{\overline{\imath}_{a_{0}}})$. Note that the null $\mathds{H}_{0}$ is equivalent to the hypothesis that all the true positives (i.e., the support of $\beta^{0}$) are among the first $a_{0}$ variables selected by LARS, namely $\{\overline{\imath}_{1},\ldots,\overline{\imath}_{a_{0}}\}$. Following the original idea of Lockhart et al., (2014), we study testing procedures of $\mathds{H}_{0}$ based on the knots of the LARS path. Note that, conditional on the selection event, the law of three consecutive knots has been studied by Lockhart et al., (2014), where it was referred to as the spacing test (ST) (Tibshirani et al.,, 2015). The article Azaïs et al., (2018) proved that the spacing test is unbiased, and introduced a Studentized version of this test. In the same direction, inference after model selection has been studied in several papers, such as Fithian et al., (2014); Taylor and Tibshirani, (2015) and respectively Tian et al., (2018) for selective inference and respectively a joint estimate of the noise level. This raises the following questions.

• [Q2]

Can we provide exact testing procedures based on knots that are not consecutive?

• [Q3]

What is the most powerful test among these spacing tests?

• [Q4]

Can we provide exact testing procedures when the noise level is not known?

$\bullet\ \mathbf{{Contribution\ }(i)}:$ First, our test is based on the conditional joint law of three, not necessarily consecutive, knots $a_{0}\leq a<b<c\leq K+1$. In this way, we extend the work of Tibshirani et al., (2016), where the spacing test works for consecutive knots. We present this framework in Section 3.2 and we refer to these new tests are Generalized Spacing tests (GSt). At the first reading of the next theorem, one can set $\widehat{m}=a_{0}$ for a fixed value $0\leq a_{0}\leq a$. The selection procedure, defining $\widehat{m}$, will be presented in Section 2.5 with the notion of an ‘admissible procedure’ ($\bf\mathcal{A}_{\mathrm{\,Stop}}$).

Theorem 1.2**.**

Let $a,b,$ and $c$ be such that $0\leq a<b<c\leq K+1$. Let $(\lambda_{1},\ldots,\lambda_{K},\lambda_{K+1})$ be the first knots and let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first variables entering along the LARS path. If $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$) and $\widehat{m}$ is chosen according to a procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), then under the null hypothesis

[TABLE]

and conditional on the selection event \big{\{}\,\widehat{m}\leq a\big{\}}, it follows that

[TABLE]

namely, it is uniformly distributed over $(0,1)$.

The proof of Theorem 1.2 is presented in Section 3.2.3. The construction of the $p$-values $\widehat{\alpha}_{abc}$ and of $\mathds{F}_{abc}$ is given in (3.19) and Section 3.2.2 respectively. We consider the following Generalized Spacing test procedures (GSt):

[TABLE]

that rejects if the $p$-value $\widehat{\alpha}_{abc}$ is less than the level $\alpha$ of the test. One may remark that

the $p$-value $\widehat{\alpha}_{abc}$ detects abnormally large values of $\lambda_{b}$ relatively to the interval $(\lambda_{a},\lambda_{c})$.

When the noise variance is unknown, we introduce the Generalized $t$-Spacing tests (GtSt) whose theoretical guarantees are given in the next theorem. The estimator of the variance $\widehat{\sigma}^{2}$ is given in Section 2.4.

Theorem 1.3**.**

Let $a,b,$ and $c$ be such that $0\leq a<b<c\leq K+1$. Let $(\lambda_{1},\ldots,\lambda_{K+1})$ be the first knots and let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first variables entering along the LARS path. If $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$) and $\widehat{m}$ is chosen according to a procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$) then under the null hypothesis

[TABLE]

and conditional on the selection event \big{\{}\widehat{m}\leq a\big{\}}, it follows that

[TABLE]

where $\Lambda_{k}:=\lambda_{k}/\widehat{\sigma}$.

The proof of Theorem 1.3 is presented in Section 3.4.1. The construction of the $p$-values $\widehat{\beta}_{abc}$, of $\widetilde{\mathds{F}}_{abc}$ and of the estimation of the noise $\widehat{\sigma}$ is given in (3.24), Section 3.23 and Section 2.4 respectively. One may remark that

the $p$-value $\widehat{\beta}_{abc}$ detects abnormally large values of $\Lambda_{b}$ relatively to the interval $(\Lambda_{a},\Lambda_{c})$.

$\bullet\ \mathbf{{Contribution\ }(ii)}:$ Working with three consecutive knots, we recover the spacing test of Tibshirani et al., (2016) and even in this framework, the present paper improves the current state of knowledge. We prove that:

* $\circ$ Under (${\mathcal{A}}_{\mathrm{Irr.}}$), the Spacing test procedure defined in (Tibshirani et al.,, 2016, Theorem 1) is exact, and is equal to the so-called ‘conservative’ spacing test defined in (Tibshirani et al.,, 2016, Theorem 2).*

The exact formulation of the spacing test of the pioneering work of Tibshirani et al., (2016) requires extra computations of the term denoted by $M^{+}$ in (Tibshirani et al.,, 2016, Lemma 5). They proved that the spacing test is asymptotically equivalent to the conservative spacing test. We remove this restriction and we prove that it suffices to check wether the Irrepresentable Check Condition (${\mathcal{A}}_{\mathrm{Irr.}}$) holds to get a non-asymptotic equivalence between the spacing test and the conservative spacing test.

$\bullet\ \mathbf{{Contribution\ }(iii)}:$ We theoretically prove that working with non-consecutive knots can allow obtaining higher power for the testing procedure.

Theorem 1.4**.**

Assume that the design $X$ is orthogonal, namely $X^{\top}X=\mathrm{Id}_{p}$. Let $a_{0}$ be an integer such that $0\leq a_{0}\leq K-1$. If $\widehat{m}$ is chosen according to a procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), then under the null hypothesis

[TABLE]

and conditional on the selection event \big{\{}\widehat{m}=a_{0}\big{\}}, it follows that the test $\mathcal{S}_{a_{0},a_{0}+1,K+1}$ is uniformly more powerful than any of the tests $\mathcal{S}_{a,b,c}$ for $a_{0}\leq a<b<c\leq K+1$.

The proof of this result is given in Appendix 7.5. This shows that the most powerful test among the set of tests $(\mathcal{S}_{a,b,c})_{a_{0}\leq a<b<c\leq K+1}$ is given by

the GSt test $\mathcal{S}_{a_{0},a_{0}+1,K+1}$ with the smallest $a$ and the largest $c$.

More precisely, in the proof of Theorem 1.4, it is shown that

[TABLE]

for orthogonal designs, where $\preccurlyeq$ denotes stochastic ordering.

1.3 A new exact testing procedure on false negatives (FN) after support selection

One specific task is to estimate the support of the target sparse vector, namely identify the true positives in the context of a multiple testing procedure. In particular, one may take the support of the LASSO (or SLOPE) solution as an estimate of the support of the solution. This strategy has been intensively studied in the literature, one may consider Wainwright, (2009); Bogdan et al., (2015); van de Geer, (2016); Bellec et al., (2018) and references therein. Support selection has been studied under the so-called ‘Irrepresentable Condition’ (IC), as presented for instance in (van de Geer,, 2016, Page 53) and (Bühlmann and van de Geer,, 2011, Sec. 7.5.1) and also referred to as the ‘Mutual Incoherence Condition’ (Wainwright,, 2009). Under the so-called ‘Beta-Min Condition’, one may prove (Bühlmann and van de Geer,, 2011; van de Geer,, 2016) that the LASSO asymptotically returns the true support. Following this line of thought, one may ask:

• [Q5]

Can we provide a false negative testing procedure with a controlled Type I error?

In this article, we build an exact non-asymptotic multiple test for false non-negatives, see Sections 3.3 and 3.4. The control of the false negatives after model selection in the case of an unknown noise level is given in Section 3.4 and the procedure is introduced in Algorithm 1. We assume ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), which assumes that the model has been selected using an ‘admissible’ procedure, which basically means that the decision to select a model of size $a$ only depends on the orthogonal projection of the observation $Y$ onto $H_{a}=\mathrm{Span}(X_{\overline{\imath}_{1}},\ldots,X_{\overline{\imath}_{a}})$. Assuming further that (${\mathcal{A}}_{\mathrm{Irr.}}$) holds, we provide an exact testing method for false negatives. In order to reach high power, the test statistic is the $p$-value of three non-consecutive knots of the LARS path. To compute this $p$-value, one needs to marginalize the joint law of the knots, leading to a numerical integration whose complexity grows exponentially with the space between the indices of the knots. We propose to use QMC techniques to compute the statistic, see Appendix 8.

1.4 False Discovery Rate control for LARS

Simultaneous controls of confidence intervals independently of the selection procedure have been studied under the concept of post-selection constants as introduced in Berk et al., (2013) and studied for instance in Bachoc et al., (2018). Asymptotic confidence intervals can be build using the de-sparsified LASSO, the reader may refer to (van de Geer,, 2016, Chapter 5) and references therein. We also point a recent study (Javanmard et al.,, 2019) of the FDR control as the sample size tends to infinity using de-biased LASSO, which has been implemented in Section 4.2. Asymptotic FDR control has been studied in Barber et al., (2015) and references therein, which has been implemented in Section 4.2. Let us point recent control of the Joint family-wise Error Rate as in Blanchard et al., (2017) and references therein. Following these lines of work, one may ask:

• [Q6]

Can we provide multiple Spacing Tests with a controlled False Discovery Rate (FDR)?

To the best of our knowledge, this paper is the first to study the joint law and an exact control of multiple spacing tests of LARS knots in a non-asymptotic frame, see Sections 3.2 and 3.6. We investigate the consecutive spacings of the knots of the LARS as test statistics and we prove an exact FDR control using a Benjamini–Hochberg procedure (Benjamini and Hochberg,, 1995) in the orthogonal design case, see Theorem 3.23 and Section 3.6. Our proof (see Appendix 7.8) is based on the Weak Positive Regression Dependency (WPRDS), the reader may consult Blanchard et al., (2008) or the survey Roquain, (2011), and Knothe-Rosemblatt transport, see for instance (Santambrogio,, 2015, Sec.2.3, Page 67) or (Villani,, 2008, Page 20), which is based on conditional quantile transforms.

1.5 Additional related works on high-dimensional statistics

Parsimonious models have become ubiquitous tools to tackle high-dimensional representations with a small budget of observations. Successful applications may be found in signal processing (see for instance the pioneering works of Chen et al., (1998); Candès et al., (2006) and references therein) and biology (see for instance Barber et al., (2015) or (Bühlmann and van de Geer,, 2011, Chapter 1.4) and references therein). These applications have shown that there are interesting almost sparse representations in some well chosen basis. Nowadays, in many practical situations, this sparsity assumption is recognized as reasonable.

These important successes have put a focus on High-Dimensional Statistics and Compressed Sensing in the past decades, which may be due to the deployment of tractable algorithms with strong theoretical guarantees. Among the large panoply of methods, one may consider $\ell_{1}$-regularization, which benefits from a remarkable tractability, empirical performance, and theoretical guarantees. Nowadays, sparse regression techniques based on $\ell_{1}$-regularization are a common and powerful tool in high-dimensional settings. Popular estimators, among which one may point to the LASSO (Tibshirani,, 1996) and SLOPE (Bogdan et al.,, 2015), are known to achieve a minimax rate of prediction and to satisfy the sharp oracle inequalities under conditions on the design, such as Restricted Eigenvalue (Bickel et al.,, 2009; Bellec et al.,, 2018) or Compatibility (Bühlmann and van de Geer,, 2011; van de Geer,, 2016). The sharp oracle inequalities show that the estimation errors, in $\ell_{1}$ and $\ell_{2}$ norm, of these estimators are optimal, see for instance (van de Geer,, 2016, Chapter 2.7).

Variable selection has also been investigated, and it has been proven, see for instance (Bühlmann and van de Geer,, 2011, Theorem 7.5), that the LASSO selects the true variables (i.e., there are no false negatives) under the Compatibility condition and the so-called beta-min condition (which assumes that the true parameters are large enough with respect to some threshold that scales linearly with the regularization parameter $\lambda$ of LASSO). Under a stronger assumption, referred to as the ‘irrepresentable condition’, one can prove, see for instance (Bühlmann and van de Geer,, 2011, Theorem 7.1), that the $\ell_{\infty}$-estimation error scales linearly with the regularization parameter $\lambda$ of LASSO. As the regularization parameter $\lambda$ tends to zero, when the number of observation goes to infinity and under some assumption on the noise, these results show that LASSO produces a consistent selection of the variables (it asymptotically finds the true support with no errors).

1.6 Outline of the paper

1.6.1 Detailed outline

Section 2 introduces the notation (see also Section 1.6.3 for a summary), assumptions (${\mathcal{A}}_{\mathrm{Irr.}}$) and ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), and variance estimate $\widehat{\sigma}^{2}$. The variance estimate is a key step in our testing procedures: we introduce new variance estimate with properties useful for deriving exact and non-asymptotic post-selection laws, see Section 2.4.

The main assumption is based on the Irrepresentable Check condition (${\mathcal{A}}_{\mathrm{Irr.}}$), which can be checked in practice, see Section 2.3. Under (${\mathcal{A}}_{\mathrm{Irr.}}$), we obtain a new characterization of the selection event in Proposition 3.3 of Section 3.1.

Section 3 gives the main results: Section 3.2 describes the joint distribution of the LARS knots as a mixture of Gaussian order statistics and the GST and GtST tests. The power in the orthogonal case is considered in Section 3.5. The control of the false negatives in a post selection inference with estimation of the variance is presented in Section 3.4 (when the variance is known, this procedure is studied in Section 3.3). A procedure to control the FDR in the orthogonal case is presented in Section 3.6.

Illustrations of our method, both on simulated data and on real data, are presented in Section 4. A Zenodo repository of the code used in all our experiments can be found at De Castro, (2021).

1.6.2 Dependency diagram

The outline can be depicted by the following dependency diagram:

Conditional hypotheses $H_{a}$ in Sec. 2.2Proposition 3.3, selection event (Sec. 3.1)Irrepresentable Check (${\mathcal{A}}_{\mathrm{Irr.}}$) in Sec. 2.3Variance estimate $\widehat{\sigma}^{2}$ in Sec. 2.4Admissible procedures ($\bf\mathcal{A}_{\mathrm{\,Stop}}$) in Sec. 2.5Orthogonal design, $X^{\top}X=\mathrm{Id}_{p}$Assumptions & toolsTheorem 1.1, conditional joint law (Sec. 3.2)Theorem 1.4, power/Type II error (Sec. 3.5)Theorem 1.2, GSt & FN testing (Sec. 3.3)Theorem 1.3, GtSt & FN testing (Sec. 3.4)Theorem 3.23, FDR control (Sec. 3.6)Results

1.6.3 Notation and commands

[TABLE]

[TABLE]

2 Assumptions, Variance Estimation and Admissible Procedures

2.1 Signed variables of LARS

We give some notation that will be useful. We denote by $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k})\in[2p]^{k}$ the ‘signed’ variables that enter the model along the LARS path with the convention that

[TABLE]

so that $\widehat{\imath}_{k}\in[2p]$ is a useful way of encoding both the variable $\overline{\imath}_{k}\in[p]$ and its sign $\varepsilon_{k}=\pm 1$ as used in Algorithm 4. We denote by $\overline{Z}:=X^{\top}Y$ the correlation vector such that $\overline{Z}_{k}$ is the scalar product between the $k^{\text{th}}$ predictor and the response variable, and we denote by $\sigma^{2}\overline{R}$ its variance-covariance matrix. For the sake of presentation, we may consider the $2p$-vector

[TABLE]

whose mean is given by

[TABLE]

and whose variance-covariance matrix is $\sigma^{2}R$ with

[TABLE]

We also denote by

1. $\circ$

$\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k}$, the first $k$ signed variables entering the LARS, 2. $\circ$

$i_{1},\ldots,i_{k}$, a generic value of the sequence above, 3. $\circ$

\overline{\imath}_{1},\ldots,\overline{\imath}_{k}\leavevmode\nobreak\, the first $k$ variables entering the LARS, 4. $\circ$

$j_{1},\ldots,j_{k}$, a generic value of the sequence above, 5. $\circ$

$\varepsilon_{1},\ldots,\varepsilon_{k}$, the first $k$ signs of the coefficients of the variables entering in the LARS, 6. $\circ$

$s_{1},\ldots,s_{k}$, a generic value of the sequence above.

The quantities above are related by (2.1) and

[TABLE]

2.2 Models, conditional hypotheses, and the notation $K$

We are interested in selecting the true support $S^{0}$ of $\beta^{0}$, where the support is defined by

[TABLE]

To estimate this support, we will consider the models that appear along the LARS path: the selected model $\widehat{S}$ would be chosen from the family of nested models

[TABLE]

where $K$ denotes the maximal model size. We denote by $\widehat{m}$ the size of the selected model $\widehat{S}$, and then

[TABLE]

Respectively, denote

[TABLE]

the corresponding family of nested subspaces of $\mathds{R}^{n}$. Once the model has been selected, we will construct tests based on the $K+1$ first knots of the LARS.

Remark 2.1**.**

The testing procedures under consideration are not standard since the $(H_{a})_{a=1}^{K}$ are random subspaces. We are interested in the framework of selective testing, namely, testing procedures conditional on the selection event $\{\widehat{m}=a,\widehat{\imath}_{1}=i_{1},\ldots,\widehat{\imath}_{K}=i_{K}\}$, for some fixed $a\in[K-1]$. Conditional on the event, note that $H_{a}$ is fixed. By convention, we may consider the case $a=0$, that is, testing the global null hypothesis.

Throughout this paper, we assume that

[TABLE]

In practice, $K$ can be considerably much smaller than $n$. Our analysis is conditional on $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ and in this spirit it can be referred to as a ‘Post-Section’ procedure, see e.g. Taylor and Tibshirani, (2015); Tibshirani et al., (2016, 2015).

2.3 Irrepresentable Check on the Active sets

We define the set of Active Sets $\mathcal{A}_{K}$ as all the sequences $i_{1},\ldots,i_{K}$ of signed variables such that $j_{1},\ldots,j_{K}$ are pairwise different, where the $j$’s are defined by (2.5), namely

[TABLE]

Sometimes it would be useful to consider $\mathcal{A}_{K+1}$, the set of active sets of size $K+1$. We introduce the notion of ‘Irrepresentable Check’, which is the only assumption on the design and the selected active set in most of our results.

Definition 2.2** (Irrepresentable Check).**

An active set $(i_{1},\ldots,i_{K})\in\mathcal{A}_{K}$ is said to satisfy the Irrepresentable Check condition if

[TABLE]

where $j_{k}$ and $s_{k}$ are defined from $i_{k}$ using (2.5). By a slight abuse of notation we will denote by (${\mathcal{A}}_{\mathrm{Irr.}}$) the set of sequences $(i_{1},\ldots,i_{K})$ that satisfy this property.

In our procedures and theoretical results, we will limit our attention to sequences $\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K}\,$ chosen by LARS that satisfy (${\mathcal{A}}_{\mathrm{Irr.}}$). A particular case is when the property is true for all possible active sets. This is equivalent to the Irrepresentable Condition that we will now recall.

Definition 2.3** (Irrepresentable Condition of order $K$).**

The design matrix $X$ satisfies the Irrepresentable Condition of order $K$ if and only if

[TABLE]

where $X_{j}$ denotes the $j^{\text{th}}$ column of $X$ and $X_{S}$ the sub-matrix of $X$ obtained by keeping the columns indexed by $S$.

Remark 2.4**.**

Note that the Irrepresentable Condition is a standard condition, as presented for instance, in (van de Geer,, 2016, Page 53) and (Bühlmann and van de Geer,, 2011, Sec. 7.5.1), and is also referred to as the Mutual Incoherence Condition (Wainwright,, 2009).

Remark 2.5**.**

This condition has been intensively studied in the literature and it is now well established that some random matrix models satisfy it with high probability. For instance, one may refer to Wainwright, (2009), where it is shown that a design matrix $X\in\mathds{R}^{n\times p}$ whose rows are drawn independently with respect to a centered Gaussian distribution with variance-covariance matrix satisfying (${\mathrm{Irrep.}}$) (for instance the identity matrix) satisfies (${\mathrm{Irrep.}}$) with high probability when $n\gtrsim K\log(p-K)$, where $\gtrsim$ denotes an inequality up to some multiplicative constant.

In practice, the Irrepresentable Condition (${\mathrm{Irrep.}}$) is a strong requirement on the design $X$ and, in addition, this condition cannot be checked in polynomial time. One important feature of our results is that we do not require the Irrepresentable Condition (${\mathrm{Irrep.}}$) but only the weaker requirement of Irrepresentable Check (${\mathcal{A}}_{\mathrm{Irr.}}$) on the selected active set. Namely, we will assume that

[TABLE]

Given $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$, note that this condition can be checked in polynomial time.

Example 2.6**.**

Taking the (signed) variables entering the model with LARS in an iid Gaussian design and as response variable a centered Gaussian vector with iid entries from $10,000$ Monte Carlo repetitions, Figure 1 illustrates the law of the maximal order $K_{\mathrm{max}}$ for which the Irrepresentable Check condition holds. For example, we found that for $p=1,000$ and $n=100$ $($with ratio $n/p=0.1)$ resp. $n=500$ $($with ratio $n/p=0.5)$, the Irrepresentable Check condition (${\mathcal{A}}_{\mathrm{Irr.}}$) of order $K_{\mathrm{max}}$ holds when $K_{\mathrm{max}}$ is about $K_{\mathrm{max}}\simeq 0.16\times n=16$ respectively $K_{\mathrm{max}}\simeq 0.12\times n=60$, see Figure 1.

2.3.1 Irrepresentable Check: An equivalent formulation

Now, we can define

[TABLE]

where $(1,\ldots,1)$ is the column vector of size $k$ whose entries are equal to one; $\sigma^{2}M_{i_{1},\ldots,i_{k}}$ is the variance-covariance matrix of the vector $(Z_{i_{1}},\cdots,Z_{i_{k}})$ and $(R_{j,i_{1}}\cdots R_{j,i_{k}})$ is a row vector of size $k$. Note that $M_{i_{1},\ldots,i_{k}}$ is the submatrix of $R$ obtained by keeping the columns and the rows indexed by $\{i_{1},\ldots,i_{k}\}$, namely

[TABLE]

Remark that

[TABLE]

when $\mathds{E}Z=0$. Then Proposition 2.7 shows that the Irrepresentable Condition (${\mathrm{Irrep.}}$) of order $K$ is equivalently given by

[TABLE]

where $\theta_{j}(i_{1},\ldots,i_{k})$ is given by (2.10).

Proposition 2.7**.**

Let $X$ and $R$ be defined by (2.4). Then, the following assumptions are equivalent:

• •

the design matrix $X$ satisfies (${\mathrm{Irrep.}}$) of order $K$,

• •

the variance-covariance matrix $R$ satisfies (2.12) of order $K$.

Furthermore, they imply that for all $(i_{1},\ldots,i_{K})\in\mathcal{A}_{K}$ one has

[TABLE]

which is an equivalent formulation of $(i_{1},\ldots,i_{K})$ satisfying (${\mathcal{A}}_{\mathrm{Irr.}}$).

Proof 2.8**.**

Let $S=\{j_{1},\ldots,j_{k}\}\subset[p]$ and $j\in[2p]\setminus S$. Let $\overline{v}=(\overline{v}_{1},\ldots,\overline{v}_{k})\in\{-1,1\}^{k}$ and define $i_{\ell}=j_{\ell}+p(1-\overline{v}_{\ell})/2$ for $\ell\in[k]$. Note that

[TABLE]

Now, observe that

[TABLE]

showing the equivalence between the two assumptions.

Remark 2.9**.**

One may require that the design be ‘normalized’ so that $R_{i,i}=1$, namely its columns have unit Euclidean norm. Under this normalization, one can check that $R$ satisfies (${\mathrm{Irrep.}}$) of order $K=1$. Hence, up to some normalization, one can always assume (${\mathrm{Irrep.}}$) of order $K=1$.

Remark 2.10**.**

When computing the LARS path, one has to compute the values

[TABLE]

see for instance Algorithm 2 or Algorithm 3, where these values are given by $\theta$, as shown by Proposition 2.7. This implies that, in practice, along the LARS path, one witnesses the maximal order $K$ for which Irrepresentable Check (${\mathcal{A}}_{\mathrm{Irr.}}$) holds.

2.4 The estimator of the variance

In our analysis, we introduce an estimate of the variance $\widehat{\sigma}^{2}$ to perform post-selection inference when the noise level $\sigma$ is unknown. The degree of freedom to the estimation of the variance is $n-K$. Let us fix, for the moment, $j_{1},\ldots,j_{K}$, the indices that are the putative indices for the selected variables. Let $P^{\perp}_{K}$ be the orthogonal projection on the orthogonal to $H_{K}^{j_{1},\ldots,j_{K}}:=\mathrm{Span}(X_{j_{1}},\ldots,X_{j_{K}})$. We define

[TABLE]

By a slight abuse of notation, we can index the estimator above by the signed indexes $i_{1},\ldots,i_{K}$. Eventually, we set

[TABLE]

the estimates of the standard deviation $\sigma$.

2.5 Admissible Selection Procedures

Note that choosing a model $\widehat{S}$ is equivalent to choosing a model size $\widehat{m}$ so that

[TABLE]

Our procedure is flexible on this point and allows any choice of $\widehat{m}$ as long as the following property ($\bf\mathcal{A}_{\mathrm{\,Stop}}$) is satisfied:

Stopping Rule: The estimated model size $\widehat{m}$ is a ‘stopping time’: $\widehat{m}\in[K-1]$ and, for all $a\in[K-1]$,

[TABLE]

In other words, the decision to select a model of size $\{\widehat{m}=a\}$ depends only on the first $a$ variables entering the LARS.

Remark 2.11**.**

We now give an example to show that ($\bf\mathcal{A}_{\mathrm{\,Stop}}$) implies some restriction. Suppose, for example, that we want to decide wether the target $\beta^{0}$ is two sparse or one sparse. A natural decision rule is to look at large values of the second knot $\lambda_{2}$, if $\text{\textquoteleft\textquoteleft}\,\lambda_{2}>\mathrm{(some\ threshold)}\text{''}$ choose $m=2$ otherwise choose $m=1$. This rule does not satisfy ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), since looking at $\lambda_{2}$ we can choose only sizes $m$ greater than or equal to 2.

Denote by $P_{k}(Y)$ (resp. $P^{\perp}_{k}(Y)$) the orthogonal projection of the observation $Y$ onto $H_{k}$ (resp. the orthogonal of $H_{k}$) for all $k\geq 1$ where $H_{k}$ are defined by (2.8). Given $h$ any measurable function,

[TABLE]

determines a class of selection procedures satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$). These procedures decide whether to stop at $\{\widehat{m}=a\}$ based on the information given by $P_{a}(Y)$. Once one has selected a model of size $\widehat{m}$, one may be willing to test if $\widehat{S}$ contains the true support $S^{0}$ by considering the null hypothesis

[TABLE]

namely there are no false negatives. Equivalently, one aims at testing the null hypothesis

[TABLE]

at an exact significance level $\alpha\in(0,1)$, where $(H_{a})_{a=0}^{K-1}$ is defined by (2.8).

3 Exact Controls using Least Angle Regression: Main Results

3.1 Key notion: the ‘frozen’ knots, their means and variances

3.1.1 Frozen knots

Given $K$ as defined in (2.9) and fixed $i_{1},\ldots,i_{K+1}\in[2p]$, one may define

[TABLE]

where

[TABLE]

and $\theta_{j}(i_{1},\ldots,i_{k})$ is given by (2.10). When $\mathds{E}Z=0$, one may remark that $\Pi_{i_{1},\ldots,i_{k}}(Z_{j})$ is the regression of $Z_{j}$ on the vector $(Z_{i_{1}},\cdots,Z_{i_{k}})$ whose variance-covariance matrix is $\sigma^{2}M_{i_{1},\ldots,i_{k}}$, namely

[TABLE]

From this point on, we introduce

[TABLE]

and we emphasize that

[TABLE]

as proven in Appendix 5.4 (Eq. (5.3)) and Proposition 3.3. We are now able to define the “frozen” values of the knots:

[TABLE]

They are the Gaussian random variables that coincide with $\lambda_{1},\lambda_{2},\cdots,\lambda_{K+1}$ when the random variables defined by the signed indices \widehat{\imath}_{1},\widehat{\imath}_{2},\ldots,\widehat{\imath}_{K+1}\leavevmode\nobreak\ take the particular values $i_{1},i_{2},\ldots,i_{K+1}$.

Remark 3.1**.**

An interesting feature of the LARS knots is that they have a simple expression in terms of the partition given by the identity

[TABLE]

As we have seen in (3.5),

[TABLE]

giving the definition of the frozen knots $\lambda^{f}_{k}$ above.

3.1.2 Mean and centering of the frozen knots

Now, write

[TABLE]

for the orthogonal projection of $y$ onto $\mathrm{Span}(X_{j_{1}},\ldots,X_{j_{k}})$. Recall that $P_{k}$ is the orthogonal projection onto $H_{k}$, which is a random subspace. Recall also that, conditional on the event $\{\widehat{\imath}_{1}=i_{1},\ldots,\widehat{\imath}_{k}=i_{k}\}$ the subspace $H_{k}$ is fixed. Note that

[TABLE]

where $i=j+p({1-s})/2$. The mean $m^{f}_{k}$ and standard deviation $\sigma\,\rho^{f}_{k}$ of $\lambda^{f}_{k}$ are important values defined for all $k\in[K]$:

[TABLE]

with

[TABLE]

and this definition is equivalent to (3.12) and (3.13), see Section 3.2. Recall that $\overline{\mu}^{0}=\overline{R}\beta^{0}$ as defined in (2.3) and note that

[TABLE]

which is true when the true support $S^{0}$ of $\beta^{0}$ is included in $\overline{S}^{k}$, defined by (2.6). This proves the next proposition.

Proposition 3.2**.**

For fixed $0\leq a\leq K-1$, conditional on the selection event $\{\widehat{\imath}_{1}=i_{1},\ldots,\widehat{\imath}_{K}=i_{K}\}$, the hypothesis

[TABLE]

implies that $m^{f}_{k}=0$ for all $a<k\leq K$, namely $(Z_{i_{a+1}}^{(i_{1},\ldots,i_{a})},\ldots,Z_{i_{K}}^{(i_{1},\ldots,i_{K-1})})$ is centered.

This proposition is important for defining the hypothesis under consideration, see also Remark 2.1.

3.1.3 A key result: The characterization of the selection event

Regarding the joint law of the frozen knots, one has the following important proposition whose proof can be found in Section 7.1.

Proposition 3.3**.**

Let $(i_{1},\ldots,i_{K},i_{K+1})\in\mathcal{A}_{K+1}$, that is, a fixed active set of size $K+1$.

• •

If $(i_{1},\ldots,i_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$) then

[TABLE]

for any $0\leq a<k\leq K$ with the convention $\lambda_{0}=\infty$.

• •

It holds that

[TABLE]

are mutually independent, for any $k\in[K]$.

Furthermore, if $X\beta^{0}\in H_{K}$ then

[TABLE]

• •

*If $(i_{1},\ldots,i_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$) then *

[TABLE]

Remark 3.4**.**

Is Proposition 3.3 a new polyhedral lemma?* The characterization of the selection event for the inference of a single testing statistic has been known as the ‘polyhedral lemma’ in the literature, see for instance (Tibshirani et al.,, 2015, Figure 6.9) and references therein. This result is the cornerstone of selective inference with sparse models. It is based on two ingredients: First, the selection event can be expressed as a polyhedra; Second, conditional on the selection event, any linear statistics is distributed according to a truncated Gaussian with independent truncation bounds.*

A first remark is that the polyhedral lemma is shown for one linear statistic and, as far as we known, there is no polyhedral lemma for multiple linear statistics. The interesting point is that our result (Proposition 3.3) can be seen as a polyhedral lemma for multiple linear statistics. Under (${\mathcal{A}}_{\mathrm{Irr.}}$), the selection event $\{\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K}\}$ corresponds to a polyhedra described by the $Z_{i_{k}}^{(i_{1},\cdots,i_{k-1})}$ variables in the third point of Proposition 3.3. Our main result shows that the joint law of these multiple linear statistics are the Gaussian distribution restricted to the polyhedra $\{\ell_{1}\geq\ldots\geq\ell_{K}\geq\lambda_{K+1}\}$, see Theorem 1.1.

Note that the selection event has to include $\lambda_{K+1}$. As discussed above, our polyhedral lemma (Theorem 1.1 and Proposition 3.3) shows that, conditional on the selection event, $\lambda_{1},\ldots,\lambda_{K}$ are distributed with respect to a Gaussian distribution restricted to the polyhedra $\{\ell_{1}\geq\ldots\geq\ell_{K}\geq\lambda_{K+1}\}$. If one does not include $\lambda_{K+1}$ in the selection event, then one has to integrate this latter conditional law with respect to the distribution of $\lambda_{K+1}$ which is not known.

Proposition 3.5**.**

Assume that the design $X$ is such that the Irrepresentable Condition (${\mathrm{Irrep.}}$) of order $K$ holds. Almost surely, one has

• •

Among all possible sets $(i_{1},\ldots,i_{K})\in\mathcal{A}_{K}$, there is one and only one such that

[TABLE]

• •

This set is the set selected by LARS, namely $\widehat{\imath}_{1}=i_{1},\ldots,\widehat{\imath}_{K}=i_{K}$,

• •

and, for all $(i_{1},\ldots,i_{K})\in\mathcal{A}_{K}$,

[TABLE]

Proof 3.6**.**

Note that (${\mathrm{Irrep.}}$) implies (${\mathcal{A}}_{\mathrm{Irr.}}$) by Proposition 2.7. Then apply the first point of Proposition 3.3 to conclude.

Finding the set $\{\widehat{\imath}_{1}=i_{1},\ldots,\widehat{\imath}_{K}=i_{K}\}$ selected by LARS may be related to a combinatorial search testing (3.11) all possible candidates $(i_{1},\ldots,i_{K})\in\mathcal{A}_{K}$. Under the Irrepresentable Condition, the support selected by LARS is given by (3.11), which can be seen as the extension of (1.2) introducing [Q1] in Section 1.1.

3.2 Main results: Joint law and construction of post-selection tests

We assume that $K$ is defined as in (2.9). Except in Section 3.4, $\sigma^{2}$ is assumed to be known. Let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first signed variables entering along the LARS path. In this section, we are interested in the joint law of the LARS knots $(\lambda_{1},\ldots,\lambda_{K})$ conditional on $\lambda_{K+1}$ and $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$. To determine this joint law, we need to make precise the centering parameters $m_{k}$, by (see also (3.7))

[TABLE]

the first standard deviation $\sigma\rho_{1}$ with $\rho_{1}:=\sqrt{R_{\widehat{\imath}_{1},\widehat{\imath}_{1}}}$, and the others $\sigma\rho_{k}$ by (see also (3.8))

[TABLE]

where

[TABLE]

is defined by (2.10) and $M_{\widehat{\imath}_{1},\ldots,\widehat{\imath}_{\ell-1}}$ is defined by (2.11).

3.2.1 Proof of Theorem 1.1

From the definition of the Gaussian random variable Z_{i_{k}}^{(i_{1},\ldots,i_{k-1})}\leavevmode\nobreak\ in (3.1) one can deduce that its mean $m_{k}$ is given by (3.12) and its standard deviation $v_{k}$ by (3.13), considering putative indices for the selected variables. By the second point of Proposition 3.3, we know that these variables are independent. We deduce that their joint density $(Z_{i_{1}},Z_{i_{2}}^{(i_{1})},\ldots,Z_{i_{K}}^{(i_{1},\ldots,i_{K-1})})$ is

[TABLE]

with respect to Lebesgue measure. For now on, we condition on $\mathcal{E}:=\{\widehat{\imath}_{1}=i_{1},\ldots,\widehat{\imath}_{K}=i_{K},\lambda_{K+1}\}$ and we assume that $(i_{1},\ldots,i_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$). By the first equality of the third point of Proposition 3.3 we known that \mathcal{E}=\big{\{}\lambda_{K+1}\leq Z_{i_{K}}^{(i_{1},\ldots,i_{K-1})}\leq\cdots\leq Z_{i_{1}}\big{\}}, and on the event $\mathcal{E}$,

[TABLE]

Conditional on $\mathcal{E}$, the joint density of $(Z_{i_{1}},Z_{i_{2}}^{(i_{1})},\ldots,Z_{i_{K}}^{(i_{1},\ldots,i_{K-1})})$ is proportional to

[TABLE]

with respect to Lebesgue measure, and by (3.14) it is the conditional density of the knots.

3.2.2 Construction of the Generalized Spacing test

A useful consequence of Theorem 1.1 is that one can explicitly describe the joint law of the LARS knots after having selected a support $\widehat{S}$ of size $\widehat{m}$ with any procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$). In the sequel, we write

[TABLE]

where $\lambda_{0}=\infty$ and $F_{0}=1$ by convention.

Proposition 3.7**.**

Let $a\in\mathds{N}$ be such that $0\leq a\leq K-1$. Let $\widehat{m}$ be a selection procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$). Under the conditions of Theorem 1.1, under the null hypothesis

[TABLE]

and conditional on the selection event \big{\{}\widehat{m}=a,F_{a},F_{K+1},\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K}\big{\}}, we have that $(F_{a+1},\ldots,F_{K})$ is uniformly distributed on

[TABLE]

where the $\mathcal{P}_{i,j}$ are described in (3.16).

A proof of this proposition can be found in Appendix 7.2.

Remark 3.8**.**

The previous statement is consistent with the case $a=0$ corresponding to the global null hypothesis $\mathds{H}_{0}\,:\,\text{`}X\beta^{0}=0\text{'}$ $($or equivalently $\mathds{E}Z=0)$. Therefore, if $Z$ is centered, then, conditional on $F_{K+1}$, one has that $(F_{1},\ldots,F_{K})$ is uniformly distributed on

[TABLE]

Remark 3.9**.**

In the orthogonal case, where $\overline{R}=\mathrm{Id}$, note that $\theta_{j}(i_{1},\ldots,i_{\ell})=0$ for all $\ell\geq 1$ and all $i_{1},\ldots,i_{\ell}\neq j$, $\rho_{j}=1$ and $\mathcal{P}_{i,j}(f)=f$. We recover that $\mathcal{D}_{1,K}$ is the set of order statistics

[TABLE]

In this case, the knots $\lambda_{i}$ are Gaussian order statistics $\lambda_{1}=Z_{\widehat{\imath}_{1}}\geq\lambda_{2}=Z_{\widehat{\imath}_{2}}\geq\ldots\geq\lambda_{K}=Z_{\widehat{\imath}_{K}}\geq\lambda_{K+1}$ for the vector $Z$.

From Theorem 1.1, we deduce several test statistics. To this end, we introduce some notation defining

[TABLE]

and also

[TABLE]

Remark 3.10**.**

On the numerical side, note that this quantity can be computed using Quasi Monte Carlo (QMC) methods as in (Genz and Bretz,, 2009, Chapter 5.1) or Appendix 8. The function $\mathds{F}_{abc}$ gives the CDF of $\lambda_{b}$ conditional on $\lambda_{a},\lambda_{c}$ and on some selection event, as shown in the next proposition.

Proposition 3.11**.**

Let $a,b,$ and $c$ be such that $0\leq a<b<c\leq K+1$. Let $(\lambda_{1},\ldots,\lambda_{K},\lambda_{K+1})$ be the first knots and let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first variables entering along the LARS path. If $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$) and $\widehat{m}$ is chosen according to a procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), then under the null hypothesis

[TABLE]

it holds that

[TABLE]

A proof of this proposition can be found in Appendix 7.3.

Remark 3.12**.**

Note that the deterministic choice $\widehat{m}=a$, for a fixed $a\in[K-1]$, is a procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$) and Proposition 3.11 holds. This shows that if $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$), then, under the null hypothesis $\mathds{H}_{0}\,:\ \text{\textquoteleft\textquoteleft}X\beta^{0}\in H_{a}\,\text{''}$,

[TABLE]

for any $0\leq a<b<c\leq K+1$.

3.2.3 Proof of Theorem 1.2

From Proposition 3.11 we know that under $\mathds{H}_{0}$ and conditional on the selection event \big{\{}\widehat{m}=a\big{\}}, Eq. (3.20) gives the conditional CDF of $\lambda_{b}$. As a consequence,

[TABLE]

Finally, considerations of the distribution under the alternative show that to obtain a $p$-value, we must consider the complement to 1 of the quantity above.

3.2.4 Monte Carlo simulations, Spacing tests, and Generalized Spacing tests

Theorem 1.2 is illustrated numerically in Figure 2. Note that we have a perfect fit with the uniform law: the conditional law of the LARS knots obtained theoretically is numerically validated333A reproducible experiment given in a Python notebook is available at https://github.com/ydecastro/lar_testing/blob/master/Law_LAR.ipynb. This test statistic generalizes previous test statistics that appeared in ‘Spacing Tests’, as presented in (Tibshirani et al.,, 2015, Chapter 5) for instance, and will be referred to as the Generalized Spacing test.

Remark 3.13**.**

If one takes $a=0$, $b=1$, and $c=2$ then

[TABLE]

Similarly, taking $b=a+1$ and $c=a+2$,

[TABLE]

which is the conservative spacing test, see (Tibshirani et al.,, 2016, Theorem 2).

3.3 Exact false negative testing after model selection

We return to the case of a general design. Given $\alpha\in(0,1)$ and using Theorem 1.2, one can consider the following exact testing procedure at level $\alpha$ on false negatives, see the pseudo-code in Algorithm 1. The theoretical guarantee of this algorithm is given by the next proposition. It shows that conditional on the event that ‘there are no false negatives’, namely $\text{`}X\beta^{0}\in H_{\widehat{m}}\,\text{'}$, the observed significance $\widehat{\alpha}$ obeys the uniform law and hence $\mathds{1}_{\{\widehat{\alpha}\leq\alpha\}}$ is a testing procedure with level exactly $\alpha$.

Corollary 3.14**.**

Let $(\lambda_{1},\ldots,\lambda_{K},\lambda_{K+1})$ be the first knots and let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first variables entering along the LARS path. If $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$) and $\widehat{m}$ is chosen according to a procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), then, conditional on the null hypothesis

[TABLE]

it holds that

[TABLE]

that is, it is uniformly distributed over $(0,1)$.

Proof 3.15**.**

By Theorem 1.2, the conditional law of $\widehat{\alpha}_{a(a+1)(K+1)}$ with respect to $\{\widehat{m}\leq a\}$ is the uniform distribution. Note that the conditional law (1.3) does not depend on $a,b=a+1,c=K+1$, hence this law is unconditional on $\widehat{m}$.

When the variance $\sigma^{2}$ is unknown, one can ‘Studentize’ this test, as presented in the next section. The reader may consult Section 3.4 for a definition and check that the quantities $\widehat{\beta}_{abc},\widetilde{\mathds{F}},\Lambda_{k}$ do not require $\sigma$ to be computed.

Corollary 3.16**.**

Let $(\lambda_{1},\ldots,\lambda_{K},\lambda_{K+1})$ be the first knots and let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first variables entering along the LARS path. If $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$) and $\widehat{m}$ is chosen according to a procedure satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), then, conditional on the null hypothesis

[TABLE]

it holds that

[TABLE]

that is, it is uniformly distributed over $(0,1)$.

Proof 3.17**.**

By Theorem 1.3, the conditional law of $\widehat{\beta}_{a(a+1)(K+1)}$ with respect to $\{\widehat{m}\leq a\}$ is the uniform distribution. Note that the conditional law (3.24) does not depend on $a,b=a+1,c=K+1$, hence this law is unconditional on $\widehat{m}$.

3.4 Exact Testing Procedure for False Negatives with Variance Estimation

From the results of Section 3.2, one can present a method to select a model and propose an exact test of false negatives in the case of a general design, when the variance is unknown. We introduce a new exact testing procedure that can be deployed when ($\bf\mathcal{A}_{\mathrm{\,Stop}}$) holds, namely an ‘admissible’ selection procedure is used to build $\widehat{S}$. We start by a preliminary result whose proof is in Appendix 7.4.

Proposition 3.18**.**

Let $(\lambda_{1},\ldots,\lambda_{K},\lambda_{K+1})$ be the first knots of LARS and let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first variables entering along the LARS path. If $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$), then

• •

conditional on $\{\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K},\lambda_{K+1}\}$, the random variables $(\lambda_{1},\ldots,\lambda_{K})$ and $\widehat{\sigma}$ are independent;

• •

conditional on $\{\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K}\}$ and under the null hypothesis $\mathds{H}_{0}\,:\ \text{`}X\beta^{0}\in H_{K}\text{'}$, the random variables $(\lambda_{K+1}/\widehat{\sigma})$ and $\widehat{\sigma}$ are independent;

• •

conditional on $\{\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K},\lambda_{K+1}\}$ and under the null hypothesis $\mathds{H}_{0}\,:\ \text{`}X\beta^{0}\in H_{K}\text{'}$, the distribution of $(\lambda_{1},\ldots,\lambda_{K})$ is given by Theorem 1.1, while the distribution of $\widehat{\sigma}/\sigma$ is the same as the random variable

[TABLE]

where $W:=(w_{1},\ldots,w_{n-K})$ is a ‘truncated’ standard Gaussian vector with the truncation given by

[TABLE]

where $U\in\mathds{R}^{n\times(n-K)}$ is any matrix such that $UU^{\top}=\mathrm{Id}_{n}-P^{(\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})}$, $\theta^{K}:=(\theta_{j}(\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K}))_{j\in[p]}$, and with the convention $0/0=0$.

Remark 3.19**.**

Under the null hypothesis $\mathds{H}_{0}\,:\ \text{`}X\beta^{0}\in H_{K}\text{'}$, the Gaussian vectors $P_{K}^{\perp}(Y)$ $($see (2.13)$)$, defining the variance estimate $\widehat{\sigma}^{2}$, is centered. This null hypothesis means that the true support is included in the set of the $K$ first indices chosen by LARS. One may choose $K$ large enough to guarantee this null hypothesis.

Recall that, up to some positive numerical constant, the probability density function of the multivariate $t$-distribution with $(n-K)$ degrees of freedom, mean $m=(m_{1},\ldots,m_{K})$ and variance-covariance matrix $\mathrm{Diag}(\rho_{1}^{2},\ldots,\rho_{K}^{2})$ is given by

[TABLE]

We have an analogue to Theorem 1.1 giving the joint law of

[TABLE]

where $\widehat{\sigma}$ is given by (2.13) has $n-K$ degrees of freedom, see Proposition 3.18.

Theorem 3.20** (Conditional Joint Law of the Studentized LARS knots).**

Let $(\lambda_{1},\ldots,\lambda_{K},\lambda_{K+1})$ be the first knots and let $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ be the first variables entering along the LARS path. If $(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$) then, under the null hypothesis

[TABLE]

and conditional on the selection event $\{\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K},\Lambda_{K+1}\}$, the vector $(\Lambda_{1},\ldots,\Lambda_{K})$ obeys a law with the density $($w.r.t. Lebesgue measure$)$

[TABLE]

at point $(t_{1},t_{2},\ldots,t_{K})$, where ${\mathrm{P}}_{(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{K},\Lambda_{K+1})}$ is a normalizing constant, $m_{k}$ and $\rho_{k}$ are as in (3.12) and (3.13).

Proof 3.21** (Proof of Theorem 3.20).**

Let us fix some values $i_{1},\ldots,i_{K}$. From the definition of the Gaussian random variable Z_{i_{k}}^{(i_{1},\ldots,i_{k-1})}\leavevmode\nobreak\ in (3.1), one can deduce that its mean $m_{k}$ is given by (3.12) and its standard deviation $v_{k}:=\sigma\rho_{k}$ by (3.13), considering putative indices for the selected variables. By the proof of Proposition 3.18, we know that these variables are independent of $\widehat{\sigma}^{i_{1},\ldots,i_{K}}$. We deduce that the vector $(Z_{i_{1}}/\widehat{\sigma}^{i_{1},\ldots,i_{K}},\ldots,Z_{i_{K}}^{(i_{1},\ldots,i_{K-1})}/\widehat{\sigma}^{i_{1},\ldots,i_{K}})$ has density a multivariate $t$-distribution with ${n-K}$ degrees of freedom, mean $m=(m_{1},\ldots,m_{K})$ and variance-covariance matrix $\mathrm{Diag}(\rho_{1},\ldots,\rho_{K})$. Furthermore, by (3.10) of Proposition 3.18, we know that this vector is independent of $({Z^{(i_{1},\ldots,i_{K})}_{j}}/{\widehat{\sigma}^{i_{1},\ldots,i_{K}}})_{j\neq i_{1},\ldots,i_{K}}$, and, in particular, independent of $\Lambda_{K+1}^{(i_{1},\ldots,i_{K})}:=\max_{j}\{{Z^{(i_{1},\ldots,i_{K})}_{j}}/{\widehat{\sigma}^{i_{1},\ldots,i_{K}}}\}$. Recall that, conditional on

[TABLE]

and assuming that $(i_{1},\ldots,i_{K})$ satisfies (${\mathcal{A}}_{\mathrm{Irr.}}$), Proposition 3.3 implies that

[TABLE]

Furthermore, on the event $\mathcal{E}$ we have

[TABLE]

Because of the independence above, this implies that the conditional distribution is the one claimed.

For $0\leq a<b\leq K+1$, we introduce

[TABLE]

with the convention $\widetilde{\mathcal{I}}_{ab}(s,t)=1$ when $b=a+1$: and also

[TABLE]

When $m_{a+1}=\cdots=m_{c-1}=0$, the function $\widetilde{\mathds{F}}_{abc}$ gives the CDF of $\Lambda_{b}$ conditional on $\Lambda_{a},\Lambda_{c}$ and on some selection event, as shown below in Theorem 1.3 and (3.25). For $0\leq a<b<c\leq K+1$, we introduce the $p$-value

[TABLE]

On the numerical side, note that this quantity can be computed using Quasi Monte Carlo (QMC) methods as in (Genz and Bretz,, 2009, Chapter 5.1).

3.4.1 Proof of Theorem 1.3

Fix $a$ such that $0\leq a\leq K-1$ and consider any selection procedure $\widehat{m}$ satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$). From Proposition 3.2, conditional on

[TABLE]

and under the null hypothesis $\mathds{H}_{0}\,:\ \text{`}X\beta^{0}\in H_{a}\,\text{'}$, we know that $m_{a+1}=\ldots=m_{K}=0$. From Theorem 3.20 we know that the density of $(\Lambda_{a+1},\Lambda_{a+2},\ldots,\Lambda_{K})$ conditional on $\mathcal{F}$ is given by

[TABLE]

From the definition of assumption ($\bf\mathcal{A}_{\mathrm{\,Stop}}$), and on the event $\mathcal{F}$, we know that the indicator $\mathds{1}_{\{\widehat{m}=a\}}$ is a measurable function of $\lambda_{1},\ldots,\lambda_{a-1}$, which are respectively equal to $Z_{i_{1}},\ldots,Z_{i_{a-1}}^{(i_{1},\ldots,i_{a-2})}$ on $\mathcal{F}$ by (3.5). By (3.10) of Proposition 3.3, we deduce that $\mathds{1}_{\{\widehat{m}=a\}}$ is independent of $(\lambda_{a+1}/\widehat{\sigma}^{i_{1},\ldots,i_{K}},\ldots,\lambda_{K}/\widehat{\sigma}^{i_{1},\ldots,i_{K}})$ and of $\Lambda_{K+1}:=\lambda_{K+1}/\widehat{\sigma}^{i_{1},\ldots,i_{K}}$ conditional on $\mathcal{F}$. We deduce that the conditional density above is also the conditional density on the event

[TABLE]

Now, a simple integration shows that

[TABLE]

As a consequence and under the same conditioning, one has

[TABLE]

Finally, considerations of the distribution under the alternative show that to obtain a $p$-value we must consider the complement to $1$ of the quantity above.

3.4.2 $t$-Spacing tests and Generalized $t$-Spacing tests

Consider the following testing procedures:

[TABLE]

that rejects if the $p$-value $\widehat{\beta}_{abc}$ is less than the level $\alpha$ of the test. This test statistic generalizes previous test statistics that appeared in $t$-Spacing Tests, as presented in Azaïs et al., (2018) for instance, and will be referred to as the Generalized $t$-Spacing test (GtSt).

Remark 3.22**.**

If one takes $a=0$, $b=1$ and $c=2$, then one gets

[TABLE]

Similarly, taking $b=a+1$ and $c=a+2$, one gets

[TABLE]

which is the $t$-spacing test as presented in Azaïs et al., (2018), where

[TABLE]

is, up to some positive numerical constant, the CDF of a centered $t$-Student distribution with variance $\rho_{k}^{2}$ and $n-K$ degrees of freedom.

3.5 Power studies

3.5.1 Power when the design is orthogonal

One may investigate the power of these tests at detecting false negatives, namely, the alternatives given by: there exists $k\in S^{0}$ such that $k\notin\{\overline{\imath}_{1},\ldots,\overline{\imath}_{a}\}$. In particular, what is the most powerful test among these latter (1.4) testing procedures? A comprehensive study for the case of orthogonal designs is given by Theorem 1.4.

3.5.2 Numerical studies on the power for the general design case

In the case of an orthogonal design, Theorem 1.4 shows that the test based on $\widehat{\alpha}_{a,a+1,K+1}$ is uniformly more powerful than tests based on $\widehat{\alpha}_{x,y,z}$ with $a\leq x<y<z\leq K+1$. Numerical experiments on the power of these tests are presented in Figure 3 and they witness the same phenomenon for Gaussian designs. It presents the CDF of the $p$-value $\widehat{\alpha}_{abc}$ under the null and under two $2$-sparse alternatives, one with low signal and one with $5$ times more signal. The numerical results show that all the tests are exact (leftmost panel) and the test $\mathcal{S}_{125}$ is the most powerful. A detailed presentation of this is given in Section 4.1.

More precisely, one has, as proved in the orthogonal case by Theorem 1.4 (and its proof), that

• •

$\widehat{\alpha}_{125}\preccurlyeq\widehat{\alpha}_{124}\preccurlyeq\widehat{\alpha}_{123}$;

• •

$\widehat{\alpha}_{125}\preccurlyeq\widehat{\alpha}_{135}\preccurlyeq\widehat{\alpha}_{235}\preccurlyeq\widehat{\alpha}_{234}$;

• •

$\widehat{\alpha}_{125}\preccurlyeq\widehat{\alpha}_{135}\preccurlyeq\widehat{\alpha}_{145}\preccurlyeq\widehat{\alpha}_{245}$.

where $\preccurlyeq$ denotes stochastic ordering. In the proof of Theorem 1.4, it was shown that

[TABLE]

for orthogonal designs.

3.6 Control of False Discovery Rate in the Orthogonal Design case

3.6.1 Presentation in the general case

For the sake of readability, we will assume, for the moment, that $\sigma$ is known. We understand that the law of test statistics are parametrized by the hypotheses $(m_{k})_{k\in[K]}$, where $m_{k}$ is given by (3.12).

We recall that we write $\overline{\mu}^{0}=X^{\top}X\beta^{0}$ and $\overline{\mu}^{0}_{i}$ for its $i$th coordinate. Assuming that the predictors are normalised, in the general case, this quantity is the sum of $\beta_{i}^{0}$ and a linear combination of the $\beta_{j}^{0}$’s whose predictors $X_{j}$ are highly correlated with the predictor $X_{i}$. Now, given the variables $\overline{\imath}_{1},\ldots,\overline{\imath}_{k}\in[p]$ and signs $\varepsilon_{1},\ldots,\varepsilon_{{k}}\in\{\pm 1\}^{k}$, we denote by $(\Pi_{\overline{\imath}_{1},\ldots,\overline{\imath}_{k-1}}^{\perp}(\overline{\mu}^{0}))_{\overline{\imath}_{k}}$ the orthogonal projection given by

[TABLE]

The tested null hypotheses are conditional on some sub-sequence of variables $(\overline{\imath}_{1},\ldots,\overline{\imath}_{K+1})\in[p]^{K+1}$ and signs $\varepsilon_{1},\ldots,\varepsilon_{{K+1}}\in\{\pm 1\}^{K+1}$ entering the model. The $p$-values under consideration are

[TABLE]

We write $I_{0}$ of the set I_{0}=\big{\{}k\in[K]\ :\ \mathds{H}_{0,k}\text{ is true}\big{\}}. Given a subset $\widehat{R}\subseteq[K]$ of hypotheses that we consider as rejected, we call false positive ($\mathrm{FP}$) and true positive ($\mathrm{TP}$) the quantities $\mathrm{FP}=\mathrm{card}(\,\widehat{R}\cap I_{0})$ and $\mathrm{TP}=\mathrm{card}(\,\widehat{R}\setminus I_{0})$. Denote by $\widehat{p}_{(1)}\leq\ldots\leq\widehat{p}_{(K)}$ the $p$-values ranked in a nondecreasing order. Let $\alpha\in(0,1)$ and consider the Benjamini–Hochberg procedure, see for instance Benjamini and Hochberg, (1995), defined by a rejection set $\widehat{R}\subseteq[K]$ such that $\widehat{R}=\emptyset$ when $\{k\in[K]\ :\ \widehat{p}_{(k)}\leq\alpha k/K\}=\emptyset$ and

[TABLE]

Recall the definition of the FDR as the mean of the False Discovery Proportion ($\mathrm{FDP}$), namely

[TABLE]

where the expectation is unconditional on the sequence of variables entering the model, while the hypotheses that are being tested are conditional on the sequence of variables entering the model. This FDR can be understood by invoking the following decomposition

[TABLE]

where \overline{\pi}_{(\imath_{1},\ldots,\imath_{K})}=\mathds{P}\big{\{}\overline{\imath}_{1}=\imath_{1},\ldots,\overline{\imath}_{K}=\imath_{K}\big{\}}.

3.6.2 Control of the FDR by the Benjamini–Hochberg procedure in the orthogonal design case

We now consider the case of an orthogonal design where $X^{\top}X=\mathrm{Id}_{p}$ and the set of $p$-values is given by (3.29). Note that $I_{0}$ is simply the set of null coordinates of $\beta^{0}$. Remark also that the Irrepresentable Condition (${\mathrm{Irrep.}}$) of order $p$ holds and so does Empirical Irrepresentable Check (${\mathcal{A}}_{\mathrm{Irr.}}$), see Proposition 2.7.

Theorem 3.23**.**

Assume that the design is orthogonal, i.e. $X^{\top}X=\mathrm{Id}_{p}$, and let $K\in[p]$. Let $(\overline{\imath}_{1},\ldots,\overline{\imath}_{K})$ be the first variables entering along the LARS path. Consider the $p$-values given by (3.29) and the set $\widehat{R}$ given by (3.30). Then

[TABLE]

and so $\mathrm{FDR}$ is bounded above by $\alpha$.

The proof of this result is given in Appendix 7.8. One interpretation of post-selection type may be given as follows: if one looks at all the experiments giving the same sequence of variables entering the model $\{\overline{\imath}_{1}=i_{1},\ldots,\overline{\imath}_{K}=i_{K}\}$ and if one considers the Benjamini–Hochberg procedure for the hypotheses described in Section 3.6.1, then the FDR is exactly controlled by $\alpha$.

4 Testing procedures: Numerical studies

4.1 Power in the non-orthogonal case

To study the power in the case of a non-orthogonal design, we built a Monte-Carlo experiment with:

• •

a model with $n=200$ observations and $p=300$ predictors,

• •

a random design matrix $X$ given by $300$ independent column vectors uniformly distributed on the Euclidean sphere $\mathds{S}^{199}$,

• •

and we ran $3,000$ Monte Carlo experiments.

• •

The results are presented in Figure 3.

The computation of the function $\mathds{F}_{abc}$ given by (3.19) requires multivariate integration tools. All our test statistics can be efficiently computed using Quasi Monte Carlo methods (QMC) for Multi-Variate Normal (MVN) and $t$ (MVT) distributions, see the book Genz and Bretz, (2009) for a comprehensive treatment of this topic or Appendix 8 for a short overview of the method we used. We compute spacings of length at most $4$, which implies that $c\leq 5$ when $a=1$ in our experimental framework.

A Python notebook and codes are given at https://github.com/ydecastro/lar_testing. The base function is

observed_significance_CBC(lars, sigma, start, end, middle)

in the file multiple_spacing_tests.py. It gives the $p$-value $\widehat{\alpha}_{(\mathtt{start})(\mathtt{middle})(\mathtt{end})}$ of the knots and indices given by lars and an estimate of (or the true) standard deviation given by sigma. We ran $3,000$ repetitions of this function to get the laws displayed in Figure 3. It presents the CDF of the $p$-value $\widehat{\alpha}_{abc}$ under the null and under two $2$-sparse alternatives, one with low signal and one with $5$ times more signal. The results show, in our particular case, that all the tests are exact and the test $\mathcal{S}_{125}$ is the most powerful, see Section 3.5.2 for further details.

4.2 A comparison of FDR control and power on simulated data

We take the experiments introduced in (Javanmard et al.,, 2019, Section 5). As in this reference, we consider a linear model with design $X$ with independent rows drawn with respect to $\mathcal{N}_{p}(0,\Sigma)$. The covariance $\Sigma\in\mathds{R}^{p\times p}$ is such that $\Sigma_{ij}=r^{|i-j|}$, for some parameter $r\in(0,1)$. We then normalize the columns of $X$ to have unit Euclidean norm. We draw a $k$-sparse vector $\beta^{0}\in\mathds{R}^{p}$ by choosing a support of size $k$ at random with values $\{\pm A\}$ uniformly at random, where $A>0$ denotes the absolute value of the amplitudes. The Gaussian noise term $\eta$ is drawn from $\mathcal{N}_{n}(0,\mathrm{Id}_{n})$.

We compare the performances of three procedures:

• •

[Knockoff] Knockoff filters for FDR control (Barber et al.,, 2015) and we use knockoff+ as implemented on https://web.stanford.edu/group/candes/knockoffs/;

• •

[FCD] False Discovery Control via Debiasing (Javanmard et al.,, 2019, Section 5) and we use the implementation of debiased lasso presented on the webpage https://web.stanford.edu/~montanar/sslasso/ with the theoretical value $\overline{\lambda}=2\sqrt{(2\log p)/n}$ for the regularizing parameter. When the sample size is larger than the number of predictors ($n\geq p$), the debiasing step in FCD is superfluous as the decorrelating matrix ($M$) can be the inverse of the sample covariance. So, in this case, we start with an unbiased estimator upfront (which is Ordinary Least Squares OLS). The FCD then becomes thresholding the test statistics $|T_{i}|$ obtained from OLS;

• •

[GtSt-BH] Generalized $t$-Spacing tests on successive entries of the LARS path combined with a Benjamini–Hochberg procedure Benjamini and Hochberg, (1995) based on the sequence of spacings $\widehat{\beta}_{012},\widehat{\beta}_{123},\ldots,\widehat{\beta}_{a(a+1)(a+2)},\ldots$ as described in Section 3.4 with nominal value $\alpha=0.1$;

we numerically investigate the effects of the level of sparsity, the magnitude of the signal, the correlation between the features, and the empirical power. In all simulations, we set the target level FDR to $\alpha=0.1$.

4.2.1 The effect of the amplitude of the signal

We chose $n=200$ (sample size), $p=100$ (predictors), $k=20$ (sparsity), $\eta=0.1$ (features correlation) and varied the amplitude within the set $A\in[25]$. We computed the FDR and power by averaging over $3,000$ realizations of the noise and generations of the coefficients of the vector $\beta^{0}$. The results are plotted in Figure 4. Recall that, in the case $n\geq p$, FCD is a thresholded OLS and it might be considered as the best test here. One may note that it presents the best features (low FDR and high power). GtSt controls the FDR below the nominal value $\alpha=0.1$ with a slightly lower power than FCD. Knockoff+ has controlled FDR and matches the power of FCD.

4.2.2 Effect of feature correlation

We test the effect of correlations between the features with $n=200$, $p=100$, $k=20$, and $A=10$. Recall that the rows of the design matrix $X$ are generated from an $\mathcal{N}_{p}(0,\Sigma)$ distribution, with $\Sigma_{ij}=\eta^{|i-j|}$, and then the columns of $X$ are normalized to have unit norm. We vary the parameter $\eta$ within the set $\{0.1,0.15,0.2,\dotsc,0.75,0.8,0.85\}$. For each value of $\eta$, we compute the FDR and power by averaging over $3,000$ realizations of the noise and design matrix $X$. The results are displayed in Figure 5. One may note that, in the case $n\geq p$, FCD is a thresholded OLS and might be considered as the best estimation here (with low FDR and high power). Knockoff+ has controlled FDR and matches the power of FCD for small feature correlations. GtSt controls the FDR below the nominal value $\alpha=0.1$ with a lower power than Knockoff+.

4.2.3 Effect of sparsity

We set $n=200$, $p=100$, $A=10$, and $\eta=0.1$, and varied the level of sparsity of the coefficients within the set $k\in(10,40)$. The power and the FDR are computed by averaging over $3,000$ trials of the noise and generations of the coefficients of the vector $\beta^{0}$. The results are displayed in Figure 6. Knockoff+ has controlled FDR and matches the power of FCD. GtSt controls the FDR below the nominal value $\alpha=0.1$ with a lower power than Knockoff+.

4.3 FDR on real data

A detailed presentation in a Python notebook is available at https://github.com/ydecastro/lar_testing/blob/master/multiple_spacing_tests.ipynb. We consider a data set about HIV drug resistance extracted from Barber et al., (2015) and Rhee et al., (2006). The experiment consists in identifying mutations of the genes of the HIV that are involved with drug resistance. The data set contains about $p=200$ and $n=700$ observations. Since some protocol was used to remove some genes and some individuals, the exact numbers depend on the considered drug.

The methods considered are [Knockoff], [FCD], [GtSt-BH], and:

• •

[Slope] Slope for FDR control, as presented in (Bogdan et al.,, 2015).

The comparison is displayed in Figure 7. It appears that GtSt-BH and FCD procedures are more conservative but they give a better control of the False Discovery Proportion (FDP). SLOPE and Knockoff are more powerful but their FDP is greater than the expected FDR $\alpha=0.2$ (in $6$ experiments out of $7$ for SLOPE, in $3$ out of $7$ for Knockoff1, and in $5$ out of $7$ for Knockoff2).

Acknowledgements

The authors would like to thank anonymous referees for their time, comments and useful remarks on preliminary versions of this paper. The authors are in debt to Quentin Duchemin for his valuable comments and remarks on this work.

Supplement to\bUnif Multiple Testing and Variable Selection along\bUnif the path of the Least Angle Regression

5 Representing the LARS knots

5.1 The equivalent formulations of the LARS algorithm

We present here three equivalent formulations of the LARS that are a consequence of the analysis provided in Appendices 5 and 6. One new formulation is given by Algorithm 4.

5.2 A new formulation of Least Angle Regression algorithm

The Least Angle Regression (LARS) algorithm has been introduced in the seminal article Efron et al., (2004). In the context of linear regression in high dimensions, the LARS algorithm can be used to identify a subset of potential covariates. The LARS outputs a piecewise affine solutions path, and the knots $\lambda_{1}\geq\lambda_{2}\geq\dots>0$ are the change points of the LARS path that are built by tracking the $\ell_{\infty}$ of the residual. At each knot, the LAR algorithm adds to the active set of variables the covariate the most correlated with the actual residual. In that way, the descent direction is always equiangular to all variables present in the current active set. . This sequence of knots is closely related to the sequence of knots of LASSO (Tibshirani,, 1996), as they differ by only one rule: “Only in the LASSO case, if a nonzero coefficient crosses zero before the next variable enters, drop it from the active set and recompute the current joint least-squares direction”, as mentioned in (Tibshirani et al.,, 2015, Page 120) or (Efron et al.,, 2004, Theorem 1) for instance.

Theorem 5.1**.**

Let $n,p$ be integers. Given a vector $Y\in\mathbb{R}^{n}$ and matrix $X\in\mathds{R}^{n\times p}$ of rank $r$ then Algorithm 2 (LARS standard formulation), Algorithm 3 (LARS projected formulation) and Algorithm 4 (LARS recursive formulation) output the same sequence $((\lambda_{k},\overline{\imath}_{k},\varepsilon_{k}))_{k=1}^{r}$ from the input given by $\overline{Z}=X^{\top}Y$ and $\overline{R}=X^{\top}X$, where $\lambda_{1}\geq\lambda_{2}\geq\lambda_{2}\ldots\geq 0$ are the knots, and $\overline{\imath}_{1},\overline{\imath}_{2},\ldots,\overline{\imath}_{r}$ are the variables that enter the model with signs $\varepsilon_{1},\varepsilon_{2},\ldots,\varepsilon_{r}$ $(\varepsilon_{k}=\pm 1)$.

The formulation of Algorithm 2 (LARS standard formulation), Algorithm 3 (LARS projected formulation) and the proof of Theorem 5.1 are given in Appendices 5 and 6.

5.3 Initialization: First Knot

The first step of the LARS algorithm (Step 1 in Algorithm 2) seeks the most correlated predictor with the observation. In our formulation, introduce the first residual $N^{(1)}:=Z$ and observe that $N^{(1)}:=(\overline{N}^{(1)},-\overline{N}^{(1)})$. We define the first knot $\lambda_{1}>0$ as

[TABLE]

One may see that this definition is consistent with $\lambda_{1}$ in Algorithm 2 and note that $\widehat{\imath}_{1}$ and $(\overline{\imath}_{1},\varepsilon_{1})$ are related as in (2.1).

The LARS algorithm is a forward algorithm that selects a new variable and maintains a residual at each step. We also define

[TABLE]

and one can check that $N^{(2)}(\lambda)=(\overline{N}^{(2)}(\lambda),-\overline{N}^{(2)}(\lambda))$ where $\overline{N}(\lambda)$ is defined in Algorithm 2. It is clear that the coordinate $\widehat{\imath}_{1}$ of $N^{(2)}(\lambda)$ is equal to $\lambda$. On the other hand $N^{(1)}=Z$ attains its maximum at the single point $\widehat{\imath}_{1}$. By continuity this last property is kept for $\lambda$ in a left neighborhood of $\lambda_{1}$. We search for the first value of $\lambda$ such that this property is not met, i.e. the largest value of $\lambda$ such that

[TABLE]

as in Step ${\bf 3}$ of Algorithm 2. We call this value $\lambda_{2}$ and one may check that this definition is consistent with $\lambda_{2}$ in Algorithm 2.

Now, we can be more explicit about the expression of $\lambda_{2}$. Indeed, we make the following discussion on the values of \theta_{j}(\,\widehat{\imath}_{1})\leavevmode\nobreak\.

• •

If $\theta_{j}(\,\widehat{\imath}_{1})\geq 1$ , since $N^{(1)}_{j}<N^{(1)}_{\widehat{\imath}_{1}}$ for $j\neq\widehat{\imath}_{1}$ there is no hope to achieve the equality between N^{(2)}_{j}(\lambda)\leavevmode\nobreak\ and $N^{(2)}_{\widehat{\imath}_{1}}(\lambda)=\lambda$ for $0<\lambda\leq\lambda_{1}$ in view of (5.1).

• •

Thus we limit our attention to the $j$’s such that $\theta_{j}(\,\widehat{\imath}_{1})<1$. We have equality $N_{j}^{(2)}(\lambda)=\lambda$ when

[TABLE]

So we can also define the second knot $\lambda_{2}$ of the LARS as

[TABLE]

where $\Pi_{i_{1}}(Z_{j}):=Z_{i_{1}}\theta_{j}(i_{1})$. Remark that $\Pi_{i_{1}}(Z_{j})=\mathds{E}(Z_{j}\ |\ Z_{i_{1}})$ is the regression of $Z_{j}$ on $Z_{i_{1}}$ when $\mathds{E}Z=0$.

5.4 Recursion: Next Knots

The loop $({\bf 2}\rightleftarrows{\bf 3})$ in Algorithm 2 builds iteratively the knots $\lambda_{1},\lambda_{2}\ldots$ of the LARS algorithm and some “residuals” $\overline{N}^{(1)},\overline{N}^{(2)},\ldots$ defined in Step 3. We will present here an equivalent formulation of these knots.

Assume that $k\geq 2$ and we have build $\lambda_{1},\ldots,\lambda_{k-1}$ and selected the “signed” variables $\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1}$. Introduce $N^{(k-1)}:=(\overline{N}^{(k-1)},-\overline{N}^{(k-1)})$ and define

[TABLE]

Check that $\theta_{j}(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1})=(\overline{\theta}_{j},-\overline{\theta}_{j})$ where we recall that we define

[TABLE]

at Step ${\bf 2}$ and it holds that $\widehat{\imath}_{\ell}$ and $(\overline{\imath}_{\ell},\varepsilon_{\ell})$ are related as in (2.1). From this equality, we deduce that it holds $N^{(k)}(\lambda)=(\overline{N}^{(k)}(\lambda),-\overline{N}^{(k)}(\lambda))$. One may also check that the coordinates $\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1}$ of $N^{(k)}(\lambda)$ are equal to $\lambda$.

Again if we want to solve $N_{j}^{(k)}(\lambda)=\lambda$ for some $j$, we have to limit our attention to $j$’s such that $\theta_{j}(\,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1})<1$. Solving this latter equality yields to

[TABLE]

This expression is consistent with $\lambda_{k}$ in Algorithm 2.

Now, we can give an other expression of $\lambda_{k}$ that will be useful in the proofs of our main theorems. Note that the residuals satisfy the relation

[TABLE]

and that $N^{(k-1)}_{j}=\lambda_{k-1}$ for $j=\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1}$. The following lemma permits a drastic simplification of the expression of the knots. Its proof is given in Appendix 7.6.

Lemma 5.2**.**

It holds

[TABLE]

where we denote $\Pi_{i_{1},\ldots,i_{k-1}}(Z)=(\Pi_{i_{1},\ldots,i_{k-1}}(Z_{1}),\ldots,\Pi_{i_{1},\ldots,i_{k-1}}(Z_{2p}))$ and note that, for all $j\in[2p]$, one has \Pi_{i_{1},\ldots,i_{k-1}}(Z_{j})=\big{(}R_{j,i_{1}}\cdots R_{j,i_{k-1}}\big{)}M^{-1}_{i_{1},\ldots,i_{k-1}}(Z_{i_{1}},\ldots,Z_{i_{k-1}}).

Using Lemma 5.2 we deduce that $\lambda_{k}$ in Algorithm 2 is consistent with

[TABLE]

where \Pi_{\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1}}(Z_{j})=\big{(}R_{j,\widehat{\imath}_{1}}\cdots R_{j,\widehat{\imath}_{k-1}}\big{)}M^{-1}_{\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1}}(Z_{\widehat{\imath}_{1}},\ldots,Z_{\widehat{\imath}_{k-1}}). When $\mathds{E}Z=0$, one may remark that $\Pi_{i_{1},\ldots,i_{k-1}}(Z_{j})$ is the regression of Z_{j}\leavevmode\nobreak\ on the vector $(Z_{i_{1}},\cdots,Z_{i_{k-1}})$ whose variance-covariance matrix is $M_{i_{1},\ldots,i_{k-1}}$. This analysis leads to an equivalent formulation of the LARS algorithm (Algorithm 2). We present this formulation in Algorithm 3.

Remark 5.3**.**

Note that Algorithm 2 implies that $\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k}$ are pairwise different, but also that they differ modulo $p$.

6 First Steps to Derive the Joint Law of the LARS knots

6.1 Law of the First Knot

One has the following lemma governing the law of $\lambda_{1}$.

Lemma 6.1**.**

It holds that

• •

$Z_{i_{1}}$* is independent of $(Z_{j}^{(i_{1})})_{j\neq i_{1}}$,*

• •

If $\theta_{j}(i_{1})<1$ for all $j\neq i_{1}$ then \big{\{}\widehat{\imath}_{1}=i_{1}\big{\}}=\big{\{}\lambda_{2}^{(i_{1})}\leq Z_{i_{1}}\big{\}},

• •

If $\theta_{j}(i_{1})<1$ for all $j\neq i_{1}$ then, conditional on $\{\widehat{\imath}_{1}=i_{1}\}$ and $\lambda_{2}$, $\lambda_{1}$ is a truncated Gaussian random variable with mean $\operatorname{\mathds{E}}(Z_{i_{1}})$ and variance $\rho_{1}^{2}:=R_{\widehat{\imath}_{1},\widehat{\imath}_{1}}$ subject to be greater than $\lambda_{2}$.

Proof 6.2**.**

The first point is a consequence or the properties of Gaussian regression. Now, observe that

[TABLE]

as claimed. The last statement is a consequence of the two previous points.

6.2 Recursive Formulation of the LARS

One has the following proposition whose proof can be found in Section 7.7. As we will see in this section, this intermediate result as a deep consequence, the LARS algorithm can be stated in a recursive way applying the same function repeatedly, as presented in Algorithm 4.

Proposition 6.3**.**

Set

[TABLE]

and observe that $\tau_{j,i_{k}}$ is the covariance between $Z^{(i_{1},\ldots,i_{k-1})}_{j}$ and $Z^{(i_{1},\ldots,i_{k-1})}_{i_{k}}$. Furthermore, it holds

[TABLE]

and

[TABLE]

Now, we present Algorithm 4. Define $R(0):=R$, $Z(0)=Z$ and $T(0)=0$. For $k\geq 1$ and $i_{1},\ldots,i_{k}\in[2p]$, introduce

[TABLE]

and note that $R(k)$ is the variance-covariance matrix of the Gaussian vector $Z(k)$. The key property is following. Let $v_{1},\ldots,v_{k},$ be $k$ linearly independent vectors of an Euclidean space and let $u$ be any vector of the space. Set

[TABLE]

the projection of $v_{k}$ orthogonally to $v_{1},\ldots,v_{k-1}$.Then

[TABLE]

Using this result we deduce that

[TABLE]

where ${\bf x}(k)=R_{i_{k}}(k-1)/R_{i_{k},i_{k}}(k-1)$. It yields that

[TABLE]

Using (6.1) (or (7.22)), remark that

[TABLE]

These relations give a recursive formulation of the LARS as presented in Algorithm 4.

7 Proofs

7.1 Proof of Proposition 3.3

**First and third points: ** The first point works by induction. The initialization of the proof is given by the second point of Lemma 6.1. We will use Proposition 6.3 to prove the first point. We have

[TABLE]

using (6.2) and that $1-{\tau_{j,i_{k}}}/{\tau_{i_{k},i_{k}}}>0$ (which is a consequence of (6.1) and (${\mathcal{A}}_{\mathrm{Irr.}}$)) in (7.1). By induction and using (7.1), it holds that

[TABLE]

Now, observe that $\widehat{i}_{k+1}$ is the (unique) arg max of $\lambda_{k+1}^{(i_{1},\ldots,i_{k})}$ on the event \big{\{}\widehat{\imath}_{1}=i_{1},\ldots,\widehat{\imath}_{k}=i_{k}\big{\}}.

It yields that

[TABLE]

as claimed. Stopping at $a$ as in ($s_{a}$) gives the second part of the statement. The third point of the proposition is a direct consequence of the first point.

**Second point: **The proof of the second point can be lead by induction. The initialization of the proof is given by the first point of Lemma 6.1. Now, observe that $Z_{i_{k}}^{(i_{1},\ldots,i_{k-1})},\cdots,Z_{i_{2}}^{(i_{1})},Z_{i_{1}}$ are measurable functions of $(Z_{i_{1}},\ldots,Z_{i_{k}})$ and one may check that $(Z_{i_{1}},\ldots,Z_{i_{k}})$ is independent of $(Z^{(i_{1},\ldots,i_{k})}_{j})_{j\neq i_{1},\ldots,i_{k}}$. We deduce that $Z_{i_{k}}^{(i_{1},\ldots,i_{k-1})},\cdots,Z_{i_{2}}^{(i_{1})},Z_{i_{1}}$ are independent of $(Z^{(i_{1},\ldots,i_{k})}_{j})_{j\neq i_{1},\ldots,i_{k}}$. One can also check that $Z_{i_{k+1}}^{(i_{1},\ldots,i_{k})}$ is independent of $(Z^{(i_{1},\ldots,i_{k})}_{j})_{j\neq i_{1},\ldots,i_{k+1}}$. Deduce that the variables $Z_{i_{K}}^{(i_{1},\ldots,i_{K-1})},\cdots,Z_{i_{2}}^{(i_{1})},Z_{i_{1}}$ are independent and independent of $(Z^{(i_{1},\ldots,i_{K})}_{j})_{j\neq i_{1},\ldots,i_{K}}$. This gives the first part of the second point of the proposition.

Now observe that $(Z^{(i_{1},\ldots,i_{K})}_{j})_{j\neq i_{1},\ldots,i_{K}}$ is a linear function of $P^{\perp}_{K}Y$. More precisely, consider the vector $\widetilde{W}:=(Z^{(i_{1},\ldots i_{K})}_{j})_{j\in[p]}$ and note that

[TABLE]

where $\theta^{K}:=(\theta_{j}(i_{1},\ldots,i_{K}))_{j\in[p]}$, $P_{K}^{\perp}=\mathrm{Id}_{n}-P^{(i_{1},\ldots,i_{K})}$, and with the convention $0/0=0$. The key remark is that the direction and the norm of a centered Gaussian vector are independent, namely $\|P_{K}^{\perp}(Y)\|$ is independent of $P_{K}^{\perp}(Y)/\|P_{K}^{\perp}(Y)\|$ when $\mathds{E}P_{K}^{\perp}(Y)=0$. Observe that $\mathds{E}P_{K}^{\perp}(Y)=P_{K}^{\perp}(X\beta^{0})$ and $\mathds{E}P_{K}^{\perp}(Y)=0$ is equivalent to $X\beta^{0}\in H_{K}$. We deduce that if $X\beta^{0}\in H_{K}$ then $\widehat{\sigma}^{i_{1},\ldots,i_{K}}$ is independent of $\widetilde{W}/\sigma^{i_{1},\ldots,i_{K}}$, as claimed.

7.2 Proof of Proposition 3.7

Fix $a$ such that $0\leq a\leq K-1$ and consider any selection procedure $\widehat{m}$ satisfying ($\bf\mathcal{A}_{\mathrm{\,Stop}}$). From Theorem 1.1 (more precisely (3.15)) we know that the density of $(\lambda_{a+1},\lambda_{a+2},\ldots,\lambda_{K})$ conditional on

[TABLE]

is given by

[TABLE]

From Proposition 3.2, conditional on $\mathcal{F}$ and under the null hypothesis described in Proposition 3.7, we know that $m_{a+1}=\ldots=m_{K}=0$. It implies that $\Phi_{k}$ is the CDF of $\lambda_{k}$ for $a<k\leq K$.

From the definition of a stopping time given by ($\bf\mathcal{A}_{\mathrm{\,Stop}}$) and on the event $\mathcal{F}$, we know that $\mathds{1}_{\{\widehat{m}=a\}}$ is a measurable function of $\lambda_{1},\ldots,\lambda_{a-1}$ which are respectively equal to $Z_{i_{1}},\ldots,Z_{i_{a-1}}^{(i_{1},\ldots,i_{a-2})}$ on $\mathcal{F}$ by (3.5) (as proven in Appendix 5.4 and Eq. (5.3)). By Proposition 3.3 (more precisely (3.10)), we also know that this function is independent of $(\lambda_{a+1},\lambda_{a+2},\ldots,\lambda_{K})$ conditional on $\mathcal{F}$. Remark that its is also independent of $\widehat{\sigma}^{i_{1},\ldots,i_{K}}$ conditional of $\mathcal{F}$ for the same reason (it would be useful later, when we will build testing procedures when the variance is unknown). We deduce that the conditional density above is also the conditional density on the event

[TABLE]

From $F_{k}=\Phi_{k}(\lambda_{k})$ (i.e., applying the CDF) we deduce by a change of variables that conditional on the selection event $\mathcal{G}$, the vector $(F_{a+1},\ldots,F_{K})$ is uniformly distributed on

[TABLE]

where $\mathcal{P}_{i,j}$ are described in (3.16).

7.3 Proof of Proposition 3.11

By Proposition 3.7, a simple integration shows that

[TABLE]

under the null hypothesis of Proposition 3.11 (which implies that $m_{a+1}=\ldots=m_{K}=0$). Then note that the function $\mathds{F}_{abc}$ is defined by $\sigma,\lambda_{a},\lambda_{c},\widehat{\imath}_{a+1},\ldots,\widehat{\imath}_{c-1}$ only. We deduce that we can de-condition on $\widehat{m}=a,\widehat{\imath}_{1},\ldots,\widehat{\imath}_{a},\widehat{\imath}_{c},\ldots,\widehat{\imath}_{K}$, which gives the result.

7.4 Proof of Proposition 3.18

Let us fix some values $i_{1},\ldots,i_{K+1}$. Recall that the frozen values of the knots

[TABLE]

are Gaussian, independent, and \widehat{\sigma}^{i_{1},\ldots,i_{K}}\mathchoice{\mathrel{\hbox to0.0pt{\displaystyle\perp\hss}\mkern 2.0mu{\displaystyle\perp}}}{\mathrel{\hbox to0.0pt{\textstyle\perp\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptstyle\perp\hss}\mkern 2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptscriptstyle\perp\hss}\mkern 2.0mu{\scriptscriptstyle\perp}}}\big{(}Z_{i_{K+1}}^{(i_{1},\ldots,i_{K})}/\widehat{\sigma}^{i_{1},\ldots,i_{K}}\big{)}\mathchoice{\mathrel{\hbox to0.0pt{\displaystyle\perp\hss}\mkern 2.0mu{\displaystyle\perp}}}{\mathrel{\hbox to0.0pt{\textstyle\perp\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptstyle\perp\hss}\mkern 2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptscriptstyle\perp\hss}\mkern 2.0mu{\scriptscriptstyle\perp}}}Z_{i_{K}}^{(i_{1},\ldots,i_{K-1})}\mathchoice{\mathrel{\hbox to0.0pt{\displaystyle\perp\hss}\mkern 2.0mu{\displaystyle\perp}}}{\mathrel{\hbox to0.0pt{\textstyle\perp\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptstyle\perp\hss}\mkern 2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptscriptstyle\perp\hss}\mkern 2.0mu{\scriptscriptstyle\perp}}}\cdots\mathchoice{\mathrel{\hbox to0.0pt{\displaystyle\perp\hss}\mkern 2.0mu{\displaystyle\perp}}}{\mathrel{\hbox to0.0pt{\textstyle\perp\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptstyle\perp\hss}\mkern 2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptscriptstyle\perp\hss}\mkern 2.0mu{\scriptscriptstyle\perp}}}Z_{i_{2}}^{(i_{1})}\mathchoice{\mathrel{\hbox to0.0pt{\displaystyle\perp\hss}\mkern 2.0mu{\displaystyle\perp}}}{\mathrel{\hbox to0.0pt{\textstyle\perp\hss}\mkern 2.0mu{\textstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptstyle\perp\hss}\mkern 2.0mu{\scriptstyle\perp}}}{\mathrel{\hbox to0.0pt{\scriptscriptstyle\perp\hss}\mkern 2.0mu{\scriptscriptstyle\perp}}}Z_{i_{1}}, see Proposition 3.3 and (3.10). Let us condition by $\{\widehat{\imath}_{1}=i_{1},\ldots,\widehat{\imath}_{K}=i_{K},\lambda_{K+1}=\ell_{K+1}\}$. Note that, on this event, $\lambda_{h}^{f}=\lambda_{h}$ , $h\in[K]$ and by Proposition 3.3 this event is equivalent to

[TABLE]

Because of the independence above, we get the conditional independence stated in the first and second point of the proposition.

For the last point, consider $\widetilde{W}:=(Z^{(i_{1},\ldots i_{K})}_{j})_{j\in[p]}$ and note that

[TABLE]

where $\theta^{K}:=(\theta_{j}(i_{1},\ldots,i_{K}))_{j\in[p]}$, $P_{K}^{\perp}=\mathrm{Id}_{n}-P^{(i_{1},\ldots,i_{K})}$, and with the convention $0/0=0$. Now, let $U\in\mathds{R}^{n\times(n-K)}$ be any matrix such that $UU^{\top}=P_{K}^{\perp}$ and define $W:=U^{\top}Y/\sigma$, then

[TABLE]

Because of the independence above and (7.2), the distribution of $\widehat{\sigma}$ is independent of the other variables and such that

[TABLE]

This implies that the conditional distribution is the one claimed.

7.5 Orthogonal Case: Proof of Theorem 1.4

Let $\mathcal{I}$ the set of admissible indexes

[TABLE]

$\circ$ **Step 1: ** We prove that, when the considered indexes such that $c+1\leq K+1$ belong to $\mathcal{I}$, $\mathcal{S}_{a,b,c+1}$ is more powerful than $\mathcal{S}_{a,b,c}$. Our proof is conditional to F_{a}=f_{a},F_{c+1}=f_{c+1}\leavevmode\nobreak\. Note that (F_{a+1},\ldots,F_{c})\leavevmode\nobreak\ has for distribution the uniform distribution on the simplex

[TABLE]

This implies by direct calculations that

[TABLE]

and that

[TABLE]

where $\bm{F}_{\beta}$ is the cumulative distribution of the Beta distribution in reference. Using monotony of this function the $\mathcal{S}_{abc}$ test has for rejection region

[TABLE]

where $z_{1}$ is some threshold, depending on $\alpha$, that belongs to $(0,1)$.

Similarly $\mathcal{S}_{ab(c+1)}$ has for rejection region

[TABLE]

where $z_{2}$ is some other threshold belonging to $(0,1)$. We use the following lemmas.

Lemma 7.1**.**

Let $c\leq K$. The density $h_{\mu}$ of $f_{1},\ldots,f_{c}$, conditional on $F_{c+1}$ with respect of the Lebesgue measure under the alternative is coordinate-wise non-decreasing and given by (7.6).

Proof 7.2**.**

Observe that it suffises to prove the result when $\sigma=1$. Note that

[TABLE]

Thus its density $p_{\mu^{0},i_{1},\ldots,i_{c}}$ does not depend on $\mu^{0}_{\overline{i}_{1}},\ldots,\mu^{0}_{\overline{i}_{c}}$. As a consequence the following variables have the same distribution ; \lambda_{c+1}^{i_{1}+\epsilon_{1}p,\dots,i_{c}+\epsilon_{c}p}\leavevmode\nobreak\ , where $\epsilon_{1},\ldots,\epsilon_{c}$ take the value 0 or 1 and indices are taken modulo $p$.

Because of the independence of the different variables, the joint density, under the alternative hypothesis, of $\lambda_{1},\ldots,\lambda_{c+1}$ taken at $\ell_{1},\ldots,\ell_{c+1},$ on the domain $\{\lambda_{1}>\cdots>\lambda_{c+1}\}$ takes the value

[TABLE]

Here the sum $\displaystyle\sum^{\prime}$ is taken over all different $j_{1},\ldots,j_{c}$ belonging to $\llbracket 1,p\rrbracket$.

Then the density, conditional on $F_{c+1}=f_{c+1}$, of  F_{1},\ldots,F_{c}\leavevmode\nobreak\ at $f_{1},\ldots,f_{c}$ takes the value

[TABLE]

implying that this density is coordinate-wise non-decreasing.

Lemma 7.3**.**

Let $\nu_{0}$ the image on the plane $(F_{b},F_{c})$ on the uniform probability on $\mathcal{S}$: it is the distribution under the null of $(F_{b},F_{c})$. The two rejection regions : $\mathcal{R}_{1}$ associated to (7.4) and $\mathcal{R}_{2}$ associated to (7.5) have of course the same probability $\alpha$ under $\nu_{0}$. Let $\eta_{\mu^{0}}$ the density w.r.t. $\nu_{0}$ of the distribution of $(F_{b},F_{c})$ under the alternative. Then $\eta_{\mu^{0}}$ is non decreasing coordinate-wise.

Proof 7.4**.**

Integration yields that density of $\nu_{0}$ w.r.t. the Lebesgue measure taken at point $(f_{b},f_{c})$ is

[TABLE]

The density of $\nu_{\mu^{0}}$ w.r.t. Lebesgue measure is

[TABLE]

Thus $\eta_{\mu^{0}}$ which is the quotient of these two quantities is just a mean value of $h_{\mu^{0}}$ on the domain of integration $\mathcal{D}_{f_{b},f_{c}}$ in (7.7).

Suppose that $f_{b}$ and $f_{c}$ increase, then all the borns of the domain \mathcal{D}_{f_{b},f_{c}}\leavevmode\nobreak\ increase also. By Lemma 7.1 the mean value increases.

We finish now the proof of Step 1: For a given level $\alpha$ let us consider the two rejection regions $R_{a,b,c}$ and R_{a,b,(c+1)}\leavevmode\nobreak\ of the two considered tests in the plane $F_{b},F_{c}$ and set

[TABLE]

see Figure 8. These two regions have the same $\nu_{0}$ measure. By elementary geometry there exist a point $K=(K_{b},K_{c})$ in the plane such that

• •

For every point of $A$, $F_{b}\leq K_{b}$, $F_{c}\leq K_{c}$,

• •

For every point of $B$, $F_{b}\geq K_{b}$, $F_{c}\geq K_{c}$,

By transport of measure there exists a transport function $\mathcal{T}$ that preserve the measure $\nu_{0}$ and that is one-to one $A\to B$. As a consequence the transport by $\mathcal{T}$ improve the probability under the alternative: the power of $\mathcal{S}_{a,b,c+1}$ is larger than that of  $\mathcal{S}_{a,b,c}$.

$\circ$ **Step 2: ** We prove that, when the considered indexes belong to $\mathcal{I}$ such that $a<b-1$, $\mathcal{S}_{a,(b-1),c}$ is more powerful than $\mathcal{S}_{a,b,c}$. Our proof is conditional on F_{a}=f_{a},F_{b}=f_{b}\leavevmode\nobreak\ and is located in the plane $(F_{b-1},F_{c})$.

The rejection region $R_{a,b,c}$ takes the form  F_{c}\leq\frac{1}{1-z_{1}}f_{b}-\frac{z_{1}}{1-z_{1}}f_{a}\leavevmode\nobreak\ for some threshold $z_{1}$ belonging to $(0,1)$.

The rejection region $R_{a,(b-1),c}$ takes the form F_{c}\leq\frac{1}{1-z_{2}}F_{b-1}-\frac{z_{2}}{1-z_{2}}f_{a}\leavevmode\nobreak\ for some other threshold $z_{2}$ belonging to $(0,1)$.

These regions as well as the regions $A$ and $B$ and the point $K$ are indicated in Figure 8.

Transport of measure and the convenient modification of Lemma 7.3 imply that the power of the test $\mathcal{S}_{a,(b-1),c}$ is greater of equal than that of $\mathcal{S}_{a,b,c}$.

$\circ$ **Step 3: ** We prove that, when the considered indexes belong to $\mathcal{I}$ such that $a+1<b$, $\mathcal{S}_{a,b,c}$ is more powerful than $\mathcal{S}_{(a+1),b,c}$. Our proof is conditional on F_{a}=f_{a},F_{c}=f_{c}\leavevmode\nobreak\ and is located in the plane $F_{a+1},F_{b}$.

The rejection region $R_{a,b,c}$ takes the form F_{b}\geq z_{1}f_{a}+(1-z_{1})f_{c}\leavevmode\nobreak\ for some threshold $z_{1}$belonging to $(0,1)$.

The rejection region $R_{a+1,b,c}$ takes the form F_{b}\geq z_{2}F_{a+1}+(1-z_{2})f_{c}\leavevmode\nobreak\ for some other threshold $z_{2}$ belonging to $(0,1)$.

These regions as well as the regions $A$ and $B$ and the point $K$ are indicated in Figure 8.

Transport of measure and the convenient modification of Lemma 7.3 imply that the power of $\mathcal{S}_{a,b,c}$ is greater of equal that that of $\mathcal{S}_{(a+1),b,c}$.

Considering the three cases above, we get the desired result.

7.6 Proof of Lemma 5.2

The proof works by induction. Let us check the relation for $k=2$, namely

[TABLE]

Now, let $k\geq 3$. First, the three perpendicular theorem implies that for every j,i_{1},\ldots,i_{k-1}\leavevmode\nobreak\,

[TABLE]

By induction, using (5.2), we get that

[TABLE]

Then, recall that $N^{(k-1)}_{j}=\lambda_{k-1}$ for $j=\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1}$ and remark that

[TABLE]

Using (7.8) at indices $j=\widehat{\imath}_{1},\ldots,\widehat{\imath}_{k-1}$, we deduce that

[TABLE]

Namely

[TABLE]

Using again (7.8) we get that

[TABLE]

as claimed.

7.7 Proof of Proposition 6.3

We denote

[TABLE]

and observe that

[TABLE]

using Schur complement of block $M$ of the matrix $\overline{M}$ and a LU decomposition. Note also that

[TABLE]

To prove (6.2), it suffices to show that the R.H.S term above is equal to the following R.H.S term

[TABLE]

We will prove that numerators are equal and that denominators are equal. For denominators,

[TABLE]

using (7.18). Furthermore, it proves (6.1). For the numerators, we use that

[TABLE]

using (7.18).

7.8 Proof of Theorem 3.23

We rely on the Weak Positive Regression Dependency (WPRDS) property to prove the result, one may consult (Giraud,, 2014, Page 173) for instance. We say that a function $g:[0,1]^{K}\to\mathds{R}^{+}$ is nondecreasing if for any $p,q\in[0,1]^{K}$ such that $p_{k}\geq q_{k}$ for every $k=1,\ldots,K$, we have $g(p)\geq g(q)$. We say that a Borel set $\Gamma\in[0,1]^{K}$ is nondecreasing if $g=\mathds{1}_{\Gamma}$ is nondecreasing. In other words if $y\in\gamma$ and if $z\geq 0$, then $y+z\in\gamma$. We say that the $p$-values $(\,\widehat{p}_{1}=\widehat{\alpha}_{0,1,K+1},\ldots,\widehat{p}_{K}=\widehat{\alpha}_{K-1,K,K+1})$ satisfy the WPRDS property if for any nondecreasing set $\Gamma$ and for all $k^{0}\in I_{0}$, the function

[TABLE]

where $\mu^{0}=\beta^{0}$ in our orthogonal design case, and we recall that

[TABLE]

To prove Theorem 3.23, note that it is sufficient (Giraud,, 2014, Chapter 8) to prove that

[TABLE]

where $\overline{\mathds{E}},\overline{\operatorname{\mathds{P}}}$ will denote that expectations and probabilities are conditional on $\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},\lambda_{K+1}\}$ and under the hypothesis that $\mu^{0}=X^{\top}X\beta^{0}$. Note that one can integrate in $\lambda_{K+1}$ to get the statement of Theorem 3.23.

$\circ$ **Step 1: **We start by giving the joint law of the LARS knots under the alternative in the orthogonal design case. Lemma 7.1 and (7.6) show that, conditional on $\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},\lambda_{K+1}\}$, $(\lambda_{1},\ldots,\lambda_{K})$ is distributed on the set $\lambda_{1}\geq\lambda_{2}\geq\cdots\geq\lambda_{K}\geq\lambda_{K+1}$ and it has a coordinate-wise nondecreasing density. Now we can assume without loss of generality that $\sigma^{2}=1$, in addition because of orthogonality  $\rho_{k}^{2}=1$ implying that $F_{k}=\Phi(\lambda_{k})$ $\mathcal{P}_{i,j}=\Phi_{i}\circ\Phi_{j}^{-1}=\mathrm{Id}$. We deduce that, conditional on $\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1}\}$, $(F_{1},\ldots,F_{K})$ is distributed on the set

[TABLE]

it has an explicit density given by (7.6), and we denote it by $h_{\mu^{0}}$. By the change of variables $G_{k}:=\frac{F_{k}-F_{K+1}}{F_{k-1}-F_{K+1}}$ one obtains that the distribution of $(G_{1},\ldots,G_{K})$ is supported on $[0,1]^{K}$. More precisely, define

[TABLE]

whose inverse Jacobian determinant is

[TABLE]

We deduce that the density of $(G_{1},\ldots,G_{K})|\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1}\}$ at point $g$ with respect to Lebesgue measure is

[TABLE]

where we have used (7.6). From (1.3) and (7.3), one has

[TABLE]

where $\bm{F}_{\beta}$ is the cumulative distribution of the Beta distribution in reference. We deduce that for any $v\in(0,1)$ and for any $\ell\in[K]$,

[TABLE]

so that

[TABLE]

where $\overline{\Gamma}$ can be proved to be a nonincreasing Borel set from (7.25).

$\circ$ **Step 2: ** Let $0<x<y<1$ and denote by $\mu_{x}$ the following conditional law

[TABLE]

Remark that if there exists a measurable $T:[0,1]^{K}\mapsto[0,1]^{K}$ such that

• •

$T$ is nondecreasing, meaning that for any $g\in[0,1]^{K}$, $T(g)\geq g$;

• •

$T$ is such that push-forward of $\mu_{x}$ by $T$ gives $\mu_{y}$, namely $T_{\#}\mu_{x}=\mu_{y}$;

then it holds

• •

$\mathds{1}_{\{T(g)\in\overline{\Gamma}\}}\leq\mathds{1}_{\{g\in\overline{\Gamma}\}}$;

• •

\mathrm{law}\big{[}T(G)|\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1},G_{k^{0}}\geq x\}\big{]}=\mathrm{law}\big{[}G|\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1},G_{k^{0}}\geq y\}\big{]} where $G=(G_{1},\ldots,G_{K})$.

In this case, we deduce that

[TABLE]

If one can prove that such function $T$ exists for any $0<x<y<1$, it proves that

[TABLE]

and, in view of (7.26), it proves (7.23). Proving that such function $T$ exists is done in the next step.

$\circ$ **Step 3: ** Let $0<x<y<1$. Consider the Knothe-Rosenblatt transport map $T$ of $\mu_{x}$ toward $\mu_{y}$ following the order

[TABLE]

It is based on a sequence of conditional quantile transforms defined following the ordering above. Its construction is presented for instance in (Santambrogio,, 2015, Sec.2.3, P.67) or (Villani,, 2008, P.20). The transport $T$ is defined as follows. Given $z,z^{\prime}\in[0,1]^{K}$ such that $z^{\prime}=T(z)$ it holds

[TABLE]

where $T^{(k^{0})},T^{(k^{0}+1)},\ldots,T^{(K)},T^{(k^{0}-1)}\ldots,T^{(1)}$ will be build in the sequel, in which we will drop their dependencies in the $z^{\prime}_{k}$’s to ease notations. It remains to prove that

• •

$T$ is nondecreasing, meaning that for any $g\in[0,1]^{K}$, $T(g)\geq g$;

• •

$T$ is such that push-forward of $\mu_{x}$ by $T$ gives $\mu_{y}$, namely $T_{\#}\mu_{x}=\mu_{y}$;

to conclude. The last point is a property of the Knothe-Rosenblatt transport map. Proving the first point will be done in the rest of the proof.

$\circ$ *Step 3.1: * We start by the first transport map $T^{(k^{0})}:[0,1]\mapsto[0,1]$. Denote $\mu_{x}^{(k^{0})}$ the following conditional law

[TABLE]

and $\mathds{F}_{x}^{(k^{0})}$ its cdf. Note that the Knothe-Rosenblatt construction gives $T^{(k^{0})}=(\mathds{F}_{y}^{(k^{0})})^{-1}\circ\mathds{F}_{x}^{(k^{0})}$. We would like to prove that $T^{(k^{0})}(t)\geq t$ for all $z\in(0,1)$. This is equivalent to prove that it holds $\mathds{F}_{x}^{(k^{0})}\geq\mathds{F}_{y}^{(k^{0})}$. For $t\leq y$, $\mathds{F}_{y}^{(k^{0})}(t)=0$ and it implies that $\mathds{F}_{x}^{(k^{0})}(t)\geq\mathds{F}_{y}^{(k^{0})}(t)$. Let $t>y$, using the conditional density $\mathbf{p}$ defined in (7.24), note that

[TABLE]

where, for example

[TABLE]

A simple calculation (see also Figure 9) gives that

[TABLE]

and it proves that $\mathds{F}_{x}^{(k^{0})}\geq\mathds{F}_{y}^{(k^{0})}$.

$\circ$ *Step 3.2: * We continue with the second transport map in Knothe-Rosenblatt construction. Let $z_{k^{0}}\in(x,1)$ and denote $\mu_{z_{k^{0}}}^{(k^{0}+1)}$ the following conditional law

[TABLE]

and $\mathds{F}_{z_{k^{0}}}^{(k^{0}+1)}$ its cdf. Let $z^{\prime}_{k^{0}}:=T^{(k^{0})}(z^{\prime}_{k^{0}})$ and denote $\mu_{z^{\prime}_{k^{0}}}^{(k^{0}+1)}$ the following conditional law

[TABLE]

and $\mathds{F}_{z^{\prime}_{k^{0}}}^{(k^{0}+1)}$ its cdf. Note that $x<z_{k^{0}}\leq z^{\prime}_{k^{0}}=T^{(k^{0})}(z_{k^{0}})\leq 1$. Again, we would like to prove that $\mathds{F}_{z_{k^{0}}}^{(k^{0}+1)}\geq\mathds{F}_{z^{\prime}_{k^{0}}}^{(k^{0}+1)}$ which implies that the transport map T^{(k^{0}+1)}:=\big{(}\mathds{F}_{z^{\prime}_{k^{0}}}^{(k^{0}+1)}\big{)}^{-1}\circ\mathds{F}_{z_{k^{0}}}^{(k^{0}+1)} satisfies $T^{(k^{0}+1)}(u)\geq u$ for all $u\in(0,1)$.

Recall that the conditional density $\mathbf{p}$ of $G|\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1}\}$ is given by (7.24) and recall that $k^{0}\in I_{0}$. Observe that $\mu^{0}_{k^{0}}=0$, so that the conditional density of $G|\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1},G_{k^{0}}=\mathbf{z}\}$ is

[TABLE]

Set $\tau:=z^{\prime}_{k^{0}}/z_{k^{0}}\geq 1$ and $G^{\prime}_{k^{0}+1}=\tau G_{k^{0}+1}$ so that

[TABLE]

Denote $G^{\prime}:=(G_{1},\ldots,G_{k^{0}-1},G^{\prime}_{k^{0}+1},G_{k^{0}+2},\ldots,G_{K})\in(0,1)^{k^{0}-1}\times(0,\tau)\times(0,1)^{K-k^{0}-1}$ and note that the conditional density of $G^{\prime}|\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1},G_{k^{0}}=\tau\mathbf{z}\}$ is

[TABLE]

which, up to some normalising constant, is the same as (7.27) up to the following change of support

[TABLE]

By an abuse of notation, we denote by $\mathbf{p}$ this function, namely

[TABLE]

We deduce that

[TABLE]

where

[TABLE]

We now present an inequality on the to conclude. Observe that we are integrating on domains depicted in Figure 10. The two boxes have same area for the uniform measure and we would like to compare their respective measure for the $\mathbf{p}\otimes\mathbf{p}$ measure. We start by the next lemma whose proof is omitted.

Lemma 7.5**.**

Let $a,b\geq 0$. The function

[TABLE]

is non-decreasing on the domain $[1,\infty)$.

Now, let $(g_{1},\ldots,g_{k^{0}-1},g_{k^{0}+2}\ldots,g_{K})\in(0,1)^{K-1}$ be fixed in the integrals (7.28). We are the looking at the weights of the domains $(h_{1},h_{2})\in(0,t)\times(0,\tau)$ and $(h_{3},h_{4})\in(0,\tau t)\times(0,1)$ for the weight function $w$ given by

[TABLE]

where the constant $\mathrm{C}_{1}$ depends on $(g_{1},\ldots,g_{k^{0}-1},g_{k^{0}+2}\ldots,g_{K})\in(0,1)^{K-1}$. By the change of variables $h^{\prime}_{1}=h_{3}/t$ and $h^{\prime}_{2}=th_{4}$, the right hand term of (7.28) is given by the integration on the domain $(h^{\prime}_{1},h^{\prime}_{2})\in(0,t)\times(0,\tau)$ of the weight function $w^{\prime}$ given by

[TABLE]

Now, invoke Lemma 7.5 with

[TABLE]

where $h=h_{1}$ or $h_{2}$, to get that $w^{\prime}\geq w$ and so

[TABLE]

which concludes this part of the proof.

$\circ$ *Step 3.3: * We continue by induction with the other transport maps in Knothe-Rosenblatt’s construction. Assume that we have built $z^{\prime}:=(z^{\prime}_{k},\ldots,z^{\prime}_{k^{0}})$ and $z:=(z_{k},\ldots,z_{k^{0}})$ for some $k>k^{0}$. Denote $\mu_{z}^{(k+1)}$ the following conditional law

[TABLE]

and $\mathds{F}_{z}^{(k+1)}$ its cdf. Denote $\mu_{z^{\prime}}^{(k+1)}$ the following conditional law

[TABLE]

and $\mathds{F}_{z^{\prime}}^{(k+1)}$ its cdf. Note that $z\leq z^{\prime}=T^{(k)}(z)\leq 1$. Again, we would prove that $\mathds{F}_{z}^{(k+1)}\geq\mathds{F}_{z^{\prime}}^{(k+1)}$ which implies that the transport map T^{(k+1)}:=\big{(}\mathds{F}_{z^{\prime}}^{(k+1)}\big{)}^{-1}\circ\mathds{F}_{z}^{(k+1)} satisfies $T^{(k+1)}(u)\geq u$ for all $u\in(0,1)$.

For $\mathbf{z}\in(0,1)^{k-k^{0}}\times(x,1)$, the conditional density of $G|\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1},G^{[k,k^{0}]}=\mathbf{z}\}$ is

[TABLE]

Set $\tau:=\prod_{n=k^{0}}^{k}z^{\prime}_{n}/\prod_{n=k^{0}}^{k}z_{n}\geq 1$ and $G^{\prime}_{k}=\tau G_{k^{0}+1}$ so that

[TABLE]

Then the proof follows the same idea as in Step 3.2 and we will not detail it here.

$\circ$ *Step 3.4: * This is the last step of the proof. Assume that we have built $z^{\prime}:=(z^{\prime}_{K},\ldots,z^{\prime}_{k^{0}})$ and $z:=(z_{K},\ldots,z_{k^{0}})$. Denote $\mu_{z}^{(k^{0}-1)}$ the following conditional law

[TABLE]

and $\mathds{F}_{z}^{(k^{0}-1)}$ its cdf. Denote $\mu_{z^{\prime}}^{(k^{0}-1)}$ the following conditional law

[TABLE]

and $\mathds{F}_{z^{\prime}}^{(k^{0}1)}$ its cdf. Note that $z\leq z^{\prime}=T^{(K)}(z)\leq 1$. Again, we would prove that $\mathds{F}_{z}^{(k^{0}-1)}\geq\mathds{F}_{z^{\prime}}^{(k^{0}-1)}$ which implies that the transport map T^{(k^{0}-1)}:=\big{(}\mathds{F}_{z^{\prime}}^{(k^{0}-1)}\big{)}^{-1}\circ\mathds{F}_{z}^{(k^{0}-1)} satisfies $T^{(k^{0}-1)}(u)\geq u$ for all $u\in(0,1)$.

For $\mathbf{z}\in(0,1)^{K-k^{0}}\times(x,1)$, the conditional density of $G|\{\overline{\imath}_{1},\ldots,\overline{\imath}_{K},F_{K+1},G^{[K,k^{0}]}=\mathbf{z}\}$ is

[TABLE]

Now, let $(g_{1},\ldots,g_{k^{0}-2})\in(0,1)^{k^{0}-2}$ be fixed and denote by

[TABLE]

and, substituting $z$ by $z^{\prime}$, define $w_{z^{\prime}}$ as well. Let $t\in(0,1)$. Following the idea of Step 3.2, one can check that it is sufficient to prove that

[TABLE]

Substituting

[TABLE]

on both parts, one is reduced to prove that

[TABLE]

Observe that $g\leq g^{\prime}$ in the last two integrals. Now, we have this lemma whose proof is omitted.

Lemma 7.6**.**

Let $0<a\leq a^{\prime}$ and $b>0$. The function

[TABLE]

is non-increasing on the domain $(0,\infty)$.

Let $g\leq g^{\prime}$. From Lemma 7.6, we deduce that $\cosh(ag+b)\cosh(a^{\prime}g^{\prime}+b)\geq\cosh(ag^{\prime}+b)\cosh(a^{\prime}g+b)$, proving that $w_{z}(g)w_{z^{\prime}}(g^{\prime})\geq w_{z^{\prime}}(g)w_{z}(g^{\prime})$. It proves that $T^{(k^{0}-1)}(u)\geq u$ for all $u\in(0,1)$.

We then proceed by induction for $k^{0}-1\to k^{0}-2\to\cdots\to 1$. The proof follows the same line as above, Step 3.4.

8 A Quasi Monte Carlo (QMC) method: Cubature by lattice rule

Our goal is to compute the integral of some function$f$ on the hypercube of dimension $d$, namely

[TABLE]

We want to approximate it by a finite sum over $n$ points

[TABLE]

A convenient way of constructing the sequence $x^{(i)},i=1,\dots,n$ is the so-called lattice rule: from the first point $x^{(1)}$ we deduce the others $x^{(i)}$ by

[TABLE]

where the $\{\}$ brackets mean that we take the fractional part coordinate by coordinate. In such a case the error given by

[TABLE]

is a function, in particular, of starting point $x^{(1)}$ .

The Fast-rank algorithm (Nuyens and Cools,, 2006) is a fast algorithm that finds, component by component and as a function of the prime $n$, the sequence of coordinates of $x^{(1)}$ that minimizes the maximal error when f\leavevmode\nobreak\ varies in a unit ball $\mathcal{E}$ of some RKHS, namely a tensorial product of Koborov spaces. In addition it gives an expression of its minimax error, namely

[TABLE]

In practice, very few properties are known on the function $f$, so the result above is not directly applicable. Nevertheless for many functions $f$, it happens that the convergence of $I_{n}$ to $I$ is “fast”: typically of the order $1/n$ while the Monte-Carlo method (choosing the $x^{(i)}$ at random) converges at rate $1/\sqrt{n}$.

A reliable estimate of the error is obtained by adding a Monte-Carlo layer as in Genz, (1992) for instance. This can be done as follows. Let $U$ a unique uniform variable on $[0,1]^{d}$, we define

[TABLE]

Classical computations show that $I_{n,U}$ is now an unbiased estimator of $I$. In a final step, we perform $N$ (in practice 15-20) independent repetitions of the experiment above an we compute usual asymptotic confidence intervals for independent observations.