Camille Charbonnier
email: camille.charbonnier@ensae.org

Nicolas Verzelen
email: nicolas.verzelen@supagro.inra.fr

Fanny Villers
email: fanny.villers@upmc.fr

A Global Homogeneity Test for High-Dimensional Linear Regression

Keywords: AMS 2000 subject classifications: Primary 62H15; secondary 62P10 Gaussian graphical model, two-sample hypothesis testing, high-dimensional statistics, multiple testing, adaptive testing, minimax hypothesis testing, detection boundary

This paper is motivated by the comparison of genetic networks based on microarray samples. The aim is to test whether the differences observed between two inferred Gaussian graphical models come from real differences or arise from estimation uncertainties. Adopting a neighborhood approach, we consider a twosample linear regression model with random design and propose a procedure to test whether these two regressions are the same. Relying on multiple testing and variable selection strategies, we develop a testing procedure that applies to high-dimensional settings where the number of covariates p is larger than the number of observations n 1 and n 2 of the two samples. Both type I and type II errors are explicitely controlled from a non-asymptotic perspective and the test is proved to be minimax adaptive to the sparsity. The performances of the test are evaluated on simulated data. Moreover, we illustrate how this procedure can be used to compare genetic networks on Hess et al breast cancer microarray dataset.

Introduction

The recent flood of high-dimensional data has motivated the development of a vast range of sparse estimators for linear regressions, in particular a large variety of derivatives from the Lasso. Although theoretical guarantees have been provided in terms of prediction, estimation and selection performances (among a lot of others [START_REF] Bickel | Simultaneous analysis of Lasso and Dantzig selector[END_REF][START_REF] Meinshausen | Lasso-type recovery of sparse representations for high-dimensional data[END_REF][START_REF] Wainwright | Sharp thresholds for high-dimensional and noisy sparsity recovery using 1 -constrained quadratic programming (lasso)[END_REF]), the research effort has only recently turned to the construction of high-dimensional confidence intervals or parametric hypothesis testing schemes [START_REF] Bühlmann | Statistical significance in high-dimensional linear models[END_REF][START_REF] Javanmard | Confidence Intervals and Hypothesis Testing for High-Dimensional Regression[END_REF][START_REF] Lockhart | A significance test for the lasso[END_REF][START_REF] Meinshausen | Assumption-free confidence intervals for groups of variables in sparse high-dimensional regression[END_REF][START_REF] Zhang | Confidence intervals for low-dimensional parameters with high-dimensional data[END_REF]]. Yet, quantifying the confidence surrounding coefficient estimates and selected covariates is essential in areas of application where these will nourrish further targeted investigations.

In this paper we consider the two-sample linear regression model with Gaussian random design.

Y (1) = X (1) β (1) + (1) (1)

Y (2) = X (2) β (2) + (2) . (2)
In this statistical model, the size p row vectors X (1) and X (2) follow Gaussian distributions N (0 p , Σ (1) ) and N (0 p , Σ (2) ) whose covariance matrices remain unknown. The noise components (1) and (2) are independent from the design matrices and follow a centered Gaussian distribution with unknown standard deviations σ (1) and σ (2) . In this formal setting, our objective is to develop a test for the equality of β (1) and β (2) which remains valid in high-dimension. Suppose that we observe an n 1 -sample of (Y (1) , X (1) ) and an n 2 -sample of (Y (2) , X (2) ) noted Y (1) , X (1) , and Y (2) , X (2) , with n 1 and n 2 remaining smaller than p. Defining analogously (1) and (2) , we obtain the decompositions Y (1) = X (1) β (1) + (1) and Y (2) = X (2) β (2) + (2) . Given these observations, we want to test whether models (1) and ( 2) are the same, that is H 0 : β (1) = β (2) , σ (1) = σ (2) , and Σ (1) = Σ (2) H 1 : β (1) = β (2) or σ (1) = σ (2) .

(3)

In the null hypothesis, we include the assumption that the population covariances of the covariates are equal (Σ (1) = Σ (2) ), while under the alternative hypothesis the population covariances are not required to be the same. This choice of assumptions is primarely motivated by our final objective to derive homogeneity tests for Gaussian graphical models (see below). A discussion of the design hypotheses is deferred to Section 7.

Connection with two-sample Gaussian graphical model testing

This testing framework is mainly motivated by the validation of differences observed between Gaussian graphical models (modelling regulation networks) inferred from transcriptomic data from two samples [START_REF] Chiquet | Inferring multiple graphical structures[END_REF][START_REF] Friedman | Sparse inverse covariance estimation with the graphical Lasso[END_REF][START_REF] Meinshausen | High-dimensional graphs and variable selection with the Lasso[END_REF] when looking for new potential drug or knock-out targets [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF]. Following the development of univariate differential analysis techniques, there is now a surging demand for the detection of differential regulations between pairs of conditions (treated vs. placebo, diseased vs. healthy, exposed vs. control, . . . ). Given two gene regulation networks inferred from two transcriptomic data samples, it is however difficult to disentangle differences in the estimated networks that are due to estimation errors from real differences in the true underlying networks.

We suggest to build upon our two-sample high-dimensional linear regression testing scheme to derive a global test for the equality of high-dimensional Gaussian graphical models inferred under pairs of conditions.

Formally speaking, the global two-sample GGM testing problem is defined as follows. Consider two Gaussian random vectors Z (1) ∼ N (0, [Ω (1) ] -1 ) and Z (2) ∼ N (0, [Ω (2) ] -1 ). The dependency graphs are characterized by the non-zero entries of the precision matrices Ω (1) and Ω (2) [START_REF] Lauritzen | Graphical Models[END_REF]. Given an n 1 -sample of Z (1) and an n 2 -sample of Z (2) , the objective is to test H G 0 : Ω (1) = Ω (2) versus H G 1 : Ω (1) = Ω (2) ,

where Ω (1) and Ω (2) are assumed to be sparse (most of their entries are zero). This testing problem is global as the objective is to assess a statistically significant difference between the two distributions. If the test is rejected, a more ambitious objective is to infer the entries where the precision matrices differ (ie Ω

i,j = Ω

i,j ). Adopting a neighborhood selection approach [START_REF] Meinshausen | High-dimensional graphs and variable selection with the Lasso[END_REF] as recalled in Section 6, high-dimensional GGM estimation can be solved by multiple high-dimensional linear regressions. As such, two-sample GGM testing (4) can be solved via multiple two-sample hypothesis testing as [START_REF] Arias-Castro | Global testing under sparse alternatives: Anova, multiple comparisons and the higher criticism[END_REF] in the usual linear regression framework. This extension of two-sample linear regression tests to GGMs is described in Section 6.

Related work

The literature on high-dimensional two-sample tests is very light. In the context of highdimensional two-sample comparison of means, [START_REF] Bai | Effect of high dimension: by an example of a two sample problem[END_REF][START_REF] Chen | A two-sample test for high-dimensional data with applications to gene-set testing[END_REF][START_REF] Lopes | A more powerful two-sample test in high dimensions using random projection[END_REF][START_REF] Srivastava | A test for the mean vector with fewer observations than the dimension[END_REF] have introduced global tests to compare the means of two high-dimensional Gaussian vectors with unknown variance. Recently, [START_REF] Cai | Two-sample covariance matrix testing and support recovery in high-dimensional and sparse settings[END_REF][START_REF] Li | Two sample tests for high-dimensional covariance matrices[END_REF] developped two-sample tests for covariance matrices of two high-dimensional vectors.

In contrast, the one-sample analog of our problem has recently attracted a lot of attention, offering as many theoretical bases for extension to the two-sample problem. In fact, the high-dimensional linear regression tests for the nullity of coefficients can be interpreted as a limit of the two-sample test in the case where β (2) is known to be zero, and the sample size n 2 is considered infinite so that we perfectly know the distribution of the second sample.

There are basically two different objectives in high-dimensional linear testing: a local and a global approach. In the local approach, one considers the p tests for the nullity of each coefficient H 0,i : β

(1) i = 0 (i = 1, . . . , p) with the purpose of controling error measures such as the false discovery rate of the resulting multiple testing procedures. In a way, one aims to assess the individual statistical significance of each of the variables. This can be achieved by providing a confidence region for β (1) [START_REF] Bühlmann | Statistical significance in high-dimensional linear models[END_REF][START_REF] Javanmard | Confidence Intervals and Hypothesis Testing for High-Dimensional Regression[END_REF][START_REF] Lockhart | A significance test for the lasso[END_REF][START_REF] Meinshausen | Assumption-free confidence intervals for groups of variables in sparse high-dimensional regression[END_REF][START_REF] Zhang | Confidence intervals for low-dimensional parameters with high-dimensional data[END_REF]. Another line of work derives p-values for the nullity of each of the coefficients. Namely, [START_REF] Wasserman | High dimensional variable selection[END_REF] suggests a screen and clean procedure based upon half-sampling. Model selection is first applied upon a random half of the sample in order to test for the significance of each coefficient using the usual combination of ordinary least squares and Student t-tests on a model of reasonnable size on the remaining second half. To reduce the dependency of the results to the splitting, [START_REF] Meinshausen | P-values for highdimensional regression[END_REF] advocate to use half-sampling B times and then aggregate the B p-values obtained for variable j in a way which controls either the family-wise error rate or false discovery rate.

In the global approach, the objective is to test the null hypothesis H 0 : β (1) = 0. Although global approaches are clearly less informative than approaches providing individual significance tests like [START_REF] Bühlmann | Statistical significance in high-dimensional linear models[END_REF][START_REF] Meinshausen | P-values for highdimensional regression[END_REF][START_REF] Zhang | Confidence intervals for low-dimensional parameters with high-dimensional data[END_REF], they can reach better performances from smaller sample sizes. Such a property is of tremendous importance when dealing with high-dimensional datasets. The idea of [START_REF] Verzelen | Goodness-of-fit tests for high-dimensonal gaussian linear models[END_REF], based upon the work of [START_REF] Baraud | Adaptive tests of linear hypotheses by model selection[END_REF], is to approximate the alternative H 1 : β (1) = 0 by a collection of tractable alternatives {H S 1 : ∃j ∈ S , β

= 0, S ∈ S} working on subsets S ⊂ {1, . . . , p} of reasonable sizes. The null hypothesis is rejected if the null hypothesis H S 0 is rejected for at least one of the subsets S ∈ S. Admittedly, the resulting procedure is computationally intensive. Nonetheless it is non-asymptotically minimax adaptive to the unknown sparsity of β (1) , that is it achieves the optimal rate of detection without any assumption on the population covariance Σ (1) of the covariates. Another series of work relies on higher-criticism. This last testing framework was originally introduced in orthonormal designs [START_REF] Donoho | Higher criticism for detecting sparse heterogeneous mixture[END_REF], but has been proved to reach optimal detection rates in high-dimensional linear regression as well [START_REF] Arias-Castro | Global testing under sparse alternatives: Anova, multiple comparisons and the higher criticism[END_REF][START_REF] Ingster | Detection boundary in sparse regression[END_REF]. In the end, highercriticism is highly competitive in terms of computing time and achieves the asymptotic rate of detection with the optimal constants. However, these nice properties require strong assumptions on the design.

While writing this paper, we came across the parallel work of Städler and Mukherjee [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF], which adopts a local approach in an elegant adaptation of the screen and clean procedure in its simple-split [START_REF] Wasserman | High dimensional variable selection[END_REF] and multi-split [START_REF] Meinshausen | P-values for highdimensional regression[END_REF] versions to the two-sample framework. Interestingly, their work also led to an extension to GGM testing [START_REF] Städler | Network-based multivariate gene-set testing[END_REF].

In contrast we build our testing strategy upon the global approach developed by [START_REF] Baraud | Adaptive tests of linear hypotheses by model selection[END_REF] and [START_REF] Verzelen | Goodness-of-fit tests for high-dimensonal gaussian linear models[END_REF]. A more detailed comparison of [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF][START_REF] Städler | Network-based multivariate gene-set testing[END_REF] with our contribution is deferred to simulations (Section 5) and discussion (Section 7).

Our contribution

Our suggested approach stems from the fundamental assumption that either the true supports of β (1) and β (2) are sparse or that their difference β (1) β (2) is sparse, so that the test can be successfully led in a subset S ⊂ {1, . . . , p} of variables with reasonnable size, compared to the sample sizes n 1 and n 2 . Of course, this low dimensional subset S is unknown. The whole objective of the testing strategy is to achieve similar rates of detection (up to a logarithmic constant) as an oracle test which would know in advance the optimal low-dimensional subset S .

Concretely, we proceed in three steps :

1. We define algorithms to select a data-driven collection of subsets S identified as most informative for our testing problem, in an attempt to circumvent the optimal subset S . 2. New parametric statistics related to the likelihood ratio statistic between the conditional distributions Y (1) |X

S and Y (2) |X

S are defined for S ∈ S. 3. We define two calibration procedures which both guarantee a control on type-I error:

• we use a Bonferroni calibration which is both computationally and conceptually simple, allowing us to prove that this procedure is minimax adaptive to the sparsity of β (1) and β (2) from a non-asymptotic point of view;

• we define a calibration procedure based upon permutations to reach a fine tuning of multiple testing calibration in practice, for an increase in empirical power.

The resulting testing procedure is completely data-driven and its type I error is explicitely controlled. Furthermore, it is computationally amenable in a large p and small n setting. Interestingly, the procedure does not require any half-sampling steps which are known to decrease the robustness when the sample size is small. The procedure is described in Section 2 while Section 3 is devoted to technical details, among which theoretical controls on Type I error, as well as some useful empirical tools for interpretation. Section 4 provides a non-asymptotic control of the power. Section 5 provides simulated experiments comparing the performances of the suggested procedures with the approach of [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF]. In Section 6, we detail the extension of the procedure to handle the comparison of Gaussian graphical models. The method is illustrated on Transcriptomic Breast Cancer Data. Finally, all the proofs are postponed to Section 8.

The R codes of our algorithms are available at [START_REF]R code[END_REF].

Notation

In the sequel, p norms are denoted | • | p , except for the l 2 norm which is referred as . to alleviate notations. For any positive definite matrix Σ, . Σ denotes the Euclidean norm associated with the scalar product induced by Σ: for every vector x, x 2 Σ = x Σx. Besides, fo every set S, |S| denote its cardinality. For any integer k, I k stands for the identity matrix of size k. For any square matrix A, ϕ max (A) and ϕ min (A) denote respectively the maximum and minimum eigenvalues of A. When the context makes it obvious, we may omit to mention A to alleviate notations and use ϕ max and ϕ min instead. Moreover, Y refers to the size n 1 + n 2 concatenation of Y (1) and Y (2) and X refers to the size (n 1 + n 2 ) × p the concatenation of X (1) and X (2) . To finish with, L refers to a positive numerical constant that may vary from line to line.

Description of the testing strategy

Likelihood ratio statistics used to test hypotheses like H 0 in the classical large n, small p setting are intractable on high-dimensional datasets for the mere reason that the maximum likelihood estimator is not itself defined under high-dimensional design proportions.

Our approach approximates the intractable high-dimensional test by a multiple testing construction, similarly to the strategy developped by [START_REF] Baraud | Adaptive tests of linear hypotheses by model selection[END_REF] in order to derive statistical tests against non-parametric alternatives and adapted to one sample tests for high-dimensional linear regression in [START_REF] Verzelen | Goodness-of-fit tests for high-dimensonal gaussian linear models[END_REF].

For any subset S of {1, . . . , p} satisfying 2|S| ≤ n 1 ∧ n 2 , denote X

(1) S and X

(2) S the restrictions of the random vectors X (1) and X (2) to covariates indexed by S. Their covariance structure is noted Σ

S β

(1)

S + (1) S Y (2) = X (2) S β (2)
S +

(2)

S ,
where the noise variables

S and

(2)

S are independent from X

(1) S and X

(2) S and follow centered Gaussian distributions with new unkwown conditional standard deviations σ

(1) S and σ

(2) S . We now state the test hypotheses in reduced dimension:

H 0,S : β (1) S = β (2) S , σ (1)
S = σ (2)

S , and Σ

(1)

S = Σ (2) S , H 1,S : β (1) S = β (2) S or σ (1) S = σ (2) S .
Of course, there is no reason in general for β S to coincide with the restrictions of β (1) and β (2) to S, even less in high-dimension since variables in S can be in all likelihood correlated with covariates in S c . Yet, as exhibited by Lemma 2.1, there is still a strong link between the collection of low dimensional hypotheses H 0,S and the global null hypothesis H 0 .

Lemma 2.1. The hypothesis H 0 implies H 0,S for any subset S ⊂ {1, . . . p}.

Proof. Under H 0 , the random vectors of size p + 1 (Y (1) , X (1) ) and (Y (2) , X (2) ) follow the same distribution. Hence, for any subset S, Y (1) follows conditionally on X

(1) S the same distribution as Y (2) conditionnally on X

(2) S . In other words, β

S = β (2) S . (1)
By contraposition, it suffices to reject at least one of the H 0,S hypotheses to reject the global null hypothesis. This fundamental observation motivates our testing procedure. As summarized in Algorithm 1, the idea is to build a well-calibrated multiple testing procedure that considers the testing problems H 0,S against H 1,S for a collection of subsets S. Obviously, it would be prohibitive in terms of algorithmic complexity to test H 0,S for every S ⊂ {1, . . . , p}, since there would be 2 p such sets. As a result, we restrain ourselves to a relatively small collection of hypotheses {H 0,S , S ∈ S}, where the collection of supports S is potentially data-driven. If the collection S is judiciously selected, then we can manage not to lose too much power compared to the exhaustive search.

We now turn to the description of the three major elements required by our overall strategy (see Algorithm 1):

1. a well-targeted data-driven collection of models S as produced by Algorithm 2; 2. a parametric statistic to test the hypotheses H 0,S for S ∈ S, we resort actually to a combination of three parametric statistics F S,V , F S,1 and F S,2 ; 3. a calibration procedure guaranteeing the control on type I error as in Algorithm 3 or 4.

Algorithm 1 Overall Adaptive Testing Strategy

Require: Data X (1) ,X (2) ,Y (1) ,Y (2) , collection S and desired level α

Step 1 -Choose a collection S of low-dimensional models (as e.g. S Lasso in Algorithm 2) procedure ModelChoice(X (1) ,X (2) ,Y (1) ,Y (2) ,S) Define the model collection S ⊂ S end procedure

Step 2 -Compute p-values for each test in low dimension procedure Test(X

S ,X

S ,Y (1) ,Y (2) , S) for each subset S in S do Compute the p-values qV,S , q1,S , q2,S associated to the statistics F S,V , F S,1 , F S,2 end for end procedure

Step 3 -Calibrate decision thresholds as in Algorithms 3 (Bonferroni) or 4 (Permutations)

procedure Calibration(X (1) ,X (2) ,Y (1) ,Y (2) , S,α) for each subset S in S and each i = V, 1, 2 do Define a threshold α i,S . end for end procedure

Step 4 -Final Decision if there is a least one model S in S such that there is at least one p-value for which qi,S < α i,S then Reject the global null hypothesis H 0 end if

Choices of Test Collections (Step 1)

The first step of our procedure (Algorithm 1) amounts to select a collection S of subsets of {1, . . . , p}, also called models. A good collection S of subsets must satisfy a tradeoff between the inclusion of the maximum number of relevant subsets S and a reasonable computing time for the whole testing procedure, which is linear in the size | S| of the collection. The construction of S proceeds in two steps: (i) One chooses a deterministic collection S of models. (ii) One defines an algorithm (called ModelChoice in Algorithm 1) mapping the raw data (X, Y) to some collection S satisfying S ⊂ S. Even though the introduction of S as an argument of the mapping could appear artificial at this point, this quantity will be used in the calibration step of the procedure. Our methodology can be applied to any fixed or data-driven collection. Still, we focus here on two particular collections. The first one is useful for undertaking the first steps of the mathematical analysis. For practical applications, we advise to use the second collection.

Deterministic Collections S ≤k . By deterministic, we mean the model choice step is trivial in the sense ModelChoice(X, Y, S) = S. Among deterministic collections, the most straightforward collections consist of all size-k subsets of {1, . . . , p}, which we denote S k . This kind of family provides collections which are independent from the data, thereby reducing the risk of overfitting. However, as we allow the model size k or the total number of candidate variables p to grow, these deterministic families can rapidly reach unreasasonble sizes. Admittedly, S 1 always remains feasible, but reducing the search to models of size 1 can be costly in terms of power. As a variation on size k models, we introduce the collection of all models of size smaller than k, denoted S ≤k = k j=1 S j , which will prove useful in theoretical developments.

Lasso-type Collection S Lasso . Among all data-driven collections, we suggest the Lassotype collection S Lasso . Before proceeding to its definition, let us informally discuss the subsets that a "good" collection S should contain. Let supp(β) denote the support of a vector β. Intuitively, under the alternative hypothesis, good candidates for the subsets are either subsets of S * ∨ := supp(β (1) ) ∪ supp(β (2) ) or subsets of S * ∆ := Supp(β (1) β (2) ).

The first model S * ∨ nicely satisfies β

(1) S * ∨ = β (1) and β

(2) S * ∨ = β (2) . The second subset has a smaller size than S * ∨ and focuses on covariates corresponding to different parameters in the full regression. However, the divergence between effects might only appear conditionally on other variables with similar effects, this is why the first subset S * ∨ is also of interest. Obviously, both subsets S * ∨ and S * ∆ are unknown. This is why we consider a Lasso methodology that amounts to estimating both S * ∨ and S * ∆ in the collection S Lasso . Details on the construction of S Lasso are postponed to Section 3.1.

Parametric Test Statistic (Step 2)

Given a subset S, we consider the three following statistics to test H 0,S against H 1,S :

F S,V := -2 + Y (1) -X (1) β (1) S 2 /n 1 Y (2) -X (2) β (2) S 2 /n 2 + Y (2) -X (2) β (2) S 2 /n 2 Y (1) -X (1) β (1) S 2 /n 1 , (5)
F S,1 := X (2) S ( β (1)
S -β (2) S ) 2 /n 2 Y (1) -X (1) β (1) S 2 /n 1 , F S,2 := X (1) S ( β (1)

Sβ

(2)

S ) 2 /n 1 Y (2) -X (2) β (2) S 2 /n 2 . ( 6
)
As explained in Section 3, these three statistics derive from the Kullback-Leibler divergence between the conditional distributions Y (1) |X

S and

Y (2) |X (2)
S . While the first term F S,V evaluates the discrepancies in terms of conditional variances, the last two terms F S,1 and F S,2 address the comparison of β (1) to β (2) .

Because the distributions under the null of the statistics F S,i , for i = V, 1, 2, depend on the size of S, the only way to calibrate the multiple testing step over a collection of models of various sizes is to convert the statistics to a unique common scale. The most natural is to convert observed F S,i 's into p-values. Under H 0,S , the conditional distributions of F S,i for i = V, 1, 2 to X S are parameter-free and explicit (see Proposition 3.1 in the next section). Consequently, one can define the exact p-values associated to F S,i , conditional on X S . However, the computation of the p-values require a function inversion, which can be computationally prohibitive. This is why we introduce explicit upper bounds qi,S (Equations [START_REF] Chung | Development of an orallyadministrative melk-targeting inhibitor that suppresses the growth of various types of human cancer[END_REF][START_REF] Friedman | Sparse inverse covariance estimation with the graphical Lasso[END_REF]) of the exact p-values.

Combining the parametric statistics (Step 3)

The objective of this subsection is to calibrate a multiple testing procedure based on the sequence of p-values {(q V,S , q1,S , q2,S ), S ∈ S}, so that the type-I error remains smaller than a chosen level α. In particular, when using a data-driven model collection, we must take good care of preventing the risk of overfitting which results from using the same dataset both for model selection and hypothesis testing.

For the sake of simplicity, we assume in the two following paragraphs that ∅ S, which merely means that we do not include in the collection of tests the raw comparison of Var(Y (1) ) to Var(Y (2) ).

Testing Procedure Given a model collection S and a sequence α = (α i,S ) i=V,1,2, S∈ S , we define the test function:

T α S = 1 if ∃S ∈ S, ∃i ∈ {V, 1, 2} qi,S ≤ α i,S . 0 otherwise. (7)
In other words, the test function rejects the global null if there exists at least one model S ∈ S such that at least one of the three p-values is below the corresponding threshold α i,S . In Section 3.3, we describe two calibration methods for choosing the thresholds (α i,S ) S∈ S . We first define a natural Bonferroni procedure, whose conceptual simplicity allows us to derive non-asymptotic type II error bounds of the corresponding tests (Section 4). However, this Bonferroni correction reveals too conservative in practice, in part paying the price for resorting to data-driven collections and upper bounds on the true p-values. This is why we introduce as a second option the permutation calibration procedure. This second procedure controls the type I error at the nominal level and therefore outperforms the Bonferroni calibration in practice. Nevertheless, the mathematical analysis of the corresponding test becomes more intricate and we are not able to provide sharp type II error bounds.

Remark: In practice, we advocate the use of the Lasso Collection S Lasso (Algorithm 2) combined with the permutation calibration method (Algorithm 4). Henceforth, the corresponding procedure is denoted T P SLasso .

Discussion of the procedure and Type I error

In this section, we provide remaining details on the three steps of the testing procedure. First, we describe the collection S Lasso and provide an informal justification of its definition. Second, we explain the ideas underlying the parametric statistics F S,i , i = V, 1, 2 and we define the corresponding p-values qi,S . Finally, the Bonferroni and permutation calibration methods are defined, which allows us to control the type I error of the corresponding testing procedures.

Collection S Lasso

We start from S ≤Dmax , where, in practice, D max = (n 1 ∧ n 2 )/2 and we consider the following reparametrized joint regression model.

Y (1) Y (2) = X (1) X (1) X (2) -X (2) θ (1) * θ (2) * + (1) (2)
.

In this new model, θ

* captures the mean effect (β (1) + β (2) )/2, while θ

(2) * captures the discrepancy between the sample-specific effect β (i) and the mean effect θ (1) * , that is to say θ

(2) * = (β (1) -β (2) )/2. Consequently, S * ∆ := Supp(β (1) -β (2) ) = supp(θ (2)
* ) and S * ∨ := supp(β (1) ) ∪ supp(β (2) ) = supp(θ

(1) * ) ∪ supp(θ (2)
* ). To simplify notations, denote by Y the concatenation of Y (1) and Y (2) , as well as by W the reparametrized design matrix of [START_REF] Bühlmann | Statistical significance in high-dimensional linear models[END_REF]. For a given λ > 0, the Lasso estimator of θ * is defined by

θ λ := θ (1) λ θ (2) λ := arg min θ∈R 2p Y -Wθ + λ θ 1 , (9)
Vλ := supp( θ λ ),

V (i) λ := supp( θ (i) λ ), i = 1, 2 . ( 10
)
For a suitable choice of the tuning parameter λ and under assumptions of the designs, it is proved [START_REF] Bickel | Simultaneous analysis of Lasso and Dantzig selector[END_REF][START_REF] Meinshausen | Lasso-type recovery of sparse representations for high-dimensional data[END_REF] that θ λ estimates well θ * and Vλ is a good estimator of supp(θ * ). The Lasso parameter λ tunes the amount of sparsity of θ λ : the larger the parameter λ, the smaller the support Vλ . As the optimal choice of λ is unknown, the collection S Lasso is built using the collection of all estimators ( Vλ ) λ>0 , also called the Lasso regularization path of θ * . Below we provide an algorithm for computing S Lasso along with some additional justifications.

Algorithm 2 Construction of the Lasso-type Collection S Lasso Require: Data X (1) ,X (2) ,Y (1) ,Y (2) , Collection S ≤Dmax

Y ← Y (1) Y (2)
W ← X (1) X (1) X (2) -X (2) Compute the function f : λ → Vλ (defined in (9,10)) using Lars-Lasso Algorithm [START_REF] Efron | Least angle regression[END_REF] Compute the decreasing sequences (λ k ) 1≤k≤q of jumps in

f k ← 1, S (1)
L ← ∅, S (2)
L ← ∅ while | V (1) λ k ∪ V (2)
λ k | < Dmax do S (1)
L ← S

(1)

L ∪ { V (1) λ k ∪ V (2) λ k } S (2)
L ← S

(2)

L ∪ { V (2) λ k } k ← k + 1 end while S Lasso ← S (1) L ∪ S (2) L ∪ S 1
It is known [START_REF] Efron | Least angle regression[END_REF] that the function f : λ → Vλ is piecewise constant. Consequently, there exist thresholds λ 1 > λ 2 > . . . such that Vλ changes on λ k 's only. The function f and the collection (λ k ) are computed efficiently using the Lars-Lasso Algorithm [START_REF] Efron | Least angle regression[END_REF]. We build two collections of models using ( V (1) λ k ) k≥1 and ( V (2) λ k ) k≥1 . Following the intuition described above, for a fixed λ k , V

λ k is an estimator of supp(β (1) β (2) ) while

V (1) λ k ∪ V (2)
λ k is an estimator of supp(β (1) ) ∪ supp(β (2) ). This is why we define

S (1) L := kmax k=1 V (1) λ k ∪ V (2) λ k , S (2)
L := kmax k=1 V (2) λ k
,

where k max is the smallest integer q such that |

V (1) λq+1 ∪ V (2) λq+1 | > D max .
In the end, we consider the following S Lasso data-driven family,

S Lasso := S (1) L ∪ S (2) L ∪ S 1 . (11)
Recall that S 1 is the collection of the p models of size 1. Recently, data-driven procedures have been proposed to tune the Lasso and find a parameter λ is such a way that θ λ is a good estimator of θ * (see e.g. [START_REF] Baraud | Estimator selection in the gaussian setting[END_REF][START_REF] Sun | Scaled sparse linear regression[END_REF]). We use the whole regularization path instead of the sole estimator θ λ , because our objective is to find subsets S such that the statistics F S,i are powerful. Consider an example where β (2) = 0 and β (1) contains one large coefficient and many small coefficients. If the sample size is large enough, a well-tuned Lasso estimator will select several variables. In contrast, the best subset S (in terms of power of F S,i ) contains only one variable. Using the whole regularization path, we hope to find the best trade-off between sparsity (small size of S) and differences between β

S and β

S . This last remark is formalized in Section 4.4. Finally, the size of the collection S Lasso is generally linear with (n 1 ∧n 2 )∨p, which makes the computation of (q i,S ) S∈ SLasso,i=V,1,2 reasonnable.

Parametric statistics and p-values

Symmetric conditional likelihood

In this subsection, we explain the intuition behind the choice of the parametric statistics (F S,V , F S,1 , F S,2 ) defined in Equations [START_REF] Baraud | Estimator selection in the gaussian setting[END_REF][START_REF] Baraud | Adaptive tests of linear hypotheses by model selection[END_REF]. Let us denote by L (1) (resp. L (2) ) the loglikelihood of the first (resp. second) sample normalized by n 1 (resp. n 2 ). Given a subset S ⊂ {1, . . . , p} of size smaller than n 1 ∧ n 2 , ( β

(1) S , σ (1)
S ) stands for the maximum likelihood estimator of (β (1) , σ (1) ) among vectors β whose supports are included in S. Similarly, we note ( β

(2) S , σ (2)
S ) for the maximum likelihood corresponding to the second sample. Statistics F S,V , F S,1 and F S,2 appear as the decomposition of a two-sample likelihoodratio, measuring the symmetrical adequacy of sample-specific estimators to the opposite sample. To do so, let us define the likelihood ratio in sample i between an arbitrary pair (β, σ) and the corresponding sample-specific estimator ( β

(i) S , σ (i) S ): D (i) ni (β, σ) := L (i) ni β (i) S , σ (i) S -L (i) ni (β, σ) .
With this definition, D

n1 ( β (2) , σ (2) ) measures how far ( β (2) , σ (2) ) is from ( β (1) , σ (1) ) in terms of likelihood within sample 1. The following symmetrized likelihood statistic can be decomposed into the sum of F S,V , F S,1 and F S,2 : 2 D (1) n1 ( β (2) , σ (2) ) + D (2) n2 ( β (1) , σ (1)

) = F S,V + F S,1 + F S,2 . (12)
Instead of the three statistics (F S,i ) i=V,1,2 , one could use the symmetric likelihood [START_REF] Chiquet | Inferring multiple graphical structures[END_REF] to build a testing procedure. However, we do not manage to obtain an explicit and sharp upper bound of the p-values associated to the statistic [START_REF] Chiquet | Inferring multiple graphical structures[END_REF], which makes the resulting procedure either computationally intensive if one estimated the p-values by a Monte-Carlo approach or less powerful if one uses a non-sharp upper bound of the p-values. In contrast, we explain below how, by considering separately F S,V , F S,1 and F S,2 , one upper bounds sharply the exact p-values.

Definition of the p-values

Denote by g(x) = -2 + x + 1/x the non-negative function defined on R + . Since the restriction of g to [1; +∞) is a bijection, we note g -1 the corresponding reciprocal function.

Proposition 3.1 (Conditional distributions of F S,V , F S,1 and F S,2 under H 0,S ).

1. Let Z denote a Fisher random variable with (n 1 -|S|, n 2 -|S|) degrees of freedom. Then, under the null hypothesis,

F S,V |X S ∼ H 0,S g Z n 2 (n 1 -|S|) n 1 (n 2 -|S|) .
2. Let Z 1 and Z 2 be two centered and independent Gaussian vectors with covariance

X (2) S (X (1)T S X (1)
S ) -1 + (X (2)T S X (2)
S ) -1 X

(2)T S and I n1-|S| . Then, under the null hypothesis,

F S,1 |X S ∼ H 0,S Z 1 2 /n 2 Z 2 2 /n 1 .
A symmetric result holds for F S,2 .

Although the distributions identified in Proposition 3.1 are not all familiar distributions with ready-to-use quantile tables, they all share the advantage that they do not depend on any unknown quantity, such as design variances Σ (1) and Σ (2) , noise variances σ (1) and σ (2) , or even true signals β (1) and β (2) . For any i = V, 1, 2, we note Q i,|S| (u|X S ) for the conditional probability that F S,i is larger than u. By Proposition 3.1, the exact p-value qV,S = Q V,|S| (F S,V |X S ) associated to F S,V is easily computed from the distribution function of a Fisher random variable:

qV,S = F n1-|S|,n2-|S| g -1 (F S,V ) n 1 (n 2 -|S|) n 2 (n 1 -|S|) +F n2-|S|,n1-|S| g -1 (F S,V ) n 2 (n 1 -|S|) n 1 (n 2 -|S|) , (13)
where F m,n (u) denotes the probability that a Fisher random variable with (m, n) degrees of freedom is larger than u.

Since the conditional distribution of F S,1 given X S only depends on |S|, n 1 , n 2 , and X S , one could compute an estimation of the p-value Q 1 (u|X S ) associated with an observed value u by Monte-Carlo simulations. However, this approach is computationally prohibitive for large collections of subsets S. This is why we use instead an explicit upper bound of Q 1,|S| (u|X S ) based on Laplace method, as given in the definition below and justified in the proof of Proposition 3.3.

Definition 3.2 (Definition of Q 1,|S| and Q 2,|S| ). Let us note a = (a 1 , . . . , a |S| ) the positive eigenvalues of n 1 n 2 (n 1 -|S|) X (2) S (X (1)T S X (1)
S ) -1 + (X (2)T S X (2) S ) -1 X (2)T S . For any u ≤ |a| 1 , define Q 1,|S| (u|X S ) := 1. For any u > |a| 1 , take Q 1,|S| (u|X S ) := exp   - 1 2 |S| i=1 log (1 -2λ * a i ) - n 1 -|S| 2 log 1 + 2λ * u n 1 -|S|   , (14)
where λ * is defined as follows. If all the components of a are equal, then λ * :=

u-|a|1 2u(|a|∞+ |a| 1 n 1 -|S| ) .
If a is not a constant vector, then we define

b := |a| 1 u |a| ∞ (n 1 -|S|) + u + a 2 |a| ∞ -|a| 1 , ∆ := b 2 - 4u (u -|a| 1 ) (n 1 -|S|)|a| ∞ |a| 1 - a 2 |a| ∞ , (15)
λ * := 1 4u n1-|S| |a| 1 -a 2 |a|∞ b - √ ∆ . (16)
Q 2,|S| is defined analogously by exchanging the role of X

S and X

S .

Proposition 3.3. For any u ≥ 0, and

for i = 1, 2, Q i,|S| (u|X S ) ≤ Q i,|S| (u|X S ).
Finally we define the approximate p-values q1,S and q2,S by

q1,S := Q 1,|S| (F S,1 |X S ) , q2,S := Q 2,|S| (F S,2 |X S ) . (17)
Although we use similar notations for qi,S with i = V, 1, 2, this must not mask the essential difference that q1,S is the exact p-value of F S,1 whereas q1,S and q2,S only are upper-bounds on F S,2 and F S,2 p-values. The consequences of this asymetry in terms of calibration of the test is discussed in the next subsection.

Comparison of the calibration procedures and Type I error

Bonferroni Calibration (B)

Recall that a data-driven model collection S is defined as the result of a fixed algorithm mapping a deterministic collection S and (X, Y) to a subcollection S. The collection of

thresholds α B = α i,S , S ∈ S is chosen such that S∈S i=V,1,2 α i,S ≤ α . ( 18
)
For the collection S ≤k , or any data-driven collection derived from S ≤k , a natural choice is

α V,S := α 2k p |S| -1 , α 1,S = α 2,S := α 4k p |S| -1 , (19)
which puts as much weight to the comparison of the conditional variances (F S,V ) and the comparison of the coefficients (F S,1 , F S,2 ). Similarly for the collection S Lasso , a natural choice is [START_REF] Giraud | Supplement to 'Highdimensional regression with unknown variance[END_REF] with k replaced by D max (which equals (n 1 ∧ n 2 )/2 in practice).

Algorithm 3 Bonferroni Calibration for a collection S ⊂ S ≤Dmax

Require:

maximum model dimension Dmax, model collection S, desired level α for each subset S in S do α V,S ← α(2Dmax) -1 p |S| -1 α 1,S ← α(4Dmax) -1 p |S| -1 , α 2,S ← α 1,S
end for

Given any data-driven collection S, denote by T B S the multiple testing procedure calibrated by Bonferroni thesholds α B (18). Proposition 3.4 (Size of T B S ). The test function T B S satisfies P H0 [T B S = 1] ≤ α. Remark 3.1 (Bonferroni correction on S and not on S). Note that even though we only compute the statistics F S,i for models S ∈ S, the Bonferroni correction (18) must be applied to the initial deterministic collection S including S. Indeed, if we replace the condition [START_REF] Giraud | Graph selection with GGMselect[END_REF] by the condition S∈ S 3 i=1 α i,S ≤ α, then the size of the corresponding is not constrained anymore to be smaller than α. This is due to the fact that we use the same data set to select S ⊂ S and to perform the multiple testing procedure. As a simple example, consider any deterministic collection S and the data-driven collection

S = arg min S∈S min i=V,1,2 qi,S ,
meaning that S only contains the subset S that minimizes the p-values of the parametric tests. Thus, computing T B S for this particular collection S is equivalent to performing a multiple testing procedure on S.

Although procedure T B

S is computationally and conceptually simple, the size of the corresponding test can be much lower than α because of three difficulties:

1. Independently from our problem, Bonferroni corrections are known to be too conservative under dependence of the test statistics.

2. As emphasized by Remark 3.1, whereas the Bonferroni correction needs to be based on the whole collection S, only the subsets S ∈ S are considered. Provided we could afford the computational cost of testing all subsets within S, this loss cannot be compensated for if we use the Bonferroni correction. 3. As underlined in the above subsection, for computational reasons we do not consider the exact p-values of F S,1 and F S,2 but only upper bounds q1,S and q2,S of them. We therefore overestimate the type I error due to F S,1 and F S,2 .

In fact, the three aforementionned issues are addressed by the permutation approach.

Calibration by permutation (P).

The collection of thresholds α P = {α i,S , S ∈ S} is chosen such that each α i,S remains inversely proportional to p |S| in order to put all subset sizes at equal footage. In other words, we choose a collection of thresholds of the form

α i,S = C i p |S| -1 , (20)
where C i 's are calibrated by permutation to control the type I error of the global test.

Given a permutation π of the set {1, . . . , n 1 +n 2 }, one gets Y π and X π by permuting the components of Y and the rows of X. This allows us to get a new sample (Y π, (1) , Y π, (2) , X π, (1) , X π, (2) ). Using this new sample, we compute a new collection S π , parametric statistics (F π S,i ) i=V,1,2 and p-values (q i,S ) i=V,1,2 . Denote P the uniform distribution over the permutations of size n 1 + n 2 .

We define C V as the α/2-quantiles with respect to P of

min S∈ S π qV,S p |S| . (21)
Similarly, C 1 = C 2 are the α/2-quantiles with respect to P of min

S∈ S π (q 1,S ∧ q2,S ) p |S| . ( 22
)
In practice, the quantiles C i are estimated by sampling a large number B of permutations. The permutation calibration procedure for the Lasso collection S Lasso is summarized in Algorithm 4.

Given any data-driven collection S, denote by T P S the multiple testing procedure calibrated by the permutation method [START_REF] Heidel | Potent sirna inhibitors of ribonucleotide reductase subunit rrm2 reduce cell proliferation in vitro and in vivo[END_REF]. Remark 3.2. Through the three constants C V , C 1 and C 2 (Eq. [START_REF] Hess | Pharmacogenomic predictor of sensitivity to preoperative chemotherapy with paclitaxel and fluorouracil, doxorubicin, and cyclophosphamide in breast cancer[END_REF][START_REF] Ingster | Detection boundary in sparse regression[END_REF]), the permutation approach corrects simultaneously for the losses mentioned earlier due to the Bonferroni correction, in particular the restriction to a data-driven class S and the approximate pvalues q1,S and q2,S . Yet, the level of T P S is not exactly α because we treat separately the the statistics F S,V and (F S,1 , F S,2 ) and apply a Bonferroni correction. It would be possible to calibrate all the statistics simultaneously in order to constrain the size of the corresponding test to be exactly α. However, this last approach would favor the statistic F S,1 too much, because we would put on the same level the exact p-value qV,S and the upper bounds q1,S and q2,S .

Algorithm 4 Calibration by Permutation for S Lasso

Require: Data X (1) ,X (2) ,Y (1) ,Y (2) , maximum model dimension Dmax, number B of permutations, desired level α for b = 1,. . . , B do Draw π a random permutation of {1, . . . ,

n 1 + n 2 } X (b) , Y (b) ← π-permutation of (X, Y)
procedure LassoModelChoice(X (1,b) ,X (2,b) ,Y (1,b) ,Y (2,b) ,S ≤Dmax ) Define S (b) Lasso (as in Algorithm 2) end procedure procedure Test(X (1,b) ,X (2,b) ,Y (1,b) ,Y (2,b) ,

S (b) Lasso ) for each subset S in S (b) Lasso do Compute the p-values q(b) i,S for i = V, 1, 2. end for C (b) V ← min S∈ S (b) Lasso q(b) V,S p |S| C (b) 1 ← min S∈ S (b) Lasso q(b) 1,S ∧ q(b) 2,S p |S| end procedure end for Define C V as the α/2-quantile of the (C (1) V , . . . , C (B) V ) distribution Define C 1 = C 2 as the α/2-quantile of the (C (1) 1 , . . . , C (B) 1 ) distribution for each subset S in S Lasso , each i = V, 1, 2, do α i,S ← C i p |S| -1
end for

Interpretation tools

Empirical p-value When using a calibration by permutations, one can derive an empirical p-value p empirical to assess the global significance of the test. In contrast with model and statistic specific p-values qi,S , this p-value provides a nominally accurate estimation of the type-I error rate associated with the global multiple testing procedure, every model in the collection and test statistic being considered. It can be directly compared to the desired level α to decide about the rejection or not of the global null hypothesis. This empirical p-value is obtained as the fraction of the permuted values of the statistic that are less than the observed test statistic. Keeping the notation of Algorithm 4, the empirical p-value for the variance and coefficient parts are given respectively by :

p empirical V = 1 B B b=1 1 C (b) V < min S∈ SLasso qV,S p |S| , p empirical 1-2 = 1 B B b=1 1 C (b) 1 < min S∈ SLasso (q 1,S ∧ q2,S ) p |S| .
The empirical p-value for the global test is then given by the following equation.

p empirical = 2 min(p empirical V , p empirical 1-2 ). ( 23
)
Rejected model Moreover, one can keep track of the model responsible for the rejection, unveiling sensible information on which particular coefficients most likely differ between samples. The rejected models for the variance and coefficient parts are given respectively by :

S R V = arg min S∈ SLasso qV,S p |S| S R 1-2 = arg min S∈ SLasso (q 1,S ∧ q2,S ) p |S|
We define the rejected model S R as model S R V or S R 1-2 according to the smallest empirical p-value p empirical V or p empirical 1-2 .

Power and Adaptation to Sparsity

Let us fix some number δ ∈ (0, 1). The objective is to investigate the set of parameters (β (1) , σ (1) , β (2) , σ (2) ) that enforce the power of the test to exceed 1δ. We focus here on the Bonferroni calibration (B) procedure because the analysis is easier. Section 5 will illustrate that the permutation calibration (P) outperforms the Bonferroni calibration (B) in practice. In the sequel, A B (resp. A B) means that for some positive constant L(α, δ) that only depends on α and δ, A ≤ L(α, δ)B (resp. A ≥ L(α, δ)B).

We first define the symmetrized Kullback-Leibler divergence as a way to measure the discrepancies between (β (1) , σ (1) ) and (β (2) , σ (2) ). Then, we consider tests with deterministic collections in Sections 4.2-4.3. We prove that the corresponding tests are minimax adaptive to the sparsity of the parameters or to the sparsity of the difference β (1)β (2) . Sections 4.4-4.5 are devoted to the analysis T B SLasso . Under stronger assumptions on the population covariances than for deterministic collections, we prove that the performances of T B SLasso are nearly optimal.

Symmetrized Kullback-Leibler divergence

Intuitively, the test T B S should reject H 0 with large probability when (β (1) , σ (1) ) is far from (β (2) , σ (2) ) in some sense. A classical way of measuring the divergence between two distributions is the Kullback-Leibler discrepancy. In the sequel, we note K P Y (1) |X ; P Y (2) |X the Kullback discrepancy between the conditional distribution of Y (1) given X (1) = X and conditional distribution of Y (2) given X (2) = X. Then, we denote K 1 the expectation of this Kullback divergence when X ∼ N (0 p , Σ (1) ). Exchanging the roles of Σ (1) and Σ (2) , we also define K 2 :

K 1 := E X (1) K P Y (1) |X ; P Y (2) |X , K 2 := E X (2) K P Y (2) |X ; P Y (1)
|X .

The sum K 1 + K 2 forms a semidistance with respect to (β (1) , σ (1) ) and (β (2) , σ (2) ) as proved by the following decomposition

2 (K 1 + K 2 ) = σ (1) σ (2) 2 + σ (2) σ (1) 2 -2 + β (2) -β (1) 2 Σ (2) (σ (1) ) 2 + β (2) -β (1) 2 Σ (1) (σ (2) ) 2 .
When Σ (1) = Σ (2) , we quantify the discrepancy between these covariance matrices by 1) .

ϕ Σ (1) ,Σ (2) := ϕ max Σ (2) (Σ (1) ) -1 Σ (2) + Σ (1) (Σ (2) ) -1 Σ (
Observe that the quantity ϕ Σ (1) ,Σ (2) can be considered as a constant if we assume that the smallest and largest eigenvalues of Σ (i) are bounded away from zero and infinity.

Power of T B S ≤k

First, we control the power of T B S for a deterministic collection S = S ≤k (with some k ≤ (n 1 ∧ n 2 )/2) and the Bonferroni calibration weights α i,S as in [START_REF] Giraud | Supplement to 'Highdimensional regression with unknown variance[END_REF]). For any β ∈ R p , |β| 0 refers to the size of its support and |β| stands for the vector (|β i |), i = 1, . . . , p. We consider the two following assumptions

A.1 : log(1/(αδ)) n 1 ∧n 2 .
A.2 :

|β (1) | 0 + |β (2) | 0 k ∧ n 1 ∧ n 2 log(p) , log(p) ≤ n 1 ∧ n 2 .
Remark 4.1. Condition A.1 requires that the type I and type II errors under consideration are not exponentially smaller than the sample size. Condition A.2 tells us that the number of non-zero components of β (1) and β (2) has to be smaller than (n 1 ∧ n 2 )/ log(p). This requirement has been shown [START_REF] Verzelen | Minimax risks for sparse regressions: Ultra-high-dimensional phenomenons[END_REF] to be minimal to obtain fast rates of testing of the form [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF] in the specific case β (2) = 0, σ (1) = σ (2) and n 2 = ∞.

K 1 + K 2 ϕ Σ (1) ,Σ (2)
|β (1)

| 0 ∨ |β (2) | 0 ∨ 1 log (p) + log 1 αδ n 1 ∧ n 2 . ( 24
)
If we further assume that

Σ (1) = Σ (2) := Σ, then P[T B S ≤k = 1] ≥ 1 -δ as long as β (1) -β (2) 2 Σ
Var[Y (1) ] ∧ Var[Y (2) ]

|β (1) -β (2) | 0 log (p) + log 1 αδ n 1 ∧ n 2 . ( 25
)
Remark 4.2. The condition Σ (1) = Σ (2) is not necessary to control the power of T B S ≤k in terms of |β (1) β (2) | 0 as in [START_REF] Laurent | Adaptive estimation of a quadratic functional by model selection[END_REF]. However, the expression (25) would become far more involved.

Remark 4.3. Before assessing the optimality of Theorem 4.1, let us briefly compare the two rates of detection [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF] and [START_REF] Laurent | Adaptive estimation of a quadratic functional by model selection[END_REF]. According to [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF], T B S ≤k is powerful as soon as the symmetrized Kullback distance is large compared to {|β (1) | 0 ∨ |β (2) | 0 } log (p) /(n 1 ∧ n 2 ). In contrast, [START_REF] Laurent | Adaptive estimation of a quadratic functional by model selection[END_REF] tells us that T B S ≤k is powerful when β (1)β (2) 2 Σ /(Var[Y (1) ] ∧ Var[Y (2) ]) is large compared to the sparsity of the difference:

|β (1) -β (2) | 0 log (p) /(n 1 ∧ n 2 ).
When β (1) and β (2) have many non-zero coefficients in common, |β (1) β (2) | 0 is much smaller than |β (1) | 0 ∨ |β (2) | 0 . Furthermore, the left-hand side of (25) is of the same order as K 1 + K 2 when Σ (1) = Σ (2) , σ (1) = σ (2) and β (i) Σ /σ (i) 1 for i = 1, 2, that is when the conditional variances are equal and when the signals β (i) Σ are at most at the same order as the noises levels σ (i) . In such a case, (25) outperforms [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF] and only the sparsity of the difference β (1)β (2) plays a role in the detection rates. Below, we prove that [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF] and (25) are both optimal from a minimax point of view but on different sets. Proposition 4.2 (Minimax lower bounds). Assume that p ≥ 5, Σ (1) = Σ (2) = I p , fix some γ > 0, and fix (α, δ) such that α + δ < 53%. There exist two constants L(α, δ, γ) and L (α, δ, γ) such that the following holds.

• For all 1 ≤ s ≤ p 1/2-γ no level-α test has a power larger than 1δ simultaneously over all s-sparse vectors (β (1) , β (2) ) satisfying A.2 and

K 1 + K 2 ≥ L(α, δ, γ) s n 1 ∧ n 2 log (p) . (26)
• For all 1 ≤ s ≤ p 1/2-γ , no level-α test has a power larger than 1δ simultaneously over all sparse vectors (β (1) , β (2) ) satisfying A.2, |β (1) β (2) | 0 ≤ s and

β (1) -β (2) 2

Ip

Var[Y (1) ] ∧ Var[Y (2) ]

≥ L (α, δ, γ) s n 1 ∧ n 2 log (p) . (27)
The proof (in Section 8) is a straightforward application of minimax lower bounds obtained for the one-sample testing problem [START_REF] Arias-Castro | Global testing under sparse alternatives: Anova, multiple comparisons and the higher criticism[END_REF][START_REF] Verzelen | Goodness-of-fit tests for high-dimensonal gaussian linear models[END_REF]. Remark 4.4. Equation [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF] together with [START_REF] Lauritzen | Graphical Models[END_REF] tell us that T B S ≤k simultaneously achieves (up to a constant) the optimal rates of detection over s-sparse vectors β (1) and β (2) for all

s k ∧ p 1/2-γ ∧ n 1 ∧ n 2 log(p) ,
for any γ > 0. Nevertheless, we only managed to prove the minimax lower bound for Σ (1) = Σ (2) = I p , implying that, even though the detection rate ( 24) is unimprovable uniformly over all (Σ (1) , Σ (2) ), some improvement is perhaps possible for specific covariance matrices. Up to our knowledge, there exist no such results of adaptation to the population covariance of the design even in the one sample problem.

Remark 4.5. Equation (25) together with [START_REF] Li | Two sample tests for high-dimensional covariance matrices[END_REF] tells us that T B S ≤k simultaneously achieves (up to a constant) the optimal rates of detection over s-sparse differences β (1) β (2) satisfying

β (1) Σ σ (1) ∨ β (2) Σ σ (2) ≤ 1 for all s k ∧ p 1/2-γ ∧ n1∧n2 log(p) .
Remark 4.6 (Informal justification of the introduction of the collection S Lasso ). If we look at the proof of Theorem 4.1, we observe that the power (24) is achieved by the statistics (F S∨,V , F S∨,1 , F S∨,2 ) where S ∨ is the union of the support of β (1) and β (2) . In contrast, (25) is achieved by the statistics (F S∆,V , F S∆,1 , F S∆,2 ) where S ∆ is the support of β (1)β (2) . Intuitively, the idea underlying the collection S

L in the definition (11) of S Lasso is to estimate S ∨ , while the idea underlying the collection S L is to estimate S ∆ .

Power of T B

S for any deterministic S Theorem 4.3 belows extends Theorem 4.1 from deterministic collections of the form S ≤k to any deterministic collection S, unveiling a bias/variance-like trade-off linked to the cardinality of subsets S of collection S. To do so, we need to consider the Kullback discrepancy between the conditional distribution of Y (1) given X

S = X S and the conditional distribution of Y (2) given X

(2) S = X S , which we denote K P Y (1) |X S ; P Y (2) |X S . For short, we respectively note K 1 (S) and K 2 (S)

K 1 (S) := E X (1) S K P Y (1) |X S ; P Y (2) |X S , K 2 (S) := E X (2) S K P Y (2) |X S ; P Y (1) |X S .
Intuitively, K 1 (S) + K 2 (S) corresponds to some distance between the regression of Y (1) given X

S and of Y (2) given X

S . Noting Σ

(1)

S (resp. Σ (2)
S ) the restriction of Σ (1) (resp. Σ (2) ) to indices in S, we define

ϕ S := ϕ max Σ (2) S (Σ (1)
S ) -1 Σ

(2)

S + Σ (1)
S (Σ (2) S ) -1 Σ (1) S . ( 28
)
Theorem 4.3 (Power of T B S for any deterministic S). For any S ∈ S, we note α S = min i=V,1,2 α i,S . The power of T B S is larger than 1δ as long as there exists S ∈ S such that |S| n 1 ∧ n 2 and

1 + log[1/(δα S )] n 1 ∧ n 2 , (29)
and

K 1 (S) + K 2 (S) ϕ S 1 n 1 + 1 n 2 |S| + log 1 α S δ . ( 30
)
Remark 4.7. Let us note ∆(S) the right hand side of [START_REF] Lopes | A more powerful two-sample test in high dimensions using random projection[END_REF]. According to Theorem 4.3, The term ∆(S) plays the role of a variance term and therefore increases with the cardinality of S. Furthermore, the term

K 1 -K 1 (S) + K 2 -K 2 (S)
plays the role of a bias. Let us note S * the subcollection of S made of sets S satisfying [START_REF] Lockhart | A significance test for the lasso[END_REF]. According to theorem 4.3, T B S is powerful as long as

K 1 + K 2 is larger (up to constants) to inf S∈S * {K 1 -K 1 (S) + K 2 -K 2 (S)} + ∆(S) (31)
Such a result is comparable to oracle inequalities obtained in estimation since the test T B S is powerful when the Kullback loss K 1 + K 2 is larger than the trade-off (31) between a bias-like term and a variance-like term without requiring the knowledge of this trade-off in advance. We refer to [START_REF] Baraud | Adaptive tests of linear hypotheses by model selection[END_REF] for a thorough comparison between oracle inequalities in model selection and second type error terms of this form.

Power of T B

SLasso

For the sake of simplicity, we restrict in this subsection to the case n 1 = n 2 := n, more general results being postponed to the next subsection. The test T B S ≤n/2 is computationally expensive (non polynomial with respect to p). The collection S Lasso has been introduced to fix this burden. We consider T B SLasso with the prescribed Bonferroni calibration weights α i,S (as in [START_REF] Giraud | Supplement to 'Highdimensional regression with unknown variance[END_REF] with k replaced by (n 1 ∧ n 2 )/2 . In the statements below, ψ

Σ (1) ,Σ (2) , ψ (1)
Σ (1) ,Σ (2) ,. . . refer to positive quantities that only depend on the largest and the smallest eigenvalues of Σ (1) and Σ (2) . Consider the additional assumptions

A.3 : |β (1) | 0 ∨|β (2) | 0 ψ (1)
Σ (1) ,Σ (2) n log(p) .

A.4 :

|β (1) | 0 ∨ |β (2) | 0 ψ (2)
Σ (1) ,Σ (2) n log(p) .

Theorem 4.4. Assuming that A.1 and A.3 hold, we have P[T B SLasso = 1] ≥ 1δ as long as

K 1 + K 2 ψ
(3) Σ (1) ,Σ (2) |β (1)

| 0 ∨ |β (2) | 0 ∨ 1 log (p) + log 1 αδ n . ( 32
)
If Σ (1) = Σ (2) = Σ and if A.1 and A.4 hold, then P[T B SLasso = 1] ≥ 1δ as long as

β (1) -β (2) 2 Σ
Var[Y (1) ] ∧ Var[Y (2) ] ψ

(4) Σ,Σ |β (1) -β (2) | 0 log (p) + log 1 αδ n . ( 33
)
Remark 4.8. The rates of detection (32) and the sparsity condition A.3 are analogous to [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF] and Condition A.2 in Theorem 4.1 for T B S ≤(n 1 ∧n 2 )/2 . The second result (33) is also similar to [START_REF] Laurent | Adaptive estimation of a quadratic functional by model selection[END_REF]. As a consequence, T B SLasso is minimax adaptive to the sparsity of (β (1) , β (2) ) and of β (1) β (2) . Remark 4.9. Dependencies of A.3, A.4, (32) and (33) on Σ (1) and Σ (2) are unavoidable because the collection S Lasso is based on the Lasso estimator which require design assumptions to work well [START_REF] Candes | Near-ideal model selection by 1 minimization[END_REF]. Nevertheless, one can improve all these dependencies using restricted eigenvalues instead of largest eigenvalues. This and other extensions are considered in next subsection.

Sharper analysis of T B

SLasso

Given a matrix X, an integer k, and a number M , one respectively defines the largest and smallest eigenvalues of order k, the compatibility constants κ[M, k, X] and η[M, k, X] (see [START_REF] Van De Geer | On the conditions used to prove oracle results for the Lasso[END_REF]) by

Φ k,+ (X) = sup θ,1≤|θ|0≤k Xθ 2 θ 2 , Φ k,-(X) = inf θ,1≤|θ|0≤k Xθ 2 θ 2 , κ[M, k, X] = min T,θ: |T |≤k, θ∈C(M,T ) Xθ θ , η[M, k, X] = min T,θ: |T |≤k, θ∈C(M,T ) √ k Xθ |θ| 1 , (34)
where

C(M, T ) = {θ : |θ T c | 1 < M |θ T | 1 }. Given an integer k, define γ Σ (1) ,Σ (2) ,k := i=1,2 κ 2 6, k * , √ Σ (i) i=1,2 Φ k * ,+ ( √ Σ (i) ) , γ Σ (1) ,Σ (2) ,k := i=1,2 Φ 2 k,+ ( √ Σ (i) ) i=1,2 Φ k,-( √ Σ (i) ) i=1,2 κ 2 [6, k, Σ (i) ] ,
that measure the closeness to orthogonality of Σ (1) and Σ (2) . Theorem 4.4 is straightforward consequence of the two following results.

Proposition 4.5. There exist four positive constants L * , L * 1 , L * 2 , and L * 3 such that following holds. Define k * as the largest integer that satisfies

(k * + 1) log(p) ≤ L * (n 1 ∧ n 2 ) , (35)
and assume that

1 + log [1/(αδ)] < L * 1 (n 1 ∧ n 2 ) . (36)
The hypothesis H 0 is rejected by T B SLasso with probability larger than 1-δ for any (β (1) , β (2) ) satisfying

|β (1) | 0 + |β (2) | 0 ≤ L * 2 γ Σ (1) ,Σ (2) ,k k * n 1 n 2 ∧ n 2 n 1 . ( 37
)
and

K 1 + K 2 ≥ L * 3 γ Σ (1) ,Σ (2) ,k * |β (1) | 0 ∨ |β (2) | 0 ∨ 1 log(p) + log{1/(αδ)} n 1 ∧ n 2 n 1 n 2 ∨ n 2 n 1 .
This proposition tells us that T B SLasso behaves nearly as well as what has been obtained in [START_REF] Jeanmougin | Defining a robust biological prior from pathway analysis to drive network inference[END_REF] for T B S ≤(n 1 ∧n 2 )/2 , at least when n 1 and n 2 are of the same order. In the next proposition, we assume that Σ (1) = Σ (2)

:= Σ. Given an integer k, define γΣ,k := κ[6, k, √ Σ]Φ 1/2 k,-( √ Σ) Φ 1,+ ( √ Σ) , γ (2)
Σ,k :=

κ 2 6, k, √ Σ Φ k,+ ( √ Σ) , γ (3)
Σ,k := Φ 2 1,+ ( √ Σ) κ 2 [6, k, √ Σ] .
(k * + 1) log(p) ≤ L * (n 1 ∧ n 2 ) , k * ≤ L * γΣ,k * n 1 ∧ n 2 |n 1 -n 2 | ∧ n 1 ∧ n 2 log(p) , (38)
with the convention x/0 = ∞. Assume that

1 + log [1/(αδ)] < L * 1 (n 1 ∧ n 2 ) .
The hypothesis H 0 is rejected by T B SLasso with probability larger than 1-δ for any (β (1) , β (2) ) satisfying

|β (1)

| 0 + |β (2) | 0 ≤ L * 2 γ(2) Σ, k * k * . ( 39
)
and

β (1) -β (2) 2 Σ Var(Y (1) ) ∧ Var(Y (2) ) ≥ L * 3 γ(3) Σ,k * |β (1) -β (2) | 0 ∨ 1 log(p) + log{1/(αδ)} .
λ of θ

(2) *

for some λ > 0 (see the proof in Section 8 for further details). Still, this is not clear to what extent this assumption is necessary.

Numerical Experiments

This section evaluates the performances of the suggested test statistics along with aforementioned test collections and calibrations on simulated linear regression datasets.

Synthetic Linear Regression Data

In order to calibrate the difficulty of the testing task, we simulate our data according to the rare and weak parametrization adopted in [START_REF] Arias-Castro | Global testing under sparse alternatives: Anova, multiple comparisons and the higher criticism[END_REF]. We choose a large but still reasonable number of variables p = 200, and restrict ourselves to cases where the number of observations n = n 1 = n 2 in each sample remains smaller than p. The sparsity of sample-specific coefficients β (1) and β (2) is parametrized by the number of non zero common coefficients p 1-η and the number of non zero coefficients p 1-η2 which are specific to β (2) . The magnitude µ r of all non zero coefficients is set to a common value of √ 2r log p, where we let the magnitude parameter range from r = 0 to r = 0.5:

β (1) = (µ r µ r . . . µ r 0 . . . 0 0 . . . 0) β (2) = (µ r µ r . . . µ r p 1-η common coefficients µ r . . . µ r p 1-η 2 sample-2-specific coefficients 0 . . . 0)
We consider three sample sizes n = 25, 50, 100, and generate two sub-samples of equal size n 1 = n 2 = n according to the following sample specific linear regression models: Y (1) = X (1) β (1) + ε (1) , Y (2) = X (2) β (2) + ε (2) .

Design matrices X (1) and X (2) are generated by multivariate Gaussian distributions, X (j) i ∼ N (0, Σ (j) ) with varying choices of Σ (j) , as detailed below. Noise components ε

(1) i and ε

(2) i are generated independantly from X (1) and X (2) according to a standard centered Gaussian distribution.

The next two paragraphs detail the different design scenarios under study as well as test statistics, collections and calibrations in competition. Each experiment is repeated 1000 times.

Design Scenarios Under Study.

Sparsity Patterns. We study six different sparsity patterns as summarized in Table 1. The first two are meant to validate type I error control. The last four allow us to compare the performances of the various test statistics, collections and calibrations under different sparsity levels and proportions of shared coefficients. In all cases, the choices of sparsity parameters η and η 2 lead to strong to very strong levels of sparsity. The last column of Table 1 illustrates the signal sparsity patterns of β (1) and β (2) associated with each scenario. In scenarios 1 and 2, sample-specific signals share little, if not none, non zero coefficient. In scenarios 3 and 4, sample-specific coefficients show some overlap. Scenario 4 is the most difficult one since the number of sample-2-specific coefficients is much smaller than the number of common non zero coefficients: the sparsity of the difference between β (1) and β (2) is much smaller than the global sparsity of β (2) . This explains why the illustration in the last column might be misleading: the two patterns are not equal but do actually differ by only one covariate.

Beyond those six varying sparsity patterns, we consider three different correlation structures Σ (1) and Σ (2) for the generation of the design matrix. In all three cases, we assume that Σ (1) = Σ (2) = Σ. On top of the basic orthogonal matrix Σ (1) = Σ (2) = I p , we investigate two randomly generated correlation structures.

Power Decay Correlation Structure. First, we consider a power decay correlation structure such that Σ i,j = ρ |i-j| . Since the sparsity pattern of β (1) and β (2) is linked to the order of the covariates, we randomly permute at each run the columns and rows of Σ in order to make sure that the correlation structure is independent from the sparsity pattern.

Gaussian Graphical Model Structure. Second, we simulate correlation structures with the R package GGMselect. The function simulateGraph generates covariance matrices corresponding to Gaussian graphical model structure made of clusters with some intracluster and extra-cluster connectivity coefficients. See Section 4 of [START_REF] Giraud | Graph selection with GGMselect[END_REF] for more details. A new structure is generated at each run.

Both random correlation structures are calibrated such that, on average, each covariate is correlated with 10 other covariates with correlations above 0.2 in absolute value. This corresponds to fixing ρ at a value of 0.75 in the power decay correlation structure and the intra-cluster connectivity coefficient to 5% in the Gaussian graphical model structure. With the default option of the function simulateGraph the extra-cluster connectivity coefficient is taken five times smaller.

Test statistics, collections and calibrations in competition

In the following, we present the results of the proposed test statistics combined with two test collections, namely a deterministic and data-driven model collection, respectively S 1 and S Lasso , as well as with a Bonferroni (B) or Permutation (P) calibration (computed with 100 random permutations). 2) specific Signals

Setting η common η 2 β (
H 00 - 0 - 0 β (2)
β (1)

H 0 5/8 7 - 0 β (2)
β (1)

1 - 0 5/8 7 β (2)
β (1)

2 7/8 1 5/8 7 β (2)
β (1)

3 5/8 7 5/8 7 β (2)
β (1)

4 5/8 7 7/8 1 β (2)
β (1) Table 1 Summary of the six different sparsity patterns under study.

Furthermore, to put those results in perspective, we compare our suggested test statistic to the usual Fisher statistic and we compare our approach with the parallel work of [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF].

Fisher statistic. For a given support |S| of reduced dimension the usual likelihood ratio statistic for the equality of β

Fi S = Y -X S β S 2 -Y (1) -X (1)
S β

(1) S 2 -Y (2) -X (2)
S β

(2) S 2 Y (1) -X (1)
S β

(1) S 2 + Y (2) -X (2)
S β

(2) S 2

n 1 + n 2 -2|S| |S| , (40)
where β S is the maximum likelihood estimator restricted to covariates in support S on the concatenated sample (X, Y). While this statistic Fi S is able to detect differences between β (1) and β (2) , it is not really suited for detecting differences between the standard deviations σ (1) and σ (2) . The Fisher statistic Fi S is adapted to the high-dimensional framework similarly as the suggested statistics (F S,V , F S,1 , F S,2 ), except that exact p-values are available. The corresponding test with a collection S and a Bonferroni (resp. permutation) calibration is denotes T B,Fisher S (T P,Fisher

S

). Procedure of Städler and Mukherjee [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF]. The DiffRegr procedure of Städler and Mukherjee performs two-sample testing between high-dimensional regression models. The procedure is based on sample-splitting: the data is split into two parts, the first one allowing to reduce dimensionality (screening step) and the second being used to compute p-values based on a restricted log-likelihood-ratio statistic (cleaning step). To increase the stability of the results the splitting step is repeated multiple times and the resulting p-values must be aggregated. We choose to use the p-value calculations based on permutations as it remains computationally reasonable for the regression case, see [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF]. The single-splitting and multi-splitting procedures are denoted respectively SS(perm) and MS(perm).

Validation of Type I Error Control

Control Under the Global Null Hypothesis H 00 . Table 2 presents estimated type I error rates, that is the percentage of simulations for which the null hypothesis is rejected, based upon 1000 simulations under the restricted null hypothesis H 00 , where β (1) = β (2) = 0 and under orthogonal correlation structure. The desired level is α = 5%, and the estimated levels are given with a 95% Gaussian confidence interval.

As expected under independence, the combination of the S 1 collection with Bonferroni correction gives accurate alpha-level when applied to the usual Fisher statistic. On the contrary when applied to the suggested statistics, the use of upper bounds on p-values leads to a strong decrease in observed type-I error. This decrease is exacerbated when using the S Lasso collection, since we are accounting for many more models than the number actually tested in order to prevent overfitting. This effect can be seen both on the Fisher statistic and our suggested statistic. Even with the usual Fisher statistic, for which we know the exact p-value, it is unthinkable to use Bonferroni calibration as soon as we adopt data-driven collections instead of deterministic ones.

On the contrary, a calibration by permutations restores a control of type-I error at the desired nominal level, whatever the test statistic or model collection.

As noted by [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF], the multi-splitting procedure yields conservative results in terms of type I error control at level 5%.

Model collection

S 1 S Lasso Calibration (B) (P) (B) (P) n= [START_REF] Laurent | Adaptive estimation of a quadratic functional by model selection[END_REF] 1 ± 0.6 6.9 ± 1.6 0 ± 0 6.9 ± 1.6 n= 50

1.8 ± 0. Control Under the Global Equality of Non Null Coefficients H 0 . Figures 1 and3 present level checks under H 0 but with non null β (1) = β (2) = 0, under respectively orthogonal and non-orthogonal correlation structures. Conclusions are perfectly similar to the case H 00 : all methods behave well, except the multi-split DiffRegr procedure and the Bonferroni calibration-based procedures T B S (for any collection S) and T B,Fisher SLasso . In particular, the Fisher statistic combined with S 1 and Bonferroni calibration is more conservative than the desired nominal level under correlated designs. Power Analysis. We do not investigate the power of the Bonferroni-based procedures T B S and T B,Fisher S as they have been shown to be too conservative in the above Type I error analysis. Figure 4 represents power performances for the test T P S and the usual likelihood ratio test T P,Fisher S combined with either S 1 or S Lasso test collections using a calibration by permutation under an orthogonal covariance matrix Σ, as well as power performance for the DiffRegr procedure. Figure 5 represents equivalent results for power decay and GGM covariance structures when n = 50.

In the absence of common coefficients (scenarios 1 and 2), the test T P S reaches 100% power from very low signal magnitudes and small sample sizes. Compared to the test based on usual likelihood ratio statistics, which does not reach more than 40% power when . The deterministic collection S 1 is drawn in empty points, while the data-driven collection S Lasso is in plain points. Green circles represent the DiffRegr procedure, respectively plain and empty for single-splitting and multi-splitting.

n = 25 given the signal magnitudes under consideration, the suggested statistics proves itself extremely efficient. Under these settings as well, any subset of size 1 containing one of the variables activated in only β (2) can suffice to reject the null, which is why collection S 1 performs actually very well when associated with (F S,V , F S,1 , F S,2 ) and not so badly when associated with Fi S . However, in more complex settings 3 and 4, where larger subsets are required to correct for strong and numerous common effects, subset collection S Lasso yields a higher power than the collection S 1 .

For small n, the test T P SLasso outperforms the procedure DiffRegr, whose limitation likely stems from the half sampling step. This limitation of sample splitting approaches has already been noticed by [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF]. However, for n = 100 the procedure DiffRegr performs better than our procedure in the highly challenging setting 4.

Figure 5 provide similar results under respectively power decay correlated designs and GGM-like correlated designs for a sample size of n = 50, leading to similar conclusions as in the uncorrelated case.

Application to GGM

The following section explicits the extension of the two-sample linear regression testing framework to the two-sample Gaussian graphical model testing framework. We describe the tools and guidelines for a correct interpretation of the results and illustrate the approach on a typical two-sample transcriptomic data-set.

How to Apply this Strategy to GGM Testing

Neighborhood Selection Approach The procedure developed in Section 2 can be adapted to the case of Gaussian graphical models as in [START_REF] Verzelen | Tests for gaussian graphical models[END_REF]. We quickly recall why estimation of the Gaussian graphical model amounts to the estimation of p independent linear regressions when adopting a neighborhood selection approach [START_REF] Meinshausen | High-dimensional graphs and variable selection with the Lasso[END_REF].

Consider two Gaussian random vectors Z (1) ∼ N (0, Ω (1) -1 ) and Z (2) ∼ N (0, Ω (2) -1 ). Their respective conditional independence structures are represented by the graphs G (1) and G (2) , which consist of a common set of nodes Γ = {1, . . . , p} and their respective sets

S

. The deterministic collection S 1 is drawn in empty points, while the data-driven collection S Lasso is in plain points. Results for the DiffRegr procedure are represented by green circles, respectively plain and empty for single-splitting and multi-splitting approaches.

Settings

Power decay GGM 1.

β (2)
β (1) 0.0

S

. The deterministic collection S 1 is drawn in empty points, while the data-driven collection S Lasso is in plain points. Results for the DiffRegr procedure are represented by green circles, respectively plain and empty for single-splitting and multi-splitting approaches.

of edges E (1) and E (2) . When speaking of gene regulation networks, each node represents a gene, and edges between genes are indicative of potential regulations. In contrast with gene co-expression networks, edges in Gaussian graphical models do not reflect correlations but partial correlations between gene expression profiles. Formally, an edge (i, j) belongs to the edge set E (1) (resp. E (2) ) if

Z (1) i (resp. Z (2) i ) is independent from Z (1) j (resp. Z (2) j ) conditional on all other variables Z (1)
\i,j (resp. Z

(2) \i,j ). When the precision matrix Ω (k) is nonsingular, the edges are characterized by its non zero entries.

The idea of neighborhood selection is to circumvent the intricate issue of estimating the precision matrix by recovering the sets of edges neighborhood by neighborhood, through the conditional distribution of Z (k) i given all remaining variables Z (k) \i . Indeed, this distribution is again a Gaussian distribution, whose mean is a linear combination of Z (k) \i while its variance is independent from

Z (k) \i . Hence, Z (k) i
can be decomposed into the following linear regression:

Z (k) i = j =i Z (k) j β (k) ij + ε (k) i = Z (k) \i β (k) i + ε (k) i , (41)
where

β (k) ij = -Ω (k) ij /Ω (k)
ii and Var[ε

(k) i ] = (Ω (k) ii ) -1 .
Given an n 1 -sample of Z (1) and an n 2 -sample of Z (2) , we recall that our objective is to test as formalized in ( 4)

H G 0 : Ω (1)
= Ω (2) versus H G 1 : Ω (1) = Ω (2) .

As a result of Equation ( 41), testing for the equality of the matrix rows Ω

(1)

i. = Ω (2) i.
is equivalent to testing for β

(1) i = β (2) i
and Var[ε

(1) i ] = Var[ε (2) i ]. Denote by Σ (k) \i the covariance of Z (k)
\i . Under the null H G 0 , we have that for any i, Σ

\i = Σ (2) (1)
\i . Consequently, we can translate the GGM hypotheses given in Equation ( 4) into a conjunction of twosample linear regression tests:

H G 0 : i β (1) i = β (2) i , Ω (1)
ii = Ω

(2)

ii , Σ

\i = Σ (2) \i (1)
H G 1 : i β (1) i = β (2) i ∪ Ω (1) ii = Ω (2) ii (42)
.

Concretely, we apply the previous two-sample linear regression model with

X (1) = Z (1) ,\i ,X (2) = Z (2) ,\i , Y (1) = Z (1)
,i , and

Y (2) = Z (2)
,i for every gene i and combine multiple neighborhood tests using a Bonferroni calibration as presented in Algorithm 5. The equality of σ (k) 's in H 0 models the equality of Ω (k)

ii 's in H G 0 while the equality of Σ (k) 's accounts for the equality of Σ (k) \i 's.

Interpretation. Because we need Ω

(1)

ii = Ω (2)
ii and Σ

(1) \i = Σ (2)
\i for every neighborhood in the two-sample GGM null hypothesis H N 0 (42), the assumptions that σ (1) = σ (2) and Σ (1) = Σ (2) in the two-sample linear regression null hypothesis H 0 (3) are crucial for each neighorbood test to be interpreted correctly. As a result, only the global test can be strictly speaking interpreted in a statistically correct sense.

However in practice, when the global null hypothesis is rejected, our construction of neighborhood tests provides helpful clues on the location of disruptive regulations. In particular, for each rejected neighorhood test i, one can keep track of the rejected model

Algorithm 5 Gaussian Graphical Model Testing Strategy

Require: Data Z (1) ,Z (2) , maximum model dimension Dmax and desired level α for each gene i = 1, . . . , p do procedure Neighborhood Test Define X (1) = Z

,\i ,

X (2) = Z (2) ,\i Define Y (1) = Z (1) ,i , Y (2) = Z (2) ,i
Apply the Adaptive Testing Strategy of Algorithm 1 to X (1) ,X (2) ,Y (1) ,Y (2) end procedure end for Reject the global null hypothesis if at least one Neighborhood Test is rejected at level α/p S i R , retaining sensible information on which particular regulations are most likely altered between samples.

Illustration on Real Transcriptomic Breast Cancer Data

We apply this strategy to the full (training and validation) breast cancer dataset studied by [START_REF] Hess | Pharmacogenomic predictor of sensitivity to preoperative chemotherapy with paclitaxel and fluorouracil, doxorubicin, and cyclophosphamide in breast cancer[END_REF] and [START_REF] Natowicz | Prediction of the outcome of preoperative chemotherapy in breast cancer using dna probes that provide information on both complete and incomplete responses[END_REF], whose training subset was originally published in [START_REF] Rouzier | Breast cancer molecular subtypes respond differently to preoperative chemotherapy[END_REF]. The full dataset consists of microarray gene expression profiles from 133 patients with stage I-III breast cancer undergoing preoperative chemotherapy. A majority of patients (n=99) presented residual disease (RD), while 34 patients demonstrated a pathologic complete response (pCR). The common objective of [START_REF] Hess | Pharmacogenomic predictor of sensitivity to preoperative chemotherapy with paclitaxel and fluorouracil, doxorubicin, and cyclophosphamide in breast cancer[END_REF] and [START_REF] Natowicz | Prediction of the outcome of preoperative chemotherapy in breast cancer using dna probes that provide information on both complete and incomplete responses[END_REF] was to develop a predictor of complete response to treatment from gene expression profiling. In particular, [START_REF] Hess | Pharmacogenomic predictor of sensitivity to preoperative chemotherapy with paclitaxel and fluorouracil, doxorubicin, and cyclophosphamide in breast cancer[END_REF] identified an optimal predictive subset of 30 probes, mapping to 26 distinct genes.

[2] inferred Gaussian graphical models among those 26 genes on each patient class using weighted neighborhood selection. The corresponding graphs of conditional dependencies for medium regularization are presented in Figure 6. Those two graphs happen to differ dramatically from one another. The question we tackle is whether those differences remain when taking into account estimation uncertainties.

Pathologic Complete Response (pCR)

Residual Disease (RD) Figure 6. Graphs of conditional dependencies among the 26 genes selected by [START_REF] Hess | Pharmacogenomic predictor of sensitivity to preoperative chemotherapy with paclitaxel and fluorouracil, doxorubicin, and cyclophosphamide in breast cancer[END_REF] on patients with pathologic complete response or residual disease with medium regularization as presented in Figure 3 of [START_REF] Ambroise | Inferring sparse Gaussian graphical models with latent structure[END_REF].

We run for each of the p = 26 genes a neighborhood test T P ŜLasso at level 0.05/26. We associate to each neighborhood test the empirical p-value defined in 23 that has to be be compared to α/p.

Most of the graph estimation methods proposed in the literature, such as the procedure of [START_REF] Ambroise | Inferring sparse Gaussian graphical models with latent structure[END_REF] leading to Figure 6, rely on the assumption that observations are i.i.d. Yet the training and validation datasets have been collected and analysed separately by two different clinical centers. We therefore start by checking whether the pooled sample can be considered as homogeneous. Within each group of patients (RD and pCR), we lead a test for the homogeneity of Gaussian graphical models between the training and validation subsets.

Within pCR patients (3), two neighborhood tests corresponding to CA12 and PDGFRA are rejected at level 0.05/26. Within RD patients (4), half of the neighborhoods happen to differ significantly between the training and validation datasets. Genes CA12 and JMJD2B are responsible for the rejection of respectively seven and six neighborhoods. Because of these surprisingly significant divergences between training and validation subsets, we restrict the subsequent analysis to the training set (n=82 patients, among which 61 RD and 21 pCR patients).

To roughly check that we got rid of the underlying heterogeneity, we create an artificial dataset under H 0 by permutation of the patients, regardless of their class. No neighborhood test is rejected at a level corrected for multiple testing. We also cut the group of patients with residual disease artificially in half. When testing for the difference between the two halves, no significant heterogeneity remains, whatever the neighborhood.

Within the training set, the comparison of Gaussial graphical structures between pCR and RD patients leads to the rejection of all neighborhood tests after Bonferroni correction for multiple testing of the 26 neighborhoods, as summarized in Table 5. RRM2, MAPT and MELK genes appear as responsible for the rejection of respectively nine, nine and four of these neighborhood tests. Quite interestingly, these three genes have all been described in clinical literature as new promising drug targets. [START_REF] Heidel | Potent sirna inhibitors of ribonucleotide reductase subunit rrm2 reduce cell proliferation in vitro and in vivo[END_REF] exhibited inhibitors of RRM2 expression, which reduced in vitro and in vivo cell proliferation. [START_REF] Rouzier | Microtubule-associated protein tau: a marker of paclitaxel sensitivity in breast cancer[END_REF] led functional biology experiments validating the relationship between MAPT expression levels and response to therapy, suggesting to inhibit its expression to increase sensivity to treatment. More recently, [START_REF] Chung | Development of an orallyadministrative melk-targeting inhibitor that suppresses the growth of various types of human cancer[END_REF] Although this setting is particularly suited to consider the two-sample GGM testing (Section 6), one may wonder whether one can circumvent these two restrictions. We doubt that this is possible without making the testing problem much more difficult. First, the formulation (3) allows the null hypothesis to be interpreted as a relevant intermediary case between two extreme fixed design settings: design equality (X (1) = X (2) ) and arbitrary different design (X (1) = X (2) ). In the first case, the two-sample problem amounts to a one-sample problem by considering Ỹ = Y (1) -Y (2) and it has therefore already been thoroughly studied. The second case is on the contrary extremely difficult as illustrated by the proposition below. Proposition 7.1. Consider the design matrices X (1) and X (2) as fixed and assume that σ (1) = σ (2) = 1. If the (n 1 + n 2 ) × p matrix formed by X (1) and X (2) has rank n 1 + n 2 , then any test T of β (1) = β (2) vs β (1) = β (2) based on the data (Y, X) satisfies:

sup β∈R p P β,β [T = 1] + inf β (1) =β (2) ∈R p P β (1) ,β (2) [T = 0] ≥ 1 ,
where P β (1) ,β (2) (.) denotes the distribution of (Y (1) , Y (2) ). In other words, any level-α test T has a type II error larger than 1-α, and this uniformly over β (1) and β (2) . Consequently, any test in this setting cannot perform better than complete random guess.

Furthermore, if Σ (1) = Σ (2) is allowed in the null (3), then the two-sample testing problem becomes much more difficult in the sense that it is impossible to reformulate the null hypothesis into a conjonction of low-dimensional hypotheses as done in Lemma 2.1. Indeed, consider the following toy example:

σ (1) = σ (2) = 1, β (1) = β (2)
= (a, 0, 0, . . .) T for some a > 0, Σ (1) = I p and Σ (2) = (ρ + 1 i=j ) 1≤i,j≤p for some ρ > 0. Then, for any subset S that does not contain the first component, the parameters β

S = β (2) (1)
S does not imply that β (1) = β (2) and one should not rule out the parameter equality hypothesis relying on some low-dimensional regressions.

Comparison with related work [START_REF] Städler | Two-sample testing in high-dimensional models[END_REF][START_REF] Städler | Network-based multivariate gene-set testing[END_REF] Städler and Mukherjee propose a very general approach to high-dimensional two-sample testing, being applicable to a wide range of models. In particular this approach allows for the direct comparison of two-sample Gaussian graphical models without adopting a neighborhood selection approach. This avoids the burden of multiple neighborhood linear regression and the multiple testing correction which follows.

Because they estimate the supports of sample-specific estimators and joint estimator separately in the screening step, they resort to an elegant estimation of the p-values for the non-nested likelihood ratio test in the cleaning step. Yet, they do not provide any theoretical controls on type I error rate or power for their overall testing strategy.

Finally, as it appears in the numerical experiments, their approach is based on halfsampling and can thus suffer from an acute reduction of power on small samples. On the bright side, the multi-split procedure shows stable results and performs well even in difficult scenarios as soon as n is sufficiently large.

Non asymptotic bounds and constants. In the spirit of [START_REF] Baraud | Adaptive tests of linear hypotheses by model selection[END_REF], our type II error analysis is completely non-asymptotic. However, the numerical constants involved in the bounds are clearly not optimal. Another line of work initiated by [START_REF] Donoho | Higher criticism for detecting sparse heterogeneous mixture[END_REF] considers an asymptotic but high-dimensional framework and aims to provide detection rates with optimal constants. For instance [START_REF] Arias-Castro | Global testing under sparse alternatives: Anova, multiple comparisons and the higher criticism[END_REF][START_REF] Ingster | Detection boundary in sparse regression[END_REF] have derived such results in the one-sample high-dimensional linear regression testing problem under strong assumptions on the design matrices. In our opinion, both analyses are complementary. While deriving sharp detection rates (under perhaps stronger assumptions on the covariance) is a stimulating open problem, it is beyond the scope of our paper.

Loss functions and interpretation. The Kullback discrepancies considered in the power analysis of the test depend on β (1) and β (2) through the prediction distances β (1)β (2) Σ i , i = 1, 2 rather than the l 2 distance β (1)β (2) . On the one hand, such a dependency on the prediction abilities is natural, as our testing procedures relies on the likelihood ratio. On the other hand, it is possible to characterize the power of our testing procedures as in Theorems 4.1 and 4.4 in terms of the distance β (1) β (2) by inverting Σ (1) and Σ (2) at β (1) β (2) . However, the inversion would lead to an additional factor of the form Φ

-1 |β (1) -β (2) |0,-( √ Σ (i)
) in the testing rates. In terms of interpretation, even though our procedure adopts a global testing approach through prediction distances, our real dataset example illustrates that identifying which subset in the collection is responsible for rejecting the null hypothesis provides clues into which specific coefficients are most likely to differ between samples.

Gene network inference. Thinking of gene network inference by Gaussian graphical modeling, the high levels of correlations encountered within transcriptomic datasets and the potential number of missing variables result in highly unstable graphical estimations. Our global testing approach provides a way to validate whether sample-specific graphs eventually share comparable predictive abilities or disclose genuine structural changes. Such a statistical validation is obviously crucial before translating any graphical analysis into further biological experiments. Interestingly, the three main genes pointed out by our testing strategy have been validated as promising therapeutic targets by functional biology experiments.

Finally, this test should also facilitate the validation of the fundamental i.i.d. assumption across multiple samples, paving the way to pooled analyses when possible. In that respect, we draw attention to the significant heterogeneity detected between the training and validation subsets of the well-known Hess et al dataset, suggesting that these samples should be used separately as originally intended. Methods which require i.i.d. observations should only be applied with caution to this dataset if considered as a single large and homogeneous sample.

Proofs

Two-sample testing for fixed and differents designs

Proofs of Proposition 7.1. Using the rank condition, we derive that for any vector (a, b) in R n1 × R n2 , there exists β ∈ R p such that X (1) β = a and X (2) β = b. Consequently, under the null hypothesis, (Y (1) , Y (2) ) follows any distributions N (a, I n1 )⊗N (b, I n2 ) with (a, b) arbritary in R n1 × R n2 . Hence, for any β (1) = β (2) ∈ R p , the distribution P β (1) ,β (2) of (Y (1) , Y (2) ) is not distinguishable from the null hypothesis. The result follows.

Upper bounds of the quantiles

Proof of Proposition 3.3. For the sake of simplicity, we note N = n 1 -|S|, (Z 1 , . . . , Z |S| ) a standard Gaussian random vector and W N a χ 2 random variable with N degrees of freedom. We apply Laplace method to upper bound P[F S,1 ≥ u]:

P[F S,1 ≥ u] = P   |S| i=1 a i Z 2 i ≥ uW N /N   ≤ inf λ>0 E exp   λ |S| i=1 a i Z 2 i -λuW N /N   ≤ inf 0<λ<|a|∞/2 exp [ψ u (λ)] ,
where

ψ u (λ) = - 1 2 |S| i=1 log(1 -2λa i ) - N 2 log 1 + 2λu N .
The sharpest upper-bound is given by the value λ * which minimizes ψ u (λ). We obtain an approximation of λ * by cancelling the second-order approximation of its derivative. Deriving ψ u gives

ψ u (λ) = |S| i=1 a i 1 -2λa i - u 1 + 2λu N ,
which admits the following second order approximation :

|a| 1 + 2λ a 2 1 -2|a| ∞ λ - u 1 + 2λu N . (43)
Cancelling this quantity amounts to solving a polynomial equation of the second degree. The smallest solution of this equation leads to the desired λ * .

Additional Notations. Given a subset S, Π

S (resp. Π (1)
S ) stands for the orthogonal projection onto the space spanned by the rows of X (1)

S (resp. X (2) S ). Moreover, Π (1)
S ⊥ denotes the orthogonal projection along the space spanned by the rows of X Let us consider the regression of Y (1) (resp. Y (2) ) with respect to X

Y (1) = X (1) S β (1) S + (1) S , Y (2) = X (2) S β (2)
S +

S .

Under the null hypothesis H 0,S , we have β

(1) S = β (2) S and σ (1) S = σ (2)
S . For the sake of simplicity, we write β S and σ S for these two quantities. Define the random variable T 1 and T 2 as

T 1 = Π (1) S ⊥ (1) S 2 (n 1 -|S|)σ 2 S , T 2 = Π (2) S ⊥ (2) S 2 (n 2 -|S|)σ 2 S . (44)
Given X, T 1 /T 2 follows a Fisher distribution with (n 1 -|S|, n 2 -|S|) degrees of freedom.

Observing that under the null hypothesis

F S,V = -2 + T 1 T 2 n 2 (n 1 -|S|) n 1 (n 2 -|S|) + T 2 T 1 n 1 (n 2 -|S|) n 2 (n 1 -|S|)
allows us to prove the first assertion of Proposition 3.1. Let us turn to the second statistic:

F S,1 = n 1 n 2 (n 1 -|S|) U T 1 ,
where

U = X (2)
S (X

S X

(2)

S ) -1 X (2) S (2)
S -X

S (X

S X

(1)

S ) -1 X (1) S (1) S 2 σ 2 S . Given X, U is independent from T 1 since T 1 is a function of Π (1) S ⊥ (1)
S while U is a function of (

S , Π (2)
S ). Furthermore, U is the squared norm of a centered Gaussian vector with covariance X

(2) S (X

S X

(1)

S ) -1 + (X (2)
S X

(2)

S ) -1 X (2)
P H0 Q i,|S| (F S,i |X S }) ≤ α i,S |X S ≤ α i,S .
Applying a union bound and integrating with respect to X allows us to control the type I error:

P H0 [T B S = 1] = E   S∈ S i=V,1,2 P Q i,|S| (F S,i |X S }) < α i,S   ≤ S∈S i=V,1,2 P Q i,|S| (F S,i |X S }) < α i,S ≤ S∈S i=V,1,2 E X S P Q i,|S| (F S,i |X S }) < α i,S ≤ S∈S α i,S ≤ α ,
where we have upper bounded the sum over the random collection S by the sum over S.

Proof of Proposition 3.5. Consider i ∈ {V, 1, 2}. Under H 0 , the distributions of min

S∈ Sπ Q V,|S| (F S,V (π)|X π S ) p |S| , min S∈ Sπ Q 1,|S| (F S,1 (π)|X π S ) Q 2,|S| (F S,2 (π)|X π S )
p |S| are invariant with respect to the permutation π. Hence, we derive

P H0 min S∈ S Q V,|S| (F S,V |X S ) p |S| ≤ C 1 X S = α/2 , P H0 min S∈ Sπ Q 1,|S| (F S,1 (π)|X π S ) Q 2,|S| (F S,2 (π)|X π S ) p |S| ≤ C 2 X S = α/2 .
Applying a union bound and integrating with respect to X allows us to conclude.

Proof of Theorem 4.3

The objective is to exhibit a subset for which the power of T B S is larger than 1δ. This subset is such that the distance between the two sample-specific distributions is large enough that we can actually reject the null hypothesis with large probability. As exposed in Theorem 4.3, we rely on the semi-distances K 1 (S) + K 2 (S) for S ∈ S:

2(K 1 (S) + K 2 (S)) = σ (1) S σ (2) S 2 + σ (2) S σ (1) S 2 -2 + β (2) S -β (1) S 2 Σ (2) (σ (2)
S ) 2 + β (2)
Sβ

(1) S 2 Σ (1) (σ (1)
S ) 2 . ( 45
)
The proof is split into five main lemmas. First, we upper bound

Q -1 V,|S| (x|X S ), Q -1
) < L • 1 (n 1 ∧ n 2 ), log(1/α S ) ≤ L • 2 (n 1 ∧ n 2 ) , |S| ≤ L • 3 (46)
where the numerical constants L V,|S| (x|X S ) ). There exists a positive universal constant L such that the following holds. Consider some 0 < x < 1 such that 16 log(2/x) ≤ n 1 ∧n 2 . For any subset S of size smaller than (n 1 ∧ n 2 )/2, we have

Q -1 V,|S| (x|X S ) ≤ L |S|(n 1 -n 2 ) n 1 n 2 2 + log(2/x) 1 n 1 + 1 n 2 . ( 47
)
We recall that a = (a 1 , . . . , a |S| ) denotes the positive eigenvalues of

n 1 n 2 (n 1 -|S|) X (2) S (X (1)
S X

(1)

S ) -1 + (X (2)
S X

(2)

S ) -1 X (2) S . Lemma 8.2 (Upper-bound of Q -1 1,|S| (x|X S )
). There exist two positive universal constants L 1 and L 2 such that the following holds. If

|a| 1 < u ≤ (n 1 -|S|)|a| ∞ and if |S| ≤ L 1 n 1 , log Q 1,|S| (u|X S ) ≤ - (u -|a| 1 ) 2 4 [|a| ∞ (u -|a| 1 ) + a 2 ] + (u -|a| 1 )u 3 2(n 1 -|S|) [|a| ∞ (u -|a| 1 ) + a 2 ] 2 .
For any 0 < x < 1, satisfying

L 2 log(1/x) ≤ n 1 -|S| , (48)
we have the following upper bound

Q -1 1,|S| (x|X S ) ≤ |a| ∞ 2|S| + 2 2|S|
|a| ∞ ≤ L 2     1 n 2 + ϕ max Σ (2)
S (Σ

S ) -1 Σ (2) S n 1     . (1)
Lemma 8.4 (Deviations of F S,V ). There exist three positive universal constants L 1 , L 2 and L 3 such that the following holds. Assume that L 1 log(1/δ) ≤ n 1 ∧ n 2 . With probability larger than 1δ, we have

F S,V ≥ L 2 (σ (1)
S ) 2 -(σ (2)
S ) 2 σ (1) S σ (2) S 2 -L 3 |S| 2 1 n 2 1 + 1 n 2 2 + log 1 δ 1 n 1 + 1 n 2 . ( 50
)
Lemma 8.5 (Deviations of F S,1 ). There exist two positive universal constants L 1 and L 2 such that the following holds. Assume that

L 1 log(12/δ) < n 1 ∧ n 2 . ( 51
)
With probability larger than 1δ/2, we have

F S,1 ≥ β (2) S -β (1) S 2 Σ (2) 8(σ (1)
S ) 2 -log (6/δ) L 2 1 n 2 (σ (2)
S ) 2 (σ (1) S ) 2 + ϕ S n 1 , (52)
where ϕ S is defined in [START_REF] Lin | Domain and functional analysis of a novel breast tumor suppressor protein, scube2[END_REF].

Consider some S ∈ S . Combining Lemmas 8.1 and 8.4, we derive that Q V,|S| (F S,V |X S ) ≤ α S holds with probability larger than 1δ if (σ

(1) S ) 2 -(σ (2) S ) 2 2 (σ (1) S ) 2 (σ (2) S ) 2 ≥ L |S| 2 1 n 2 1 + 1 n 2 2 + log[1/(α S δ)] 1 n 1 + 1 n 2 .
Similarly, combining Lemmas 8.2, 8.3, and 8.5, we derive that Q 1,|S| (F S,1 |X S ) ≤ α S with probability larger than 1δ if

β (2) S -β (1) S 2 Σ (2) (σ (1)
S ) 2 ≥L 1 (ϕ S + 1) 1 n 1 + 1 n 2 |S| + log 6 δα S + L 2 n 2 σ (2) S σ (1) S 2 log 6 δ . A symmetric result holds for Q 2,|S| (F S,2 |X S ). Consequently, Q V,|S| (F S,V |X S ) ∧ Q 1,|S| (F S,1 |X S ) ∧ Q 2,|S| (F S,2 |X S ) ≤ α S with proba- bility larger than 1 -δ if K 1 (S) + K 2 (S) ≥ L * 1 ϕ S 1 n 1 + 1 n 2 |S| + log 6 α S δ +L * 2 log(6/δ) 1 n 1 + 1 n 2   σ (2) S σ (1) S 2 + σ (1) S σ (2) S 2   . (53)
Since we assume that 4L * 2 log(6/δ) ≤ n 1 ∧ n 2 in (46), the last condition is fulfilled if

K 1 (S) + K 2 (S) ≥ L * ϕ S 1 n 1 + 1 n 2 [|S| + log{6/(α S δ)}] .
We now proceed to the proof of the five previous lemmas.

Proof of Lemma 8.1. Let u ∈ (0, 1) and F -1 D,N (u) be the 1u quantile of a Fisher random variable with D and N degrees of freedom. According to [START_REF] Baraud | Adaptive tests of linear hypotheses by model selection[END_REF], we have

F -1 D,N (u) ≤ 1 + 2 1 D + 1 N log 1 u + N 2D + 1 exp 4 N log 1 u -1 .
Let us assume that 8/N log(1/u) ≤ 1. By convexity of the exponential function it holds that

F -1 D,N (u) ≤ 1 + 2 1 D + 1 N log 1 u + 4 D + 8 N log 1 u .
Recall T 1 and T 2 defined in [START_REF] Verzelen | Minimax risks for sparse regressions: Ultra-high-dimensional phenomenons[END_REF]. Under hypothesis H 0 ,

T 1 T 2 ∼ Fisher(n 1 -|S|, n 2 -|S|) .
Consider some x > 0 such that [8/(n 1 -|S|) ∨ 8/(n 2 -|S|)] log(2/x) ≤ 1. Then, with probability larger than 1x/2 we have,

T 1 T 2 n 2 (n 1 -|S|) n 1 (n 2 -|S|) ≤ 1 + |S|(n 1 -n 2 ) n 1 (n 2 -|S|) 1 + 8 log(2/x) n 1 -|S| + 8 log(2/x) n 2 -|S| ≤ 1 + |S|(n 1 -n 2 ) n 1 (n 2 -|S|)   1 + 12 log(2/x) n 1 + 12 log(2/x) n 2   ≤ L , since |S| ≤ (n 1 ∧ n 2 )/2.
Similarly, with probability at least 1x/2, we have

T 2 T 1 n 1 (n 2 -|S|) n 2 (n 1 -|S|) ≤   1 + |S|(n 2 -n 1 ) n 2 (n 1 -|S|)   1 + 12 log(2/x) n 1 + 12 log(2/x) n 2     ∧ L . ( 54
)
Depending on the sign of T1 T2 n2(n1-|S|) n1(n2-|S|) -1, we apply one the two following identities:

T 1 T 2 n 2 (n 1 -|S|) n 1 (n 2 -|S|) + T 2 T 1 n 1 (n 2 -|S|) n 2 (n 1 -|S|) -2 = T 1 T 2 n 2 (n 1 -|S|) n 1 (n 2 -|S|) -1 2 T 2 T 1 n 1 (n 2 -|S|) n 2 (n 1 -|S|) , T 1 T 2 n 2 (n 1 -|S|) n 1 (n 2 -|S|) + T 2 T 1 n 1 (n 2 -|S|) n 2 (n 1 -|S|) -2 = T 2 T 1 n 1 (n 2 -|S|) n 2 (n 1 -|S|) -1 2 T 1 T 2 n 2 (n 1 -|S|) n 1 (n 2 -|S|) .
Combining the different bounds, we conclude that with probability larger than 1x,

F S,V := T 1 T 2 n 2 (n 1 -|S|) n 1 (n 2 -|S|) + T 2 T 1 n 1 (n 2 -|S|) n 2 (n 1 -|S|) -2 ≤ L |S|(n 1 -n 2 ) n 1 n 2 2 + log(2/x) n 1 + n 2 n 1 n 2 .
Proof of Lemma 8.2. As in the proof of Proposition 3.3, we note N = n 1 -|S|. Recall that Q 1,|S| (u|X S ) is defined as exp ψ u (λ ) (see Definition 3.2). We start by upperbounding ψ u (λ ), which proves the first upper-bound of the logarithm of the tail probability log Q 1,|S| (u|X S ). We then exhibit a value u x such that ψ ux (λ ) ≤ log x.

Upper-bound of the tail probability. Since Equation ( 43) is increasing with respect to λ and with respect to N , λ * decreases with N . Consequently,

λ * ≤ λ + := u -|a| 1 2 [|a| ∞ (u -|a| 1 ) + a 2 ]
.

By convexity, 1 -√ 1x ≥ x/2 for any 0 ≤ x ≤ 1. Applying this inequality, we upper bound √ ∆ and derive that

λ * ≥ λ -:= u -|a| 1 2 |a| ∞ (u -|a| 1 ) + a 2 + |a|1u N . Since u ≤ N |a| ∞ , 2λ * u ≤ N . Observing that -log(1 -2x)/2 ≤ x + x 2 /(1 -2x
) for any 0 < x < 1/2 and that log(1 + x) ≥ xx 2 for any x > 0, we derive

ψ u (λ * ) ≤ |a| 1 λ + + λ 2 + a 2 1 -2|a| ∞ λ + -λ * u + 2 (λ * ) 2 u 2 N ≤ - (u -|a| 1 ) 2 4 [|a| ∞ (u -|a| 1 ) + a 2 ] + 2λ 2 + u 2 N + (λ + -λ -)u ≤ - (u -|a| 1 ) 2 4 [|a| ∞ (u -|a| 1 ) + a 2 ] + (u -|a| 1 )u 3 2N [|a| ∞ (u -|a| 1 ) + a 2 ] 2 . (55)
Upper-bound of the quantile. Let us turn to the upper bound of Q -1 1,|S| (x|X S ). Consider u x the solution larger than |a| 1 of the equation

(u -|a| 1 ) 2 4 [|a| ∞ (u -|a| 1 ) + a 2 ] = 2 log(1/x) ,
and observe that

2 a log(1/x) ≤ u x -|a| 1 ≤ 2 √ 2 a log(1/x) + 8|a| ∞ log(1/x) .
Choosing L 1 and L 2 large enough in the condition |S| ≤ L 1 n 1 and in condition [START_REF] Wasserman | High dimensional variable selection[END_REF] leads us to u x ≤ N |a| ∞ . We now prove that ψ ux∨2|a|1 (λ * ) ≤ log x. If u x ≥ 2|a| 1 , then u 3

x ≤ 8(u x -|a| 1 ) 3 and it follows from (55) that

ψ ux (λ * ) ≤ log(1/x) -2 + 2 8 log(1/x) N ≤ -log(1/x)
if we take L 2 large enough in Condition [START_REF] Wasserman | High dimensional variable selection[END_REF]. If

u x ≤ 2|a| 1 , then |a| 2 1 /(|a| ∞ |a| 1 + a 2 ) ≥ 8 log(1/x) and ψ ux∨2|a|1 (λ * ) ≤ - |a| 2 1 4 [|a| ∞ |a| 1 + a 2 ] 1 - 2 4 |a| 2 1 N [|a| ∞ |a| 1 + a 2 ] ≤ -log(1/x) ,
if we take L 1 and L 2 large enough in the two aforementionned condition. since |S| ≤ 2 -6 n 1 . Thus, we conclude that

Q -1 1,|S| (x|X S ) ≤ u x ∨ 2|a| 1 ≤ |a| 1 + 2 √ 2 a log(1/x) + 8|a| ∞ log(1/x) ∨ |a| 1 .
Proof of Lemma 8.3. Upon defining Z

(1)

S = X (1) S Σ (1) S -1/2
and Z

(2)

S = X (2) S Σ (2) S -1/2
, it follows that Z

S and Z

S follow standard Gaussian distributions.

|a| ∞ ≤ n 1 n 2 (n 1 -|S|) 1 + ϕ max Z (2) S Σ (2)
S (Σ

S ) -1 Z (1)
S Z

(1) S -1

(Σ

S ) -1 Σ (2) S Z (2) S ≤ 2 n 2 + 2 ϕ max [Z (1)
S Z

(2)

S ] n 2 ϕ max [Z (1) S Z (1) S ] ϕ max Σ (2)
S (Σ

S ) -1 Σ

(2) S

.

In order to conclude, we control the largest and the smallest eigenvalues of Standard Wishart matrices applying Lemma 8.12.

Proof of Lemma 8.4. By symmetry, we can assume that σ

S /σ

S ≥ 1. Recall the definition of T 1 and T 2 in the proof of Proposition 3.1

CASE 1. Suppose that T 1 /T 2 ≥ 1. -2 + (σ (1)
S ) 2 (σ (2) S ) 2 T 1 T 2 + (σ (2)
S ) 2 (σ (1) S ) 2 T 2 T 1 ≥ [(σ (1)
S ) 2 -(σ (2)
S ) 2 ] 2 (σ (1) S ) 2 (σ (2) S ) 2 + (σ (1)
S ) 2 (σ (2) S ) 2 T 1 T 2 -1 + (σ (2)
S ) 2 (σ (1) S ) 2 T 2 T 1 -1 ≥ [(σ (1)
S ) 2 -(σ (2)
S ) 2 ] 2 (σ (1) S ) 2 (σ (2) S ) 2 . ( 56
) CASE 2. Suppose that T 1 /T 2 ≤ 1. -2 + (σ (1)
S ) 2 (σ (2) S ) 2 T 1 T 2 + (σ (2)
S ) 2 (σ (1) S ) 2 T 2 T 1 = (σ (1)
S ) 2 (σ (2) S ) 2 - T 2 T 1 2 (σ (2)
S ) 2 (σ (1) S ) 2 T 1 T 2 ≥ T 1 T 2 [(σ (1)
S ) 2 -(σ (2)
S ) 2 ] 2 4(σ (1) S ) 2 (σ (2) S ) 2 1 (σ (1) S ) 2 (σ (2) S ) 2 -1≥2 T 2 T 1 -1
.

We need to control the deviations of T 2 /T 1 . Using bound (54), we get

T 2 T 1 ≤ 1 + |S|(n 2 -n 1 ) n 2 (n 1 -|S|)   1 + 12 log(1/δ) n 1 + 12 log(1/δ) n 2   ,
with probability larger than 1δ. Since |S| ≤ (n 1 ∧ n 2 )/2, we derive that

T 2 T 1 -1 ≤ 2|S| n 1 + 24 log(1/δ) n 1 + 24 log(1/δ) n 2 ≤ 3 ,
for L 1 large enough in the statement of the lemma. In conclusion, we have

-2 + (σ (1)
S ) 2 (σ (2) S ) 2 T 1 T 2 + (σ (2) S ) 2 (σ (1) S ) 2 T 2 T 1 ≥ [(σ (1)
S ) 2 -(σ (2) S ) 2 ] 2 16(σ (1) S ) 2 (σ (2) S ) 2 , ( 57
)
with probability larger than 1δ, as long as

[(σ (1)
S ) 2 -(σ (2)
S ) 2 ] 2 (σ (1) S ) 2 (σ (2) S ) 2 ≥ L |S| 2 n 2 1 + |S| 2 n 2 2 + log(1/δ) 1 n 1 + 1 n 2 . (58)
Combining ( 56), (57), and (58), we derive

-2 + (σ (1)
S ) 2 (σ (2) S ) 2 T 1 T 2 + (σ (2) S ) 2 (σ (1) S ) 2 T 2 T 1 ≥ [(σ (1)
S ) 2 -(σ (2) S ) 2 ] 2 16(σ (1) S ) 2 (σ (2) S ) 2 -L |S| 2 n 2 1 + |S| 2 n 2 2 + log(1/δ) 1 n 1 + 1 n 2 ,
with probability larger than 1δ.

Proof of Lemma 8.5. We want to lower bound the random variable F S,1 =

Rn1 (σ (1) S ) 2 T1(n1-|S|)
where R is defined by

R := X (2) S (β (2) S -β (1) S ) + Π (2) S (2) S -X (2) S (X (1) S X (1) S ) (-1) X (1) S (1) S 2 /n 2 .
Let us first work conditionally to X

S ) 2 Π

(2)

S + (σ (1) S ) 2 X (2)
S (X

S X

(1)

S ) (-1) X (2) S , we get R = X (2)
S (β

Sβ

S ) + W 2 /n 2 . We have the following lower bound:

R ≥ X

(2)

S (β (2)
Sβ

S ) + W, X

S (β

Sβ

S ) X

(2)

S (β (2) S -β (1) S ) 2 /n 2 ≥ X (2)
S (β

Sβ

S ) 2 2n 2 - 1 n 2 W, X (1)
S (β

Sβ

S ) X

(2)

S (β (2) S -β (1) S ) 2 The random variable X (2) S (β (2) S -β (1) S ) 2 / β (2) S -β (1) S 2 Σ (2) follows a χ 2 distribution with n 2 degrees of freedom. Given (X (1)
S , X

S ), W, X

S (β

S -β

S ) X

(2)

S (β (2) S -β (1) S )
2 is proportional to a χ 2 distributed random variable with one degree of freedom and its variance is smaller than (σ

(2) S ) 2 + ϕ max [X (2) S (X (1) S X (1) S ) (-1) X (2) S ](σ (1)
S ) 2 . Applying Lemma 8.11, we derive that with probability larger than 1x/6, R ≥ β

(2)

S -β (1) S 2 Σ (2) 2   1 -2 log(12/x) n 2   -4 log (12/x) n 2 (σ (2)
S ) 2 + (σ (1)
S ) 2 ϕ max {X (2)
S (X

S X

S ) (-1) X

S } .

Using the upper bound |S| ≤ (n 1 ∧ n 2 )/2 and Lemma 8.12, we control the last term

ϕ max X (2)
S (X

S X

(1)

S ) (-1) X (2) S ≤ Lϕ S n 2 n 1 , with probability larger than 1 -2 exp[-(n 1 ∧ n 2 )L ].
If we take the constant L 1 large enough in condition (51), then we get

R ≥ β (2) S -β (1) S 2 Σ (2) 4 -log (12/δ) L (σ (2)
S ) 2 n 2 + (σ (1)
S ) 2 n 1 ϕ S , (59)
with probability larger than 1δ/3.

Let us now upper bound the random variable T 1 (n 1 -|S|)/n 1 . Since (n 1 -S)T 1 follows a χ 2 distribution with n 1 -|S| degrees of freedom, we derive from Lemma 8.11 that

T 1 (n 1 -|S|)/n 1 ≤ 1 + 2 log(6/δ) n 1 + 2 n 1 log(6/δ) ≤ 2 , (60)
with probability larger than 1δ/6. Gathering (59) and (60), we conclude that

F S,1 ≥ β (2) S -β (1) S 2 Σ (2) 8(σ (1)
S ) 2 -log (6/δ) L   σ (2) S 1 n2 σ (1) S 2 + ϕ S n 1   ,
with probability larger than 1δ/2.

Proof of Theorem 4.1: Power of T B S ≤k

This proposition is a straightforward corollary of Theorem 4.3. Consider the subsets S ∨ and S ∆ of {1, . . . , p} such that S ∨ is the union of the support of β (1) and β (2) and S ∆ is the supports of β (2)β (1) . 2) | 0 allows to prove the first result. Let us turn to the second result. According to Theorem 4.3, T B S ≤k rejects H 0 with probability larger than 1δ when 2) | 0 , the second result follows.

-δ when K 1 (S ∨ ) + K 2 (S ∨ ) ≥ ϕ S∪ 1 n 1 + 1 n 2 |S ∨ | + log 1 α S∪ . Observing that ϕ S∪ ≤ ϕ Σ (1) ,Σ (2) , K 1 (S ∨ ) = K 1 , K 2 (S ∨ ) = K 2 and that |S ∨ | ≤ |β (1) | 0 + |β (
K 1 (S ∆ ) + K 2 (S ∆ ) ≥ ϕ S∆ 1 n 1 + 1 n 2 |S ∆ | + log 1 α S∆ . Since K 1 (S ∆ )+K 2 (S ∆ ) ≥ β (1) -β (2) 2 Σ 2[Var(Y (1) )∧Var(Y (2) )] and since |S ∆ | = |β (1) -β (
If S ∨ = ∅, then we can consider any subset of size 1 to prove the first result. If S ∆ = ∅, then β (1) = β (2) and the second result does not tell us anything.

Proof of Proposition 4.5

For simplicity, we assume in the sequel that β (1) = 0 or β (2) = 0, the case β (1) = β (2) = 0 being handled by any set S ∈ S 1 ⊂ S Lasso .

This proof is divided into two main steps. First, we prove that with large probability the collection S Lasso contains some set S λ close to the union S ∨ of the supports of β (1) and β (2) . Then, we show that the statistics (F S λ ,V , F S λ ,1 , F S λ ,2 ) allow to reject H 0 with large probability.

Recall that the collection S Lasso is based on the Lasso regularization path of the following heteroscedastic Gaussian linear model,

Y (1) Y (2) = X (1) X (1) X (2) -X (2) θ (1) * θ (2) * + (1) (2)
(61) which we denote for short Y = Wθ * + . Given a tuning parameter λ, θ λ refers to the Lasso estimator of θ:

θ λ = arg inf θ∈R 2p Y -Wθ 2 + λ|θ| 1 .
In order to analyze the Lasso solution θ λ , we need to control how W acts on sparse vectors.

Lemma 8.6 (Control of the design W). If we take the constants L * , L * 1 , and L * 2 in Proposition 4.5 small enough then the following holds. The event

A := ∀θ s.t. |θ| 0 ≤ k * , 1/2 ≤ X (1) θ 2 n 1 θ 2 Σ (1) ≤ 2 and 1/2 ≤ X (2) θ 2 n 2 θ 2 Σ (2) ≤ 2    κ 6, |β (1) | 0 + |β (2) | 0 , X (1) / √ n 1 κ 6, |β (1) | 0 + |β (2) | 0 , √ Σ (1) κ 6, |θ * | 0 , X (2) / √ n 1 κ 6, |θ * | 0 , √ Σ (2) ≥ 2 -3
   has probability larger than 1δ/4. Furthermore, on the event A,

Φ k,+ (W) ≤ 4(n 1 + n 2 ) Φ k,+ ( Σ (1) ) ∨ Φ k,+ ( Σ (2) ) , Φ k,-(W) ≥ (n 1 ∧ n 2 ) Φ k,-( Σ (1) ) ∧ Φ k,-( Σ (2) ) , for any k ≤ k * .
The following property is a slight variation of Lemma 11.2 in [START_REF] Van De Geer | On the conditions used to prove oracle results for the Lasso[END_REF] and Lemma 3.2 in [START_REF] Giraud | Supplement to 'Highdimensional regression with unknown variance[END_REF]. Lemma 8.7 (Behavior of the Lasso estimator θ λ ). If we take L * 2 in Proposition 4.5 small enough then the following holds. The event

B = |W T | ∞ ≤ 2(σ (1) ∨ σ (2) ) 2Φ 1,+ (W) log(p)
occurs with probability larger than 1 -1/p. Assume that

λ ≥ 8(σ (1) ∨ σ (2) ) 2Φ 1,+ (W) log(p) .
Then, on the event A ∩ B we have

W( θ λ -θ * ) 2 ≤ L 1 λ 2 /(n 1 ∧ n 2 ) κ 2 [6, |θ * | 0 , √ Σ (1) ] ∧ κ 2 [6, |θ * | 0 , √ Σ (2) ] |θ * | 0 , (62)
| θ λ | 0 ≤ L 2 n 1 ∨ n 2 n 1 ∧ n 2 Φ k * ,+ ( √ Σ (1) ) ∨ Φ k * ,+ ( √ Σ (2) ) κ 2 [6, |θ * | 0 , √ Σ (1) ] ∧ κ 2 [6, |θ * | 0 , √ Σ (2) ] |θ * | 0 ≤ k * /2 . (63)
In the sequel, we fix

λ = 16(σ (1) ∨ σ (2) ) 2(n 1 + n 2 ) Φ 1,+ ( Σ (1) ) ∨ Φ 1,+ ( Σ (2) ) log(p) .
Consider an integer k ≤ k * and a k-sparse vector θ = θ (1) θ (2) in R 2p . Under event A, we have

Wθ 2 = X (1) (θ (1) + θ (2) ) 2 + X (2) (θ (1) -θ (2) ) 2 ≤ 2n 1 θ (1) + θ (2) 2 Σ (1) + 2n 2 θ (1) -θ (2) 2 Σ (2) ≤ 4(n 1 + n 2 ) Φ k,+ ( Σ (1) ) ∨ Φ k,+ ( Σ (2) ) θ 2 Wθ 2 ≥ 1 2 n 1 (θ (1) + θ (2) 2 Σ (1) + n 2 θ (1) -θ (2) 2 Σ (2) ≥ (n 1 ∧ n 2 ) Φ k,-( Σ (1) ) ∧ Φ k,-( Σ (2) ) θ 2 .
Proof of Lemma 8.7. Observe that the variance of [W ] i given W is smaller than Φ 1,+ (W)(σ (1) ∨ σ (2) ) 2 . Using a union bound and the deviations of the Gaussian distribution, it follows that P(B) ≥ 1 -1/p. Recall the definition of η[., .] in [START_REF] Meinshausen | Lasso-type recovery of sparse representations for high-dimensional data[END_REF]. A slight variation of Lemma 11.2 in [START_REF] Van De Geer | On the conditions used to prove oracle results for the Lasso[END_REF] ensures that

W( θ λ -θ * ) 2 ≤ L λ 2 η 2 [3, |θ * | 0 , W] |θ * | 0 (64) on event B. Fix k = |θ * | 0 and consider some θ = θ (1) θ (2) ∈ C (3, T ) with |T | = k. Define T ⊂ {1, . . . , p} by i ∈ T if i ∈ T or i + p ∈ T . We have |(θ (1) + θ (2) ) T c | 1 ∨ |(θ (1) -θ (2) ) T c | 1 ≤ |θ (1)
T c | 1 + |θ (2)
T c | 1 ≤ |θ T c | 1 ≤ 3|θ T | 1 ≤ 3 |θ (1) T | 1 + |θ (2) T | 1 ≤ 6 |(θ (1) + θ (2) ) T | 1 ∨ |(θ (1) -θ (2) ) T | 1
It follows that θ (1) + θ (2) ∈ C(6, T ) or θ (1) θ (2) ∈ C(6, T ). By symmetry, we assume that |(θ (1) + θ (2) ) T | 1 ≥ |(θ (1) θ (2) ) T | 1 . Let us lower bound the l 1 norm of θ (1) + θ (2) in terms of θ.

2|θ (1) + θ (2) | 1 ≥ θ (1) + θ (2) T 1 + θ (1) -θ (2) T 1 ≥ |θ T | 1 ≥ |θ| 1 4 ,
since θ belongs to C(3, T ). Thus, we derive the lower bound

k Wθ 2 |θ| 2 1 ≥ k X (1) (θ (2) + θ (1) ) 2 |θ| 2 1 + k X (2) (θ (2) -θ (1) ) 2 |θ| 2 1 ≥ (n 1 ∧ n 2 )|θ (2) + θ (1) | 2 1 |θ| 2 1   i=1,2 η 2 6, k, X (i) / √ n i   ≥ L(n 1 ∧ n 2 )   i=1,2 κ 2 6, k, X (i) / √ n i   ≥ L(n 1 ∧ n 2 ) κ 2 6, k, Σ (1)
∧ κ 2 6, k, Σ (2) ,

where the last inequality proceeds from Lemma 8.6. We conclude that L κ 2 [3, |θ * | 0 , W] ≥ (n 1 ∧ n 2 ) κ 2 6, k, Σ (1) ∧ κ 2 6, k, Σ (2) .

Gathering this bound with (64), it follows that 1) ] ∧ κ 2 [6, |θ * | 0 , Σ (2) ] |θ * | 0 , which allows us to prove (62). Lemma 3.1 in [START_REF] Giraud | Supplement to 'Highdimensional regression with unknown variance[END_REF] tells us that on event B,

W( θ λ -θ * ) 2 ≤ L λ 2 /(n 1 ∧ n 2 ) κ 2 [6, |θ * | 0 , √ Σ (
λ 2 | θ λ | 0 ≤ 16Φ | θ λ |0,+ (W) W( θ λ -θ * ) 2 .
Gathering the last two bounds and Lemma 8.6, we obtain

| θ λ | 0 ≤ L Φ | θ λ |0,+ (
| θ λ | 0 ≤ L n 1 ∨ n 2 n 1 ∧ n 2 Φ k * ,+ ( √ Σ (1) ) ∨ Φ k * ,+ ( √ Σ (2) ) κ 2 [6, |θ * | 0 , √ Σ (1) ] ∧ κ 2 [6, |θ * | 0 , √ Σ (2) ] |θ * | 0 1 + | θ λ | 0 k * ≤ k * + | θ λ | 0 /2 ,
where the last inequality holds if we take L * 2 in (37) small enough. Hence, | θ λ | 0 ≤ k * . Coming back to (65), we prove (63).

Proof of Lemma 8.8. Given the Lasso estimator θ λ of θ * in model (61), we define β

λ = θ (2)
λθ

λ .

On event A ∩ B, we upper bound the difference between (β (1) , β (2) ) and ( β

λ , β

λ ).

β (1) β

(1) λ 2 Σ (1) + β (2) -β (2) λ 2 Σ (2)
≤ 2 X (1) √ n 1 (β (1) β

(1)

λ ) 2 + X (2) √ n 2 (β (2) -β (2) λ ) 2 ≤ 2 n 1 ∧ n 2 W(θ * -θ λ ) 2 ≤ L i=1,2 Φ 1,+ ( √ Σ (i) ) i=1,2 κ 2 [6, |θ * | 0 , √ Σ (i) ] |S ∨ |(n 1 ∨ n 2 ) (n 1 ∧ n 2 ) 2 log(p)(σ (1) ∨ σ (2) ) 2 ,
where the last inequality follows from Lemma 8. (σ

(1) S λ ) 2 -2 ≥ (σ (1) ) 2 (σ (2) ) 2 + (σ (2) ) 2 (σ (1) ) 2 /(1 + A) -2 ≥ (σ (1) ) 2 (σ (2) ) 2 + (σ (2) ) 2 (σ (1) ) 2 -2 - A 1 + A (σ (1) ) 2 (σ (2) ) 2 + (σ (2) ) 2 (σ (1) ) 2 ≥ (σ (1) ) 2 (σ (2) ) 2 + (σ (2) ) 2 (σ (1) ) 2 -2 -3A . ( 67
)
Let us now lower bound the remaining part of K 1 ( S λ )+K 2 ( S λ ). For i = 1, 2, |β (i)β β (1) β (2) 2 Σ (2) (σ (1) ) 2 + β (1)β (2) 2 Σ (1) (σ (2) ) 2 ≤ 3 (σ (1)

) 2 ∧ (σ (2) ) 2 2 i=1 β (1) -β (1) S λ 2 Σ (i) + β (2) -β (2) S λ 2 Σ (i) + β (1) S λ -β (2) S λ 2 Σ (i) ≤ L 1   β (1) S λ -β (2) S λ 2 Σ (1) (σ (2) ) 2 + β (1) S λ -β (2) S λ 2 Σ (2) (σ (1) ) 2   + L 2 (σ (1) ∧ σ (2) ) 2 i=1,2 Φ k * ,+ ( √ Σ (i) ) i=1,2 Φ k * ,-( √ Σ (i) ) 2 i=1 β (i) -β (i) λ 2 Σ (i) ≤ L 1   β (1) S λ -β (2) S λ 2 Σ (1) (σ (2) ) 2 + β (1) S λ -β (2) S λ 2 Σ (2) (σ (1) ) 2   + L 2 i=1,2 Φ k * ,+ ( √ Σ (i) ) i=1,2 Φ k * ,-( √ Σ (i) ) A
Gathering the last inequality with (67) yields

K 1 ( S λ ) + K 2 ( S λ ) ≥ L 1 [K 1 + K 2 ] -L 2 i=1,2 Φ k * ,+ ( √ Σ (i) ) i=1,2 Φ k * ,-( √ Σ (i) ) A .
Proof of Lemma 8.9. For any non empty set S of size smaller or equal to k * , define δ S = δ 2 |S| p k * -1

. If we take L * and L * 1 in (35-36) small enough, then 1+log[1/(α S δ S )]/(n 1 ∧ n 2 ) is smaller than some constant L small enough so that we can apply Theorem 4. [START_REF] Arias-Castro | Global testing under sparse alternatives: Anova, multiple comparisons and the higher criticism[END_REF] This proof follows the same steps as above. Taking L * small enough, we can assume that n 1 ∨ n 2 ≤ 2(n 1 ∧ n 2 ). Rewrite the linear regression model Y = Wθ * + as follows: 1) θ

Y = W (
(1) * + W (2) θ

(2) * + .

From the definition of the Lasso estimator θ λ = θ

We fix λ = 16(σ (1) ∨ σ (2) ) 2(n 1 + n 2 )Φ 1,+ ( √ Σ) log(p) .

and we suppose that event A ∩ B (defined in the last proof) holds. Recall that P[A ∩ B] ≥ 1δ/4 -1/p. We consider the set S

(2) λ defined as the support of θ

λ . Note that S L ⊂ S Lasso . Lemma 8.10. If we take constants L * and L * 2 in Proposition 4.6 small enough, then the following holds. There exists an C of probability larger than 1 -1/p such that, under A ∩ B ∩ C, we have |W (2) W (1) θ

(1) *θ

(1) λ | ∞ ≤ λ/8 (69)
It follows that on A ∩ B ∩ C:

W (2)
+ W (1) θ

(1) * -θ (1) λ ∞ ≤ λ/4
Arguing as in the proof of Lemma 8.7 and taking L * 2 small enough, we derive that under A ∩ B,

W (2) (θ (2) * -θ (2) λ ) 2 ≤ L 1 λ 2 /(n 1 ∧ n 2 ) κ 2 [6, k * , √ Σ] |θ (2) * | 0 , (70) | θ
(2)

λ | 0 ≤ L 2 Φ k * ,+ ( √ Σ) κ 2 [6, k * , √ Σ] |θ (2) * | 0 ≤ k * /2 ≤ k * /2 . (71)
This allows us to upper bound θ

(2) *θ

Σ under event A.

θ (2) * -θ (2) λ 2 Σ ≤ L n 1 ∧ n 2 X (1) (θ (2) * -θ (2)
λ ) 2 + X (2) (θ

(2) * -θ (2) λ ) 2 ≤ L n 1 ∧ n 2 W (2) (θ (2) * -θ (2) λ ) 2 .
Pythagorean inequality then gives

β (1) -β (2) 2 Σ = β (1)
S (2) λ -β (2)

S

(2) λ 2 Σ + β (1) β (2)β

S

(2) λ + β

(2)

S (2) λ 2 Σ ≤ β (1)
S (2) λ -β (2)

S

(2)

λ 2 Σ + θ (2) * -θ (2) λ 2 Σ ≤ β (1)
S (2) λ -β (2)

S

(2)

λ 2 Σ + L |θ (2) * | 0 log(p) n 1 ∧ n 2 Φ 1,+ ( √ Σ) κ 2 [6, k * , √ Σ] (σ (1) ∨ σ (2) ) 2 ,
where we use the two previous upper bounds in the last line. Consequently, we obtain

K 1 ( S (2)
λ ) + K 2 ( S (2)
λ ) ≥ L β (1) -β (2) 2 Σ
Var(Y (1) ) ∨ Var(Y (2) ) -L |θ

(2) * | 0 log(p) n 1 ∧ n 2 Φ 1,+ ( √ Σ) κ 2 [6, k * , √ Σ] .
Applying Lemma 8.9 to S

λ , using (71) and taking L * 3 large enough then allows us to conclude.

Proof of Lemma 8.10. Given any matrix A, we define the norm A ∞ = max i,j |A i,j |. Suppose that we are under events A ∩ B defined previously. Arguing as in the proof of Lemma 8. (1) *θ

(1) λ ∞ = X (1) X (1) -X (2) X (2) θ

(1) *θ

(1) λ ∞ ≤ θ

(1) *θ

(1) λ k * X (1) X (1) -X (2) X (2)

∞ ≤ W(θ * -θ) Φ k * ,-(W) k * X (1) X (1) -X (2) X (2) ∞ ≤ L λ k * X (1) X (1) -X (2) X (2) ∞ √ n 1 ∧ n 2 κ[6, |θ * | 0 , √ Σ] Φ k * ,-(W) , (73)
where we used (72) in the last line.

Combining deviations inequality for χ 2 distributions (Lemma 8.11) and for Gaussian distributions and a union bound, we derive that

X (1) X (1) -X (2) X (2) ∞ ≤ Φ 1,+ ( √ Σ) |n 1 -n 2 | + L (n 1 ∨ n 2 ) log(p) , (74)
with probability larger than 1 -1/p. Consider some θ with |θ| 0 ≤ k * . When event A defined in Lemma 8.6 holds, we have By symmetry, we can assume that n 1 ≤ n 2 . Let us fix β (2) = 0 and σ (2) = 1. Fix some positive integer s ≤ p 1/2-γ and fix r ∈ (0, 1/ √ 2). We consider the test of hypotheses H 0 : β (1) = 0, σ (1) = 1 against H 1 : |β (1) | 0 = s, β (1) = r 2 , and σ (1) = √ 1r 2 . Note that for this problem, the data (Y (2) , X (2) ) do not bring any information on the hypotheses. This one-sample testing problem is a specific case of the two-sample testing problem considered in the proposition. Thus, a minimax lower bound for the one-sample problem provides us a minimax lower bound for the two-sample problem.

Wθ 2 θ 2 = X (1) (θ (1) + θ (2) ) 2 θ 2 + X (2) (θ (1) -θ (2) ) 2 θ 2 ≥ Φ k * ,-( √ Σ) 2 n 1 θ (1) + θ (2) 2 + n 2 θ (1) -θ (2) 2 θ 2 ≥ Φ k * ,-( √ Σ)(n 1 ∧ n 2 ) .
According to Theorem 4.3 in [START_REF] Verzelen | Goodness-of-fit tests for high-dimensonal gaussian linear models[END_REF], no level α test has power larger than 1δ if

r 2 1 -r 2 ≤ s 2n 1 log 1 + p s 2 + 2p s 2

S

  ). Consider the linear regression of Y(1) by X

Proposition 3 . 5 (

 35 Size of T P S ). The test function T P S satisfies α/2 ≤ P H0 T P S = 1 ≤ α .

Theorem 4 . 1 (

 41 Power of T B S ≤k ). Assuming that A.1 and A.2 hold, P[T B S ≤k = 1] ≥ 1δ as long as

Proposition 4 . 6 .

 46 Let us assume that Σ(1) = Σ(2) := Σ. There exist five positive constants L * , L * , L * 1 , L * 2 , and L * 3 such that following holds. Define k * and k * as the largest positive integers that satisfy

Remark 4 . 10 .

 410 The definition[START_REF] Rouzier | Microtubule-associated protein tau: a marker of paclitaxel sensitivity in breast cancer[END_REF] of k * together with Condition (39) restrict the number of non-zero components |β(1) | 0 + |β(2) | 0 to be small in front of (n 1 ∧ n 2 )/|n 1n 2 |. This technical assumption enforces the design matrix in the reparametrized model (8) to be almost block-diagonal and allows us to control efficiently the Lasso estimator θ

S

  follows a Fisher distribution with |S| and n 1 + n 2 -2|S| degrees of freedom:

5 ( 2

 52 ± 1.4 6.8 ± 1.6 0.5 ± 0.4 6.5 ± 1.5 n= 50 4.5 ± 1.3 5.5 ± 1.4 0.1 ± 0.2 5.3 ± 1.4 n= 100 4.8 ± 1.3 6.6 ± 1.5 0.1 ± 0.2 6.5 ± 1.Estimated test levels in percentage along with 95% Gaussian confidence interval (in percentage) under H 00 based upon 1000 simulations.

Figure 1 .Figure 2 .

 12 Figure 1. Estimated test levels in percentage under H 0 for varying magnitudes of common non null coefficients, based upon 1000 simulations. Bonferroni calibration in dotted lines, calibration by permutation in plain lines. Blue squares represent the suggested test T * S , red triangles stand for the Fisher

Figure 3 .

 3 Figure 3. Estimated test levels in percentage under H 0 for varying magnitudes of common non null coefficients, based upon 1000 simulations, under power decay and GGM correlation structures when n = 50. Bonferroni calibration in dotted lines, calibration by permutation in plain lines. Blue squares represent the suggested test T * S , red triangles stand for the Fisher test T * ,Fisher S

Figure 4 .

 4 Figure 4. Power (in percentage) as a function of signal magnitude parameter r for various sparsity pattern under the assumption of uncorrelated designs Σ (1) = Σ (2) = Ip. Results for the suggested test T P S

Figure 5 .

 5 Figure 5. Power (in percentage) as a function of signal magnitude parameter r for various sparsity patterns under power decay and GGM correlated designs, at n = 50 observations. Results for the suggested test T P S and the test T P,FisherS

3 .

 3 Distributions of F S,V , F S,1 and F S,2 (Proposition 3.1)

4 .

 4 By definition of the p-values Q i,|S| , we have under H 0 for each S ∈ S and each i ∈ {V, 1, 2}

• 1 , L • 2 ,

 2 and L • 3 only depend on L * 2 in (53) and on the constants introduced in Lemmas 8.1-8.5. These conditions allow us to fix the constants in the statement (29) of Theorem 4.3. Lemma 8.1 (Upper-bound of Q -1

S

  . Upon defining the Gaussian vector W by W ∼ N 0, (σ

W) κ 2 [ 6 ,

 26 |θ * | 0 , √ Σ (1) ] ∧ κ 2 [6, |θ * | 0 , Σ (2) ] |θ * | 0 . (65) Recall that |θ * | 0 ≤ |β (1) | 0 +|β (2) | 0 . The upper-bound Φ | θ λ |0,+ (W) ≤ (1+| θ λ | 0 /k * )Φ k * ,+ (W)and Lemma 8.6 enforce

  |θ * | 0 + | θ λ | 0 ≤ k * byLemma 8.7 and Condition (37).

+

  W(2) θ (2) * + W (1) θ (1) *θ (1) λ -W (2) θ + λ|θ | 1 .

λ

  7, we derive that |θ* | 0 + | θ λ | 0 ≤ k * and W( θ λθ * ) 2 ≤ L 1 λ 2 κ 2 [6, |θ * | 0 , √ Σ](n 1 ∧ n 2 ) k * . | 0 ≤ k * and we derive W (2) W (1) θ

|n 1 -n 2 | n 1 ∧ n 2 + log(p) n 1 ∧ n 2  T Σ , since n 1 ∨ 8 . 9 .

 1222189 n 2 ≤ 2(n 1 ∧ n 2 ). Taking L * small enough in definition (38) of k * allows us to conclude. Proof of Proposition 4.2

Table 3

 3 Homogeneity test between training and test samples among pCR patients. Summary of test decisions after Bonferroni multiple testing correction and empirical p-values for each neighborhood test as defined in Section 3.4.

		AMFR	BB S4	BECNI	BTG3	CA12	CTNND2	E2F3
	decision	0	0	0	0	1	0	0
	p empirical	0.0492	0.0072	0.1972	1	0.0018	0.0100	0.1080
		ERBB4	FGFRIOP FLJ10916	FLJI2650	GAMT	GFRAI	IGFBP4
	decision 0	0	0	0	0	0	0
	p empirical	0.5610	0.0242	0.2542	0.0312	0.1158	0.5318	0.0458
		JMJD2B	KIA1467	MAPT	MBTP SI	MELK	METRN	PDGFRA
	decision	0	0	0	0	0	0	1
	p empirical	0.0128	0.0272	0.0178	0.0062	0.5602	1	0.0012
		RAMPI	RRM2	SCUBE2	THRAP2	ZNF552
	decision	0	0	0	0	0
	p empirical	0.0444	0.0022	0.2372	0.0228	0.0028
		AMFR	BB S4	BECNI	BTG3	CA12	CTNND2	E2F3
	decision	0	1	1	0	1	0	0
	p empirical	0.0046	<0.0001	<0.0001	0.0202	<0.0001	0.0684	0.0428
		ERBB4	FGFRIOP FLJ10916	FLJI2650	GAMT	GFRAI	IGFBP4
	decision	0	1	1	1	1	0	0
	p empirical	0.26	<0.0001	<0.0001	0.002	<0.0001	0.3606	0.389
		JMJD2B	KIA1467	MAPT	MBTP SI	MELK	METRN	PDGFRA
	decision	1	1	0	1	0	0	1
	p empirical	<0.0001	2e-04	0.006	6e-04	0.1556	0.1054	<0.0001
		RAMPI	RRM2	SCUBE2	THRAP2	ZNF552
	decision	0	0	0	1	1
	p empirical	0.2288	0.2988	0.3552	<0.0001	<0.0001

Table 4

 4 Homogeneity test between training and test samples among RD patients. Summary of test decisions after Bonferroni multiple testing correction and empirical p-values for each neighborhood test as defined in Section 3.4

Table 6

 6 developed a therapeutic candidate inhibiting MELK expression that was proved to suppress the growth of tumour-initiating cells in mice with various cancer types, including breast cancer. Summary of neighborhood tests between RD and pCR patients within the validation set (n=51).Decision is made at level 0.05/26 to correct for multiple testing. The empirical p-value and the rejected model are defined in Section 3.4.

		AMFR	BB S4	BECNI	BTG3	CA12	CTNND2	E2F3
	decision	0	0	0	1	0	0	1
	p empirical	0.0024	0.0028	0.0048	0.0018	0.0028	0.0082	<0.0001
	rejected model	-	-	-	SCUBE2	-	-	METRN
		ERBB4	FGFRIOP FLJ10916	FLJI2650	GAMT	GFRAI	IGFBP4
	decision	1	0	1	0	0	1	1
	p empirical	0.0014	0.0072	8e-04	0.0142	0.0046	8e-04	2e-04
	rejected model SCUBE2	-	SCUBE2	-	-	E2F3	SCUBE2
		JMJD2B	KIA1467	MAPT	MBTP SI	MELK	METRN	PDGFRA
	decision	0	1	0	0	0	1	0
	p empirical	0.0054	0.0018	0.0032	0.0078	0.0036	4e-04	0.0104
	rejected model	-	SCUBE2	-	-	-	E2F3	-
		RAMPI	RRM2	SCUBE2	THRAP2	ZNF552
	decision	0	0	1	0	0
	p empirical	0.0056	0.0034	2e-04	0.0024	0.006
	rejected model	-	-	FLJ10916	-	-

  1,|S| (x|X S ), and Q -1 2,|S| (x|X S ) in Lemmas 8.1, 8.2 and 8.3. Then, we control the deviations of F S,V , F S,1 , and F S,2 under H 1,S in Lemmas 8.4 and 8.5. In the sequel, we call S the subcollection of S made of subsets S satisfying |S| ≤ (n 1 ∧ n 2 )/2 and log(12/δ

  Lemma 8.3 (Upper-bound of |a| ∞). There exist two positive universal constants L 1 and L 2 such that the following holds. Consider δ a positive number sastifying L 1 log(4/δ) < n 1 ∧ n 2 . With probability larger than 1δ/2, we have

	log(1/x) + 8 log(1/x) .	(49)

  Assume first that S ∨ and S ∆ are non empty.

			By Definition
	(19) of the weights, we have
	log	1 α i,S∪	≤ log(4k) + log(1/α) + |S

∨ | log(p) ≤ 2|S| ∪ log(p) + log(1/α) .

A similar upper bound holds for log(1/α i,S∆ ). If we choose the numerical constants large enough in Conditions A.1 and A.2, then the sets S ∨ and S ∆ follow the conditions of Theorem 4.3.

Applying Theorem 4.3, we derive that T B S ≤k rejects H 0 with probability larger than 1

  7. Let us now lower bound the Kullbackdiscrepancy 2[K 1 ( S λ ) + K 2 ( S λ )] which equals CASE 1: σ (1) ∨σ (2) σ (1) ∧σ (2) ≥ √ 2.By symmetry, we can assume that σ (1) > σ(2) .) 2 = (σ (1) ) 2 + β (1)β ,2 κ 2 [6, |θ * | 0 , √ Σ (i) ] |S ∨ |(n 1 ∨ n 2 ) (n 1 ∧ n 2 ) 2 log(p)(σ (1) ) 2 (66) ≤ (σ (2) ) 2 + (σ (1) ) 2 4 ,where we used conditions (35) and[START_REF] Rouzier | Breast cancer molecular subtypes respond differently to preoperative chemotherapy[END_REF] in the last inequality assuming that we have taken L * and L * 2 small enough in these two conditions. This enforces2 K 1 S λ + K 2 S λ ≥ 1 12 .Let us lower bound K 1 ( S λ ) + K 2 ( S λ ) in terms of K 1 + K 2 . First, we consider the ratio of the variances

	(σ S λ (1)											(1) S λ	2 Σ (1) ≥ (σ (1) ) 2
	(σ S λ (2)	) 2 = (σ (2) ) 2 + β (2) -β ≤ (σ (2) ) 2 + L i=1,2 Φ 1,+ ( (2) S λ 2 Σ (2) ≤ (σ (2) ) 2 + β (2) -β √ Σ (i) )	(2) λ	2 Σ (2)
	CASE 2: σ (1) ∨σ (2) σ (1) ∧σ (2) ≤ A = 2L √	2. Let us note i=1,2 Φ 1,+ ( √	Σ (i) )
	(1) S (σ (σ (1) S λ ) 2 (σ (2) S λ ) 2 + (σ (2) S λ ) 2
	 	σ σ	(1) S λ (2) S λ	  2	+	 	σ S λ (1) σ (2) S λ	 	2	-2 +	β	(2) S λ (σ -β (1) (1) S λ S λ ) 2	2 Σ (2)	+	β	(2) S λ (σ -β (2) (1) S λ S λ ) 2	2 Σ (1)	.

i=1i=1,2 κ 2 [6, |θ * | 0 , √ Σ (i) ] |S ∨ |(n 1 ∨ n 2 ) (n 1 ∧ n 2 ) 2 log(p) ,

with L as in (66). Arguing as in Case 1, we derive that (σ

(1) S λ ) 2 ≤ (σ (1) ) 2 [1 + A] ≤ 2(σ (1) ) 2 , λ ) 2 ≤ (σ (2) ) 2 [1 + A] ≤ 2(σ (2) ) 2 .

  . Arguing as in the proof of this Theorem, we derive that Applying a union bound over all sets S of size smaller or equal to k * allows us to prove i, S λ (F S λ ,i |X S λ ) < α S λ ≥ 1δ . 8.8. Proof of Proposition 4.6

		P	min i∈{V,1,2}
			1 n 1	+	1 n 2	|S| log(p) + log	1 αδ	+ log(p) .
	P	min i∈{V,1,2}

Q i,S (F S,i |X S ) < α S ≥ 1δ S if K 1 (S) + K 2 (S) ≥ Lϕ S Q

Acknowledgements

We are grateful to Christophe Giraud, two anonynous reviewers, the Associate Editor and the Editor for careful suggestions. We also thank Nicolas Städler for providing us the code of DiffRegr. The research of N. Verzelen is partly supported by the french Agence Nationale de la Recherche (ANR 2011 BS01 010 01 projet Calibration).

For comprehensiveness, we add that similar analysis of the validation set (n=51 patients, among which 38 RD and 13 pCR patients) leads to the identification of only 9 significantly altered neighborhoods between pCR and RD patients 6. This difference in the number of significantly altered neighborhoods can be explained by the reduced size of the sample. Yet, genes responsible for the rejection of the tests differ from those identified on the training set. In particular, five of the significant tests are rejected because of SCUBE2, which has been recently recognised as a novel tumor suppressor gene [START_REF] Lin | Domain and functional analysis of a novel breast tumor suppressor protein, scube2[END_REF].

Discussion

Design hypotheses. In this work, we have made two main assumptions on the design matrices:

(i) The design matrices X (1) and X (2) are random. (ii) Under the null hypothesis (3), we further suppose that the population covariances Σ (1) and Σ (2) are equal. and we consider the set S λ defined by the union of the support of θ

λ and θ

(2)

λ . On the event A ∩ B, Lemma 8.7 tells us that | S λ | ≤ k * . Thus, S λ belongs to the collection S Lasso . We shall prove that min i∈{V,1,2}

with probability larger than 1-δ/2. In the following lemma, we compare

.

Lemma 8.8. On the event A ∩ B, we have

Then, we closely follow the arguments of Theorem 4.3 to state that T B SLasso rejects H 0 with large probability as long as K 1 ( S λ ) + K 2 ( S λ ) is large enough. Lemma 8.9. If on the event A ∩ B, we have

We derive from (63) that on the event A ∩ B,

Since |S ∨ | ≤ |β (1) | 0 + |β (2) | 0 , it follows from Condition (37) that | S λ | ≤ k * . Gathering Lemmas 8.8 and 8.9 allows us to conclude if we take L * 3 in Proposition 4.5 large enough. Proof of Lemma 8.6. In order to bound P(A), we apply Lemma 8.12 to simultaneously control ϕ max (X

S ), and ϕ min (X

S ) for all sets S of size k * . Combining a union bound with Conditions [START_REF] Natowicz | Prediction of the outcome of preoperative chemotherapy in breast cancer using dna probes that provide information on both complete and incomplete responses[END_REF] and [START_REF] Raskutti | Restricted eigenvalue properties for correlated Gaussian designs[END_REF] allows us to prove that

Applying Corollary 1 in [START_REF] Raskutti | Restricted eigenvalue properties for correlated Gaussian designs[END_REF], we derive that there exist three positive constant c 1 , c 2 and c 3 such that the following holds. With probability larger than 1 -

. Hence, we conclude that P [A] ≥ 1δ/4.

Since s ≤ p 1/2-γ , no level α test has power larger than 1δ if

By Assumption (A.2), one may assume that that the right-hand side term is smaller than 1/2. Observe that

The result follows.

Technical lemmas

In this section, some useful deviation inequalities for χ 2 random variables [START_REF] Laurent | Adaptive estimation of a quadratic functional by model selection[END_REF] and for Wishart matrices [START_REF] Davidson | Local operator theory, random matrices and Banach spaces[END_REF] are reminded.