A Conditional Randomization Test for Sparse Logistic Regression in High-Dimension
Authors: Binh Nguyen
Bertrand Thirion
Sylvain Arlot
Venue: NeurIPS 2022 Accept (2022)
Keywords: high-dimension statistics
sparse logistic regression
variable selection
variable importance
brain-imaging analysis
statistical inference

Abstract: Identifying the relevant variables for a classification model with correct confidence levels is a central but difficult task in high-dimension. Despite the core role of sparse logistic regression in statistics and machine learning, it still lacks a good solution for accurate inference in the regime where the number of features $p$ is as large as or larger than the number of samples $n$. Here we tackle this problem by improving the Conditional Randomization Test (CRT). The original CRT algorithm shows promise as a way to output p-values while making few assumptions on the distribution of the test statistics. As it comes with a prohibitive computational cost even in mildly high-dimensional problems, faster solutions based on distillation have been proposed. Yet, they rely on unrealistic hypotheses and result in low-power solutions. To improve this, we propose \emph{CRT-logit}, an algorithm that combines a variable-distillation step and a decorrelation step that takes into account the geometry of $\ell_1$-penalized logistic regression problem. We provide a theoretical analysis of this procedure, and demonstrate its effectiveness on simulations, along with experiments on large-scale brain-imaging and genomics datasets.

A Conditional Randomization Test for Sparse Logistic Regression in High-Dimension
Binh T. Nguyen
Bertrand Thirion
LTCI, Télecom Paris, IP Paris tuanbinhs@gmail.c com
Université Paris-Saclay, Inria, CEA, Palaiseau 91120, France bertrand. thirion@inria.fr
Sylvain Arlot
Université Paris-Saclay, CNRS, Inria, Laboratoire de mathématiques d' d'Orsay, 91405, Orsay, France sylvain. cellottuniversitte wewprnneeeeeeeeii--ieee-eeeeeeee -a--p--rsssaaaay
Abstract
Identifying the relevant variables for a classification model with correct confidence levels is central but difficult task in high-dimension. Despite the role of sparse logistic regression in statistics and machine learning, it still lacks a good solution for accurate inference in the regime where the number of features is as large as or larger than the number of samples n. Here we tackle this problem by improving the Conditional Randomization Test (CRT). The original CRT algorithm shows promise as a way to output p-values while making few assumptions on the distribution of the test statistics. As it comes with prohibitive computational cost even in mildly nigh-dimensional problems, faster solutions based on distillation have been proposed. Yet, they rely on unrealistic hypotheses and result in low-power solutions. To improve this, we propose CRT-logit, an algorithm that combines variable-distillttion step and a decorrelation step that takes into account the geometry of the é -penaalieed logistic regression problem. We provide theoretical analysis of this procedure, and demonstrate its effectiveness simulations, along with experiments on large-scale brain-imaging and genomics datasets.
Introduction
Logistic regression is one of the most popular tools in modern applications of statistics and machine learning, partly due to its relative algorithmic simplicity. The method belongs to the class of generalized linear models that handle discrete outcomes, i.e. classification problems. Here, focus on the binary classification problem, where one observation of the responses {0, 1} and the data vectors [R follows the relationship:
Pl = x) = g(x 1 + exp(-xx
(I)
where g(x: 1/(1 + exp(-x)) is the sigmoid function, and 3 the vector of true regression coefficients. In the classical setting, in which the number of samples is greater than the number of features p. estimate & of the true signals @ can be obtained using maximum likelihood estimation (MLE). The asymptotic behaviour and derivation of the test statistic, confidence intervals and p-values of the MLE have been well studied, e.g. in [33. The availability of pavuuu for the test statistics makes it possible to rely on multiple hypothesis testing, where one wants to test which variables have non-zero effect the outcome, conditioned the remaining variables. Unfortunately, this line of analysis cannot be applied to the nigh-dimensional regime, where p is larger than n, as argued in (25 32 34 These works show that in the regime lim,, n/p H, the MLE estimator
36th Conference on Neural Information Processing Systems (NcurIPS 2022).

---

exists only when & > 2. However, we note that this type of analysis is done without the addition of & érregularization to the likelihood function, i.e. without using a penalized estimator to enforce sparsity.
Motivation Our focus in this paper is to do inference with statistical guarantees on high-dimensional sparse logistic regression, where p is larger much larger than n. This setting is typical in modern applications of pattern recognition, egg. in brain-imaging or genomics [3, with as large hundreds of thousands -compressible to thousands- but stays most few thousand. The family of methods we consider is the Con- ditional Randomization Test (CRT) [10]. CRT
dCRT
CRT-Logit
relies on generating multiple noisy copies of original variables to output empirical p-values in high- high-dimensional inference problems, However, prohibitive computational cost makes CRT im-
Sample Quantiles
Bamull Quentille
4
practical, as discussed at length in [10 26 7 19
There have been several lines of research at- tempting to fix this problem, most notably the distilled Conditional Randomization Test (dCRT) 191 This work introduced a distilla- tion step as a replacement for the randomized sampling step to compute the importance statis- tics (see Section for more details). It provides way to output p- p-values for multiple types of re- gression and classification problems, assuming convergence Gaussian distribution of the test statistic in large-sample regime. Yet, as shown in the left panel of Figure [] the originally pro- posed dCRT test-statistic for logistic regression does not behave well as intended. In particular, its null distribution deviates markedly from standard normal in high-dimension whenever n/p S 1.
-2 Theoretical Quantiles
-2
2
Theoretical Quantiles
Figure I: QQ-Plot for 1000 samples of test-statistic of null index for logistic regression, with simulated data, n = 200, p 400, Left: Statistics obtained from running Distilled-CRT, and Right: from our proposed algorithm. The empirical distribution of the dCRT null- statistic strays far from theoretical distribution, which standard normal, while empirical distribution of CRT- logit's 's null test score is much closer.
We propose correction for the dCRT, inspired by the decorrelation method presented in 21 The decorrelation step makes the null-distribution of the test statistics much closer to standard normal, as shown on the right panel of Figure [ and thus increases the statistical power of the method. We provide asymptotic analysis of this method, which shows that CRT-logit produces standard normal stst-statistics in the large-sample regime. In addition, we validate the high performance of CRT-logit on large-scale brain-imaginn and genetics datasets, thus showing its usefulness in practical applications.
Contribution
Related works The closest cousin of the Conditional Randomization Test Knockoff Filter recent breakthrough in the False Discovery Rate (FDR) control literature. It relies on the creation of additional noisy features, called knockoffs, to calculate variable- iimpottanceee statisties. Another extension of vanilla CRT is the Holdout Randomization Test (HRT) [26]. While still requiring multiple samplings of noisy variables. HRT solves the computational issue of original CRT by doing heavy model fitting only once on one part of the dataset, and test statistics calculation on the other part, without refitting the model. However, this method relies on sample-splittin, and hence inherently suffers from loss of statistical power. A parallel line of work has introduced the Conditional Permutation Test (CPT) 77. a non-parametric alternative to CRT that relies on random shufiling mechanism applied to original variables, instead of multiple sampling of new variables. This potentially makes CPT more robust to model mis-specification. 332 recently proposed method called SLOE, which adapts the analysis of 34], but in regime different from what we considering, where lim,,p+oo p- n/p - x E 2). and more importantly without pparsity-induodi penalty. On a separate note, we notice the similarity of dCRT [19] with debiased Lasso [16 28 33). This line of work proposed a debiasing formula for the estimator, which makes the asymptotic distribution of (GLASSO - 3r) standard normal, so that one can compute the test statistic and p-value associated with each variable.

---

2 Background
Notation We denote matrices, vectors, scalars and sets by bold uppercase, bold lowercase, script lowercase and calligraphic letters, respectively, e.g. X, X, 2, and X. The i-th row of matrix X will be denoted X.,. the j-th column X., and the (i.j)-th element Xij. For any natural number p. we denote the set [p] def. {21, ..1,,.........,,,,,,... ...i.....,,....
= {l..... p}. For each x E I' and j E [p]. we denote x-j def. p dimensional observation after removing the j-th variable. Correspondingly. x_, is the data matrix X E R=xP with column X removed. The cumulative distribution function (CDF) of the standard Gaussian distribution will be denoted ((-). The indicator function of a random event A will be denoted 14. For any two positive sequences In and Uns we write Fn X Yn if cyn. S In I Cyn for all n, for some positive constants с and C. For vector x. ||*||p denotes its Ep norm. For a function f IRP + R. v denotes its gradient w.r.l. the j-th variable, for [p].
Problem setting We consider exclusively binary classification, where the response vector is denoted y E (0, 13"h and the data matrix X E RRx consists of n p-dimensional samples. Throughout the paper, we assume the data {X..) }...)== are i..d and follow distribution with zero mean and population covariance matrix E. Moreover, we assume that ,,, and y, follow the logistic relationship Eq. ①
We denote the support set S def. E E [p] B8 * 0} and assume that it is sparse, i.e. card(S)
< p, where card denotes the cardinality of a set. Furthermore, $ def. {j E [p] : -8 * 0} indicates
an estimation of S, where 3, is an estimate of the true signal 3g. We try to obtain it through penaliize logistic estimator:
3x = argmin 1(3) + A||B111 with ((3) = 1 Σ (X,,, - log [1 + exp(Xx..
(2)
BERP
2=1
( the Fisher information matrix, and ز the partial Fisher in- formation, defined by I راا B (v,
We denote I def.
L, - -,- - --j , where ,,.,, is the row-vector made with the jth-row and the columns corresponding to -_
the sub-matris of made with the rows and columns corresponding to 3_). This quantity, defined following [13 pp. 323], plays an important role in our proposed method, detailed in Section 3
Statistical control with False Discovery Rate To quantify statistical errors, consider the False
Discovery Rate, introduced in [5]. Given estimate of the support S, the false discovery proportion (FDP) is the ratio of the number of selected features that do not belong to the true support S, divided by the total number of selected features. The False Discovery Rate is the expectation of the FDP:
FDP(S) card(lj : e S.j $ S))
and
FDR($) ((FD((S))
card($) V1
Conditional Randomization Test (CRT) and Distillation CRT (dCRT) The concept of Condi- tional Randomization Test was originally proposed in the model-X knockoff paper (10) way to output valid empirical p-values using knockoff variables. The principle of the knockoff filter is first to sample noisy copies X., of variable X..j. given a known sampling mechanism Pj I-j- One advantage of the knockoff filter is that specific assumption is placed the distribution of the inferred test statistic. However, this means that, in general, there is no mechanism to derive p-values from the knockoff statistic. This motivates the introduction of CRT, which requires running nigh-dimensional inference for each variable B times. However, the computation cost of CRT is prohibitive when p grows large: assuming that we use the Lasso program with coordinate descent to compute TCRT, its runtime would be O(Bp*) [15] pp. 93]. Moreover, CRT requires decently large B to make the empirical distribution of the p-values smooth enough. Reducing the computational cost of CRT is the main motivation of several works [7] 19 26). One of them is the introduction of distillation-CRT (dCRT) 19). The main appeal of this method is that it can output p- p-values analytically, therefore bypassing the multiple knockoffs sampling steps used to infer on each variable, and leads reasonable reduction of the computation cost.
3

---

Distillation operation The key addition of dCRT is the distillation operation: for each variable j. we want to distill all the conditional information of the remaining variables X-j to x, and y via least-squares minimization with & fe-egularrzation to enforce sparsity. For each variable j. we first solve the lasso problem by regressing X., on X
argmin BERP-1 2 Xx., X- ,x, 2 Aaz 11111
For distillation of variable j and the binary response y with logistic relationship, [9 briefly suggested to solve a penalized estimation problem, similar to Eq. (
B = argmin BER
Σ (X (X,T_ ;)yy log [1 1 exp(X 3]] -x.
111/111-
(4)
s=1
Then, test statistic is calculated for each j :.....p:
Buni (y X x- X., X - x)
T = Vn
(5)
y -- ||al>x, X -,, 112
Intuitively, Eq. (5) is the correlation of the regression residuals, calculated from Eq. 3 and 4 scaled by factor of yn. Under the null hypothesis, and more importantly, assuming linear relationship between X,, and y. this quantity follows standard normal distribution asymptotically, conditional to y and X y. It then follows that we can output a p-value for each variable j by taking P5 [ 4 (T))].
However, the formulation of test statistics in Eq. (5) not truly satisfactory in the setting of sparse logistic regression. More specifically, both the calculation of regression residuals y X and test statistics T3 do not take into account the non-linear relationship between X and the binary response y. The first row of Figure 6 plots the qq-plot of the test statistics Tj for logistic regression, which shows that even in the classical regime where >> p, its distribution is far from standard normal.
3 Decorrelating Test-Statistics for Hligh-Dimensional Logistic Regression
As we have elaborated, the formulation of dCRT is not well-suited for problems other than penalized lcast-squares regression. We therefore propose an adaptation of dCRT in the of logistic regression, inspired by the classical work of (13) and by (2 First, note that when testing H.J 599 0 under the case where n > we have the classical Rao's test statistic, defined by
ThRe vn 1/2
(6)
where v)((0, def. 13, 3 ),)
Here B_g def.
is the Fisher score.
argming B_y) is the constrained likeeiooo estimator of [_, with fixed
8j. and І is consistent estimator of the partial Fisher information L The appearance of the
term i 1/2 due to the fact that under the null hypothesis Hi. we have, by 13 Chapter 9], (4 11ーj
V VnV,e(0,,8) ( (d)
(O.
which makes the asymptotic distribution of 7 Rao standard normal. However, in the high-dimension case, where n < we do not reach this convergence in distribution. To see this, consider the Taylor expansion of the Fisher score of variable j around any given estimator B_ of the true 3Dz
(0, V ,00,32,,
-Bn,) + + ( 32, )2)
On the right-hand side, the first term converges weakly to normal distribution due to the Central Limit Theorem, the remainder term becomes negligible using é regularization to induce sparsity, but the second term does not, due to estimation bias and sparsity effect of B_,
Adattigg distillation operation for sparse logistic regression Fortunately, Eq. 7 suggests that for each variable j. we can debias the Fisher score by correcting the impact of other terms. In particular, for each variable j. replace the Fisher score by
(,((0,B)))
-
(8)

---

The inversion of the large matrix ر- E ( is computationally prohibitive, but we can estimate the term
straightforwardly by solving
W argmin WERP-1 2n
V2 - - --, + A//w|l1.
(9)
Moreover, since we have the closed-form of the derivatives of the logistic loss E(33), a simple derivation from Eq. ( suggests the following -r-diitillation operation, adapted for logistic regression:
for each variable j, where 3 given with Eq.
Bux_. argmin
E o (Xx XX.- + Ad
(10)
BERr TO
exp (X,.. ...
where the extra term (second-order derivative of the sigmoid function) g (X.,. +exp (X.
appears from differentiating twice the loss function E(3), and و 33 is defined in Eq. (2). On the other hand, we can obtain 3 from 3 by simply omitting the j-th coefficient from it, i.e.
3d def. (2,, ... 8. .p)
Finally, the equation for decorrelated test score, adapted from both Eq. (5) and 66 reads
Tykcou where the formula for the empirical partial Fisher information
1/2 Σ y 9(X, (X..,134u.) X X
(11)
X, (dx.,) X..- A summary of the full procedure, which we call CRT-logit, can found Algorithm Notice that the runtime of CRT-logit the same as dCRT, which means in general slower than KO and HRT. To speedup inference time, we introduce variable-screenigg step that eliminates potentially unimportant variables before distillation, similar to dCRT. We provide empirical benchmark of the runtime of each method in Section 4.5
Ens
Setting l-rreuularization parameter À and Adæ In general, we advise to use cross-validation for obtaining B3 in Eq. (2) and for x..-diisillltiion operator, as defined by Eq. (TO). This is in line also observe
with the theoretical argument for dCRT [9 Lemma and Theorem 3]. However, empirically that choosing the E f--regularization parameters of the distillation step can strongly affect how variables are selected by CRT-logit. We provide more details in the Supplementary Material, and leave further theoretical investigations of this phenomenon for future work.
Algorithm 1: CRT-logit
INPUT design matrix X Roxp, reponses y RR OUTPUT vector p-values (r, }/=
Ax + solve solv__pprre_______ SSCREENING + [j [p), 83 * 0}
// Using Eq.
for j SSCREENING do
gdx..; + solve_ scaled lasso_ cvass__c(X.. (x. X.,-,) ( // Using Eq.
10
Aduu + (31, . B 3p)
Tykevvrr + decorrelated test feast_scoreei X, X,y, Bdxx., ., (Br // Using Eq
p 2[1 (1, (e
10 end
11 for j & SSCREENING
12
p; =
13 end
Asymptotic analysis of the Decorrelated Test Statistic We now provide theoretical analysis of CRT-logit in large-sample regime. All the proofs can be found in the Supplementary Material. First, we introduce the following assumption.
5

---

Assumption 3.1 (Regularity conditions). Assume that
(AI) Amin (T) 2 K for some constant > O.
(A2) Sparsity of 3 and whi,, with w the ground truth weights for the distillation of x3 in Eq. TO |S| 8" and |w 1Io = s with 8* =onn>// (onn///\lo(p)) and s' = /
(A3) For all E [m]. X,,- and (- + (X.,,))))) are sub-exponential random variables, and , wod| S K' almost surely, for some constant K' > O.
We then have the following result, that states that the asymptotic distribution of the decorrelated test scores is standard normal.
Theorem 3.1. Let j E [p], and let Tyfcor be defined as in Eq. [1]. with x Ade Vn log(p)
Then, if Assumption 3.1 holds true, and we consider p p(ni), under the null hypothesis Hd. 139
O, we have
Vt R, lim [P (T S 20))=
where (.)) is the CDF of the standard Gaussian distribution. Moreover, for each j [p), if we def. 2 + (T₃)l i.e. Dj is the output of Algorithm then, under the null hypothesis
define Pj
Hly 39- = 0, we have
lim sup sup Pgo (Tg 11) st for all L e [0,1],
1-+00
that the p-values output by Algorithm 7 Mre valid asymptotically.
Remark 3.1. We also show in the proof of Theorem 3.1 (in Suppementary Material) that the rate of convergence is ((1 Vn). Compared with some of the related works (e.g., knockoffs) that come with finite sample guarantees, our theoretical analysis only works in the asymprotic regime. We leave the finite sample analysis as one of the directions for future works.
FDR control with CRT-logit As consequence of Theorem 3.1 have the following result, Yekutieli procedure 6 with the p-values output from Algorithm [ assuming that the number of tests p is fixed.
which establishes that the FDR of the test is controlled when using the Benjamini-Yekuttlli
Corollary 3.1. Under Assumptions 3.1 and logistic model defined in Eq ①. with X x Adz X
Vn log(p), assume n (/2 V s*) log(p) o(1), and assume the number of tests p is fixed. Let E 10 and SEy-CRT be given by applying following the Benjamini-Yekuse Yekutieli FDR-controlling procedure the CRT-logit p-values {e))e output from Algo Then, we have
card(Suv.ccr ( (S))
lim sup E 71-000
aa
card(Sav.cRr) V1
Remark 3.2. Assumption 3.1 is also assumed in (21 28), which also provide a detailed discussion of this regularity assumption in generalized linear models. This assumption, turn, is built on the regularity assumption in the classic work (13 Chapter 9] to establish asymptotic normality of Rao
test statistic. Theorem 3.7 is an adaptation of [21 Theorem 3.1], specialized for the case of sparse logistic regression and the p-values output from CRT-logit.
Empirical Results
We provide benchmarks of the proposed CRT-logii algorithm along with most other methods men- tioned in the introduction, in particular: model-X Knockoff (KO) [0), Debiased Lasso (dLasso)
33
[6] original CRT with 1000 samplings (10), Holdout Randomization Test with 5000 samplings 26 and Lasso-Distillation CRT (dCRT) [19]. We did not include SLOE and CPT as the provided
open-source implementation are particularly unstable and do not fit the regre--eerrssion setting (for SLOE), or implementation is not available (for CPT), For the lack of space, we leave the extra experiment with a genome-wiee wide association study in the Supplementary Material.
Remark 4.1. As a slight caveat, in the simulated and comi-realistic experiment sections (Sections 4.1 4.2 and we introduce an additional noise term to the logistic relationship of Eq. 1
P(y: = 1 xx) = g(x] 130 +o())
(12)

---

where &. ~ (0,0 1) is Gaussian noise and 0 >> O the noise magnitude. The formula in Eq. T2) has been used in previous works, e.g. [8]. There is clear justification to this: in most of the applications we consider, data are collected with measurement errors. In the case of brain-imaging, for example, recording the brain signal of the human subjects by scanners often includes noise caused either from the machine, or from the movement of the subjects f18y. Moreover, in general, this setting corresponds to a model mis spiisspecifiaaiin, which the CRT-logit is robust to under Assumption 3.1 following the same argument as in [21] Section 5].
Remark 4.2. We use Benjamini-Hochee all the empirical experiments in Section 4.2tionn422uue and App. 4.4 as we observe that the FDR is empirically controlled with this procedure, without compromising power with the conservative BY bound.
step-up procedure (5) to control FDR with the p-values in
4.1 Effectiveness of the decorrelation step
To show how decorrelating the test statistics helps, we set up a simulation with matrix X of = 400 features and vary the number of samples
{200, 400, 800, 4000]. The binary response vector y is created following Eq. (12), and the design matrix X is sampled from a multivariate normal distribution with mean, while the covariance matrix Σ E R is symmetric Toeplitz matrix, where the parameter p E (0. 1) controls the correlation between neighboring features: correlation decreases quickly when the distance between feature indices increases. The true signal 33 is picked with sparsity parameter K. = 8" /p that controls the proportion of non-zero elements with magnitude 2.0, i.e. B, 2.0 for all j E S. For the specific purpose of this experiment, non-zero indices of S kept fixed. The noise & is i.i.d. normal (0, Id,) with magnitude controlled by the SNR parameter. In short, the three main parameters controlling this simulation correlation p. sparsity degree and signal-to-noise ratio SNR. We generate randomly 1000 datasets, and run dCRT and CRT-logit algorithm to obtain one sample of test statistics (T, \ and
|| (/25 SNR),
1P (Trdecomm j=1- Then, pick 1000 samples of one null test statistic T3 and TylecoT defined in Eq. (5) and respectively, and plot the qq-plot of their empirical quantile versus the standard normal quantile. From the results in Figure 6
observe that the empirical null distribution of the test
dCRT
dCRT
dCRT
dCRT
Sample Quaztiles
Sample SamspeeQune Quantiles
Sample Quantiiee
Sample Quantiies
Theoretical Qoantiles
Theoretical Quantiles
Theoretioal Quaniiles
Theoretical Quactiles CRT-Logit
CRT -ogit
CRT -ooit
CRT -Ro-t
Sample Ouanties
Sarppte Quantile
Sample Quanttles
Sample Quannties
-25 0.0 25
-25 0.5 25 Theoretical Quantiles
-25 0.0 25
-25 00 26 Theoretical Quactties
Thecretical Quantiles
Theoretical Quantiles
(a) 200
(b)
400
(c) n 800
(d)
4000
Figure 2: QQ-Plot for one null CRT statistic for logistic regression, with varying number of samples and fixed number of variables p = 400. The theoretical quantiles are obtained from standard Gaussian distribution. The decorrelation step makes the empirical null distribution of the null statisties much closer to standard Gaussian. Parameters: SNR = 3.0, p = 0.4, sparsity = 0.06. Upper row: Distilled-CRe CRT statistic defined by Eq. 5) Bottom CKT-logit, with decorreleated test score defined by Eq. [I (ours).
statistic is much closer to a standard normal when adding the decorrelation step. In particular, when the sample size increases to 400, the decorrelated test statistic has empirical quantiles almost inline with the theoretical quantiles of the standard normal distribution, while dCRT test score strays away from the 45-degree line. Again, emphasize that the normality of T, is crucial for the p-values calculation. This outlines the importance of the decorrelating step on T,.
4.2 High-dimensional scenario with varying simulation parameters
To see how each algorithm performs under different settings, we follow the same simulation scenario as in Sec. 4.1 but vary each of the three simulation parameters, while keeping the others unchanged
7

---

at default value of SNR = 2.0, p 0.5, = 0.04. We target a control of FDR at level 0.1, using the Senjamini-Hochbbee procedure. The results in Figure that CRT-logit the most powerful method while still controlling the FDR. Moreover, in the presence of higher correlations between nearby variables (p > 0.6), other methods suffer drop in average power, but this not as severe for CRT-logii. The original CRT, in general, is conservative. We believe that this is due to using only B 500 samplings to generate empirical p-values for the two methods, due to prohibitive average runtime of the algorithm with larger B (which we provide in Section 4.5). For HRT, the conservativeness is expected, due to the usage of only half of the sample for stst-statistics calculation even though the number of samplings is bigger than original CRT (B 5000). We note that, perhaps surprisingly, the debiased lasso (cdlasso) is the most conservative. It controls FDR well in all settings. This might be due to the fact that dlasso also relies on the choice of the f--reuularization in the nodewise Lasso operation, similar to the --y-distillationn of dCRT, as noted in Section []] What makes the difference is that instead of using cross-validation for setting À for each variable J. a fixed value of x = Xmmz is used in the implementation of dlasso. We strongly suspect this fixed value is not optimal, which makes the procedure powerless.
4.3 Application: large-scale analysis on brain-imaging dataset
Description The Human Connectome Project dataset (HCP) is a collection of brain imaging data healthy young adult subjects with age ranging from 22 to 35. More specifically. the in- put X is a set of 2mm statistical maps of 400 sub- jects across 8 cognitive tasks. These are called z-maps, the data follow standard normal distribution under the null hypothesis. Each task in turn features 2 different contrasts, which ef- fectively form binary responses y E (0. 1)". In short, the goal of this fMRT data analysis is to identify voxels with task-related levels of activ- ity by fitting y through distributed brain signals. The setting is nigh-diimessionll with n = 800 samples, corresponding 400 subjects, while the total number of variables is p 200, 000 brain voxels. Following [I] 20), we use a hier-
dlasso
CRT Jogit
Avg. power
1.0
0.5
loght
FDR
0.4
0.3
SCRT
EIRT
0.6
0.8
0.2
0.1
0.2
0.0
1.0
2.0
30
5.0 6000
1.0
2.0
30
A
5.0 600
snr
snr
Avg. power
1.0
0.5
FDR
0.4
0.8
0.3
0.2
0.4
0.1
0.2
0.0
0.5
G.7
0.8
0
0.5
0.6
0.7
0.8
rho
rho
archical clustering scheme to group the variables into = 1000 spatially connected clusters. We provide details of the pre-processing step in Sup-
0.5
Avg. power
1.0
0.4
FDR
0.3
0.6
0.2
0.4
plementary Material.
0.0
0.1
0.0
0.2
0.99 0..00 8...8888.
8..0
0.0g
0.0.00 10000.o 8.888% 08
Creating semi-realistic ground-truth and re- sponse labels Since there is no ground truth for this dataset, we create synthetic true signals by fitting the data X and response y with an f1- penalized logistic classifier. In other words, the estimator giogreg will serve as true regression coefficients for each task. Then, to avoid bias in simulating label y. the z-maps matrix X of one task are used in conjunction with the discrimi- native pattern map flogros from the next task in the following order: relational, gambling. emotion, social. For instance, we use 3 flogreg of gambling ing with z-maps data matrix of relational, i.e. for all i 1,..., n, given relttiiainn
sparsity
sparsity
Figure 3: FDR/Average Power of 100 runs of simula- tions across varying parameters in high-dimensional settings. Default parameter:
= 400, P
600, SNR 2.0, 0.5, = 0.04. FDR is controlled level a 0.1. Methods: Debiased Lasso (dlaaso), model-X Knockoff (K0-1ogit)) original CRT (CRT), HRT (HRT), dCRT (dCRT), and our version of CRT (dark green line- CRT- -Rg-iogit).
(13) where Bern(a) is a Bernoulli probability mass function that takes a value 1 with probability a, a is noise magnitude and &, is a standard normal noise. Finally, we apply all inference algorithms on the semi-synthetic data (X,)). and we evaluate their performance using the ground-truth Glogre, This simulation setting is similar to [I]] except that here we consider classification and not a regression
Bern [g(
[
T
Agambing + oε.)
8

---

problem. It allows us to calculate the False Discovery Rate and average power with multiple runs of the inference procedure (across tasks).
Remark 4.3. The i...d. assumption is formally violated in this experiment, where for each subject we analyze two sample images that are not independent. Yet, this remains short-range correlation structure, and is thus not a strong challenge to the ii.d. assumption.
The results in Figure show that CRT-logit achieves better recovery compared to KO original CRT/CCRT/HRT, which results in higher statistical power. This gain comes with good control of the FDR under desired level a = 0.1. On a related note, the only analysis where dCRT makes more discoveries than CRT--Logit is in the emotion task, but at the cost of failing control FDR at the nominal level.
Results
dlasso
KO-Logit
FDR
dCRT
CRT- Logit CRT
HRT
FDR control
Avg. Power
dlasso
KO-Logit dCRT
CRT-Logit CRT
HRT
0.0
0.2 RELATIONAL
0.4
0.0
0.2
0.4
0.0
0.2
0.4
0,0
0.2
0.4
GAMBLING
EMOTION
SOCIAL
Figure FDR/Average Power of 50 runs semi-realistic experiments on four tasks of Human Connectome Project dataset. Parameters: и (taken from 400 subjects), SNR 2.0. Methods (clustering versions): Debiased Lasso (cdlasso), model-X Knockoff (cKO-logit), original CRT (CRT), HRT (HRT), dCRT (dCRT), and our version of CRT (dark green line CRT-logit).
4.4 Application: genome-wide association study with Human Brain Cancer Dataset
Description The last in our benchmark is a Genome-wide Association Study (GWAS) on the The Cancer Genome Atlas (TCGA) dataset (30 31]. We choose analyze the Glioma cohort, which consists of n = 1026 patients across wide age range, diagnosed with this type of brain tumor, with total of p 24776 genes in the data matrix recorded copy number variations (CNVs) at the gene level in log ratio format. As with the brain-imaging inference in Section 4.3 we use clustering to reduce the dimension
1000 clusters. However, we use different criterion to merge variables (genes) to clusters of variables, which is the pairwise Linkage Disequilibrium, following []] Section 4] (with available R library). For the response, long-term survivor (LTS) is defined as a patient who survived more than five years after diagnosis and would be labeled O, and otherwise would be identify significant genes that contribute
short-term survivor (STS), labeled y The objective to classification of the LTS/STS status. Similar to the Human Connectome Project dataset, there is no real ground-truth for the TCGA Glioma. However, we have the list of mutations and the frequency of those detected in the diagnosed patients. We therefore select the 1000 most frequent gene mutations that appeared in this list, iee. the ground truth list consists of 1000 genes (variables).
Table 1: List of detected genes associated with Glioma Cancer from the TCGA dataset. n 1026, p = 24776 (clustered to C' 1000). Empty line (-) signifies no detection. Methods listed in the table are the clustering version. Commonly detected genes between methods are put in bold text. Most detected genes are listed in the mutant list database that can he found in the recorded patients 30
Methods
Detected Genes
dLasso KO
ABCCIO, ANK3, CDH23, PTEN, SPEN, SVIL, ZMIZI
dCRT
ANK3, ANKRD30A. CDH23, PTEN, RET. SPEN ZZIZI ABCC10, ANKRD30A, BCOR, EPHA3, PPL, SPAG17, SPEN, SVIL, USP9X
CRT-logit
Original CRT
ABCC10, BCOR, EPHA3, SPEN, SVIL.
HRT
ABCC10, SPEN

---

Result The result from Table shows that CRT-logit finds the largest number of genes. Moreover, most of selected genes this table are detected the list of mutated genes found on recorded patients. Some genes are detected by all the benchmarked methods, most prominently SPEN, which is found on over 10 % of patients in the cohort. Furthermore, this gene is known to be associated not only with brain cancer, but also with other types of cancer in The Human Protein Atlas project [17]. Note that, in the absence of a ground-truth, this does not guarantee all genes found are associated with glioma, but this experiment demonstrates the capability of CRT-logit in GWAS studies.
4.5 Average runtime of benchmarked methods
Table 2: Average runtime of benchmarked methods for one simulation (in seconds). Standard error reported in parentheses.
Methods
Simulated (Sec. 4.2)
HCP-semi-real (Sec.
Debiased Lasso
33
16
61.83 (5.2)
154.27 (8.79)
Knockoff Filter
4
10
1.62 (0.02)
8.12 (0.62)
CRT (500 samplings)
110
2312.91 (38 21)
7069.96 (109.09)
HRT (5000 samplings)
26
14.84 (2.01)
52.17 (4.11)
dCK[[screening=ITru]
19
16.83 (1.89)
65.18 (3.91)
dCK[[[screniing=False]
370. 12 (8.18)
962. 10 (20.63)
-KT--ogitlscrveennn=Truee (this work)
14.16 (0.35)
61.26 (3.55)
CKK--logidlscreennnn--aaass (this work)
367.91 (4.11)
983.78 (17.26)
Besides statistical performance, it is equally important to assess the computational cost of inference procedures. The average runtime in Table from the two experiments shows that the original CRT is not suitable for large-scale inference: it is over 2000 times slower than the fastest method (Knockoff Filter), and over 150 times slower than dCRT//CRTlogit. The empirical runtime also confirms the effectiveness of the screening step before doing distilattio/ddecorrlatiio of the test-statistics: the step makes CRT-logit and dCRT 20 times faster than without. On related note, although in theory Debiased Lasso, dCRT and CRT-logii (both without screening) share the same runtime complexity, the latter two are slower due to the use of cross-validttion estimate the sparsity hyperparameter À and Ade (detailed in Section 3).
Discussion
We proposed an adaptation of the Conditional Randomization Test (CRT) for sparse logistic regression in the high-dimensional regime. A major improvement of CRT-logit, our proposed algorithm, compared to original CRT, comes from the decorrelation of test statisties to make their distribution closer to standard normal. Indeed, results from synthetic experiments in Figure show that in high- dimension (when 0.5 < n/p S 1.0), the empirical null distribution of CRT-logit's ' test statistic 7 deccn is much more similar to standard normal compared to the original CRT test statistic. Moreover, empirical benchmarks in Section Irationnddmman demonttrate that CRT-logit performs better than related statistical inference methods, such as the Debiased Lasso or Model-X Knockoffs. In particular, CRT-logit is the most powerful method in our synthetic experiment with sigh--iimnnsionll datasets in Section 4.2 while still keeping FDR controlled under predefined level a 0.1 (with slight caveat of using BH instead of BY procedure, as elaborated in Remark 4.2) We note that there exist some limitations to CRT-logit. The computational cost of CRT-logit, while lower than vanilla CRT, is still larger than alternative methods such as Knockoff Filter and Holdout Randomization Test. Moreover, tuning the eerrguuaainattoon Adz parameter by cross-validation, the computational cost of CRT-logit (and dCRT). Despite this, our empirical benchmarks on both simulated and real data show real promises of CRT-logit. Henceforth, believe CRT-logit is competitive for practical settings that involve structured data, such as grain-imaging and genomics applications.
is often done, can further increase
Acknowledgements
BN, BT and SA acknowledged the support of the French "Agence Nationale de la Recherche" under the project ANR-17-CE:23-00 (FastBig) and TER-00-CHIA-00255 (KARAIB AI chair). SA was also supported by Institut Universitaire de France (IUF), and BN was also supported by Research Chair ASAIDIS (Data Science and Artificial Intelligence for Digitalized Indutry and Services) of Télecom Paris.
10

---

References
[1] Christophe Ambroise, Alia Dehman, Pierre Neuvial, Guillem Rigaill, and Nathalic Vialancix. Adjacency-conntrranee hierarchical clustering of band similarity matrix with application to genomics. Algorithms for Molecular Biology, 14(1):1-14, 2019.
[2] Francis Bach. Self-concordant analysis for logistic regression. Electronic Journal of Statistics, 4(none), January 2010.
[3] Francis Bach, Rodolphe Jenatton, Julien Mairal, Guillaume Obozinski, et Optimization with sparsity-inducing penalties. Foundations and Trends@ in Machine Learning, 4(1):1-106, 2012.
[4] Rina Foygel Barber and Emmanuel J. Candès. Controlling the false discovery rate via knockoffs. The Annals of Statistics, 43(5):2055-2085, October 2015. arXiv: 1404.5609.
[5] Yoav Benjamini and Yosef Hochberg. Controlling the False Discovery Rate: Practical and Powerful Approach to Multiple Testing. Journal of the Royal Statistical Society. Series B (Methodological), 57(1):289- 300, 1995.
[6] Yoav Benjamini and Daniel Yekutieli. The control of the false discovery rate in multiple testing under dependency. Ann. Statist., 1655):16651118,, 08 2001.
[7] Thomas B. Berrett, Yi Wang. Rina Foygel Barber, and Richard Samworth. The condi- tional permutation test for independence while controlling for confounders. Journal of the Royal Statistical Sociery: Series (Statistical Methodology). 82(1):175-197, 2020. _eprint:
https///ss..ollieeiibrrryy wttp://awww..iilleooooo/o//////a
1/rssb.. 12340.
[8] Danilo Bzdok, Michael Eickenberg. Olivier Grisel, Bertrand Thirion, and Gael Varoquaux. Suni-iuper Factored Logistic Regression for High--Dimensionl Neuroimaging Data. In C. Cortes, N. Lawrence, D. Lee, M. Sugiyama, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 28, Curran Associates, Inc., 2015.
[9] Peter Bühlmann, Philipp Rütimann, Sara van de Geer, and Cun-Hui Zhang. Correlated variables in regression: clustering and sparse estimation. Journal of Statistical Planning and Inference, 143(11):1835-185) November 2013. arXiv: 1209.5908.
[10] Emmanuel Candès, Yingying Fan, Lucas Janson, and Jinchi Lv. Panning for gold: 'model-x' knockoffs for high dimensional controlled variable selection. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 80(3):551-577, 2018.
[11] Jérôme-Alexis Chevalier, Tuan-Binh Nguyen, Joseph Salmon, Gaël Varoquaux, and Bertrand Thirion. Decoding with confidence: Statistical control on decoder maps. Neurolmage, 234: 17921, 2021.
[12] dérôme-Alexis Chevalier, Tuan-Binh Nguyen, Bertrand Thirion, and Joseph Salmon. Spatially relaxed inference on high-dimensional linear models. arXiv:2106.02590 [math, stat], June 2021. arXiv: 2106.02590.
[13] David Roxbee Cox and David Victor Hinkley. Theoretical statistics. CRC Press, 1979.
[14] Wenjiang Fu and Keith Knight. Asymptotics for lasso-type estimators. The Annals of statistics, 1855)1335-13378 2000.
[15] Trevor Hastie, Robert Tibshirani, and Jerome Friedman. The elements of statistical learning: data mining, inference and prediction. 2009. OCLC: 1058138445.
[16] Adel Javanmard and Andrea Montanari. Confidence intervals and hypothesis testing for high- dimensional regression. The Journal of Machine Learning Research, (5(1):2869-290)) 2014.
[17] Stéphanie Légaré, Luca Cavallone, Aline Mamo, Catherine Chabot, Isabelle Sirois, Anthony Magliocco, Alexander Klimowicz, Patricia N Tonin, Marguerite Buchanan, Dana Keilty, et al. The estrogen receptor cofactor spen functions as tumor suppressor and candidate biomarker of drug responsiveness in bormone-dependen breast cancers. Cancer research, 75(20):4351-4333 -3663, 2015.
11

---

[18] Martin A. Lindquist. The Statistical Analysis of fMRI Data. Statistical Science, 23(4), November 2008.
[19] Molei Liu, Eugene Katsevich, Lucas Janson, and Aaditya Ramdas. Fast and powerful conditional randomization testing via distillation. Biometrika, 07 2021.
[20] Tuan-Binh Nguyen, Jérôme-Alexe Chevalier, and Bertrand Thirion. Ecko: Ensemble of clustered knockoffs for robust multivariate inference on fmni data. In International Conference on Information Processing in Medical Imaging, pages 4544666, Springer, 2019.
[21] Yang Ning and Han Liu, A general theory of hypothesis tests and confidence regions for sparse high dimensional models. Annals of Statistics, 45(1):158-195, February 2017. Publisher: Institute of Mathematical Statistics.
[22] Dmitrüi M Ostrovski and Francis Bach. Finite-sample analysis of m-estimators using self- concordance. Electronic Journal of Statistics, 1511:3326-391, 2021.
[23] Fabian Pedregosa, Gaël Varoquaux, Alexandre Gramfort, Vincent Michel, Bertrand Thirion, Olivier Grisel, Mathieu Blondel, Peter Prettenhofer, Ron Weiss, Vincent Dubourg, et al. Scikit- learn: Machine learning in python. Journal of machine learning research, 12(0ct):2825-2830 2011.
[24] Radhakrishna Rao. Large sample tests of statistical hypotheses concerning several parameters with applications to problems of estimation. In Mathematical Proceedings of the Cambridge Philosophical Society, volume 44, pages 50-57. Cambridge University Press, 1948.
[25] Pragya Sur and Emmanuel Candès. A modern maximum-likelihood theory for high- dimensional logistic regression. Proceedings of the National Academy of Sciences, 116629):1511-11455,, July 2019. Publisher: Proceedings the National Academy of Sciences.
[26] Wesley Tansey, Victor Veitch, Haoran Zhang, Raul Rabadan, and David M. Blei. The Hold- out Randomization Test for Feature Selection in Black Box Models. Journal of Computa- tional and Graphical Statistics, 31(1):151-162, 2022. Publisher: Taylor & Francis _eprint: https://doi.org/II 10. 101100000000 10011006160000 01000/00061600022211133333 102.0.9222..
[27] Bertrand Thirion, Gaël Varoquaux, Elvis Dohmatob, and Jean- Caatist Poline. Which fMRI clustering gives good brain parcellations? Frontiers in Neuroscience, 8, 2014.
[28] Sara van de Geer, Peter Bühlmann, Ya' acov Ritov, and Ruben Dezeure. On asymptotically optimal confidence regions and tests for high-dimensional models. The Annals of Statistics, (613):16661120 June 2014. arXiv: 1303.0518.
[29] Gaël Varoquaux, Alexandre Gramfort, and Bertrand Thirion. Small-sample brain mapping: sparse recovery spatially correlated designs with randomization and clustering. Proceedings of the 29th International Conference on Machine Learning, Edinburgh, Scotland, UK, 2012, page 8, 2012.
[30] Suhas V Vasaikar, Peter Straub, Jing Wang, and Bing Zhang, LinkedOmics: analyzing multi- omics data within and across 32 cancer types. Nucleic Acids Research, 46(D1):D956-D963, January 2018
[31] John N Weinstein, Eric A Collisson, Gordon B Mills, Kenna R Mills Shaw, Brad A Ozenberger, Kyle Ellrott, Ilya Shmulevich, Chris Sander, and Joshua M Stuart. The cancer genome atlas pan-cancer analysis project. Nature genetics, 45(10):1113-1120 2013.
[32] Steve Yadlowsky, Taedong Yun, Cory Y McLean, and Alexander D' Amour. SLOE: A Faster Method for Statistical Inference in ffih--Dimnsiones Logistic Regression. In M. Ranzato, A. Beygelzimer, Y. Dauphin, S. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, volume 34, pages 29517-29528. Curran Associates, Inc., 2021.
[33] Cun-Hui Zhang and Stephanie S. Zhang. Confidence intervals for low dimen- sional parameters in high dimensional linear models. Journal of the Royal Sta- ristical Society: Series B (Statistical Methodology). 76(1):217-242, 2014. _eprint: wwtps//wwiiiiiiaaaiiicceooo/oo/a http://
12

---

[34] Qian Zhao, Pragya Sur, and Emmanuel J. Candès. The asymptotic distribution of the MLE in high-dimensional logistic models: Arbitrary covariance. Bernoulli, 28(3):1835 1861, 2022. Publisher: Bernoulli Society for Mathematical Statistics and Probability.
Checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's
contributions and scope? [Yes]
(b) Did you describe the limitations of your work? [Yes]
(c) Did you discuss any potential negative societal impacts of your work? [No] We believe our work provides no potential negative societal impacts.
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yss]
2. If you are including theoretical results...
(a) Did you state the full set of assumptions of all theoretical results? [Yes]
(b) Did you include complete proofs of all theoretical results? [Yes] The proofs can be found in Supplementary Material
3. If you ran experiments...
(a) Did you include the code, data, and instructions needed to reproduce the main exper- imental results (either the supplemental material or as a URL)? [Yes] Except for datasets of realistic applications, which are not included, but references are provided
(b) Did you specify all the training details (c.g.. data splits, hyperparameters, how they were chosen)? [Yes] Both in main text and in Supplementary Material
(c) Did you report error bars (c.g.. with respect the random seed after running experi- ments multiple times)? [Yes] With the simulated and semi-rcalistic experiment
(d) Did you include the total amount of compute and the type of resources used (e.g., type of GPUs, internal cluster, cloud provider)? [No]
4. If you using existing assets (e.g., code, data, models) or curating/releasing new assets...
(a) If your work uses existing assets, did you cite the creators? [Yes]
(b) Did you mention the license of the assets? [N/A]
(c) Did you include any assets either in the supplemental material URL? [No]
(d) Did you discuss whether and how consent was obtained from people whose data you't 're using/curating? [N/A]
(e) Did you discuss whether the data you are using/curating contains personally identifiable information or offensive content? [N/A]
5. If you used crowdsourcing or conducted research with human subjects...
(a) Did you include the full text of instructions given participants and screenshots, if applicable? [N/A]
(b) Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, applicable? [N/A]
(c) Did you include the estimated hourly wage paid to participants and the total amount spent participant compensation? [N/A]
13


Reviews:
[{"review_id": "aoGxfP2HC_b", "reviewer": "Reviewer_np9G", "summary": "This paper presents a method for performing hypothesis testing for lasso-penalized logistic regression in the high dimensional (but sparse) setting.  The main result relies on the asymptotic normality of a test statistic that is essentially how correlated a given feature is ", "questions": "Comments / questions:\n\n* Figure 2 is interesting, but the asymptotics presented here -- and the really interesting case -- is considering what happens as both n and p get large.  A very small simulation study to show that the asymptotic results hold even for p >> n would, in my opinion, strengthen the work.  Furthermore, in both the brain processing application and the TCGA application the data is preprocessed so that p ~ n.  Is this because the p^3 scaling makes it infeasible to run on the full dataset, or are there additional issues with p >> n?\n\n* There is no discussion of finite sample guarantees.  e.g., Knockoffs control the FDR for finite sample sizes, whereas the present method only has such guarantees in the asymptotic regime.  This is not a big deal, but it would be good to mention this somewhere (maybe in the related works section).\n\n\n\nMinor Comments / questions :\n* The switch from $\\hat{w}^j$ to $\\hat{\\beta}^{d_{x_{*,j}}$ from equations (9) to (10) is a bit notationally confusing.  Adding a sentence beforehand would help with the transition.\n\n* In lines 171-174 the mention of a variable-screening step makes it sound like a somewhat heuristic procedure, when, in fact, it just leverages that it is trivial to compute p-values for j's where $\\hat{beta}_j$ is zero.\n\n* It would be nice to be notationally consistent between Assumption 3.1 and Assumption A.1. (e.g., $\\kappa^2$ becomes $K$ and $K$ becomes $K'$)\n\t\n* In Theorem 3.2, is p overloaded as both the number of features in the design matrix and also the number of tests? If separate they should not both be labeled $p$.  If they are the same, then isn't the theorem trivial by the assumption of $p$ being fixed (eventually $n \\gg p$ and it is essentially just standard logistic regression)\n\n* Remark 4.2 feels important and should be emphasized a bit more in my opinion -- the procedure used in the empirical tests is not proven to be valid by the analysis, but it seems to work in practice.\n\n* In Table 1, the boldface should be removed.  The standard usage of boldface in these types of tables is to highlight the fastest method.  The parenthetical \"(this work)\" is already sufficient to draw the reader's attention to CRT-logit\n\n* The notation in appendix A was a bit confusing to me.  In particular, in several places (e.g., the definition of $v^*$) a vector is defined which only makes sense later under the assumption that $j=1$.  For example $v^*$ is defined as ($1, -w^{0, j})$, but then ${v^*}^\\top \\nabla \\ell(\\beta)$ only makes sense if the 1 in $v^*$ lines up with $j$.  This is obviously not a big problem and doing something more \"correct\" might even be notationally clunky enough to not be worth it, but at the very least a quick explanation of the abuse of notation would be helpful.\n\n* $\\hat{v}$ has not yet been defined by Lemma A.4.\n\n* $v^*$ as defined on line 404 and defined (again) on line 414 have different signs for the $w^{0, j}$ part.\n\n* I may be missing something, but is some term missing from the first line in the display after line 415 (same for the first line in the display after line 417). Both of these rely on lemma A.4 and so I believe they should also pick up a term that involves $s^* \\vee s' \\log p / n$.  This obviously does not affect the result (it gets absorbed into the constant hidden by the $\\preceq$ notation at the bottom of each display) but it did make it difficult to follow the argument.  Along these lines, is Lemma A.4 missing an absolute value sign on the right hand side of both equations in the displace below line 408?\n\n* The results in Figure 5 suggest that it may be beneficial to continue increasing \\lambda beyond the range of values considered in the figure.  In most rows, the power continues to increase as \\lambda/\\lambda_{univ} increases up to the final value displayed; additionally, the FDR remains within estimation error of being well-controlled.  That is, with 100 simulations and a true FDR of 0.1, one might expect a 95% confidence interval to be estimated_FDR +/- 1.96 * sqrt(0.1 * 0.9 / 100), which is ~0.06 so any FDR below 0.106 should be tolerable.  The authors should consider continuing to increase lambda to see if further power can be gained or if there is a point at which the regularization is too strong.\n\n* In appendix G it's stated that the same clustering as 4.3 is used to reduce the dimension.  In F.1 its stated that the clustering preserves the spatial structure of the data.  Is that also the case for the genes in appendix G?  There is some clustering of genes by function in the genome, but overall not really so I'm not sure that it makes sense to include any spatial component here.  Furthermore, if the genes are clustered then why does Table 3 say n=1026 and p=24776, doesn't this contradict that the genes were clustered to 1000 clusters (p=1000)?  Finally, if the genes were clustered, then how are individual genes pulled out in Table 3?  Do all genes in a significant cluster get put in the table?\n\n\t\n\n\t\n\nTypos:\n* Line 13 \"the geometry of ell_1\" --> \"the geometry of the ell_1\"\n* Line 29 \"conditionally to\" sounds off to me. Perhaps \"conditioned on\"?\n* Lines 37-39 are difficult to parse grammatically\n* Lines 54-55 \"convergence to Gaussian\"\n* Lines 64-65: \"standard normal test-statistic\" --> either \"a standard normal test-statistic\" or \"standard normal test-statistics\"\n* Line 65 \"in large-sample regime\" --> \"in the large-sample regime\"\n* Figure 1 caption: \"The empirical distribution of dCRT null-statistic\" --> \"The empirical distribution of the dCRT null-statistic\"\n* Line 71: \"multiple sampling\" --> \"multiple samplings\"\n* Line 74: \"hence inherently\" --> \"and hence inherently\"\n* Line 169: \"where the formula for empirical\" --> \"where the formula for the empirical\"\n* Line 173-174: \"We provide empirical benchmark of runtime\" --> \"We provide an empirical benchmark of the runtime\"\n* Line 176: \"inline\" --> \"in line\"\n* Line 198: \"using Benjamini-Yekutieli\" --> \"using the Benjamini-Yekutieli\"\n* Line 215: \"implementation are\" --> \"implementation is\"\n* Line 253: \"at default value\" --> \"at a default value\"\n* Lines 253-254: \"using Benjamini-Hochberg procedure\" --> \"using the Benhamini-Hochberg procedure\"\n* Line 254: \"Results in Figure 3\" --> \"The results in Figure 3\"\n* Line 310: \"analyze two sample image\" --> \"analyze two sample images\"\n* Line 315: \"is in emotion task\" --> \"is in the emotion task\"\n* Line 315: \"failing to control FDR\" --> \"failing to control the FDR\"\n* Line 316: \"at nominal level\" --> \"at the nominal level\"\n* Line 338: \"We note that there exists\" --> \"We note that there exist\"\n* Line 404: \"under logistic model\" --> \"under the logistic model\"\n* Line 405: \"under logistic model\" --> \"under the logistic model\"\n* Line 413: \"written in more general\" --> \"written in a more general\"\n* Line 415: \"where we use triangle inequality\" --> \"where we use the triangle inequality\"\n* Line 416: \"last inequality us due\" --> \"last inequality is due\"\n* Line 456: \"from inference algorithm\" --> \"from an inference algorithm\"\n* Display below line 458: \\hat{k}_{BY} should be \\hat{k}_{BH}", "limitations": "I believe the authors have adequately addressed the limitations of their work, and I am not aware of any potential negative societal impact.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 4, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The proposed approach is interesting and theoretically well-motivated.  That it gives (asymptotically) valid p-values is a compelling advantage over the knockoff framework that only allows for FDR control.  The asymptotics are also nice because then one can avoid needing to do any resampling.  One aspect that I found missing from the paper is a discussion of (the lack of) finite-sample guarantees.  The lack of finite sample guarantees is not a concern to me, but given that that's one of the compelling aspects of some of the related works (e.g., knockoffs) it would be good to at least discuss.\n\nOne minor weakness is that the clarity of the manuscript could be improved in some places.  It is quite dense (which is understandable due to space constraints).  While some parts of the paper have really nice intuitive explanations, other parts do not (e.g., equations 4 and 5).  Furthermore, the notation seems to be shifting throughout (see some of my comments in the \"Questions\" section) which can make it somewhat difficult to follow.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "TV-tavsEhpU", "reviewer": "Reviewer_VbNR", "summary": "The authors study the problem of doing inference on a logistic regression model with an L1 penalty in high dimensions, where the number of features in the dataset is at least as large as the size of the dataset. For this problem, they develop a variant (that they call CRT-logit) of the distilled conditional randomization test (dCRT) with higher power than the latter. Their innovation is in introducing a decorrelation step that brings the null distribution of the test statistic closer its assumed distribution - a standard normal. An asymptotic analysis of the performance of CRT-logit is given.", "questions": "CRT-logit is proposed for problems in which $p>n$, but the authors give an asymptotic ($n\\rightarrow\\infty$) analysis of its performance in Theorem 3.1. At a high level, I am confused by this meta structure in the paper. Do the authors have a sense for how quickly (in $n$) the p-values supplied by the algorithm become ‘good’ in some sense?\n\nWhat is the utility of adding the extra assumption of sub-exponential tail behavior only in the appendix (I am comparing Assumption A.1, (3) to Assumption 3.1, (3))?\n\nThere is a typo in line 93. $x$ should be $x_n$.\n\nThe ‘description paragraph’ in section 4.3 could be made clearer. In particular, I suggest the authors simply state that they are working with fMRI data and that the goal of the analysis is to identify voxels with task-related levels of activity. As written, the type of data and the goal of the analysis are unclear.", "limitations": "Yes.", "rating": 7, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "I think this is a very nice paper. It focuses on an important and ubiquitous inference problem, and is generally of a high quality, largely due to its clarity and thoroughness. The central technical innovation, namely, the decorrelation procedure discussed in Eqs. 8-11, is well-motivated both analytically and in Figure 1. The new algorithm that is introduced and tested, CRT-logit, will be of interest to the broader machine learning community. Its relationship to existing algorithms is discussed; a number of experiments comparing it with existing algorithms are also presented.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "TDLn0erASxA", "reviewer": "Reviewer_XeSh", "summary": "The authors study the case of high-dimensional logistic regression when the number of features p is much greater than the number of sample n.  They propose the CRT-logit algorithm that combines a variable-distillation step and decorrelation step to keep the sparsity in l1-penalized logistic regression.  The authors provide theoretical analysis of their approach and show it is effectiveness on simulations and experiments on real-world brain-imaging and genomics data. The main contribution is the proposal of the CRT-logit method which is shown to be effective and inference cost not prohibitively high.", "questions": "1. Why do we need the noise term in Eqn. 12?\n2. Can this method be generalized to work on a family class of models in high dimensions, e.g. any generalized linear model (GLM)?\n3. What happens in the other case when n >> p?", "limitations": "Yes the authors have addressed the limitations of their work.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths\n- The authors propose CRTlogit and provide thorough theoretical and experimental results\n- Highly relevant problem since logistic regression is still one of the most commonly used methods as well as l1-penalty in machine and deep learning\n- Valid and thorough theoretical results\n- Show empirical validation of theoretical results\n- Also show experiments on average inference runtime in addition to performance results\n- Tests on real-world brain and genomics data.\n- Very nice qualitative plots e.g. Figure 1 showing the improvement on the theoretical quantiles of the proposed CRT-Logit over dCRT.\n\nWeaknesses\n- Some of the derivations in Eqn, 8-11 can be moved to the appendix.\n- Limited to only looking at logistic regression - can this extend to other models?\n- The authors should emphasize their novelty and contribution as more than just an extension of CRT.\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "8O1MMZ4dv7c", "reviewer": "Reviewer_7y9r", "summary": "This paper aims to extend the growing literature on identifying relevant features in a machine learning model. The authors focus on the setting of classification in a high-dimensional setting where the number of relevant features is sparse (i.e., less than n^{1/2}), where n is the number of measurements and p is the number of features. To do so, they extend the algorithm of distilled conditional randomization test (dCRT), which itself is an extension of the conditional randomization test (CRT) to make CRT more computationally feasible. \n\nA key step of the dCRT algorithm is to take residuals from two regressions (one regression to see if feature j lies in the linear span of the other features, and the second regression to see if the response variable lies in the linear space of the other features) and use that to compute the test statistic. In essence the key innovation of this paper is to extend the way the residuals are taken to fit a logistic model rather than a linear model. This estimator is called CRT-logic.\n\nGiven this new test statistic and an assumption of sparsity along with other regularity conditions, they prove pointwise gaussian approximation and asymptotic validity of the CRT-logit estimator. They then corroborate their results with simulations and real-world data.\n", "questions": "1.\tWhat does the notation  y_i(X^T_{I, -j} \\beta) mean? Is that an argument to y_i() or multipled by y_i? For example see (4)\n2.\tFor the comparison with dCRT in Figure 1, how does dCRT perform when not using the original features, but a dictionary mapping of the original features (e.g. low-degree polynomial)? This might be a fairer comparison as it is quite obvious that just using the residuals from doing linear regression on the original features will lead to a poor for fit for a logistic function.\n3.\tIs the result in Theorem 3.1 assuming that the solve_scaled_lasso_cv step works perfectly? The result in Theorem 3.2 is for the expected FDR. Can that result be used as is in Theorem 3.1? I don’t see how.\n", "limitations": "N/A", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths \n\n-\tThe paper tackles an important problem of understanding feature relevance in a high-dimensional classification setting with sparse features. \n-\tThe empirical results support their claim\n-\tThe paper is relatively straightforward to parse through.\n\n\nWeaknesses\n\n-\tFrom a theoretical standpoint, more work can be done in explaining the technical novelty in extending the method of dCRT to the setting of classification. Right now the revised estimator, which used the second derivative of the logistic function instead, has minimal intuition for how it is derived, and what makes this problem a technical challenge. \n-\tIt is unclear which part of the estimator / results are to do with increasing computational efficiency vs. novel statistical results for the high-dimensional classification setting. This also makes it hard to tease out how to evaluate the estimator and the results. \n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}]

Author rebuttal:

Decision: Accept