Ex1: detection of active brain regions in neuroimaging

Find voxels whose average activity differs between two groups of samples

• “Word vs face” contrast: find regions with higher signal for word than face images
• \(m=52,000\) voxels, \(n0=n1=25\) subjects

Ex2: differential expression in genomics

Find genes whose average activity differs between two groups of samples

• Leukemia gene expression data (Chiaretti et. al., Clin. cancer res, 2005)
• \(m=12,625\) genes, \(n_0 = 42\) NEG patients vs \(n_1=37\) BCR/ABL patients

FDR control for large-scale multiple testing

Strategy: one test for each feature (gene/voxel) + choose significance threshold

Post hoc inference¹

• \(\mathcal{H}=\{1, \dots m\}\) null hypotheses
• \(\mathcal{H}_0 \subset \mathcal{H}\): true null hypotheses
• For \(S \subset \mathcal{H}\), \(|S \cap \mathcal{H}_0|\) : number of “false positives” within \(S\)

Neuroimaging data: glass brain plot (Simes/ARI)

True Discovery Proportion (TDP): TDP = 1-FDP

Genomic data: volcano plot (Simes/ARI)

Post hoc inference via Joint Error Rate control

\[\textrm{Goal: } V_\alpha \quad s.t. \quad\quad \mathbb{P} \left( \textcolor{red}{\forall S \subset \mathcal{H}},\:\:\: |S \cap \mathcal{H}_0| \leq V_\alpha(S) \right) \geq 1-\alpha\]

JER calibration¹: highway to powerful post hoc inference

Find \(\mathbf{t} = (t_k)_{k}\) such that \(\mathbb{P}\left( \exists k, p_{(k:\mathcal{H}_0)} \leq t_k \right) \leq \alpha\)

Genomic data: volcano plot (Parametric Simes¹)

Genomic data: volcano plot (Semiparametric Simes¹)

IIDEA: R/shiny app for interactive differential analysis¹

Controlled variable selection

• \(X_j \in \mathbb{R}^n\): brain activity data for each voxel \(j=1 \dots p\)
• \(Y \in \mathbb{R}^n\): outcome of interest

Goal: select a subset of variables significantly associated with \(Y\).

Knockoffs in a nutshell

Knockoffs viewed by Candès

Source: https://web.stanford.edu/group/candes/knockoffs/

JER control in the knockoff framework

Associated JER: \(\quad \quad JER(\mathbf{t}) = \mathbb{P} \big(\exists k \in\{1,\dots,p_0\} \::\: \pi_{(k:\mathcal{H}_0)} < t_k \big)\)

Numerical results on semi-simulated data

• neuroimaging data from Human Connectome Project (HCP): 7 contrasts/tasks
  • \(n=778\) subjects, \(p= 156,374\) reduced to \(1,000\) by Ward parcellation
• \(\mathbf{Y} = \mathbf{X}\boldsymbol{\beta}^* + \sigma\boldsymbol{\epsilon}\), with \(X\) = dataset 1 and \(\boldsymbol{\beta}^*\) estimated from dataset 2

Empirical FDP and power for 42 contrast pairs (only CT aims for posthoc FDP control¹

Results on HCP data

Detections for the SOCIAL contrast (Social motion vs random motion)

only classical KO (FDR \(< 0.2\)) and KOPI (FDP \(< 0.2\) with proba 0.9) yield discoveries

Construction of valid knockooffs

Assuming that \(X \sim \mathcal{N}(0, \Sigma)\), Candès et al. (2018) provide a valid construction of knockoffs such that:

\[ [X, \tilde{X}] \sim \mathcal{N}(0, \mathbf{G}), \quad \text { where } \mathbf{G}=\left(\begin{array}{cc} \boldsymbol{\Sigma} & \boldsymbol{\Sigma}-\operatorname{diag}(s) \\ \boldsymbol{\Sigma}-\operatorname{diag}(s) & \boldsymbol{\Sigma} \end{array}\right) \]

Failure of Gaussian knockoffs

Simulation setup:

• \(n=500\), \(p=500\), \(X\) gaussian with spatial structure
• \(Y = X \beta^* + \sigma \epsilon\) with \(\beta \in \{0,1\}^p\) and \(\Vert \beta^* \Vert_0=0.1 p\)

Proposed approach for knockoff generation

Goals

• avoid covariance matrix estimation
• maintain low computational cost

Proposed algorithm (adapted from Candès et al. (2018))

Improved knockoff generation

No obvious violation of exchangeability for nonparametric knockoffs!

Diagnosing Knockoffs exchangeability

Tool: Classifer Two-Sample Test ¹

Diagnosing Knockoffs exchangeability: results

Conclusion

• Post hoc inference = powerful & intepretable alternative to FDR control
  • can be obtained by JER control…
  • … which can be obtained when joint null distribution of test statistics is known or can be sampled from
• sanssouci: available R and python implementation + shiny app

Some other developments on JER control

• multivariate linear models (Davenport, Thirion, and Neuvial 2025)
• conformal inference (Blain, Thirion, and Neuvial 2025)

Ongoing applications

• gene-set differential analysis of scRNA-seq (PhD thesis of Sara Fallet)
• Hi-C data (PhD thesis of Elise Jorge)

Computing the interpolation bound

\[V_\alpha(S) = \min_{1\leq k \leq |S|} (k-1) + \sum_{i \in S} \mathbb{1}\{ p_i \geq t_k \}\]

How to obtain sharper JER control?

Step 1: Semiparametric JER control

Given a family of functions \((t_k)_k\) (e.g. \(t_k(\lambda) = \lambda k/m\)), estimate from the data the largest \(\lambda\) such that

\[\mathbb{P}\left( \exists k, p_{(k:\mathcal{H}_0)} \leq t_k(\lambda) \right) \leq \alpha\]

Step 2: Nonparametric JER control

Estimate from the data the “largest” \(\mathbf{t} = (t_k)_{k}\) such that

\[\mathbb{P}\left( \exists k, p_{(k:\mathcal{H}_0)} \leq t_k \right) \leq \alpha\]