Rank tests for outlier detection

Chiara G. Magnani[Department of Economics, Management and Statistics, University of Milano-Bicocca.] Aldo Solari[Department of Economics, Management and Statistics, University of Milano-Bicocca.]

In novelty detection, the objective is to determine whether the test sample contains any outliers, using a sample of controls (inliers). This involves many-to-one comparisons of individual test points against the control sample. A recent approach applies the Benjamini-Hochberg procedure to the conformal $p$-values resulting from these comparisons, ensuring false discovery rate control.

In this paper, we suggest using Wilcoxon-Mann-Whitney tests for the comparisons and subsequently applying the closed testing principle to derive post-hoc confidence bounds for the number of outliers in any subset of the test sample. We revisit an elegant result that under a nonparametric alternative known as Lehmann's alternative, Wilcoxon-Mann-Whitney is locally most powerful among rank tests. By combining this result with a simple observation, we demonstrate that the proposed procedure is more powerful for the null hypothesis of no outliers than the Benjamini-Hochberg procedure applied to conformal $p$-values.

Keywords: Conformal p-values; Closed testing; False discovery rate; Lehmann's Alternative; Wilcoxon-Mann-Whitney test.

§ INTRODUCTION

In novelty detection, the main goal is to determine whether the test sample contains any outliers, using a sample of inliers. A portion of the inliers (the training sample) is used to fit an one-class classification model, which assigns conformity scores to the observations in the test sample. These conformity scores are then compared to the conformity scores of the remaining inliers (the control sample), leading to many-to-one comparisons of individual test points with the control sample.

A recent line of research successfully combines conformal inference [56] with the Benjamini-Hochberg (BH) [3] procedure for false discovery rate (FDR) control, as discussed in [1, 43, 15]. Since this approach requires the split of the inliers sample into a training sample and a control sample, we refer to it as “split-conformal outlier detection”.
See [41] for an extension using control-unlabeled learning and [36] for an application of conformal calibration to standard early stopping training techniques.

In this paper, we propose using Wilcoxon-Mann-Whitney (WMW) tests [58, 39] for the comparisons and then applying the closed testing principle [42] to derive simultaneous $(1-\alpha)$ confidence lower bounds for the number of outliers in any subset of the test sample, following the methodology detailed in [19]. These bounds are often referred to as post-hoc confidence bounds.

In novelty detection, we can differentiate between two objectives: outlier discovery (ODI), which aims to locate specific outliers, and outlier detection (ODE), which focuses on determining the presence of outliers without precise localization. ODI is usually preferred because it implies ODE.
However, under the rare/weak effects scenario [11], which combines the effect sparsity hypothesis (few outliers in the test set) and the effect weakness hypothesis (outliers are not individually detectable), ODI may become very challenging or even impossible. In such situations with few and elusive outliers, ODE may be more relevant. Focusing solely on discovery might overlook the simpler task of detection, risking a non-informative inference.

In § 2, we review the necessary statistical background and the notation we are going to use throughout this work.
In § 3 we will revisit an elegant result in [33], which shows that the WMW test is the optimal rank test for outlier detection under a nonparametric alternative known as Lehmann's alternative.
We make a simple observation: the Benjamini and Hochberg applied to conformal $p$-values,
is a rank test for the global null hypothesis of no outliers. Consequently, the probability of not obtaining any informative inference is greater for the Benjamini-Hochberg procedure using conformal $p$-values compared to the WMW test.
In § 4, we compare the methods in some real data experiments. Final remarks and future directions are discussed in § 5.
Let us introduce the necessary notation.

§.§ Setup and existing approaches

Suppose we fit an one-class classification model $\hat{s}: \mathbb{R}^q \mapsto \mathbb{R}$ on a training sample with $l$ independent and identically distributed (i.i.d.) observations
$tr_1,\ldots,tr_l$ from some unknown null distribution. Observations drawn from this null distribution are the inliers.
Measurements are real vectors, i.e. $tr_{i} \in \mathbb{R}^q$, whose dimension $q$ is potentially large with $q > l$.

Then, we use the trained model to assign scores on two different samples: the control sample $co_1,\ldots,co_m$ of $m$ i.i.d. inliers, providing the score $X_i=\hat{s}(co_i; tr)$, and the test sample $te_1,\ldots,te_n$ of $n$ independent observations, providing the score $Y_i=\hat{s}(te_i; tr)$; see Figure <ref>. Some of the observations of the test sample may be inliers, while others are drawn from different distributions and are called outliers. Large values of $\hat{s}(te_i;tr)$ provide evidence that $te_i$ may be an outlier.

Visualization of the training, control and test samples and the corresponding scores $\hat{s}(\cdot;tr)$.

Let $(X_1,\ldots,X_m)$ and $(Y_1,\ldots,Y_n)$ be the score vectors of the control and test samples, respectively, with combined sample of size $N=m+n$. We assume that $X_1,\ldots,X_m$ are i.i.d.
with cumulative distribution function (c.d.f.) $F \in \mathcal{F}$, where $\mathcal{F}$ denotes the family of continuous distributions on the real line. $Y_1,\ldots,Y_n$ are independent, each $Y_i$ with c.d.f. $G_i$, where $G_i$ has the same support of $F$ and is continuous so that there are no equal values among $X_1,\ldots,X_m,Y_i$. We also assume that $G_i$ is stochastically larger than $F$, i.e. $G(z) \leq F(z)$ for all $z$.

We are interested in testing the null hypothesis that the $i$th observation of the test sample is an inlier, i.e.

\begin{eqnarray}\label{eq-Hi}
H_i: G_i = F
\end{eqnarray}

versus the alternative hypothesis that $te_i$ is an oulier, i.e. $K_i: G_i < F$ , the latter interpreted as $G_i(z) \leq F(z)$ for all $z$ with strict inequality for at least one value of $z$.

Let $I_0 = \{i \in [n]: G_i = F\}$ be the index set of the true null hypotheses, where $[n]$ denotes $\{1,\ldots,n\}$. The assumptions we will require are as follows:

$(X_1,\ldots,X_m, (Y_i)_{i \in I_0})$ is exchangeable.

Assumption A<ref> is satisfied since $(X_1,\ldots,X_m, (Y_i)_{i \in I_0})$ are i.i.d. with c.d.f. $F$.
The pooled vector of scores $(X_1,\ldots,X_{m},Y_{1},\ldots,Y_n)$ has no ties almost surely.

Assumption A<ref> is satisfied since we assumed that $F$ and $G_i$ are continuous.

Assumptions A<ref> and A<ref> imply the validity of the WMW two-sample test, also known as the Wilcoxon rank sum test [58] or the Mann-Whitney $U$ test [40].
The Wilcoxon rank sum two-sample test, with the 1st sample being $X_1,\dots,X_n$ and the 2nd “sample” being just $Y_i$, rejects $H_i$ if the rank of $Y_i$ in the sequence $(X_1,\ldots,X_n,Y_i)$ is large or, equivalently, if the permutation $p$-value $p^*_i$ is small:

\begin{eqnarray}\label{eq-permp}
p^*_i &=& \frac{1}{(m+1)}\left(1+ \sum_{k = 1}^{m} \mathds{1}\{ X_k \geq Y_{i} \}\right).
\end{eqnarray}

Permutation $p$-values are also referred to as “conformal $p$-values” [1] or “empirical $p$-values” [41].
We will use the terminology “permutation $p$-values,” since rank tests are a special case of permutation tests;  for an historical overview of permutation test, refer to [23]; for similarities between conformal prediction and rank tests, see [32].

The multiple hypotheses $H_1, \ldots, H_n$ represent a well-known problem in multiple testing, known as many-to-one comparisons to a control [12]. However, since the $n$ comparisons of $(Y_i)$ to $(X_1,\ldots,X_m)$ all use the same control sample, the above operation induces a special dependence structure among the permutation $p$-values $(p^*_1,\ldots,p^*_n)$.

[1] proved that the permutation $p$-values $(p^*_1,\ldots,p^*_n)$ are positive regression dependent on a subset (PRDS) [4].
On the one hand, the PRDS property implies that Fisher's method [13] applied to permutation $p$-values $(p^*_1, \ldots, p^*_n)$, may fail due to its requirement of independence among $p$-values. Rejecting the global null hypothesis

\begin{eqnarray}\label{eq-Hg}
\end{eqnarray}

if the Fisher's test statistic $-2\sum_{i=1}^{n}\log(p^*_i)$ exceeds the $(1-\alpha)$-quantile of the $\chi^2$ distribution with $2n$ degrees of freedom does not effectively control the type I error at level $\alpha$.
On the other hand, the PRDS property ensures that the BH procedure, utilizing permutation $p$-values $(p^*_1, \ldots, p^*_n)$, called BH$^*$ procedure hereafter,
successfully controls the FDR at level $\alpha$. [43] showed that the FDR of the BH$^*$ procedure is exactly $\alpha\pi_0$ if and only if $\alpha(m+1)/n$ is an integer, otherwise it is strictly less than $\alpha\pi_0$, where $\pi_0=|I_0|/n$ is the unknown proportion of inliers in the test sample.

If $\pi_0$ is known, a simple approach is to replace the target FDR level with $\alpha/\pi_0$. However, since $\pi_0$ is not known in practice, it may be estimated as follows:

\begin{eqnarray}\label{eq-pi0Hat}
\hat{\pi}_0 &=& \frac{1+\sum_{i=1}^{n} \mathds{1}\{p^*_i > \lambda\} }{n(1-\lambda)}
\end{eqnarray}

where $\lambda = j/(n+1)$ for any pre-specified integer $j$, known as the Schweder-Spjøtvoll or Storey's estimator [47, 51, 52].
Then, the BH$^*$ procedure is applied at level $\alpha/\hat{\pi}_0$. [1] demonstrated that this $\pi_0$-adaptive BH$^*$ procedure effectively controls the FDR at level $\alpha$.

[43] in § 8.3 and [41] in § 1.1 and § 3.2 advocate utilizing the PRDS property of the permutation $p$-values to derive post-hoc bounds for the false discovery proportion (FDP). The PRDS property implies that the Simes test [49] applied to permutation $p$-values is a valid test, which can be used as local test in a closed testing procedure, resulting in upper confidence bound for the FDP, simultaneously over all subsets of hypotheses [18]. Furthermore, [17] proved that if the local tests are of exact size $\alpha$, then the  post-hoc bounds are admissible.

In the model discussed above, the unknown distributions $F$ and $G_i$ are not specified, and under $H_i$, the common distribution $F = G_i$ remains unspecified as well. Neither the general model nor the null hypothesis can be described in terms of finite-dimensional parameters, which is why these models and hypotheses are termed nonparametric.

It is worthwhile to compare the nonparametric model with the classical parametric model, where unknown distributions have a given parametric form with only a finite-dimensional parameter left unspecified. For the many-to-one comparisons problem in the one-way ANOVA setup, we assume $X_1, \ldots, X_m$ are i.i.d. $N(\mu, \sigma^2)$, and $Y_1, \ldots, Y_m$ are independent with $Y_i$ drawn from $N(\mu_i, \sigma^2)$. Testing of $H_i: \mu_i = \mu$ against $K_i: \mu_i > \mu$ is performed using test statistics $T_i = (Y_i - \bar{X})/(\hat{\sigma}\sqrt{1 + 1/m})$, following a Student $t$ distribution with $m-1$ degrees of freedom, where $\hat{\sigma}^2$ is the pooled variance estimator of the common $\sigma^2$.
Dunnett's procedure finds the  $(1-\alpha)$-quantile $c_{\alpha}$ of the maximum of the $n$-variate $t$ distribution $(T_1, \ldots, T_n)$ with a correlation matrix that has off-diagonal elements equal to $1/(1+m)$. It rejects $H_i$ if $T_i > c_{\alpha}$ while controlling the familywise error rate (FWER) at level $\alpha$. An improvement of Dunnett's procedure, known as single-step Dunnett, applies the closed testing principle [42]
(referred to as step-down Dunnett; see [7]).  A multiple comparison rank sum test, for comparing treatments with a control in a one-way classification with equal numbers of observations, is presented in [50].

Criterion     Model     2cMethod

Type I error     nonparametric     2cFisher's combination     § 2.3

FWER     parametric     2cDunnett's max$t$     § 1.1

2*FDP      2*nonparametric      2*Closed testing with      1c(adaptive) Simes tests      § 2.2

1cWMW rank tests      § 3

FDR     nonparametric     2c(adaptive) Benjamini-Hochberg      § 2.1

Overview of existing and proposed methods for sample-split outlier detection for error criterion and model considered.

Table <ref> provides an overview of existing methods (rows 1, 2 and 4) and proposed methods (row 3) for error criterion and the model considered.

§.§ Contributions

In § 2.1, we review the (adaptive) BH algorithm with permutation $p$-values $p^*_1,\ldots,p^*_n$ and the theoretical results from [1, 43]. § 2.2 is dedicated to the post-hoc confidence bounds obtained by closed testing with (adaptive) Simes combination of permutation $p$-values $p^*_1,\ldots,p^*_n$. An algorithm for computing the bounds in linear time is provided, and sufficient conditions for the admissibility of the procedure are discussed.
Our main contributions are the following:

* In § 2.3, we highlight the difference between Fisher’s combination of marginal permutation $p$-values $p^*_1,\ldots,p^*_n$ and the permutation test using Fisher’s combination as test statistic, and draw similar conclusions for Simes' test.

* In § 3 we focus on rank tests, and in § 3.1 we provide an algorithm for computing post-hoc bounds with WMW tests having the same computational complexity of the BH$^*$ procedure.
The central result of this paper is given in § 3.2: for Lehmann's alternative with $k=2$, the proposed procedure is locally most powerful for outlier detection among all procedures that are invariant under the group of continuous and strictly increasing
transformations of the conformity scores.
In particular, we demonstrate that
the probability of not making any discovery is higher for the BH$^*$ procedure than for the WMW test.

In § 4, we compare the existing and proposed methods, applying them to real datasets. R code reproducing the examples is available at https://github.com/chiaragaiamagnani/nout/tree/master/Examples. Finally, in § 5 we discuss potential future directions.

§ FDR CONTROL AND FDP POST-HOC BOUNDS WITH PERMUTATION $P$-VALUES

§.§ Benjamini-Hochberg procedure with permutation $p$-values

The BH algorithm applied to generic $p$-values $p_1,\ldots,p_n$ at level $\alpha$ returns the index set of the discoveries:

\begin{eqnarray}\label{eq-DBH}
D_{\mathrm{BH}} &=&  \{i \in [n]:  p_i \leq \alpha d_{\mathrm{BH}}/n\}
\end{eqnarray}

\begin{eqnarray}\label{eq-dBH}
d_{\mathrm{BH}} &=&  \max\left\{k \in \{0,\ldots,n\}: \sum_{i=1}^{n} \mathds{1}\{p_i \leq \alpha k/n \} \geq k \right\}.
\end{eqnarray}

is the number of discoveries made by the BH procedure, with $d_{\mathrm{BH}}=|D_{\mathrm{BH}}|$, where $|\cdot|$ denotes the cardinality of a set.

All hypotheses $H_i$ with index $i\in D_{\mathrm{BH}}$ are rejected by the BH procedure, indicating that observations $te_i$ are classified as potential outliers. Some of these discoveries may be true outliers (true discoveries) or true inliers (false discoveries), with the expected proportion of inliers among the discoveries bounded by $\alpha\pi_0$.

When $\pi_0$ is expected to be not close to 1, $\alpha\pi_0$ falls below the target level $\alpha$, making the BH procedure slightly too conservative. This motivates adjusting the level into $\alpha/\hat{\pi}_0$ with $\hat{\pi}_0$ defined in (<ref>), resulting in a $\pi_0$-adaptive version of the BH algorithm.

The properties of the (adaptive) BH$^*$ procedure are provided in the following proposition.

Under assumptions A<ref> and A<ref>, the BH$^*$ procedure  controls the
FDR at level $\pi_0\alpha$. In particular, the FDR of the BH$^*$ procedure is exactly $\pi_0\alpha$ if and only if $\alpha(m+1)/n$ is an integer.

The $\pi_0$-adaptive version of the BH$^*$ procedure and using the $\pi_0$ estimator  (<ref>) controls the FDR at level $\alpha$.

§.§ Closed testing with Simes' combination of permutation $p$-values

The closed testing procedure was initially proposed for FWER control. [19] observed that it can also be used to derive simultaneous FDP bounds.
The procedure consists of three steps:

i. Test all intersection hypotheses

\begin{eqnarray}\label{eq-HI}
H_I : \bigcap_{i \in I} H_i
\end{eqnarray}

for $I \subseteq [n]$ at level $\alpha$. Throughout the paper, we refer to such tests as local tests. We use $\phi_I \in \{0,1\}$ to denote the test result of $H_I$, with 1 indicating rejection;

ii. Adjust the local tests for multiple testing by $\bar{\phi}_I = \min\{\phi_J: J\supseteq I\}$;

iii. For any set $S \subseteq [n]$, the lower bound for the number of true discoveries in $S$ is given by

\begin{eqnarray}\label{eq-dS}
d(S) &=& \min_{K \subseteq S}\{|S \setminus K|:\bar{\phi}_K = 0\}.
\end{eqnarray}

Let  $I_1 =[n] \setminus I_0$ the index set of the false hypotheses (outliers). Then $d$ defined in (<ref>) satisfies

\begin{eqnarray}\label{eq-TDguarantee}
\mathrm{pr}(d(S) \leq |I_1 \cap S| \mathrm{\,\,for\,\,all\,\,}S  ) &\geq& 1 -\alpha.
\end{eqnarray}

We will write $d$ for $d([n])$, the lower bound on the overall number of true discoveries (cfr. the number of true outliers in the test sample).

In general, computing $d(S)$ in (<ref>) involves the evaluation of exponentially many tests, which hinders its practical application. For specific local tests, however, there exist shortcuts to derive the closed testing results in polynomial time [18, 17, 53].

Simes' test [49] for $H_I$ combines the $p$-values $(p_i,i\in I)$:

\begin{eqnarray}\label{eq-Simeslocaltest}
\phi^{\mathrm{Simes}}_I &=&  \mathds{1}\left\{ \min_{k \in \{1,\ldots,|I|\}} \{|I|p_{(k:I)}/k\} \leq \alpha \right\}
\end{eqnarray}

where $p_{(k:I)}$ is the $k$th smallest $p$-value in $\{p_i: i \in I \}$.
With Simes' local tests, the shortcut described in [18] and given in Algorithm <ref> allows calculation of $d_{\mathrm{Simes}}(S)$ in linear time, after an initial preparatory step that takes
linearithmic time (i.e. $O(n \log n)$ time in big $O$ notation).

$p$-values $p_1,\ldots,p_n$; index set $S$; significance level $\alpha$.
$(1-\alpha)$-confidence lower bound $d_{\mathrm{Simes}}(S)$ for the number of true discoveries in $S$

Sort the $p$-values $p_{(1)} \leq \ldots \leq p_{(n)}$;

Compute $h=\max\{0\leq k \leq n: p_{(n-k+i)} > i\alpha/k \,\,\mathrm{for\,\,} i=1,\ldots,k\};$

Compute  $d_{\mathrm{Simes}}(S) =  \max_{1\leq u\leq |S|} 1-u+|\{i\in S: p_{i}\leq u\alpha/h\}|$;

Shortcut for closed testing with Simes local tests.

The properties of the procedure implemented in Algorithm <ref> applied to permutation $p$-values $p^*_1,\ldots,p^*_n$ (called CT$^{*}_{\mathrm{Simes}}$ hereafter) are provided in the following proposition.

Under assumptions A<ref> and A<ref>, Simes' test $\phi^{\mathrm{Simes}}_{I_0}$ in (<ref>) combining permutation $p$-values $(p^*_i, i \in I_0)$ controls the Type I error at level $\alpha$, i.e.,
$\displaystyle \mathrm{pr}(\phi^{\mathrm{Simes}}_{I_0} = 1; H_{I_0}) \leq \alpha$, with equality if and only if $\alpha (m+1)/|I_0|$ is an integer. Consequently, CT$^{*}_{\mathrm{Simes}}$ provides the true discovery guarantee (<ref>).

The proof follows by noting that Simes' test rejects $H_{I_0}$ if and only if the BH$^{*}$ procedure applied to permutation $p$-values $(p_i^*, i \in I_0)$ rejects at least one hypothesis, which happens with probability bounded by $\alpha$ by Proposition 1. This probability is exactly $\alpha$ if and only if $\alpha(m+1)/|I_0|$ is an integer by Theorem 3.1 in [43].

If $\alpha(m+1)/i$ is an integer for all $i=1,\ldots,n$, then all local tests are of exact size $\alpha$. Then, the local tests are admissible, and admissibility of the CT$^{*}_{\mathrm{Simes}}$  follows from Theorem 3 in [17]; see Appendix C for further details. If $\alpha(m+1)/i$ is not an integer for some $i$, it is still possible to use a randomized critical region to exhaust the $\alpha$-level.

Line 2 of Algorithm <ref> gives $h$, which is a $(1-\alpha)$-confidence upper bound for the number of true hypotheses $n_0 = |I_0|$. Then  $d_{\mathrm{Simes}} = n - h$ is a $(1-\alpha)$-confidence lower bound for the number of false hypotheses (cfr. number of true outliers in the test sample) $n_1 = |I_1|$. The index set of the discoveries (cfr. localized outliers) is given by

\begin{eqnarray*}
D_{\mathrm{Simes}} &=&  \{ i \in [n]: d_{{\mathrm{Simes}}}(\{i\}) = 1\} \\
&=& \{ i \in [n]: p_{i} \leq \alpha/h\}.
\end{eqnarray*}

All hypotheses $H_i$ with index $i \in D_{\mathrm{Simes}}$ are rejected by CT$_{\mathrm{Simes}}$ - referred to as Hommel's procedure [29] - while controlling the FWER at level $\alpha$:

\begin{eqnarray*}
\mathrm{pr}( |D_{\mathrm{Simes}} \cap I_0 |>0 ) &\leq& \alpha.
\end{eqnarray*}

Let $\phi_g=\phi_{[n]}$ denotes the test for the global null hypothesis $H_g$ in (<ref>). Closed testing provides the following quantities:

* $\mathds{1}\{\phi_g=1\}$ to indicate outlier detection (ODE);

* $d$ to provide outlier enumeration (OEN);

* $|D|$ to quantify outlier discovery (ODI);

We observe that between ODE and ODI, there exists a third objective known as outlier enumeration (OEN). This involves counting the number of true outliers in the test sample without necessarily pinpointing their exact locations. OEN is implied by ODI (by using $|D|$ for OEN) and, in turn, implies ODE (by using $\mathds{1}\{d>0\}$ for ODE).
These quantities can be applied to any subset $S$ of the test sample, obtaining $\mathds{1}\{\phi_S=1\}$ for ODE in $S$, $d(S)$ for OEN in $S$, and $|D \cap S |$ for ODI in $S$.

Post-hoc confidence bounds $d(S)$ obtained by closed testing satisfy the following relations:

\begin{eqnarray}\label{CTrelations}
\mathds{1}\{\phi_S=1\} \geq \mathds{1}\{d(S)>0\}, \qquad d(S) \geq |D \cap S |
\end{eqnarray}

These relations imply that post-hoc confidence bounds may have less power than the test $\phi_S$ for $H_S$ in the context of ODE in $S$, but they can be more powerful than counting the discoveries for ODI in $S$. However,

\begin{eqnarray}\label{eq-gatekeeper1}
\mathds{1}\{\phi_g=1\} = \mathds{1}\{d>0\}.
\end{eqnarray}

This means that the global test $\phi_g$ serves as the first gatekeeper of the closed testing procedure, opening the door to exponentially many tests only if $\phi_g$ rejects the global null $H_g$.

To distinguish between discoveries with FWER and FDR guarantees, we will refer to discoveries made by closed testing as FWER discoveries, while those made by an FDR controlling procedure will be called FDR discoveries.

The Simes test for $H_g$ and the BH procedure share the same critical values, and we can expect some relationships between the two methods.
A simple observation is that the event where the CT$_{\mathrm{Simes}}$ and the BH procedure make at least one non-trivial statement is identical:

\begin{eqnarray}\label{eq-weakFWER}
\mathds{1}\{d_{\mathrm{Simes}} > 0\}=\mathds{1}\{d_{\mathrm{BH}} > 0\}.
\end{eqnarray}

Consequently, the two methods have the same weak FWER control and the same power for Outlier Detection (ODE). Furthermore, [18] showed that the confidence bound $d_{\mathrm{Simes}}$ is between the number of Hommel-FWER discoveries and the number of BH-FDR discoveries, i.e.

\begin{eqnarray}\label{eq-dSimesdBH}
| D_{\mathrm{Simes}} | \leq d_{\mathrm{Simes}} \leq | D_{\mathrm{BH}} |.
\end{eqnarray}

Finally, note that the BH discoveries do
not provide the true discovery guarantee, i.e.

\begin{eqnarray}\label{eq-BHnoguarantee}
\mathrm{pr}(| D_{\mathrm{BH}} \cap S| \leq |I_1 \cap S| \mathrm{\,\,for\,\,all\,\,}S  )
\end{eqnarray}

may be less than $1-\alpha$.

In Appendix <ref> we show that the same idea as the adaptive BH procedure applies to Simes test, resulting in a $\pi_0$-adaptive version of Simes' test. The corresponding closed testing procedure producing post-hoc confidence bounds $d_{\mathrm{ASimes}}$ will be denoted as CT$^*_{\mathrm{ASimes}}$.

§.§ Global testing with Fisher’ combination of permutation $p$-values

For testing the global null $H_g: \bigcap_{i=1}^{n}H_i$, Fisher's method combining $p$-values $p_1,\ldots,p_n$ is the test $\phi_{\mathrm{Fisher}} = \mathds{1}\{T_{\mathrm{Fisher}} \geq k_{\alpha}\}$ where

\begin{eqnarray}\label{fisher}
T_{\mathrm{Fisher}} &=& - 2 \sum_{i=1}^{n}\log(p_i)
\end{eqnarray}

is the Fisher's combination statistic and $k_{\alpha}$ is the $(1-\alpha)$-quantile of the $\chi^2$ distribution with $2n$ degrees of freedom. The test controls the probability of Type I error at level $\alpha$ if the $p$-values  $p_1,\ldots,p_n$ are independent and valid, i.e. the distribution of $p_i$ is stochastically larger than uniform when $H_i$ is true: $\mathrm{pr}(p_i \leq u ; H_i) \leq u$.

[1] showed that rejecting $H_g$ if $T^*_{\mathrm{Fisher}} = - 2 \sum_{i=1}^{n}\log(p^*_i) \geq k_{\alpha}$ does not control the Type I error because of the PRDS dependence among permutation $p$-values, and they proposed the adjusted test statistic

\begin{eqnarray}\label{adjustedFisher}
\tilde{T}^*_{\mathrm{Fisher}} &=& \frac{- 2 \sum_{i=1}^{n}\log(p^*_i) + 2(\sqrt{1+\gamma} -1)n}{\sqrt{1+\gamma}}
\end{eqnarray}

where $\gamma$ solves $m=\lfloor n\gamma \rfloor$.

For testing $H_g$, there is a distinction between

i. the test based on adjusted Fisher's combination of marginal permutation $p$-values $p_1^*,\ldots,p_n^*$; i.e. $\phi_g^{\mathrm{Fisher}}=\mathds{1}\{ \tilde{T}^*_{\mathrm{Fisher}} \geq k_\alpha\}$;

ii. the permutation test using Fisher's test statistic; i.e. $\psi_g^{\mathrm{Fisher}}=\mathds{1}\{ T^*_{\mathrm{Fisher}} \geq k^*_\alpha\}$, where $k^*_\alpha$ is the $(1-\alpha)$-quantile $k^*_\alpha$ of the permutation null distribution of $T^*_{\mathrm{Fisher}}$.

The former approach relies on the exchangeability of $(X_1,\ldots,X_m,Y_i)$ for the validity of the marginal permutation $p$-value $p^*_i$ for $i = 1,\ldots,n$ (and on the correction of Fisher’s statistic for the dependence among $p$-values $p_1^*,\ldots,p_n^*$), while the latter approach directly employs the exchangeability of the full vector $(X_1,\ldots,X_m,Y_1,\ldots,Y_n)$ in order to compute the $(1-\alpha)$-quantile $k^*_\alpha$ of the permutation distribution of $T^*_{\mathrm{Fisher}}$.
The permutation test also automatically adapts to the dependence among $p$-values $p_1^*,\ldots,p_n^*$.

The theory of permutation tests and related references can be found in the books [9],
[35],  [45], as well as in the works of [24, 25].
Permutation-based procedures for controlling the FWER, as well as bounding the FDP, are discussed in [57],  [26], [5] and [55].

Let $(Z_1,\ldots,Z_N) = (X_1,\ldots,X_m,Y_1,\ldots,Y_n)$ be the combined sample. Note that we can rewrite the permutation $p$-values in (<ref>) as

\begin{eqnarray*}\label{permpZ}
p^*_i &=& \frac{1}{(m+1)}\left(1+ \sum_{k = 1}^{m} \mathds{1}\{ Z_k \geq Z_{m+i} \}\right).
\end{eqnarray*}

Consider the permutation group $\Pi = \{\pi_1,\ldots,\pi_w\}$ whose $w=N!$ elements are permutations of the set of integers $[N]$, that is
$\pi:[N]\mapsto [N]$ with  $(\pi(1),\ldots,\pi(N))$ denoting a permutation of $[N]$.
For any $\pi \in \Pi$, the permutation $p$-values computed for the permuted pooled sample $(Z_{\pi(1)},\ldots, Z_{\pi(N)})$ are given by

\begin{eqnarray*}\label{permppermZ}
p^\pi_i &=& \frac{1}{(m+1)}\left(1+ \sum_{k = 1}^{m} \mathds{1}\{ Z_{\pi(k)} \geq Z_{\pi(m+i)} \}\right).
\end{eqnarray*}

with $p_i^* = p_i^{\pi^*}$ for the identity permutation $\pi^*$.
Then, calculate the $(1-\alpha)$-quantile $k^*_\alpha$ of the set of values
$\{ -2 \sum_{i=1}^{n} \log(p^\pi_i), \pi \in \Pi\}$ and reject $H_g$ if $T^*_{\mathrm{Fisher}}\geq k^*_\alpha$.

Under assumption A<ref> with $I_0 = [n]$ and assumption A<ref>, the permutation test  $\psi^{\mathrm{Fisher}}_{g}=\mathds{1}\{ T^*_{\mathrm{Fisher}} \geq k^*_\alpha\}$ controls the Type I error at level $\alpha$, i.e.,
$\displaystyle \mathrm{pr}(\psi^{\mathrm{Fisher}}_{g} = 1; H_{g}) \leq \alpha$.

The computational cost of this procedure can become a concern, given that computing the test statistic for $N!$ permutations may be computationally prohibitive, particularly for moderately large $N$. However, for the permutation test using Fisher's statistic, it is possible reduce the group of permutations $\Pi$ to the subgroup $\tilde\Pi$ containing only permutations
allocating $m$ observations to the control sample and the remaining $n$ observations to the test sample, reducing the number of permutations to  $|\tilde{\Pi}|= { N \choose m }$. See Algorithm <ref> for the computation of the permutation $p$-value $p^*_{\mathrm{Fisher}}$.
If $|\tilde{\Pi}|$
is still too large, use a random sample of permutations from $\tilde{\Pi}$ [24].

Control scores $X_1,\ldots,X_m$ and test scores $Y_1,\ldots,Y_n$.
$p$-value of the permutation test using Fisher's test statistic

Let $(Z_1,\ldots,Z_N) = (X_1,\ldots,X_m,Y_1,\ldots,Y_n)$ be the combined sample

$\pi \in \tilde\Pi$
Compute $p_1^{\pi}, \ldots, p_n^{\pi}$ for $(Z_{\pi(1)}, \ldots, Z_{\pi(N)}  )$

Compute $T^\pi_{\mathrm{Fisher}} = -2 \sum_{i=1}^{n} \log(p^\pi_i)$

Compute $\displaystyle p_{\mathrm{Fisher}}^{*} = \frac{1}{|\tilde{\Pi}|} \sum_{\pi \in \tilde\Pi} \mathds{1}\{T^\pi_{\mathrm{Fisher}} \geq T^*_{\mathrm{Fisher}}\}$

Permutation test using Fisher's test statistic.

The two tests $\phi_g^{\mathrm{Fisher}}$ and $\psi_g^{\mathrm{Fisher}}$ may show discrepancies in small samples due to the distinct discreteness of their permutation construction and the varied application of the exchangeability hypothesis. To illustrate the potential differences, consider the scores:

\begin{eqnarray*} (x_1,\ldots,x_7) = (5.42,5.86,6.16,6.55,6.8,7,7.11), \quad (y_1,\ldots,y_5) = (6.51,7.56,7.61,7.84,11.5)
\end{eqnarray*}

These scores are taken from an example in [38].

The permutation $p$-values are

$$p^*_1=0.625, \quad p^*_2=\ldots=p^*_5=0.125.$$
and Fisher’s combination of permutation $p$-values gives $T^*_{\mathrm{Fisher}} =  17.576$ and $\tilde{T^*}_{\mathrm{Fisher}} = 15.786$ for $\gamma = n/m$. This gives $p$-values of 0.0626 and 0.106 using the $\chi^2_{2n}$ distribution. In contrast, the permutation test using Fisher’s test statistic gives $p^*_{\mathrm{Fisher}}=0.008$. Figure <ref> displays the permutation distribution.
One reason for the discrepancy between these $p$-values - 0.0626 and 0.106 versus 0.008 - is the discreteness of the permutation $p$-values $p^*_i$, which take values in the set $\{1/7, 2/7, \ldots, 1\}$.

Permutation distribution of the Fisher's statistic $T^*_{\mathrm{Fisher}} = -2 \sum_{i=1}^{n}  \log(p^*_i)$, Simes' statistic $T^*_{\mathrm{Simes}} = \min_{ k \in \{1,\ldots,n\} }  \{ n p^*_{(k)} / k \}$, Wilcoxon rank sum statistic $T_{\mathrm{W}} = \sum_{i=1}^{n} R_i$ and of the LMPI statistic for Lehmann's alternative with $k=7$, $T_{7} = \sum_{i=1}^{n} R_i(R_i + 1) \cdots (R_i + 5)$, where $R_i$  the rank of $Y_i$ in the sequence $Z_1,\ldots,Z_N$.

Similarly, there exists a distinction between the test based on Simes' combination of marginal permutation $p$-values, i.e. $  \phi^{\mathrm{Simes}}_g =  \mathds{1}\left\{ T^*_{\mathrm{Simes}} \leq \alpha \right\}$ where $
T^*_{\mathrm{Simes}} = \min_{k \in \{1,\ldots,n\}} \{ n p^*_{(k)}/k \}$ is the Simes's test statistic and the permutation test $\psi^{\mathrm{Simes}}_g =  \mathds{1}\left\{ T^{*}_{\mathrm{Simes}} \leq a_\alpha^{*}  \right\}$ where $a_\alpha^{*}$ is the $\alpha$-quantile of the permutation distribution of  $T^{*}_{\mathrm{Simes}}$.
The Simes's combination of permutation $p$-values is $T^*_{\mathrm{Simes} } = 0.156$, while the $p$-value for the permutation test using Simes test statistic is $p^*_{\mathrm{Simes}}=0.01$; refer to Figure <ref>.

Permutation tests $\psi_I^{\mathrm{Fisher}}$ testing $H_I$ can be used as local tests in a closed testing procedure to obtain post-hoc confidence bounds. However, a naive application of closed testing involves the evaluation of exponentially many permutation tests, which limits its practical application.
For a detailed theoretical treatment and the implementation of a branch-and-bound type of shortcut for the closed testing procedure using permutation sum-based local tests, including Fisher's, refer to [55].

In the forthcoming section, we will utilize another permutation sum-based local test, specifically the WMW test $\psi_I^{\mathrm{WMW}}$ for $H_I$, using  full exchangeability of $(X_1,\ldots,X_m,(Y_i)_{i \in I})$ under $H_{I}$ and efficiently compute the post-hoc confidence bound $d_{\mathrm{WMW}}(S)$ in linear time.

§ RANK TESTS

§.§ Closed testing with WMW tests

We propose to perform closed testing using the WMW test as local test.
We reject $H_{I}$ for large values of the Mann-Whitney $U$ statistic:

\begin{eqnarray}\label{eq-UI}
U_I &=& \sum_{i \in I}\sum_{k = 1}^{m}  \mathds{1}\{ X_{k} < Y_{i}\}
\end{eqnarray}

We reject $H_{I}$ if $U_I$ is larger than or equal to $(1-\alpha)$-quantile $c_{\alpha}(m,|I|)$ of the permutation null distribution of $U_I$, which depends on $\alpha$, $m$ and $|I|$.
The WMW rank test is given by

\begin{eqnarray}\label{eq-WMWlocaltest}
\psi^{\mathrm{WMW}}_I &=& \mathds{1}\{U_I \geq c_\alpha(m,|I|)\}.
\end{eqnarray}

The permutation null distribution of $U_I$
can be derived through recursion [39] and it is implemented in several software packages. For example, the R function $\alpha$$|I|$ computes $c_{\alpha}(m,|I|)$.

Note that for $I = \{i\}$, the Mann-Whitney $U$ statistic reduces to

\begin{eqnarray}\label{eq-Ui}
U_i &=& \sum_{k=1}^m\mathds{1}\{X_k < Y_i \}
\end{eqnarray}

and its $p$-value is the permutation $p$-value $p^*_i$ in (<ref>).
The test statistic $U_I$ in (<ref>) satisfies a special
property called separability [53], that is, the statistic is a function of summation of test scores for the
individual hypotheses, i.e. $U_I = \sum_{i \in I} U_i$.

We emphasize that $\psi^{\mathrm{WMW}}_I$ is a rank test for $H_I$ that relies on the full exchangeability of $(X_1,\ldots,X_m,(Y_i)_{i \in I})$ under $H_I$, as discussed in § 2.3.

The properties of the closed testing procedure with WMW rank tests (<ref>)
(called CT$^*_{\mathrm{WMW}}$ hereafter) are provided in the following proposition.

Under assumptions A<ref> and A<ref>, the WMW test $\psi^{\mathrm{WMW}}_{I_0}$ in (<ref>) controls the Type I error at level $\alpha$, i.e.,
$\displaystyle \mathrm{pr}(\psi^{\mathrm{WMW}}_{I_0} = 1; H_{I_0}) \leq \alpha$. Consequently,  CT$^*_{\mathrm{WMW}}$ provides the true discovery guarantee (<ref>).

CT$^*_{\mathrm{WMW}}$ requires the same computational time as the BH$^{*}$ procedure.

WMW test $\psi^{\mathrm{WMW}}_{I_0}$ is a permutation test, and control of Type I error follows from Theorem 1 in [24].

The BH$^*$ procedure takes linear time after the permutation $p$-values
$p^*_1,\ldots, p^*_n$ are computed and then sorted.
The same computational time is required by the closed testing shortcut given by Algorithm 1 in [53], which
computes $d_{\mathrm{WMW}}(S)$ in  linear time after the Mann-Whitney $U$ statistics $U_1,\ldots, U_n$ are computed and sorted.

Algorithm 1 in [53] can be used since the WMW statistic $U_I$ satisfies the monotonicity, symmetry and separability requirements;
we refer the interested readers to Appendix <ref> for formal definitions of the aforementioned terms.

The test $\psi^{\mathrm{WMW}}_{I}$ is of exact size $\alpha$ for the following values of $\alpha$: $\alpha_x^{|I|} = \mathrm{pr}(U_I \geq x; H_I)$ for $x \in \{0,1,\ldots,m|I|\}$, where the null distribution of $U_I$ is defined by using a recurrence formula [39].

For large values of $m$ and $|I|$, several approximations to the exact critical values $c_{\alpha}(m,|I|)$ of the $U_I$ test statistic have been proposed. These approximations include the Normal approximation suggested in [27], the Edgeworth approximation proposed in [14], and the Uniform approximation introduced in [8]. For example, we can use $k_\alpha(m,|I|)$, the $(1-\alpha)$-quantile of the

\begin{eqnarray*}
N\left(\frac{|I|m-1}{2}, \frac{m|I|(m + |I|+ 1)}{12}\right)
\end{eqnarray*}

applying the continuity correction factor to the mean of the Normal distribution, and obtain $k_\alpha(m,|I|) \approx c_\alpha(m,|I|)$ for large values of $m$ and $|I|$.
Admissibility of the closed testing procedure with WMW tests is discussed in Appendix <ref>.

To illustrate the closed testing procedure with WMW tests, we analyze the data from [38] and compare it with the other methods.
Table <ref> reports results for different methods (CT$_{\mathrm{Dunnett}}$, CT$^*_{\mathrm{WMW}}$, CT$^*_{\mathrm{Simes}}$, BH$^*$, BH with $p$-values from $t$ tests) in terms of (i) $d$: lower 90% confidence bound for the number of true discoveries (for FDR methods, $d=|D|$); (ii) $|D|$: number of FWER/FDR discoveries at $\alpha=0.1$; (iii) $\bar{p}_i^*$: FWER/FDR adjusted $p$-value. Note that closed testing with Dunnett's tests (Dunnett step-down) and BH are valid under one-way ANOVA Gaussian setup, and utilize “parametric” $p$-values from $t$ tests, i.e.

$$p_1 = 0.446, \quad p_2=0.069, \quad p_3= 0.062, \quad p_4 = 0.039, \quad p_5 < 0.001.$$

We observe that parametric assumptions yield stronger results. However, among distribution-free tests, the closed testing procedure with WMW tests is the only one that rejects the global null $H_g$ at the $\alpha=0.1$ level and concludes that there are at least three discoveries, although without pinpointing which ones.

3cClosed testing     2cFDR

Dunnett     WMW     Simes     BH$^*$     BH

$d$      1     3     0     0     4

$|D|$      1     0     0     0     4

$\bar{p}_1^*$     0.446     0.625     0.625     0.625     0.446

$\bar{p}_2^*$     0.160     0.250     0.250     0.156     0.086

$\bar{p}_3^*$     0.160     0.250     0.250     0.156     0.086

$\bar{p}_4^*$     0.129     0.250     0.250     0.156     0.086

$\bar{p}_5^*$     0.001     0.250     0.250     0.156     0.001

[38] example: results for different methods (CT$_{\mathrm{Dunnett}}$, CT$^*_{\mathrm{WMW}}$, CT$^*_{\mathrm{Simes}}$, BH$^*$, BH) in terms of (i) $d$: lower 90% confidence bound for the number of true discoveries (for FDR methods, $d=|D|$); (ii) $|D|$: number of FWER/FDR discoveries at $\alpha=0.1$; (iii) $\bar{p}_i^*$: FWER/FDR adjusted $p$-value.

§.§ The power of rank tests

We consider a model in which the test scores are sampled independently from a mixture distribution:

\begin{eqnarray}\label{eq-mixture}
Y_i \sim (1-\theta)F + \theta G_i
\end{eqnarray}

where $\theta$ is the true proportion of outliers in the test population. The number of outliers in the test sample $n_1$ (cf. the number of false hypotheses in $H_1,\ldots,H_n$) is a realization from the Binomial distribution with $n$ trials and a probability of success $\theta$.
Note that the population global null hypothesis

\begin{eqnarray}\label{eq-populationglobalnull}
H_0:\theta = 0
\end{eqnarray}

implies the sample global null hypothesis $H_g: \bigcap_{i=1}^{n}H_i$, but not the other way around: it may be that there are outliers in the test population, i.e., $H_0$ is false, although we observe no outliers in the test sample, i.e., $H_g$ is true, which happens with probability $(1-\theta)^n$. A parallel framework is the sparse Gaussian mixtures model  with $F=N(0,1)$ and $G_i = N(\tau,1)$ discussed in [10, 11], where the Higher Criticism (HC) test was introduced to determine whether there are any nonzero effects.

From the rare/weak effects viewpoint [11], outliers discovery (ODI) may be very difficult or even impossible to achieve. [43]  discussed how large $m$ should be relative to $n$ and provided a power boundary for BH$^*$, where the power criterion is the true discorey rate (TDR), i.e,. the expected true discovery proportion
$$\mathrm{Power }_{\mathrm{TDR}} = \mathbb{E}[ |D \cap I_1 | / |I_1| ].$$

This measure serves as a natural power criterion for FDR controlling procedures, representing the expected percentage of true discoveries.

We will direct our attention to the objective of outlier detection (ODE). Among the different generalizations of power [2, 46, 21], ODE can be evaluated by the probability of rejection of the global null hypothesis, i.e.

We will refer to this criterion as Power$_g$ henceforth. Maximizing Power$_{g}$ is equivalent to minimize the probability of not making any informative inference, i.e., $\mathrm{pr}(\phi_g=0)=\mathrm{pr}(d=0)=\mathrm{pr}(D=\emptyset)$.

Permutation tests have become increasingly popular and successful in practice, motivating researchers to delve into their theoretical properties, particularly concerning power, even in more complex scenarios [30].

We now revisit the following elegant finding
on a mathematically convenient alternative hypothesis.
[33] investigated the power of rank tests with $G_i = F^k$ for some integer $k>1$, which is referred to as Lehmann's alternative [33]. Here, $F^k$ represents the distribution of the maximum of $k$ independent random variables with distribution $F$. If $X \sim F$ and $Y \sim F^k$, then $\mathrm{pr}(X<Y)=k/(k+1)$. It is worth noting that the distance of $F^k$ from $F$ remains the same for any specification of $F$.

[33] proved that in the mixture model (<ref>) with Lehmann's alternative  $G_i = F^k$ with $k=2$, the WMW test is the locally most powerful invariant (LMPI) test of $H_0:\theta = 0$,  i.e. it is invariant under the group of continuous and strictly
increasing transformations of the scores.
More generally, the following Theorem provides the test statistic of the LMPI test of $H_0:\theta = 0$ for integers $k>1$.

[33]
Assume $X_1,\ldots,X_m$ i.i.d. with continuous c.d.f. $F$, and $Y_1,\ldots,Y_n$ i.i.d. with c.d.f. $(1-\theta)F + \theta F^k$ for some $\theta \in [0,1]$ and some integer $k>1$. Write $(Z_1,\ldots,Z_N)$ for the combined sample $(X_1,\ldots,X_m,Y_1,\ldots,Y_n)$. Let $Z_{(1)} < \ldots < Z_{(N)}$ denote the ordered sequence of $Z_1,\ldots,Z_N$, and  $R_i$  the rank of $Y_i$ in this sequence, thus $Y_i = Z_{(R_i)}$.
Then, the LMPI test rejects  $H_0: \theta = 0$ when the test statistic

\begin{eqnarray}
T_k=\sum_{i=1}^n R_i(R_{i}+1)\cdots(R_{i}+k-2)
\end{eqnarray}

is large.

The proof is given in Appendix <ref>, following [33], which relies on a result proved by [28].
Note that for $k=2$, the test statistic is the Wilcoxon rank sum statistic $T_{\mathrm{W}} = \sum_{i=1}^n R_i$. The Mann-Whitney $U$ statistic can be expressed as $\sum_{i=1}^n R_i - n(n+1)/2$, differing only by a constant term.
For the data from § 2.3, Figure <ref> displays the permutation distribution of the Wilcoxon rank sum test statistic $T_{\mathrm{W}}$ and $T_{7}$, the LMPI rank test for Lehmann's alternative with $k=7$, with corresponding permutation $p$-values of 0.015 and 0.006. The value $k=7$ was chosen because it is an integer close to $-\mathrm{pr}(X < Y)/(\mathrm{pr}(X < Y)-1)$ with $\mathrm{pr}(X < Y) \approx 88.7\%$ estimated by resampling from the data.

Note that the WMW test has power $\alpha$ if $H_0$ is false but the test sample contains only inliers (i.e. $H_g$ is true) because exchangeability of $(X_1,\ldots,X_m,Y_1,\ldots,Y_n)$ holds under $H_g$. However, on average over the test samples, the power of WMW is stricly larger than $\alpha$ since the WMW test is unbiased for stochastically ordered alternatives such as Lehmann's alternative [35].

An intuitive but crucial observation is the following Lemma:

The BH$^*$ global test $\phi_{\mathrm{BH}^*}=\mathds{1}\{d_{\mathrm{BH}^{*}} > 0\}$ for $H_g$ is a rank test, i.e. it is invariant under the group of continuous and strictly increasing transformations of the scores.

Refer to Appendix <ref> for a formal proof.

We are now able to state our main result in the following Corollary.

Under the assumptions of Theorem <ref> with $k=2$ and for small values of $\theta$,
the probability of not making any informative inference is larger for the BH$^*$ procedure than for CT$^*_{\mathrm{WMW}}$:

\begin{eqnarray}
\mathrm{pr}(D_{\mathrm{BH}^*} = \emptyset; \theta, F) >  \mathrm{pr}(d_{\mathrm{WMW}} = 0; \theta, F)
\end{eqnarray}

for small values of $\theta$ and any $F \in \mathcal{F}$.

The proof of the Corollary follows from Theorem <ref> and Lemma <ref>.
Similarly, under Lehmann's alternative with $k>2$, the LMPI test using the test statistic $T_k$ outperforms BH$^*$ and WMW itself. Lemma 1 and Corollary 1 also hold with BH$^*$ replaced by ABH$^*$.

The closed testing procedure with WMW tests is optimal for outlier detection, as it maximizes $\mathrm{Power}_g$ for small values of $\theta$, and it provides the richer inference $d_{\mathrm{WMW}}(S)$ beyond merely concluding that $H_g$ is false.

§ EXPERIMENTS

In this section, we study the performance of the existing FDR-controlling methods (BH$^{*}$ and ABH$^{*}$) and the proposed FDP post-hoc bounds methods (CT$^*_{\mathrm{Simes}}$, CT$^*_{\mathrm{ASimes}}$ and CT$^*_{\mathrm{WMW}}$) for the [http://odds.cs.stonybrook.edu/pendigits-dataset] dataset under Lehmann's alternative with $k=2$ (§ 4.1) and on the , [https://www.kaggle.com/mlg-ulb/creditcardfraud] and  [http://odds.cs.stonybrook.edu/forestcovercovertype-dataset] datasets with two classes (§ 4.2). These datasets are widely used for outlier detection; see  [1] for additional details.

Further analyses on other datasets ([http://odds.cs.stonybrook.edu/shuttle-dataset],  [http://odds.cs.stonybrook.edu/mammography-dataset/] and [https://www.dbs.ifi.lmu.de/research/outlier-evaluation/DAMI/literature/ALOI]) can be found in Appendix <ref>.

FDR-controlling methods and post-hoc confidence bounds cannot be directly compared due to their control over distinct error rates. FDR methods tend to be more powerful in terms of discoveries (ODI), while post-hoc confidence bounds offer greater flexibility. A common basis for comparison lies in Power$_g$ (ODE) because under the global null hypothesis $H_g$ FDR and FWER are identical.
To assess the capability of counting outliers (OEN), we will also compare the number of FDR discoveries $|D|$ with the lower confidence bound $d$ for the number of true discoveries.

§.§ Lehmann's alternative with $k=2$

For the  dataset, we fitted a one-class model on a random sample of $l=199$ training observations from the majority class (inliers).

We trained the isolation forest algorithm [37] using the $\mathrm{R}$ library isotree.

In each replication, we drew a random sample of $m=199$ control observations (inliers), and a test sample with $n=20$ observations and varying number of outliers from 0 to $n$.
The outliers score $Y$ under Lehmann's alternative with $k=2$ is obtained as $Y = \max(X',X'')$, where $X'$ and $X''$ are  two randomly drawn inliers scores. Scores on both control and test samples are obtained by using the function predict.isolation_forest function from the isotree library.
We performed $10^5$ replications at level $\alpha=0.1$ in order to estimate $\mathrm{pr}(d>0 \, | \, n_1 = x)$ and $\mathbb{E}[ d\, | \, n_1 = x]$.

The left plot in Figure <ref>
displays the probability of rejection of $H_g$ as a function of $\theta$, i.e.,

\begin{eqnarray*}
\mathrm{Power}_{g}(\theta) &=& \mathrm{pr}( d>0 \,; \theta) \\
&=& \sum_{x=0}^{n} { n \choose x} \theta^x (1-\theta)^{(n-x)} \cdot \mathrm{pr}(d>0 \, | \, n_1 = x)
\end{eqnarray*}

while the right plot in Figure <ref> displays the average lower bound $d$ as a function of $\theta$, i.e.,

\begin{eqnarray*}
\mathrm{Average\,}d(\theta)&=& \mathbb{E}_\theta[d]\\  &=&\sum_{x=0}^{n} { n \choose x} \theta^x (1-\theta)^{(n-x)} \cdot \mathbb{E}[ d\, | \, n_1 = x].
\end{eqnarray*}

The left plot in Figure <ref> demonstrates that CT$^*_{\mathrm{WMW}}$ outperforms the other FDP post-hoc bounds methods (CT$^*_{\mathrm{Simes}}$ and CT$^*_{\mathrm{ASimes}}$), as well as the FDR-controlling methods (BH$^{*}$ and ABH$^{*}$). This is as expected from Corollary 1 for small values of $\theta$. However, for the  dataset, CT$^*_{\mathrm{WMW}}$ dominates the other procedures for all values of $\theta$, not just for small values, and  the right plot in Figure <ref> shows that CT$^*_{\mathrm{WMW}}$ also exhibits a superior average lower bound $d$.

As expected, adaptive methods tend to perform better compared to non-adaptive methods when the estimate $\hat{\pi}_0$ is less than 1, which occurs as the proportion of outliers in the test set increases. However, according to the effect sparsity hypothesis, we expect few outliers in the
test set.

xmin = 0,
xmax = 1,
xtick =0,...,1,
ymin = 0,
ymax = 1,
xlabel=  Proportion $\theta$ of outliers in the population,
ylabel = Power to reject the global null,
legend style=legend pos=north west
[color=black, no marks, line width = 1.5] coordinates  ( 0,0.0964 ) ( 0.02,0.0985 ) ( 0.04,0.1006 ) ( 0.06,0.1027 ) ( 0.08,0.1048 ) ( 0.1,0.1069 ) ( 0.12,0.1089 ) ( 0.14,0.1108 ) ( 0.16,0.1127 ) ( 0.18,0.1146 ) ( 0.2,0.1165 ) ( 0.22,0.1184 ) ( 0.24,0.1202 ) ( 0.26,0.1221 ) ( 0.28,0.124 ) ( 0.3,0.1259 ) ( 0.32,0.1278 ) ( 0.34,0.1297 ) ( 0.36,0.1317 ) ( 0.38,0.1336 ) ( 0.4,0.1356 ) ( 0.42,0.1376 ) ( 0.44,0.1395 ) ( 0.46,0.1415 ) ( 0.48,0.1434 ) ( 0.5,0.1454 ) ( 0.52,0.1473 ) ( 0.54,0.1492 ) ( 0.56,0.1511 ) ( 0.58,0.153 ) ( 0.6,0.1549 ) ( 0.62,0.1567 ) ( 0.64,0.1586 ) ( 0.66,0.1605 ) ( 0.68,0.1623 ) ( 0.7,0.1642 ) ( 0.72,0.1661 ) ( 0.74,0.168 ) ( 0.76,0.1699 ) ( 0.78,0.1718 ) ( 0.8,0.1738 ) ( 0.82,0.1757 ) ( 0.84,0.1776 ) ( 0.86,0.1795 ) ( 0.88,0.1813 ) ( 0.9,0.1832 ) ( 0.92,0.185 ) ( 0.94,0.1868 ) ( 0.96,0.1886 ) ( 0.98,0.1902 ) ( 1,0.1916 )
BH$^{*}$ / Simes$^{*}$
[color=blue, no marks, line width = 1.5] coordinates  ( 0,0.0544 ) ( 0.02,0.0576 ) ( 0.04,0.061 ) ( 0.06,0.0643 ) ( 0.08,0.0677 ) ( 0.1,0.0711 ) ( 0.12,0.0746 ) ( 0.14,0.0782 ) ( 0.16,0.0818 ) ( 0.18,0.0854 ) ( 0.2,0.0891 ) ( 0.22,0.0929 ) ( 0.24,0.0968 ) ( 0.26,0.1008 ) ( 0.28,0.1048 ) ( 0.3,0.109 ) ( 0.32,0.1132 ) ( 0.34,0.1175 ) ( 0.36,0.122 ) ( 0.38,0.1265 ) ( 0.4,0.1312 ) ( 0.42,0.1359 ) ( 0.44,0.1408 ) ( 0.46,0.1458 ) ( 0.48,0.1509 ) ( 0.5,0.1561 ) ( 0.52,0.1615 ) ( 0.54,0.167 ) ( 0.56,0.1726 ) ( 0.58,0.1783 ) ( 0.6,0.1842 ) ( 0.62,0.1903 ) ( 0.64,0.1965 ) ( 0.66,0.2028 ) ( 0.68,0.2094 ) ( 0.7,0.2161 ) ( 0.72,0.223 ) ( 0.74,0.23 ) ( 0.76,0.2373 ) ( 0.78,0.2447 ) ( 0.8,0.2524 ) ( 0.82,0.2602 ) ( 0.84,0.2683 ) ( 0.86,0.2766 ) ( 0.88,0.2851 ) ( 0.9,0.2939 ) ( 0.92,0.303 ) ( 0.94,0.3125 ) ( 0.96,0.3221 ) ( 0.98,0.332 ) ( 1,0.3418 ) ;
ABH$^{*}$ / ASimes$^{*}$
[color=red, no marks, line width = 1.5] coordinates  ( 0,0.1008 ) ( 0.02,0.1098 ) ( 0.04,0.1193 ) ( 0.06,0.1295 ) ( 0.08,0.1402 ) ( 0.1,0.1516 ) ( 0.12,0.1635 ) ( 0.14,0.176 ) ( 0.16,0.1891 ) ( 0.18,0.2027 ) ( 0.2,0.2168 ) ( 0.22,0.2314 ) ( 0.24,0.2466 ) ( 0.26,0.2622 ) ( 0.28,0.2784 ) ( 0.3,0.2951 ) ( 0.32,0.3122 ) ( 0.34,0.3298 ) ( 0.36,0.3478 ) ( 0.38,0.3662 ) ( 0.4,0.385 ) ( 0.42,0.4041 ) ( 0.44,0.4235 ) ( 0.46,0.4433 ) ( 0.48,0.4632 ) ( 0.5,0.4834 ) ( 0.52,0.5037 ) ( 0.54,0.5242 ) ( 0.56,0.5447 ) ( 0.58,0.5652 ) ( 0.6,0.5857 ) ( 0.62,0.6061 ) ( 0.64,0.6264 ) ( 0.66,0.6465 ) ( 0.68,0.6663 ) ( 0.7,0.6859 ) ( 0.72,0.7051 ) ( 0.74,0.7239 ) ( 0.76,0.7423 ) ( 0.78,0.7602 ) ( 0.8,0.7776 ) ( 0.82,0.7944 ) ( 0.84,0.8107 ) ( 0.86,0.8263 ) ( 0.88,0.8412 ) ( 0.9,0.8554 ) ( 0.92,0.8689 ) ( 0.94,0.8817 ) ( 0.96,0.8937 ) ( 0.98,0.905 ) ( 1,0.9154 ) ;

xmin = 0,
xmax = 1,
xtick =0,...,1,
ymin = 0,
ymax = 5,
xlabel= Proportion $\theta$ of outliers in the population,
ylabel= Average lower bound $d$,
legend style=legend pos=north west
[no marks, color = black, line width=1.5] coordinates  ( 0,0.1303 ) ( 0.02,0.1332 ) ( 0.04,0.1365 ) ( 0.06,0.1401 ) ( 0.08,0.1437 ) ( 0.1,0.1474 ) ( 0.12,0.1511 ) ( 0.14,0.1547 ) ( 0.16,0.1583 ) ( 0.18,0.1619 ) ( 0.2,0.1654 ) ( 0.22,0.169 ) ( 0.24,0.1727 ) ( 0.26,0.1763 ) ( 0.28,0.1801 ) ( 0.3,0.1839 ) ( 0.32,0.1878 ) ( 0.34,0.1917 ) ( 0.36,0.1957 ) ( 0.38,0.1997 ) ( 0.4,0.2038 ) ( 0.42,0.2079 ) ( 0.44,0.2121 ) ( 0.46,0.2163 ) ( 0.48,0.2205 ) ( 0.5,0.2248 ) ( 0.52,0.229 ) ( 0.54,0.2333 ) ( 0.56,0.2376 ) ( 0.58,0.2419 ) ( 0.6,0.2462 ) ( 0.62,0.2506 ) ( 0.64,0.255 ) ( 0.66,0.2594 ) ( 0.68,0.2639 ) ( 0.7,0.2684 ) ( 0.72,0.273 ) ( 0.74,0.2776 ) ( 0.76,0.2823 ) ( 0.78,0.2871 ) ( 0.8,0.2919 ) ( 0.82,0.2967 ) ( 0.84,0.3016 ) ( 0.86,0.3065 ) ( 0.88,0.3114 ) ( 0.9,0.3164 ) ( 0.92,0.3214 ) ( 0.94,0.3264 ) ( 0.96,0.3315 ) ( 0.98,0.3364 ) ( 1,0.3409 ) ;

[no marks, color = blue, line width=1.5] coordinates
( 0,0.0884 ) ( 0.02,0.0935 ) ( 0.04,0.0992 ) ( 0.06,0.1055 ) ( 0.08,0.1121 ) ( 0.1,0.1191 ) ( 0.12,0.1263 ) ( 0.14,0.1338 ) ( 0.16,0.1416 ) ( 0.18,0.1498 ) ( 0.2,0.1583 ) ( 0.22,0.1673 ) ( 0.24,0.1766 ) ( 0.26,0.1864 ) ( 0.28,0.1968 ) ( 0.3,0.2077 ) ( 0.32,0.2192 ) ( 0.34,0.2313 ) ( 0.36,0.2442 ) ( 0.38,0.2579 ) ( 0.4,0.2724 ) ( 0.42,0.2878 ) ( 0.44,0.3042 ) ( 0.46,0.3216 ) ( 0.48,0.3401 ) ( 0.5,0.3599 ) ( 0.52,0.381 ) ( 0.54,0.4034 ) ( 0.56,0.4274 ) ( 0.58,0.4529 ) ( 0.6,0.4801 ) ( 0.62,0.5092 ) ( 0.64,0.5402 ) ( 0.66,0.5733 ) ( 0.68,0.6087 ) ( 0.7,0.6465 ) ( 0.72,0.6868 ) ( 0.74,0.7299 ) ( 0.76,0.7759 ) ( 0.78,0.8251 ) ( 0.8,0.8776 ) ( 0.82,0.9336 ) ( 0.84,0.9935 ) ( 0.86,1.0575 ) ( 0.88,1.1259 ) ( 0.9,1.1989 ) ( 0.92,1.277 ) ( 0.94,1.3604 ) ( 0.96,1.4495 ) ( 0.98,1.5445 ) ( 1,1.6453 )

[no marks, color = gray, line width=1.5] coordinates
( 0,0.1068 ) ( 0.02,0.1091 ) ( 0.04,0.1117 ) ( 0.06,0.1143 ) ( 0.08,0.1169 ) ( 0.1,0.1194 ) ( 0.12,0.1219 ) ( 0.14,0.1244 ) ( 0.16,0.1267 ) ( 0.18,0.1291 ) ( 0.2,0.1314 ) ( 0.22,0.1338 ) ( 0.24,0.1361 ) ( 0.26,0.1385 ) ( 0.28,0.141 ) ( 0.3,0.1434 ) ( 0.32,0.1459 ) ( 0.34,0.1484 ) ( 0.36,0.151 ) ( 0.38,0.1535 ) ( 0.4,0.1561 ) ( 0.42,0.1587 ) ( 0.44,0.1613 ) ( 0.46,0.1639 ) ( 0.48,0.1665 ) ( 0.5,0.1691 ) ( 0.52,0.1717 ) ( 0.54,0.1743 ) ( 0.56,0.1769 ) ( 0.58,0.1794 ) ( 0.6,0.182 ) ( 0.62,0.1845 ) ( 0.64,0.1871 ) ( 0.66,0.1897 ) ( 0.68,0.1923 ) ( 0.7,0.1949 ) ( 0.72,0.1976 ) ( 0.74,0.2003 ) ( 0.76,0.2031 ) ( 0.78,0.2058 ) ( 0.8,0.2086 ) ( 0.82,0.2114 ) ( 0.84,0.2141 ) ( 0.86,0.2169 ) ( 0.88,0.2196 ) ( 0.9,0.2222 ) ( 0.92,0.2248 ) ( 0.94,0.2274 ) ( 0.96,0.2298 ) ( 0.98,0.2322 ) ( 1,0.2345 )

[no marks, color = cyan, line width=1.5] coordinates
( 0,0.0626 ) ( 0.02,0.0663 ) ( 0.04,0.0702 ) ( 0.06,0.0743 ) ( 0.08,0.0784 ) ( 0.1,0.0827 ) ( 0.12,0.087 ) ( 0.14,0.0914 ) ( 0.16,0.096 ) ( 0.18,0.1006 ) ( 0.2,0.1054 ) ( 0.22,0.1103 ) ( 0.24,0.1154 ) ( 0.26,0.1206 ) ( 0.28,0.126 ) ( 0.3,0.1317 ) ( 0.32,0.1375 ) ( 0.34,0.1435 ) ( 0.36,0.1498 ) ( 0.38,0.1563 ) ( 0.4,0.1631 ) ( 0.42,0.1702 ) ( 0.44,0.1775 ) ( 0.46,0.1852 ) ( 0.48,0.1932 ) ( 0.5,0.2015 ) ( 0.52,0.2102 ) ( 0.54,0.2193 ) ( 0.56,0.2288 ) ( 0.58,0.2388 ) ( 0.6,0.2492 ) ( 0.62,0.2601 ) ( 0.64,0.2715 ) ( 0.66,0.2836 ) ( 0.68,0.2962 ) ( 0.7,0.3096 ) ( 0.72,0.3236 ) ( 0.74,0.3384 ) ( 0.76,0.3541 ) ( 0.78,0.3706 ) ( 0.8,0.3881 ) ( 0.82,0.4066 ) ( 0.84,0.4261 ) ( 0.86,0.4469 ) ( 0.88,0.4688 ) ( 0.9,0.4922 ) ( 0.92,0.5171 ) ( 0.94,0.5436 ) ( 0.96,0.5718 ) ( 0.98,0.6021 ) ( 1,0.6344 )

[no marks, color = red, line width=1.5] coordinates
( 0,0.2115 ) ( 0.02,0.2336 ) ( 0.04,0.2574 ) ( 0.06,0.2831 ) ( 0.08,0.3106 ) ( 0.1,0.3402 ) ( 0.12,0.3718 ) ( 0.14,0.4055 ) ( 0.16,0.4414 ) ( 0.18,0.4796 ) ( 0.2,0.5202 ) ( 0.22,0.5634 ) ( 0.24,0.6091 ) ( 0.26,0.6576 ) ( 0.28,0.7088 ) ( 0.3,0.763 ) ( 0.32,0.8201 ) ( 0.34,0.8803 ) ( 0.36,0.9436 ) ( 0.38,1.0101 ) ( 0.4,1.08 ) ( 0.42,1.1532 ) ( 0.44,1.2299 ) ( 0.46,1.3101 ) ( 0.48,1.3938 ) ( 0.5,1.4811 ) ( 0.52,1.572 ) ( 0.54,1.6666 ) ( 0.56,1.7648 ) ( 0.58,1.8667 ) ( 0.6,1.9722 ) ( 0.62,2.0814 ) ( 0.64,2.1943 ) ( 0.66,2.3109 ) ( 0.68,2.4312 ) ( 0.7,2.5552 ) ( 0.72,2.6829 ) ( 0.74,2.8144 ) ( 0.76,2.9496 ) ( 0.78,3.0886 ) ( 0.8,3.2313 ) ( 0.82,3.3776 ) ( 0.84,3.5274 ) ( 0.86,3.6806 ) ( 0.88,3.8369 ) ( 0.9,3.9961 ) ( 0.92,4.1579 ) ( 0.94,4.3219 ) ( 0.96,4.4878 ) ( 0.98,4.6553 ) ( 1,4.8243 )

Power to reject the global null of no outliers (left plot) and average lower bound $d$ on the number of outliers (right plot) for different methods (BH$^*$, ABH$^*$, CT$^*_{\mathrm{Simes}}$, CT$^*_{\mathrm{ASimes}}$, CT$^*_{\mathrm{WMW}}$), as a function of the proportion $\theta$ of outliers in the population. Samples are drawn from the  dataset with $l=199$, $m=199$, $n=20$ and under Lehmann's alternative with $k=2$. The significance level is set at $\alpha=0.1$ and all the quantities are estimated with $10^5$ replications.

§.§ Real datasets with two classes

For each dataset, we fitted the same one-class model as for Lehmann's alternative with $k=2$ (see § 4.1) using a random sample of $l=199$ training observations from the majority class (inliers).
In each replication, we drew a random sample of $m=199$ control observations (inliers), and a random sample with $n=20$ test observations and varying number of outliers from 0 to $n$. In  and  datasets outliers are randomly drawn from the minority class in the dataset; in  dataset outliers are randomly picked from the
class of fraudulent credit card transactions.
We performed $10^4$ replications at level $\alpha=0.1$ in order to estimate $\mathrm{pr}(\ d>0 \, | \, n_1 = x)$ and $\mathbb{E}[d\, | \, n_1 = x]$.

Figure <ref> compares the performance of the different methods on the ,  and   datasets.
We see that at $n_1=0$, CT$^*_{\mathrm{Simes}}$ is of exact size $\alpha$ since we chose $\alpha, m, n$ such that the condition in Proposition <ref> is satisfied, while the adaptive methods are conservative because $\hat{\pi}_0$ tend to be larger than 1.
CT$^*_{\mathrm{WMW}}$ performs well for the  and  datasets, but for the  dataset, it is outperformed in terms of discoveries by the FDR-controlling methods.
Please refer to Appendix <ref> for additional results on other datasets.

title = Pendigits,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 1,
xlabel= Number of outlier $n_1$,
ylabel= Power to reject the global null,
legend style=legend pos=south east
[mark size=1.5pt, color=black, mark=o, line width=1.5] coordinates  ( 0,0.0945 ) ( 1,0.1881 ) ( 2,0.267 ) ( 3,0.3382 ) ( 4,0.4133 ) ( 5,0.4669 ) ( 6,0.5201 ) ( 7,0.5743 ) ( 8,0.6223 ) ( 9,0.6691 ) ( 10,0.7083 ) ( 11,0.7399 ) ( 12,0.7779 ) ( 13,0.7979 ) ( 14,0.8215 ) ( 15,0.8524 ) ( 16,0.8718 ) ( 17,0.8895 ) ( 18,0.9073 ) ( 19,0.9186 ) ( 20,0.9273 ) ;
BH$^{*}$ / Simes$^*$

[mark size=1.5pt, color=blue, mark=o, line width=1.5] coordinates  ( 0,0.0548 ) ( 1,0.1118 ) ( 2,0.183 ) ( 3,0.26 ) ( 4,0.3565 ) ( 5,0.4336 ) ( 6,0.518 ) ( 7,0.6042 ) ( 8,0.6806 ) ( 9,0.7544 ) ( 10,0.8085 ) ( 11,0.8607 ) ( 12,0.9074 ) ( 13,0.9411 ) ( 14,0.9628 ) ( 15,0.9817 ) ( 16,0.9938 ) ( 17,0.997 ) ( 18,0.9994 ) ( 19,0.9999 ) ( 20,1 )
ABH$^{*}$ / ASimes$^*$

[mark size=1.5pt, color=red, mark=o, line width=1.5] coordinates
( 0,0.1051 ) ( 1,0.1548 ) ( 2,0.2393 ) ( 3,0.3213 ) ( 4,0.4494 ) ( 5,0.5733 ) ( 6,0.7 ) ( 7,0.8112 ) ( 8,0.8972 ) ( 9,0.949 ) ( 10,0.9799 ) ( 11,0.9956 ) ( 12,0.9991 ) ( 13,0.9998 ) ( 14,1 ) ( 15,1 ) ( 16,1 ) ( 17,1 ) ( 18,1 ) ( 19,1 ) ( 20,1 )

title = Pendigits,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 20,
xlabel= Number of outlier $n_1$,
ylabel= Average lower bound $d$,
legend style=legend pos=north west
[mark size=1.5pt, color=black, mark=o, line width=1.5] coordinates  ( 0,0.1254 ) ( 1,0.2545 ) ( 2,0.3874 ) ( 3,0.5285 ) ( 4,0.6865 ) ( 5,0.8585 ) ( 6,1.0653 ) ( 7,1.2746 ) ( 8,1.5504 ) ( 9,1.8283 ) ( 10,2.0863 ) ( 11,2.4376 ) ( 12,2.8785 ) ( 13,3.1556 ) ( 14,3.6267 ) ( 15,4.1499 ) ( 16,4.6183 ) ( 17,5.1649 ) ( 18,5.7136 ) ( 19,6.324 ) ( 20,6.8801 )

[mark size=1.5pt, color=blue, mark=o, line width=1.5] coordinates  ( 0,0.0877 ) ( 1,0.1837 ) ( 2,0.3198 ) ( 3,0.5033 ) ( 4,0.758 ) ( 5,1.0453 ) ( 6,1.4532 ) ( 7,1.983 ) ( 8,2.6775 ) ( 9,3.5194 ) ( 10,4.4477 ) ( 11,5.7764 ) ( 12,7.304 ) ( 13,9.0032 ) ( 14,10.7471 ) ( 15,12.8319 ) ( 16,14.6059 ) ( 17,16.3966 ) ( 18,17.7686 ) ( 19,18.8183 ) ( 20,19.5768 )

[mark size=1.5pt, color=gray, mark=o, line width=1.5] coordinates
( 0,0.1025 ) ( 1,0.2096 ) ( 2,0.3103 ) ( 3,0.4116 ) ( 4,0.5195 ) ( 5,0.6211 ) ( 6,0.7264 ) ( 7,0.8516 ) ( 8,0.9808 ) ( 9,1.1078 ) ( 10,1.2347 ) ( 11,1.3691 ) ( 12,1.5507 ) ( 13,1.6719 ) ( 14,1.8393 ) ( 15,2.0577 ) ( 16,2.2464 ) ( 17,2.4734 ) ( 18,2.7144 ) ( 19,2.9991 ) ( 20,3.2433 )

[mark size=1.5pt, color=cyan, mark=o, line width=1.5] coordinates
( 0,0.0614 ) ( 1,0.1287 ) ( 2,0.2209 ) ( 3,0.331 ) ( 4,0.4721 ) ( 5,0.616 ) ( 6,0.7977 ) ( 7,1.0167 ) ( 8,1.2855 ) ( 9,1.5929 ) ( 10,1.9239 ) ( 11,2.376 ) ( 12,2.9622 ) ( 13,3.5645 ) ( 14,4.3337 ) ( 15,5.3542 ) ( 16,6.4093 ) ( 17,7.8078 ) ( 18,9.381 ) ( 19,11.2371 ) ( 20,13.4056 )

[mark size=1.5pt, color=red, mark=o, line width=1.5] coordinates
( 0,0.2171 ) ( 1,0.3425 ) ( 2,0.5584 ) ( 3,0.8331 ) ( 4,1.2222 ) ( 5,1.7174 ) ( 6,2.3667 ) ( 7,3.1176 ) ( 8,3.9548 ) ( 9,4.8416 ) ( 10,5.78 ) ( 11,6.8059 ) ( 12,7.8135 ) ( 13,8.821 ) ( 14,9.8341 ) ( 15,10.9518 ) ( 16,11.9925 ) ( 17,13.1507 ) ( 18,14.312 ) ( 19,15.5413 ) ( 20,16.861 )

title = Creditcard,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 1,
xlabel= Number of outlier $n_1$,
ylabel= Power to reject the global null,
legend style=legend pos=south east
[mark size=1.5pt, color=black, mark=o, line width=1.5] coordinates
( 0,0.1001 ) ( 1,0.395 ) ( 2,0.5918 ) ( 3,0.7378 ) ( 4,0.8407 ) ( 5,0.9065 ) ( 6,0.939 ) ( 7,0.962 ) ( 8,0.9787 ) ( 9,0.9894 ) ( 10,0.9932 ) ( 11,0.9967 ) ( 12,0.9982 ) ( 13,0.9991 ) ( 14,0.9993 ) ( 15,0.9998 ) ( 16,1 ) ( 17,0.9998 ) ( 18,1 ) ( 19,1 ) ( 20,0.9999 )
BH$^{*}$ / Simes$^*$

[mark size=1.5pt, color=blue, mark=o, line width=1.5] coordinates
( 0,0.0585 ) ( 1,0.2244 ) ( 2,0.4251 ) ( 3,0.6155 ) ( 4,0.7623 ) ( 5,0.8681 ) ( 6,0.9282 ) ( 7,0.9629 ) ( 8,0.9838 ) ( 9,0.9926 ) ( 10,0.996 ) ( 11,0.9988 ) ( 12,0.9996 ) ( 13,0.9999 ) ( 14,0.9997 ) ( 15,1 ) ( 16,1 ) ( 17,1 ) ( 18,1 ) ( 19,1 ) ( 20,1 )
ABH$^{*}$ / ASimes$^*$

[mark size=1.5pt, color=red, mark=o, line width=1.5] coordinates
( 0,0.1006 ) ( 1,0.1609 ) ( 2,0.2576 ) ( 3,0.367 ) ( 4,0.5056 ) ( 5,0.6447 ) ( 6,0.7698 ) ( 7,0.8679 ) ( 8,0.9399 ) ( 9,0.9726 ) ( 10,0.9906 ) ( 11,0.9969 ) ( 12,0.9993 ) ( 13,1 ) ( 14,1 ) ( 15,1 ) ( 16,1 ) ( 17,1 ) ( 18,1 ) ( 19,1 ) ( 20,1 )

title = Creditcard,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 20,
xlabel= Number of outlier $n_1$,
ylabel= Average lower bound $d$,
legend style=legend pos=north west
[mark size=1.5pt, color=black, mark=o, line width=1.5] coordinates
( 0,0.1323 ) ( 1,0.5376 ) ( 2,1.0587 ) ( 3,1.7281 ) ( 4,2.4999 ) ( 5,3.3952 ) ( 6,4.244 ) ( 7,5.1521 ) ( 8,6.1098 ) ( 9,7.0877 ) ( 10,8.0393 ) ( 11,9.0346 ) ( 12,9.9631 ) ( 13,10.8765 ) ( 14,11.8574 ) ( 15,12.7151 ) ( 16,13.6349 ) ( 17,14.531 ) ( 18,15.3832 ) ( 19,16.2455 ) ( 20,17.0928 )

[mark size=1.5pt, color=blue, mark=o, line width=1.5] coordinates
( 0,0.0928 ) ( 1,0.3624 ) ( 2,0.8889 ) ( 3,1.628 ) ( 4,2.5456 ) ( 5,3.6375 ) ( 6,4.7018 ) ( 7,5.8329 ) ( 8,7.0393 ) ( 9,8.1998 ) ( 10,9.3456 ) ( 11,10.5078 ) ( 12,11.6105 ) ( 13,12.7118 ) ( 14,13.8324 ) ( 15,14.8159 ) ( 16,15.8737 ) ( 17,16.8025 ) ( 18,17.6761 ) ( 19,18.4247 ) ( 20,19.049 )

[mark size=1.5pt, color=gray, mark=o, line width=1.5] coordinates
( 0,0.1093 ) ( 1,0.4463 ) ( 2,0.8037 ) ( 3,1.2137 ) ( 4,1.6524 ) ( 5,2.1251 ) ( 6,2.5908 ) ( 7,3.0706 ) ( 8,3.6353 ) ( 9,4.2569 ) ( 10,4.9042 ) ( 11,5.6325 ) ( 12,6.3644 ) ( 13,7.1723 ) ( 14,8.0829 ) ( 15,8.9012 ) ( 16,9.8759 ) ( 17,10.894 ) ( 18,11.9571 ) ( 19,13.1179 ) ( 20,14.4142 )

[mark size=1.5pt, color=cyan, mark=o, line width=1.5] coordinates
( 0,0.0665 ) ( 1,0.2633 ) ( 2,0.589 ) ( 3,1.0099 ) ( 4,1.504 ) ( 5,2.0972 ) ( 6,2.6767 ) ( 7,3.302 ) ( 8,4.015 ) ( 9,4.7564 ) ( 10,5.5414 ) ( 11,6.3929 ) ( 12,7.2709 ) ( 13,8.1863 ) ( 14,9.2096 ) ( 15,10.1421 ) ( 16,11.2234 ) ( 17,12.2953 ) ( 18,13.4191 ) ( 19,14.6094 ) ( 20,15.8884 )

[mark size=1.5pt, color=red, mark=o, line width=1.5] coordinates
( 0,0.2151 ) ( 1,0.347 ) ( 2,0.607 ) ( 3,0.9226 ) ( 4,1.4203 ) ( 5,2.0283 ) ( 6,2.738 ) ( 7,3.5794 ) ( 8,4.5264 ) ( 9,5.4556 ) ( 10,6.4529 ) ( 11,7.4432 ) ( 12,8.4716 ) ( 13,9.4676 ) ( 14,10.5334 ) ( 15,11.5067 ) ( 16,12.591 ) ( 17,13.6323 ) ( 18,14.7299 ) ( 19,15.8267 ) ( 20,16.9574 )

title = Covertype,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 1,
xlabel= Number of outlier $n_1$,
ylabel= Power to reject the global null,
legend style=legend pos=south east
[mark size=1.5pt, color=black, mark=o, line width=1.5] coordinates
( 0,0.1007 ) ( 1,0.1308 ) ( 2,0.1674 ) ( 3,0.2034 ) ( 4,0.2385 ) ( 5,0.2795 ) ( 6,0.313 ) ( 7,0.3443 ) ( 8,0.393 ) ( 9,0.4197 ) ( 10,0.462 ) ( 11,0.4976 ) ( 12,0.5232 ) ( 13,0.5741 ) ( 14,0.6021 ) ( 15,0.6308 ) ( 16,0.6634 ) ( 17,0.6983 ) ( 18,0.7206 ) ( 19,0.7556 ) ( 20,0.7786 )
BH$^{*}$ / Simes$^*$

[mark size=1.5pt, color=blue, mark=o, line width=1.5] coordinates
( 0,0.057 ) ( 1,0.0862 ) ( 2,0.1238 ) ( 3,0.1716 ) ( 4,0.2325 ) ( 5,0.29 ) ( 6,0.3633 ) ( 7,0.4346 ) ( 8,0.5207 ) ( 9,0.6077 ) ( 10,0.6962 ) ( 11,0.7732 ) ( 12,0.8515 ) ( 13,0.9113 ) ( 14,0.9558 ) ( 15,0.9777 ) ( 16,0.9935 ) ( 17,0.9985 ) ( 18,1 ) ( 19,1 ) ( 20,1 )
ABH$^{*}$ / ASimes$^*$

[mark size=1.5pt, color=red, mark=o, line width=1.5] coordinates
( 0,0.0971 ) ( 1,0.1601 ) ( 2,0.2473 ) ( 3,0.3554 ) ( 4,0.4758 ) ( 5,0.609 ) ( 6,0.7286 ) ( 7,0.8476 ) ( 8,0.9253 ) ( 9,0.9726 ) ( 10,0.9935 ) ( 11,0.9986 ) ( 12,0.9998 ) ( 13,1 ) ( 14,1 ) ( 15,1 ) ( 16,1 ) ( 17,1 ) ( 18,1 ) ( 19,1 ) ( 20,1 )

title = Covertype,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 20,
xlabel= Number of outlier $n_1$,
ylabel= Average lower bound $d$,
legend style=legend pos=north west
[mark size=1.5pt, color=black, mark=o, line width=1.5] coordinates
( 0,0.1363 ) ( 1,0.1866 ) ( 2,0.2742 ) ( 3,0.358 ) ( 4,0.4768 ) ( 5,0.6153 ) ( 6,0.7756 ) ( 7,0.9236 ) ( 8,1.202 ) ( 9,1.3927 ) ( 10,1.7234 ) ( 11,2.0704 ) ( 12,2.3886 ) ( 13,2.8652 ) ( 14,3.2828 ) ( 15,3.7763 ) ( 16,4.3163 ) ( 17,5.023 ) ( 18,5.4988 ) ( 19,6.3504 ) ( 20,7.1003 )

[mark size=1.5pt, color=blue, mark=o, line width=1.5] coordinates
( 0,0.0921 ) ( 1,0.1491 ) ( 2,0.2482 ) ( 3,0.379 ) ( 4,0.6034 ) ( 5,0.8722 ) ( 6,1.2854 ) ( 7,1.8232 ) ( 8,2.5769 ) ( 9,3.5685 ) ( 10,4.8262 ) ( 11,6.4304 ) ( 12,8.2893 ) ( 13,10.2533 ) ( 14,12.4638 ) ( 15,14.4523 ) ( 16,16.2905 ) ( 17,17.7822 ) ( 18,18.786 ) ( 19,19.4357 ) ( 20,19.9437 )

[mark size=1.5pt, color=gray, mark=o, line width=1.5] coordinates
( 0,0.112 ) ( 1,0.1481 ) ( 2,0.1975 ) ( 3,0.2451 ) ( 4,0.3047 ) ( 5,0.3753 ) ( 6,0.4333 ) ( 7,0.4994 ) ( 8,0.6031 ) ( 9,0.6764 ) ( 10,0.8037 ) ( 11,0.9208 ) ( 12,1.0274 ) ( 13,1.226 ) ( 14,1.3454 ) ( 15,1.5455 ) ( 16,1.7529 ) ( 17,2.0489 ) ( 18,2.2162 ) ( 19,2.648 ) ( 20,3.011 )

[mark size=1.5pt, color=cyan, mark=o, line width=1.5] coordinates
( 0,0.0644 ) ( 1,0.1006 ) ( 2,0.1539 ) ( 3,0.2178 ) ( 4,0.3168 ) ( 5,0.4295 ) ( 6,0.5762 ) ( 7,0.7624 ) ( 8,1.0387 ) ( 9,1.3576 ) ( 10,1.7962 ) ( 11,2.3411 ) ( 12,3.0289 ) ( 13,3.8674 ) ( 14,4.8982 ) ( 15,6.0838 ) ( 16,7.5122 ) ( 17,9.2163 ) ( 18,11.102 ) ( 19,13.2934 ) ( 20,15.6342 )

[mark size=1.5pt, color=red, mark=o, line width=1.5] coordinates
( 0,0.2049 ) ( 1,0.3574 ) ( 2,0.5976 ) ( 3,0.9051 ) ( 4,1.3431 ) ( 5,1.9079 ) ( 6,2.5663 ) ( 7,3.4292 ) ( 8,4.3501 ) ( 9,5.3292 ) ( 10,6.3995 ) ( 11,7.4419 ) ( 12,8.4958 ) ( 13,9.597 ) ( 14,10.7437 ) ( 15,11.8513 ) ( 16,12.9994 ) ( 17,14.2438 ) ( 18,15.4926 ) ( 19,16.8466 ) ( 20,18.2577 )

Power to reject the global null (left column) and average lower bound $d$ (right column) for , (first row),  (second row) and  (third row) datasets
for different methods, as a function of the number $n_1$ of outliers in the test sample.
Samples are drawn from each dataset with training size of $l=199$, control size of $m=199$ and test size of $n=20$, and varying the number of outliers from $0$ to $n$. The level is set at $\alpha=0.1$ and all the quantities are estimated with $10^4$ replications.

Figure <ref> display the (expected) lower confidence bound $d$ for the number of true discoveries and the (expected) number of FWER discoveries $|D|$ for the closed testing methods in  the  and  datasets.
For bounding the number of true discoveries (cfr. enumerating outliers), CT$^*_{\mathrm{WMW}}$ performs better than the other two, while for making discoveries (cfr. localizing outliers), CT$^*_{\mathrm{Simes}}$ and CT$^*_{\mathrm{ASimes}}$ are both much better than CT$^*_{\mathrm{WMW}}$. This is the typical behavior of post-hoc bounds based on sum-based tests such as WMW and Fisher [53].

title = Pendigits,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 18,
xlabel= Number of outlier $n_1$,
ylabel= Average $d$ and average $|D|$,
legend style=legend pos=north west

[no marks, color=gray, line width = 1.5] coordinates
( 0,0.1025 ) ( 1,0.2096 ) ( 2,0.3103 ) ( 3,0.4116 ) ( 4,0.5195 ) ( 5,0.6211 ) ( 6,0.7264 ) ( 7,0.8516 ) ( 8,0.9808 ) ( 9,1.1078 ) ( 10,1.2347 ) ( 11,1.3691 ) ( 12,1.5507 ) ( 13,1.6719 ) ( 14,1.8393 ) ( 15,2.0577 ) ( 16,2.2464 ) ( 17,2.4734 ) ( 18,2.7144 ) ( 19,2.9991 ) ( 20,3.2433 )
$d$ Simes$^*$

[no marks, color=cyan, line width = 1.5] coordinates
( 0,0.0614 ) ( 1,0.1287 ) ( 2,0.2209 ) ( 3,0.331 ) ( 4,0.4721 ) ( 5,0.616 ) ( 6,0.7977 ) ( 7,1.0167 ) ( 8,1.2855 ) ( 9,1.5929 ) ( 10,1.9239 ) ( 11,2.376 ) ( 12,2.9622 ) ( 13,3.5645 ) ( 14,4.3337 ) ( 15,5.3542 ) ( 16,6.4093 ) ( 17,7.8078 ) ( 18,9.381 ) ( 19,11.2371 ) ( 20,13.4056 )
$d$ ASimes$^*$

[no marks, color=red, line width = 1.5] coordinates
( 0,0.2171 ) ( 1,0.3425 ) ( 2,0.5584 ) ( 3,0.8331 ) ( 4,1.2222 ) ( 5,1.7174 ) ( 6,2.3667 ) ( 7,3.1176 ) ( 8,3.9548 ) ( 9,4.8416 ) ( 10,5.78 ) ( 11,6.8059 ) ( 12,7.8135 ) ( 13,8.821 ) ( 14,9.8341 ) ( 15,10.9518 ) ( 16,11.9925 ) ( 17,13.1507 ) ( 18,14.312 ) ( 19,15.5413 ) ( 20,16.861 )
$d$ WMW

[no marks, color=gray, line width = 1.5, dashed] coordinates  ( 0,0 ) ( 1,0.1925 ) ( 2,0.2846 ) ( 3,0.3774 ) ( 4,0.4704 ) ( 5,0.5571 ) ( 6,0.6395 ) ( 7,0.75 ) ( 8,0.8429 ) ( 9,0.9376 ) ( 10,1.0387 ) ( 11,1.1211 ) ( 12,1.2378 ) ( 13,1.2996 ) ( 14,1.3815 ) ( 15,1.4852 ) ( 16,1.5812 ) ( 17,1.6971 ) ( 18,1.7922 ) ( 19,1.9019 ) ( 20,1.9703 )
$|D|$ Simes$^{*}$

[no marks, color=cyan, line width = 1.5, dashed] coordinates  ( 0,0 ) ( 1,0.1008 ) ( 2,0.1751 ) ( 3,0.2631 ) ( 4,0.3667 ) ( 5,0.4773 ) ( 6,0.5969 ) ( 7,0.7436 ) ( 8,0.889 ) ( 9,1.0362 ) ( 10,1.1901 ) ( 11,1.3482 ) ( 12,1.5963 ) ( 13,1.7791 ) ( 14,2.0276 ) ( 15,2.4335 ) ( 16,2.8627 ) ( 17,3.4768 ) ( 18,4.3056 ) ( 19,5.4591 ) ( 20,7.0978 )
$|D|$ ASimes$^{*}$

[no marks, color=red, line width = 1.5, dashed] coordinates
( 0,0 ) ( 1,0 ) ( 2,0 ) ( 3,0 ) ( 4,0 ) ( 5,0 ) ( 6,0 ) ( 7,0 ) ( 8,8e-04 ) ( 9,2e-04 ) ( 10,0.0021 ) ( 11,0.0016 ) ( 12,0.009 ) ( 13,0.0132 ) ( 14,0.0348 ) ( 15,0.0802 ) ( 16,0.1624 ) ( 17,0.3975 ) ( 18,0.9642 ) ( 19,2.1279 ) ( 20,5.1422 )
$|D|$ WMW

title = Creditcard,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 18,
xlabel= Number of outlier $n_1$,
ylabel= Average $d$ and average $|D|$,
legend style=legend pos=north west

[no marks, color=gray, line width = 1.5] coordinates
( 0,0.1093 ) ( 1,0.4463 ) ( 2,0.8037 ) ( 3,1.2137 ) ( 4,1.6524 ) ( 5,2.1251 ) ( 6,2.5908 ) ( 7,3.0706 ) ( 8,3.6353 ) ( 9,4.2569 ) ( 10,4.9042 ) ( 11,5.6325 ) ( 12,6.3644 ) ( 13,7.1723 ) ( 14,8.0829 ) ( 15,8.9012 ) ( 16,9.8759 ) ( 17,10.894 ) ( 18,11.9571 ) ( 19,13.1179 ) ( 20,14.4142 )
$d$ Simes$^*$

[no marks, color=cyan, line width = 1.5] coordinates
( 0,0.0665 ) ( 1,0.2633 ) ( 2,0.589 ) ( 3,1.0099 ) ( 4,1.504 ) ( 5,2.0972 ) ( 6,2.6767 ) ( 7,3.302 ) ( 8,4.015 ) ( 9,4.7564 ) ( 10,5.5414 ) ( 11,6.3929 ) ( 12,7.2709 ) ( 13,8.1863 ) ( 14,9.2096 ) ( 15,10.1421 ) ( 16,11.2234 ) ( 17,12.2953 ) ( 18,13.4191 ) ( 19,14.6094 ) ( 20,15.8884 )
$d$ ASimes$^*$

[no marks, color=red, line width = 1.5] coordinates
( 0,0.2151 ) ( 1,0.347 ) ( 2,0.607 ) ( 3,0.9226 ) ( 4,1.4203 ) ( 5,2.0283 ) ( 6,2.738 ) ( 7,3.5794 ) ( 8,4.5264 ) ( 9,5.4556 ) ( 10,6.4529 ) ( 11,7.4432 ) ( 12,8.4716 ) ( 13,9.4676 ) ( 14,10.5334 ) ( 15,11.5067 ) ( 16,12.591 ) ( 17,13.6323 ) ( 18,14.7299 ) ( 19,15.8267 ) ( 20,16.9574 )
$d$ WMW

[no marks, color=gray, line width = 1.5, dashed] coordinates
( 0,0 ) ( 1,0.4143 ) ( 2,0.7185 ) ( 3,1.0473 ) ( 4,1.3723 ) ( 5,1.7057 ) ( 6,2.0193 ) ( 7,2.3179 ) ( 8,2.6173 ) ( 9,2.9578 ) ( 10,3.2742 ) ( 11,3.6592 ) ( 12,4.0554 ) ( 13,4.5656 ) ( 14,5.268 ) ( 15,5.9531 ) ( 16,6.7779 ) ( 17,7.7115 ) ( 18,8.8631 ) ( 19,10.2387 ) ( 20,11.9294 )
$|D|$ Simes$^{*}$

[no marks, color=cyan, line width = 1.5, dashed] coordinates
( 0,0 ) ( 1,0.2079 ) ( 2,0.4344 ) ( 3,0.7077 ) ( 4,1.0339 ) ( 5,1.4327 ) ( 6,1.8518 ) ( 7,2.2699 ) ( 8,2.7223 ) ( 9,3.2197 ) ( 10,3.7298 ) ( 11,4.3004 ) ( 12,4.906 ) ( 13,5.6108 ) ( 14,6.5082 ) ( 15,7.3408 ) ( 16,8.4123 ) ( 17,9.4907 ) ( 18,10.8208 ) ( 19,12.3968 ) ( 20,14.0441 )
$|D|$ ASimes$^{*}$

[no marks, color=red, line width = 1.5, dashed] coordinates
( 0,0 ) ( 1,0 ) ( 2,0 ) ( 3,0 ) ( 4,0 ) ( 5,0 ) ( 6,0 ) ( 7,6e-04 ) ( 8,6e-04 ) ( 9,0 ) ( 10,0.0012 ) ( 11,0.004 ) ( 12,0.0117 ) ( 13,0.0189 ) ( 14,0.063 ) ( 15,0.1292 ) ( 16,0.2747 ) ( 17,0.5681 ) ( 18,1.2943 ) ( 19,3.0587 ) ( 20,6.8191 )
$|D|$ WMW

Average lower confidence bound $d$ for the number of true discoveries and average number of FWER discoveries $|D|$ for the  dataset (left plot) and the  dataset (right plot) for the three closed testing methods.

Dataset     $q$     # inliers     # outliers     $\theta$     $\hat k$     2|c|Power$_g(\theta) \, | \, n_1 > 0 $      2|c|$\mathrm{Average} \hspace{1mm}d(\theta) \, | \, n_1 > 0 $

2*BH$^*$ / Simes$^*$     2*$\mathbf{0.3982}$     BH$^*$     $\mathbf{0.5462}$

Simes$^*$     $0.4522$

2*ABH$^*$ / ASimes$^*$     2*$0.2277$
ABH$^*$     $0.3711 $

ASimes$^*$     $0.2687$

$\mathrm{WMW}$     $0.1625$     $\mathrm{WMW}$     $0.3513$

2*BH$^*$ / Simes$^*$     2*$0.1342$     BH$^*$     $0.1948$

Simes$^*$     $0.1528$

2*ABH$^*$ / ASimes$^*$     2*$0.0898$
ABH$^*$     $0.1586$

ASimes$^*$     $0.1057$

$\mathrm{WMW}$     $\mathbf{0.1684}$     $\mathrm{WMW}$     $\mathbf{0.3804}$

5*16     5*6714     5*156
2*BH$^*$ / Simes$^*$     2*$\mathbf{0.2063}$     BH$^*$     $0.2858$

Simes$^*$     $0.2331$

2*ABH$^*$ / ASimes$^*$     2*$0.1286$
ABH$^*$     $0.2170$

ASimes$^*$     $0.1508$

$\mathrm{WMW}$     $0.1745$     $\mathrm{WMW}$      $\mathbf{0.3949}$

Summary of results on real datasets with dimension $q$, # of inliers and # of outliers. Power$_g$ and $\mathrm{Average\,} d$ are computed conditional to $n_1 > 0$, with level $\alpha=0.1$. The proportion $\theta$ is the proportion of outliers in the dataset.  $\mathrm{pr}(Y > X)$ is estimated from the dataset, providing $\hat{k}$. For each dataset, the best result is emphasized in bold.

Finally, Table <ref> summarizes the results when sampling from the population represented by the whole dataset. $\theta$ represents the actual proportion of outliers in the dataset. Then, by using the estimated $\mathrm{pr}( d>0 \, | \, n_1 = x)$ and $\mathbb{E}[ d\, | \, n_1 = x]$, we computed Power$_g(\theta)$ and $\mathrm{Average \hspace{1mm} d}(\theta)$ conditional on $n_1>0$, i.e. that the test sample contains at least one outlier.  The quantity $\hat k$ is obtained by estimating $\mathrm{pr}(X<Y)$
from the dataset and solving for $k$ in $\widehat{\mathrm{pr}}(X<Y)=k/(k+1)$ (see § 3.2).

Across different datasets, the performance of the methods varies. The BH$^*$ procedure demonstrates superiority on one dataset (), while CT$^*_\mathrm{WMW}$ proves to be better on another (). Interestingly, on the  dataset, both methods exhibit comparable performance: BH$^*$ demonstrates the highest Power$_g$, while CT$^*_\mathrm{WMW}$ achieves the highest Average $d$.

The heterogeneity of results confirms Lehmann's observation in [35]: “it seems that where nonparametric methods are appropriate, one usually does not have very precise knowledge of the alternatives.”
Adaptive methods exhibit inferior performance due to the low proportion of outliers.

The simulation study, while not exhaustive, provides valuable insights into the power properties of different

§ DISCUSSION

In this work, we have introduced post-hoc confidence bounds for novelty detection. These bounds involve a trade-off, offering reduced power for outlier discovery compared to FDR methods, but providing heightened protection against the variability of FDP around FDR. Additionally, they enable simultaneous statements about subsets of the test set. An implementation of this approach is available in the  R package https://github.com/chiaragaiamagnani/nout.

We revisited an elegant result from  [33], demonstrating that under a nonparametric alternative, the Wilcoxon-Mann-Whitney test is optimal for outlier detection among rank tests.
We made a simple observation: the Benjamini and Hochberg global test based on permutation $p$-values, being a rank test, is less powerful than the Wilcoxon-Mann-Whitney test.

However, the procedure of Benjamini and Hochberg using permutation $p$-values does not fully leverage the benefits of the permutation approach.
A FDR-controlling method based on permutations that possesses the power, simplicity, and reliability of permutation methods for FWER and FDP control remains to be found [20].

Our work inspires several avenues for future research. One such direction is to develop scalable algorithms for calculating post-hoc confidence bounds using closed testing with LMPI tests under Lehmann's alternatives with $k>2$. In the Appendix, we demonstrate that these statistics are $U$-statistics; as a result, their distribution can be approximated using the Central Limit Theorem for $U$-statistics. Additionally, further exploration could involve a comparison with the permutation approach presented by [55].

When more than one hypothesis is tested, one can replace Type I error control with FWER, FDP or FDR control. Among the many generalizations of power, power for outlier detection, i.e. the probability of rejecting the global null hypothesis of no outliers, is the simplest generalization because it reduces to power for a single test.
However, this simplicity is of great importance, as it addresses a crucial question: whether outliers exist or not.
We leave other power criteria to be considered in future research.

An intriguing question arises about the implications of employing an optimal local test within a closed testing procedure, as it may not directly lead to an optimal closed testing procedure. Nevertheless, it's worth noting that closed testing procedures hold the potential for optimality. As highlighted by [17], admissible methods are required to be closed testing procedures, and optimal methods must adhere to admissibility.

Erich Lehmann significantly influenced the development of rank tests. For an elegant summary of his main contributions, refer to [54].
According to [34], “In simple problems in which Wilcoxon type tests apply,
they may be recommended as the methods of choice”.

§ SUPPLEMENTARY MATERIAL

§ CLOSED TESTING WITH $\PI_0$-ADAPTIVE SIMES TESTS

The same idea as the adaptive BH procedure applies to Simes test, resulting in a $\pi_0$-adaptive version of Simes' test for $H_I$ that utilizes the proportion of true null hypotheses in $I$:

\begin{equation}\label{eq-ASimeslocaltest}
\phi_I^\mathrm{ASimes} = \mathds{1}\left\{\min_{k\in \{1,\ldots,|I|\}} \{|I|\hat\pi_0^Ip_{(k:I)}/k\}\leq \alpha\right\},
\end{equation}

\begin{eqnarray}\label{eq-pi0HatI}
\hat{\pi}^{I}_0 &=& \frac{1+\sum_{i\in I} \mathds{1}\{p^*_i > \lambda\} }{|I|(1-\lambda)}
\end{eqnarray}

is the estimator in (<ref>) computed considering only permutation $p$-values in subset $I$.
The test $\phi_I^\mathrm{ASimes}$ was considered in [6] and [22].
We will refer to the closed testing procedure using adaptive Simes' combination of permutation $p$-values as CT$^*_\mathrm{ASimes}$.

To reduce computational time we will derive a shortcut based on Lemma 6 in [17], which uses a conservative version of $\hat\pi_0^I$ given by:

\begin{eqnarray*}\label{eq-pi0tildeI}
\tilde\pi_0^{|I|} = \frac{1+\sum_{i=n-|I|+1}^n\mathds{1}\{p^*_{(i)}\geq\lambda\}}{|I|(1-\lambda)}.
\end{eqnarray*}

The idea is to consider the $|I|$ largest permutation $p$-values, making the estimator depend only on the cardinality of $I$.
This modification is conservative since $\tilde\pi_0^{|I|}$ is always greater than or equal to $\hat\pi_0^I$, but it allows to use Lemma 6 of [17] in order to get the post-hoc bound $\tilde d_{\mathrm{ASimes}}(I)$. The post-hoc bounds $\tilde d_{\mathrm{ASimes}}(I)$ are conservative since $\tilde d_{\mathrm{ASimes}}(I) \leq d_{\mathrm{ASimes}}(I)$ holds for every $I$, where $d_{\mathrm{ASimes}}(I)$ are the post-hoc bounds of CT$^*_\mathrm{ASimes}$.
The shortcut given in Algorithm <ref> allows calculation of $\tilde d_\mathrm{Simes}(S)$ in quadratic time.

$p$-values $p_1,\ldots,p_n$; index set $S$; significance level $\alpha$.
$(1-\alpha)$-confidence lower bound $\tilde d_{\mathrm{ASimes}}(S)$ for the number of true discoveries in $S$

Sort the $p$-values $p_{(1)} \leq \ldots \leq p_{(n)}$;

Compute $\tilde{h}=\max\{k\in\{0,\ldots,n\} : p_{(n-k+i)}>i\alpha/(k\tilde\pi_0^{k}) \,\,\mathrm{for\,\,} i=1,\ldots,k\}$

Compute  $\tilde d_{\mathrm{ASimes}}(S) = \max_{1\leq u\leq |S|} 1-u+|\{i\in S:p_{i}\leq u\alpha/(\tilde{h}\tilde\pi_0^{\tilde h})\}|$;

$\tilde d_\mathrm{ASimes}(S)$
Shortcut for closed testing with adaptive Simes local tests.

The properties of CT$^{*}_{\mathrm{ASimes}}$ are provided in the following proposition.

Under assumptions A<ref> and A<ref>, the adaptive Simes' test $\phi^{\mathrm{ASimes}}_{I_0}$ in (<ref>) combining permutation $p$-values $(p^*_i, i \in I_0)$ controls the Type I error at level $\alpha$, i.e.,
$\displaystyle \mathrm{pr}(\phi^{\mathrm{ASimes}}_{I_0} = 1; H_{I_0}) \leq \alpha$.  Consequently, CT$^{*}_{\mathrm{ASimes}}$ provides the true discovery guarantee (<ref>).

Algorithm <ref> offers conservative post-hoc bounds $\tilde d_{\mathrm{ASimes}}(I) \leq d_{\mathrm{ASimes}}(I)$ in quadratic time, with $\tilde d_{\mathrm{ASimes}} = d_{\mathrm{ASimes}} = n - \tilde{h}$ where $\tilde h$ is defined at line 2 of  Algorithm <ref>.

The proof follows by noting that adaptive Simes' test rejects $H_{I_0}$ if and only if the ABH$^{*}$ procedure applied to permutation $p$-values $(p_i^*, i \in I_0)$ rejects at least one hypothesis, which happens with probability bounded by $\alpha$  by Proposition 1.
The second part of the proposition follows from Lemma 6 in [17].

§ MONOTONICITY, SYMMETRY AND SEPARABILITY

We now recall some definitions from paper by [53] and show that WMW local test enjoys these properties. For a given $I\subseteq [n]$ they consider a local test which rejects $H_I$ according to the following indicator function

\begin{eqnarray}\label{form:localtestgeneral}
\mathds{1}\{f(|I|, (Y_i)_{i\in I})\leq \alpha\}.
\end{eqnarray}

A local test of form (<ref>) is called monotonic if for any $v\geq 1$ and any two sets of scores $(Y_1,\ldots,Y_v)$ and $(Y'_1,\ldots,Y'_v)$ such that $Y_i\leq Y'_i$ for all $i=1,\ldots,v$ we have

\begin{eqnarray*}
f(v, (Y_1,\ldots, Y_v)) \leq f(v, (Y'_1,\ldots,Y'_v)).
\end{eqnarray*}

Note that inequality among the scores is reversed with respect to definition in [53] since they consider scores whose smaller values capture evidence against the null hypothesis.
A local test of form (<ref>) is called symmetric if for any $v\geq 1$, any set of scores $(Y_1,\ldots,Y_v)$ and any permutation $\pi(1),\ldots,\pi(v)$ of $1,\ldots,v$ we have

\begin{eqnarray*}
f(v, (Y_1,\ldots, Y_v)) \leq f(v, (Y_{\pi(1)},\ldots,Y_{\pi(v)})).
\end{eqnarray*}

A local test of form (<ref>) is called separable if for any $v\geq 1$ and any set of scores $(Y_1,\ldots,Y_v)$ it exists a series of transformation functions $(h_i)_{i=1}^v$ on the scores and a function $g$ such that the local test in (<ref>) is equivalent to

\begin{eqnarray*}
\mathds{1}\left\{\sum_{i=1}^{v} h_i(Y_i)\leq g(v, \alpha)\right\}.
\end{eqnarray*}

[53] observe that a local test as in (<ref>) that is both symmetric and separable must satisfies condition of separability with transformation functions which are the same for each units in the test set, i.e.,

\begin{eqnarray} \label{form:localtestsummation}
\mathds{1}\left\{\sum_{i=1}^{v} h(Y_i)\leq g(v, \alpha)\right\}.
\end{eqnarray}

When considering the Lehmann's alternative of order $k>1$, we may choose to use the test statistics identified by Theorem <ref>

\begin{eqnarray*}
\mathds{1}\{T_{k, I}\geq c^{k}_{\alpha}(m, |I|)\},
\end{eqnarray*}

where $T_{k, I}=\sum_{i\in I}R_i(R_i+1)\cdot\ldots\cdot(R_i+k-2)$ and $c^{k}_{\alpha}(m, |I|)$ is the $(1-\alpha)$-quantile of the test statistic $T_{k,I}$ based on a control sample of size $m$. These tests have the properties of monotonicity, symmetry and separability.
Indeed, the test statistics identified by Theorem <ref> are of form (<ref>) with $h(Y_i)=R_i(R_i+1)\cdot\ldots\cdot(R_i+k-2)$, $i=1,\ldots,v$ and $g(v, \alpha)$ is the $(1-\alpha)$-quantile of $T_{k, I}$ which depends only on the significance level and on $v$. Then, the tests are separable. Since the rank $R_i$ is a non-decreasing function in $Y_i$ taking on only integer values for all $i=1,\ldots,v$, $h(Y_i)$ is a monotone non-decreasing function in $Y_i$ for all $i=1,\ldots,v$ and the tests are monotone. Finally, symmetry follows from the commutative property of summation.

§ ADMISSIBILITY

We start with recollecting some concepts and definitions. A statistical test $\delta$ is said to be admissible if it cannot be uniformly improved, i.e., if there does not exist another statistical test $\tilde\delta$ on the same hypothesis $H:P\in\mathcal{F}_H$, $\mathcal{F}_H\in\mathcal{F}$ such that $\tilde\delta \geq \delta$ and $P(\tilde\delta > \delta)>0$ for some $P\in\mathcal{F}$. Analogously a true discovery guarantee procedure which is monotone is admissible if it cannot be uniformly improved. Lemma 5 by [17] ensures that all closed testing procedures are monotone.
Moreover, a monotone procedure is admissible if all its local members are admissible. Therefore, to check the admissibility of the closed testing procedure with WMW local tests we only need to check the admissibility of the single WMW tests. A useful criterion is provided by [17] in Lemma 7 of the Supplementary Material, which we recall hereafter.

Let us suppose that for every measurable set $E$ and for every $P_1, P_2\in \mathcal{F}$ we have that $P_1(E)=0$ implies $P_2(E)=0$. Then, a statistical test $\bm{\delta}$ of hypothesis $\mathcal{F}_H\subseteq\mathcal{F}$ is admissible if $P\in \mathcal{F}_H$ exists such that $P(\bm{\delta}=1)=\alpha$.

To prove the first assumption of Lemma <ref> for the set $\mathcal{F}$ of
continuous c.d.f.'s stochastically larger than the null distribution $F$ and with same support, we consider the Lebesgue measure $\mu$ and the Borelian $\sigma$-field $\mathcal{B}(\mathds{R})$ on the real line. Let $F_1, F_2\in\mathcal{F}$ distributions whose corresponding measures on $\mathds{R}$ are denoted by $P_1$ and $P_2$. Let $E\in\mathcal{B}(\mathds{R})$ such that $P_1(E)=0$. If the intersection between the event $E$ and the support of $F$ is empty, then $P_2(E)=0$ since $P_1$ and $P_2$ have the same support. If not, let $B$ be the non-empty intersection. Thanks to the monotonicity property of the probability, since $B\subseteq E$ we have $P_1(B)=0$. In addition, also $\mu(B)=0$. Denoting by $f$ the density function of $F_1$ with respect to the Lebesgue measure, we have

\begin{equation}\label{eq:appendixAdm}
0=P_1(B) = \int_{B}dP_1 = \int_{B}f d\mu
\end{equation}

and, since $f>0$ on the support of $F$, in order to satisfy (<ref>) necessarily $\mu(B)=0$, which implies $P_2(B)=0$ because $P_2$ is
continuous with respect to Lebesgue measure. The assumption of common null sets is, then, proved.

Under the null hypothesis $H_I$ that $G_i=F$ for every $i \in I\subseteq [n]$, the Wilcoxon-Mann-Whitney test statistic $U_I$ follows the distribution defined recursively by [39]. Then,  the attainable levels are the following ones:

\begin{eqnarray*}
%pr(T_{W,\{i\}} \geq c_{\alpha}(m,|I|)), \qquad  c_{\alpha}(m,|I|) \in \{0,1,\ldots,m|I|\}, \hspace{2mm} I\subseteq [n].
\alpha^{|I|}_x = \mathrm{pr}(U_I\geq x), \qquad  x \in \{0,1,\ldots,m|I|\}.
\end{eqnarray*}

This shows that also WMW test is of exact size $\alpha$ only for particular values of $\alpha$, $\{ \alpha^{|I|}_x, x = 0,1,\ldots,m|I| \}$.

When considering integers $k>2$, LMPI test statistics are $U$-statistics as well. Indeed, $T_{k,I}$ can be expressed as the weighted sum of the test ranks $(R_i)_{i=1}^n$ to the powers from $1$ to $k-1$, e.g., for $I=[n]$ and $k=3$ and $k=5$ we have

\begin{eqnarray*}
T_3 = \sum_{i=1}^n R_i^2+R_i, \qquad T_5 = \sum_{i=1}^n R_i^4+6R_i^3+11R_i^2+6R_i.
\end{eqnarray*}

The $h$th power of rank $R_i$ can be expressed as a multiple of the sum of indicator functions involving events of the form $\{X_{s_1}< Y_i,\ldots,X_{s_h}<Y_i\}$. For example, for $h=2$

\begin{eqnarray}
R_i^2 = \left(\sum_{s=1}^m \mathds{1}\{X_s<Y_i\}\right)^2 = \sum_{s_1,s_2=1}^m \mathds{1}\{X_{s_1}<Y_i, X_{s_2}<Y_i\}
\end{eqnarray}

where the weight of the sum is $1$ in this case, and in the general case depends on the power $h$ and are linked to the weights of the expansion of the power of the sum.
Thus, for any $I\subseteq [n]$ $T_{k,I}$ can be written as the weighted sum (on $i\in I$) of the test ranks to the power of $h$ ranging from $1$ to $k-1$ and the $h$th power of $R_i$ is a multiple of an indicator function involving events of the form $\{X_{s_1}< Y_i,\ldots,X_{s_h}< Y_i\}$.
Therefore, we have a representation of $T_{k,I}$ as a $U$-statistic whose kernel function is the weighted summation of indicator functions, which is a function of all the control variables $X_1,\ldots,X_m$ and of one test variable $Y_i$, satisfying symmetry requirement (see [27]).
For example, for $k=3$ the test statistic can be written as

\begin{eqnarray*}
T_{3,I} = \sum_{i\in I}\left( \sum_{s=1}^n \mathds{1}\{X_s<Y_i\}+\sum_{s_1,s_2=1}^m\mathds{1}\{X_{s_1}<Y_i, X_{s_2}<Y_i\}\right),
\end{eqnarray*}
where the kernel function is the summation inside brakets and is a symmetric function in $X_1,\ldots, X_m$ and in $Y$.
As a consequence, their distribution can be approximated by the Central Limit Theorem for $U$-statistics <cit.>, [27]. Hence, for large values of the control size and the test size, under the null hypothesis that $H_I$, the attainable levels of the LMPI test statistics for integers $k>2$ are

\begin{eqnarray}\label{TkI.alphas}
\alpha^{k,|I|}_x = \mathrm{pr}(T_{k,I}\geq x)
\end{eqnarray}
with $x$ varying in the range of $T_{k,I}$,
proving that also WMW tests are of exact size $\alpha$ only for particular values of $\alpha$ defined in equation (<ref>).

§ INVARIANCE AND RANKS TESTS

In order to recall basic definitions and concepts about invariance we refer to [35].
Given a random variable $X\sim F \in\mathcal{F}$ where $\mathcal{F}$ is the set of
continuous c.d.f's, we consider a transformation $\rho$ of the sample space $\mathcal{X}$ of $X$. By definition a transformation is a bijective function of some space onto itself. We indicate with $\rho X$ the random variable which takes on the value $\rho x$ when $X=x$ and with $\bar\rho F$ distribution of $\rho X$, i.e. $\bar\rho$ is the transformation of the model space $\mathcal{F}$ corresponding to transformation $\rho$ in the sample space. Let us now consider the problem of testing

\begin{eqnarray}\label{appendixInv:testingproblem}
H:F\in\mathcal{F}_H \qquad \text{vs} \qquad K:F\in\mathcal{F}_K
\end{eqnarray}

with $\mathcal{F}=\mathcal{F}_H \cup \mathcal{F}_K$ and $\mathcal{F}_H\cap\mathcal{F}_K=\emptyset$.

The testing problem (<ref>) remains invariant under a transformation $\rho$ if $\bar\rho$ preserves both $\mathcal{F}_H$ and $\mathcal{F}_K$, i.e.

\begin{equation}\label{eq:probinv}
\bar\rho\mathcal{F} = \mathcal{F} \qquad \text{ and } \qquad \bar\rho\mathcal{F}_H = \mathcal{F}_H.
\end{equation}

The class of transformations $\rho$ preserving a testing problem can always be extended to an algebraic group, which will be denoted as $A$ and called transformation group. Also the class of transformations $\bar\rho$ on the model space can be extended to an algebraic group $\bar A$. We will use the classical notation $\delta:\mathcal{X}\to [0,1]$ to indicate a statistical test, with $1$ indicating rejection.

A test $\delta$ for an invariant problem (<ref>) is said to be invariant under a transformation $\rho$ if

\begin{equation} \label{eq:testinv}
\delta(\rho x) = \delta(x) \qquad \forall x\in\mathcal{X}.
\end{equation}

Moreover, we will say that the testing problem (<ref>) remains invariant under $A$ if (<ref>) holds for all $\rho\in A$ and a test $\delta$ is said to be invariant under $A$ if (<ref>) holds for all $\rho\in A$. A test invariant under the group of continuous and strictly increasing transformations is also called rank test, since it depends only on the ranks of the observations. To further explore the concept of invariance, we need to define the orbits of a transformation group. For all $x\in \mathcal{X}$ its A-orbit is defined as

\begin{equation*}
A^x = \{\rho x : \rho\in A \}.
\end{equation*}

It is known that the orbits of $A$ constitute a partition of the sample space $\mathcal{X}$.
By definition <ref> it follows that a function $\eta$ defined on the sample space $\mathcal{X}$ is invariant under A if and only if it is constant on each orbit. In addition a function $\eta$ is said to be maximal invariant with respect to A if it is invariant and if

\begin{equation*}
\eta(x_1)=\eta(x_2) \Rightarrow \exists \rho\in A : x_2=\rho x_1
\end{equation*}

that is, it is constant on the orbits and it takes on a different value for each orbit.

§ PROOFS

§.§ Proof of Theorem <ref>

To simplify notation, we write $H:\theta=0$ instead of $H_0:\theta=0$.
In the setting of Theorem <ref> we are testing the null hypothesis $H$ that the population parameter $\theta$ is zero, against the alternative $K$ that it is greater than zero.
Both hypotheses can be reformulated in terms of ranks, i.e. the null hypothesis is equivalent to require that $Y_{1},\ldots, Y_{n}$ are i.i.d. with c.d.f. $F$ and the alternative is equivalent to take $Y_{1},\ldots, Y_{n}$ i.i.d. from c.d.f. $(1-\theta)F+\theta F^k$, that is testing a simple hypothesis against a simple alternative. Thus, its solution is given by the lemma of Neyman and Pearson [44]. So, we consider the likelihood ratio

\begin{eqnarray*}
\lambda({\bf y}) =\frac{L_K({\bf y})}{L_H({\bf y})}
\end{eqnarray*}

where $L_H$ and $L_K$ denote respectively the likelihood function under the null and under the alternative and ${\bf y}$ is a realization of the scores vector of the test observations. We reject the null hypothesis at the significance level $\alpha$ when the likelihood ratio is too large and the rejection region

\begin{eqnarray*}
\mathcal{R}=\left\{{\bf y}\in\mathcal{Y} : \lambda({\bf y})>v\right\}
\end{eqnarray*}

is such that $\mathrm{pr}(\lambda({\bf y})>v;H)=\alpha$, $v>0$.
We restrict our attention to rank tests, then the likelihood ratio statistic can be expressed as a function of only the rank vector of ${\bf y}$, i.e.,

\begin{eqnarray*}
&&\lambda({\bf r}) =\frac{L_K({\bf r})}{L_H({\bf r})}, \hspace{2mm} {\bf r}\in [N]^n \\&&
\mathcal{R}=\left\{{\bf r}\in [N]^n : \lambda({\bf r})>v\right\} \text{ s.t. } \mathrm{pr}(\lambda({\bf r})>v;H)=\alpha, \hspace{2mm} v>0,
\end{eqnarray*}

where $L_H({\bf r})=:\mathrm{pr}(R_1=r_1,\ldots,R_n=r_n;H)$ and $L_K({\bf s})=:\mathrm{pr}(R_1=r_1,\ldots,R_n=r_n;K)$.
Then, after identifying the rejection region $\mathcal{R}$ we are able to compute the power of the test and maximize it in a neighborhood of the population null, i.e., $\theta \approx 0$.
The proof can be divided into two steps. First, we compute $L_H$ and $L_K$ obtaining

\begin{eqnarray*}
L_H({\bf r}) &=&
\frac{1}{\binom{m+n}{m}}, \hspace{2mm} {\bf r}\in[N]^n\\
L_K({\bf r}) &=&
\frac{1}{\binom{m+n}{m}}\mathds{E}[g'_p(U_{r_1})\cdot\ldots\cdot g'_p(U_{r_n})], \hspace{2mm} {\bf r}\in[N]^n
\end{eqnarray*}

where $U_{r_i}$ is the $r_i^{th}$ order statistics in a sample of $N$ variables distributed uniformly over $[0,1]$ and $g_{\theta}(u) = (1-\theta)+\theta u^{k}$ so that $g_{\theta}'(u) = (1-\theta)+\theta ku^{k-1}$ for integers $k>1$.
Second, we apply Neyman-Pearson lemma [44] and find the critical region of the locally most powerful rank test.

Step 1

Let ${\bf S}=(S_1,\ldots,S_m)$ be the rank vector of ${\bf X}$ in the pooled sample ${\bf Z}$ and ${\bf T}=(S_1,\ldots,S_m, R_1,\ldots,R_n)\in[N]^N$ be the pooled rank vector. Defining $\mathcal{Y}$ as the score vectors space of test observations and

\begin{eqnarray*}
\mathcal{Y}({\bf R}) &=& \{ {\bf y}\in\mathcal{Y}: (y_{(1)},\ldots, y_{(n)}) = (y_{R_1},\ldots, y_{R_n})%\text{the rank vector of }{\bf y} \text{ is } {\bf R}
\},
\end{eqnarray*}

as the subset of the score vectors with rank vector $R$, we have that $\mathrm{pr}(R_1=r_1,\ldots,R_n=r_n) = pr({\bf Y}\in\mathcal{Y}({\bf r}))$ both under the null and under the alternative and the number of disjoint sets $\mathcal{Y}({\bf R})$ in the score vectors space is equal to the number $\binom{m+n}{m}$ of permutations of the last $n$ elements of ${\bf T}$ and $\mathcal{Y}$ is the union of the $\binom{m+n}{m}$ disjoint sets ${\mathcal{Y}({\bf R})}$. Then,

\begin{eqnarray*}
Pr({\bf Y}\in \mathcal{Y})=\sum_{{\bf R}} \mathrm{pr}({\bf Y}\in\mathcal{Y}({\bf R})) = 1
\end{eqnarray*}

Under the null hypothesis $H$ each permutation of the score vector ${\bf Y}$ is equally likely so that

\begin{eqnarray}\label{eq:probSNull}
L_H({\bf r}) = Pr({\bf Y}\in\mathcal{Y}({\bf r});H) = \frac{1}{\binom{m+n}{m}}, \qquad {\bf r}\in[N]^n
\end{eqnarray}

and it does not depend on ${\bf r}$.

Shifting the attention to the alternative hypothesis, we denote by $f$ and $g_{\theta}'(F)$ the density functions of $F$ and $g_{\theta}(F)$ which are continuous c.d.f.'s and by $\mathcal{Z}$ the score vectors space of the pooled observations of control and test samples and we define

\begin{eqnarray*}
&&\mathcal{Z}({\bf T}) = \{{\bf z}\in\mathcal{Z}: (z_{(1)},\ldots, z_{(N)}) = (z_{T_1},\ldots, z_{T_N})\}, \\&&
\mathcal{Z}({\bf id}) = \{{\bf z}\in\mathcal{Z}: (z_{(1)},\ldots, z_{(N)}) = (z_{1},\ldots, z_{N})\}.
\end{eqnarray*}

By definition, $\mathcal{Z}({\bf id})$ is the product sample space given by $\mathcal{X}({\bf id}) \times \mathcal{Y}({\bf id})$, where $\mathcal{X}$ is the score vectors space of control samples and $\mathcal{X}({\bf id})$ and $ \mathcal{Y}({\bf id})$ are defined similarly to $\mathcal{Z}(id)$.
Then, under the population alternative the joint likelihood

\begin{eqnarray*}
\end{eqnarray*}

is equal to

\begin{eqnarray} \label{jointLikAlt}
pr({\bf Z}\in\mathcal{Z}(T)) &=& \int_{\mathcal{Z}({\bf T})}\prod_{j=1}^m f(x_j) \prod_{i=1}^n g_{\theta}'(F(y_i)) d{\bf z} = \nonumber \\&=&
\int_{\mathcal{Z}({\bf id})}\prod_{h=1}^{N}f(z_{t_h}) \prod_{i=1}^n \frac{g_{\theta}'(F(y_{r_i}))}{f(y_{r_i})} d{\bf z} =  \nonumber \\&=&
\frac{1}{N!}\int_{\mathcal{Z}({\bf id})}\prod_{i=1}^n \frac{g_{\theta}'(F(y_{r_i}))}{f(y_{r_i})} \underbrace{N!\prod_{h=1}^{N} f(z_{t_h})}_{f^*({\bf z})}  d{\bf z} = \\&=&
\frac{1}{N!}\mathds{E}_{f^*}\left[\prod_{i=1}^n \frac{g_{\theta}'(F(Y_{r_i}))}{f(Y_{r_i})}\right]  \nonumber
\end{eqnarray}

where ${\bf t}=(t_1,\ldots,t_N)$ is a realization of the rank vector ${\bf T}$ and

\begin{eqnarray*}
f^*({\bf z}) = \begin{cases}
N!\prod_{h=1}^{N} f(z_h), &  {\bf z}\in \mathcal{Z}(id) \\
0, & \text{otherwise}
\end{cases}
\end{eqnarray*}
is the joint density of the ordered statistic of ${\bf Z}$.
In order to find the marginal probability $\mathrm{pr}(R_1=r_1,\ldots,R_n=r_n; K)$ we need to integrate ${\bf S}$ out from the joint likelihood (<ref>). Using independence between the score vectors of the control and test samples we obtain that the marginal probability of the rank vector ${\bf R}$ under the population alternative is

\begin{eqnarray}\label{LikAlt}
&&\frac{m!n!}{N!} \underbrace{\int_{\mathcal{X}} \prod_{j=1}^m f(x_{s_j}) d{\bf x}}_{=1}  \int_{\mathcal{Y}({\bf id})} \prod_{i=1}^n \frac{g_{\theta}'(F(y_{r_i}))}{f(y_{r_i})} \frac{N!}{m!n!}\prod_{i=1}^n f(y_{r_i})d{\bf y}  = \nonumber \\
&&\frac{1}{\binom{m+n}{m}}\int_{\mathcal{Y}({\bf id})} \prod_{i=1}^n g_{\theta}'(F(y_{r_i})) \binom{m+n}{m} \underbrace{\prod_{i=1}^n \mathds{1}\{F(y_{r_i})\in[0,1]\}}_{=1} d{\bf y}  = \\
&&\frac{1}{\binom{m+n}{m}}\mathds{E}\left[\prod_{i=1}^n g_{\theta}'(U_{r_i})\right] \nonumber
\end{eqnarray}

where $\binom{m+n}{m} \prod_{i=1}^n \mathds{1}\{F(y_{r_i})\in[0,1]\}$ is the joint density functions of the ordered $(F(Y_{r_1}), \ldots, F(Y_{r_n}))$ (out of the pooled vector) which are distributed as $F$ under the population null so that $F(Y)$ is uniformly distributed over $[0,1]$.

Step 2

Since $L_H({\bf r})$ does not depend on ${\bf r}$, the rejection rule becomes

\begin{eqnarray*}
\mathcal{R} = \left\{{\bf r} \in [N]^n : \lambda({\bf r}) > \tilde v, \hspace{2mm} \tilde v = vL_H, \hspace{2mm} v>0 \right\}
\end{eqnarray*}

and the power function of the test is

\begin{eqnarray}
\beta(\theta) = \mathrm{pr}({\bf r}\in\mathcal{R}; K) = \sum_{{\bf r}\in\mathcal{R}} L_{K, \theta}({\bf r})
\end{eqnarray}

where we used the notation $L_{K,\theta}$ to emphasize the dependence of the likelihood on $\theta$. We know that for any $\theta\in[0,1]$ the power function is $\beta(0)=\mathrm{pr}({\bf r}\in\mathcal{R};H)=\alpha$.
So, to maximize it in a neighborhood of the population null we need to maximize the first derivative

\begin{eqnarray*}
\beta'(0)=\frac{d}{d\theta}\beta(\theta)\Bigr|_{\theta=0} = \sum_{{\bf s}\in\mathcal{R}}\frac{d}{d\theta}L_{K,\theta}\Bigr|_{\theta=0}.
\end{eqnarray*}

Deriving equation (<ref>) we obtain

\begin{eqnarray*}
\frac{d}{d\theta}L_{K,\theta}({\bf r}) &=& \frac{1}{\binom{m+n}{m}}\mathds{E}\left[\frac{d}{d\theta}\prod_{i=1}^n\left((1-\theta)+\theta kU_{r_i}^{k-1}\right)\right] =\\&=&
\frac{1}{\binom{m+n}{m}}\mathds{E}\left\{\sum_{i=1}^n \left[(kU_{r_i}^{k-1}-1)\prod_{\substack{h=1 \\ h\neq i}}^n (1-\theta+k\theta U_{r_h}^{k-1})\right]\right\},
\end{eqnarray*}

which evaluated in $\theta=0$ gives

\begin{eqnarray*}
\frac{d}{d\theta}L_{K,\theta}({\bf r}) &=& \frac{1}{\binom{m+n}{m}}\mathds{E}\left[\sum_{i=1}^n (kU_{r_i}^{k-1}-1)\right] =\\&=&
\frac{1}{\binom{m+n}{m}}\left[k\sum_{i=1}^n \mathds{E}\left[U_{r_i}^{k-1}\right] - n\right] =\\&=&
\frac{1}{\binom{m+n}{m}}\left[\frac{k}{(m+n+1)\cdot\ldots\cdot(m+n+k-1)}\sum_{i=1}^n r_i(r_i+1)\cdot\ldots\cdot (r_i+k-2) - n\right],
\end{eqnarray*}

where in the last step we computed the $(k-1)$th moment of the Beta distribution with parameters $r_i$ and $n+m-r_i+1$ which is the distribution of the $r_i$th ordered statistic of the Uniform distribution on $[0,1]$. Then, for all integers $k>1$ we have that

\begin{eqnarray*}
\beta'(0) = \sum_{{\bf r}\in\mathcal{R}}\frac{1}{\binom{m+n}{m}}\left[\frac{k}{(m+n+1)\cdot\ldots\cdot(m+n+k-1)}\sum_{i=1}^n r_i(r_i+1)\cdot\ldots\cdot (r_i+k-2) - n\right]
\end{eqnarray*}

is maximized by rejecting the population null when

\begin{eqnarray} \label{LMPIrejreg}
\sum_{i=1}^n r_i(r_i+1)\cdot\ldots\cdot (r_i+k-2)
\end{eqnarray}

is large. It follows that when testing $H:\theta=0$ the LMPI test is the one whose rejection region is identified by the test statistic in (<ref>).

§.§ Proof of Lemma <ref>

Given the control distribution $X \sim F$ and denoting by $G$ the distribution of the test score $Y$, according to Definition <ref> in Appendix <ref> we need to prove that problem

\begin{eqnarray} \label{appendixProofs:testingProblem}
H:G=F \qquad \text{vs} \qquad K:G<F
\end{eqnarray}

is invariant under the group $A$ of continuous and strictly increasing transformations onto the score vectors space and that BH procedure is invariant under $A$ as well. We recall that both $F$ and $G$ belong to the model space $\mathcal{F}$ of
continuous c.d.f.'s.

Step 1: invariance of the testing problem

Considering the transformed scores $\rho X\sim \bar \rho F$ and $\rho Y\sim \bar \rho G$, $\rho \in A$ we obtain a transformation of the model space and its partition through $\bar\rho$.

Referring to Definition <ref> we need to show that $\bar\rho\mathcal{F} = \mathcal{F}$ and $\bar\rho\mathcal{F}_H = \mathcal{F}_H$ for any $\rho\in A$.
The invariance of the testing problem under the group $A$ of continuous and strictly increasing transformations follows from the fact that these transformations preserve both continuity of a distribution and the property of two variables being either identically distributed or stochastically larger than the other.
It is important to observe that since $\bar\rho$ is bijective, it is sufficient to prove only one inclusion, e.g. $\bar\rho\mathcal{F}\subseteq\mathcal{F}$, in order to show the equality $\bar\rho\mathcal{F}=\mathcal{F}$.

Step 2: invariance of the test

We now have to prove that the result of BH procedure on a score vector $({\bf x}, {\bf y})$ does not change if we transform the scores with a continuous and strictly increasing transformation $\rho\in A$. Let us consider the control score vector ${\bf x}\in \mathcal{X}$ and the test score vector ${\bf y}\in \mathcal{Y}$ and their transformations $\rho{\bf x}$ and $\rho{\bf y}$ with $\rho\in A$. The conformal $p$-values computed on the original scores $p^*_i = (1+\sum_{k=1}^m \mathds{1}\{x_k\leq y_i\})/(m+1)$ and on the transformed scores $p'^{*}_i = (1+\sum_{k=1}^m \mathds{1}\{\rho x_k\leq \rho y_i\})/(m+1)$ coincide since $\rho$ is strictly increasing. Thus, the result of BH procedure on a score vector $({\bf x}, {\bf y})$ does not change if we transform the scores with a continuous and strictly increasing transformation $\rho\in A$.

§ FURTHER NUMERICAL RESULTS

We report additional analysis on ,   and  datasets. We fitted an Isolation Forest model as described in § 4.1 on a random sample of $l=199$ training observations from the class of inliers. At each replication we used a random sample of $m=199$ units as control sample and a random sample of $n=20$ units as test sample. For each value of the number of outliers $n_1\in\{0,\ldots,n\}$, we performed $10^4$ replications at level $\alpha=0.1$ in order to estimate the power to reject the global null $H_g$ that all test observations are inliers and the average lower bound $d$ on the number of discoveries, i.e., true outliers. In  dataset outliers are identified as the minority class, while in  dataset the smallest five classes are combined to form the outliers class and in  dataset the outliers class is the class of rare objects [16, 31].

CT$^*_\mathrm{WMW}$ exhausts the $\alpha$ level in the  and  datasets when $n_1=0$, but it exhibits conservatism in the  dataset due to the presence of ties in the scores of the combined sample.

Power to reject the global null (left column) and average lower bound $d$ (right column) for  , (first row),  (second row) and  (third row) datasets for different methods, as a function of the number $n_1$ of outliers in the test sample.
Samples are drawn from each dataset with training size of $l=199$, control size of $m=199$ and test size of $n=20$, and varying the number of outliers from $0$ to $n$. The level is set at $\alpha=0.1$ and all the quantities are estimated with $10^4$ replications.

title = Shuttle,
xmin = 0,
xmax = 20,
xtick =0,...,20,
xticklabels =0,,,,,,,,,,10,,,,,,,,,,20,
ymin = 0,
ymax = 20,
xlabel= Number of outlier $n_1$,
ylabel= Average $d$ and average $|D|$,
legend style=legend pos=north west
[mark size=1.5pt, color=gray, mark=o, line width=1.5] coordinates
( 0,0.1093 ) ( 1,0.9271 ) ( 2,1.7959 ) ( 3,2.6793 ) ( 4,3.567 ) ( 5,4.4812 ) ( 6,5.4048 ) ( 7,6.3282 ) ( 8,7.2621 ) ( 9,8.2547 ) ( 10,9.2338 ) ( 11,10.3598 ) ( 12,11.4111 ) ( 13,12.4024 ) ( 14,13.4029 ) ( 15,14.4897 ) ( 16,15.6145 ) ( 17,16.7268 ) ( 18,17.8592 ) ( 19,19.04 ) ( 20,19.9999 )
$d$ Simes$^*$

[mark size=1.5pt, color=red, mark=o, line width=1.5] coordinates
( 0,0.2169 ) ( 1,0.3679 ) ( 2,0.6577 ) ( 3,1.0923 ) ( 4,1.714 ) ( 5,2.4788 ) ( 6,3.4237 ) ( 7,4.4359 ) ( 8,5.5106 ) ( 9,6.647 ) ( 10,7.777 ) ( 11,8.879 ) ( 12,10.0075 ) ( 13,11.1732 ) ( 14,12.3071 ) ( 15,13.4712 ) ( 16,14.6203 ) ( 17,15.8181 ) ( 18,17.0543 ) ( 19,18.3505 ) ( 20,20 )
$d$ WMW

[no marks, color=gray, line width = 1.5, dashed] coordinates
( 0,0 ) ( 1,0.9048 ) ( 2,1.7094 ) ( 3,2.5087 ) ( 4,3.3075 ) ( 5,4.1167 ) ( 6,4.9271 ) ( 7,5.7233 ) ( 8,6.5511 ) ( 9,7.3894 ) ( 10,8.2797 ) ( 11,9.8158 ) ( 12,11.1074 ) ( 13,12.1431 ) ( 14,13.1162 ) ( 15,14.3076 ) ( 16,15.5122 ) ( 17,16.6555 ) ( 18,17.8137 ) ( 19,19.0301 ) ( 20,19.9999 )
$|D|$ Simes$^{*}$

[no marks, color=red, line width = 1.5, dashed] coordinates
( 0,0 ) ( 1,0 ) ( 2,0 ) ( 3,0 ) ( 4,0 ) ( 5,0 ) ( 6,0 ) ( 7,0 ) ( 8,0.001 ) ( 9,0.001 ) ( 10,0.009 ) ( 11,0.0152 ) ( 12,0.0254 ) ( 13,0.0398 ) ( 14,0.1163 ) ( 15,0.281 ) ( 16,0.6694 ) ( 17,1.4759 ) ( 18,3.6 ) ( 19,8.6672 ) ( 20,19.9999 )
$|D|$ WMW

Average lower confidence bound $d$ for the number of discoveries and average number of FWER discoveries $|\mathcal{D}|$ for the  for the
CT$^*_\mathrm{Simes}$ and CT$^*_\mathrm{WMW}$.

[1]
Bates, S., E. Candès, L. Lei, Y. Romano, and M. Sesia (2023).
Testing for outliers with conformal p-values.
Ann. Stat. 51(1), 149–178.

[2]
Bauer, P. (1991).
Multiple testing in clinical trials.
Statistics in medicine 10(6), 871–890.

[3]
Benjamini, Y. and Y. Hochberg (1995).
Controlling the false discovery rate: a practical and powerful
approach to multiple testing.
J. R. Statist. Soc. B 57(1), 289–300.

[4]
Benjamini, Y. and D. Yekutieli (2001).
The control of the false discovery rate in multiple testing under
Ann. Stat. 29(4), 1165–1188.

[5]
Blanchard, G., P. Neuvial, and E. Roquain (2020).
Post hoc confidence bounds on false positives using reference
Ann. Stat. 48(3), 1281–1303.

[6]
Bogomolov, M. (2023).
Testing partial conjunction hypotheses under dependency, with
applications to meta-analysis.
Electron. J. Stat. 17(1), 102–155.

[7]
Bretz, F., T. Hothorn, and P. Westfall (2016).
Multiple comparisons using R.
CRC press.

[8]
Buckle, N., H. Kraft, Charles, and C. van Eeden (1969).
An approximation to the wilcoxon-mann-whitney distribution.
J. Am. Stat. Assoc. 64(326), 225–251.

[9]
Cox, D. R. and D. V. Hinkley (1979).
Theoretical statistics.
CRC Press.

[10]
Donoho, D. and J. Jin (2004).
Higher criticism for detecting sparse heterogeneous mixtures.
Ann. Stat. 32(3), 962–994.

[11]
Donoho, D. and J. Jin (2015).
Higher criticism for large-scale inference, especially for rare and
weak effects.
Statist. Sci. 30(1), 1–25.

[12]
Dunnett, C. W. (1955).
A multiple comparison procedure for comparing several treatments with
a control.
J. Am. Stat. Assoc. 50(272), 1096–1121.

[13]
Fisher, R. A. (1970).
Statistical methods for research workers.
In Breakthroughs in statistics: Methodology and distribution,
pp.  66–70. Springer.

[14]
Fix, E. and L. J. Hodges, Joseph (1955).
Significance probabilities of the wilcoxon test.
Ann. Math. Stat. 26(2), 301–312.

[15]
Gao, Z. and Q. Zhao (2023).
Simultaneous hypothesis testing using internal negative controls with
an application to proteomics.

[16]
Geusebroek, J. M., G. J. Burghouts, and A. W. M. Smeulders (2005).
The amsterdam library of object images.
International Journal of Computer Vision 61(1),

[17]
Goeman, J. J., J. Hemerik, and A. Solari (2021).
Only closed testing procedures are admissible for controlling false
discovery proportions.
Ann. Stat. 49(2), 1218–1238.

[18]
Goeman, J. J., R. J. Meijer, T. J. Krebs, and A. Solari (2019).
Simultaneous control of all false discovery proportions in
large-scale multiple hypothesis testing.
Biometrika 106(4), 841–856.

[19]
Goeman, J. J. and A. Solari (2011).
Multiple testing for exploratory research.
Stat. Sci. 26(4), 584–597.

[20]
Goeman, J. J. and A. Solari (2014).
Multiple hypothesis testing in genomics.
Statistics in medicine 33(11), 1946–1978.

[21]
Heller, R., A. Krieger, and S. Rosset (2022).
Optimal multiple testing and design in clinical trials.
Biometrics, 1–12.

[22]
Heller, R. and A. Solari (2023).
Simultaneous directional inference.

[23]
Hemerik, J. (2023).
On the term "randomization test".

[24]
Hemerik, J. and J. Goeman (2018).
Exact testing with random permutations.
Test 27(4), 811–825.

[25]
Hemerik, J. and J. J. Goeman (2021).
Another look at the lady tasting tea and differences between
permutation tests and randomisation tests.
Int. Stat. Rev. 89(2), 367–381.

[26]
Hemerik, J., A. Solari, and J. J. Goeman (2019).
Permutation-based simultaneous confidence bounds for the false
discovery proportion.
Biometrika 106(3), 635–649.

[27]
Hoeffding, W. (1948).
A class of statistics with asymptotically normal distribution.
Ann. Math. Stat. 19(3), 293–325.

[28]
Hoeffding, W. (1951).
"optimum" nonparametric tests.
Proceedings of the Second Berkeley Symposium on Mathematical
Statistics and Probability 2, 83–92.

[29]
Hommel, G. (1988).
A stagewise rejective multiple test procedure based on a modified
bonferroni test.
Biometrika 75(2), 383–386.

[30]
Kim, I., S. Balakrishnan, and L. Wasserman (2022).
Minimax optimality of permutation tests.
Ann. Stat. 50(1), 225–251.

[31]
Kriegel, H.-P., P. Kroger, E. Schubert, and A. Zimek (2011).
Interpreting and Unifying Outlier Scores, pp.  13–24.
Society for Industrial and Applied Mathematics.

[32]
Kuchibhotla, A. K. (2021).
Exchangeability, conformal prediction, and rank tests.

[33]
Lehmann, E. L. (1953).
The power of rank tests.
Ann. Math. Stat. 24, 23–42.

[34]
Lehmann, E. L. (2012).
Parametric versus nonparametrics: two alternative
methodologies, pp.  437–445.
Springer US.

[35]
Lehmann, E. L. and J. P. Romano (2005).
Testing Statistical Hypotheses (3 ed.).
Springer Texts in Statistics. Springer New York, NY.

[36]
Liang, Z., Y. Zhou, and M. Sesia (2023).
Conformal inference is (almost) free for neural networks trained with
early stopping.

[37]
Liu, F. T., K. M. Ting, and Z.-H. Zhou (2008).
Isolation forest.
In 2008 Eighth IEEE International Conference on Data Mining,
pp.  413–422.

[38]
Ludbrook, J. and H. Dudley (1998).
Why permutation tests are superior to t and f tests in biomedical
Am. Stat. 52(2), 127–132.

[39]
Mann, H. B. and D. R. Whitney (1947a).
On a test of whether one of two random variables is stochastically
larger than the other.
Ann. Math. Stat. 18(1), 50 – 60.

[40]
Mann, H. B. and D. R. Whitney (1947b).
On a test of whether one of two random variables is stochastically
larger than the other.
Ann. Math. Stat. 18(1), 50–60.

[41]
Marandon, A., L. Lei, D. Mary, and E. Roquain (2022).
Machine learning meets false discovery rate.

[42]
Marcus, R., E. Peritz, and K. R. Gabriel (1976).
Closed testing procedures with special reference to ordered analysis
of variance.
Biometrika 1(63), 655–660.

[43]
Mary, D. and E. Roquain (2022).
Semi-supervised multiple testing.
Electron. J. Stat. 16(2), 4926–4981.

[44]
Neyman, J. and E. S. Pearson (1933).
On the problem of the most efficient tests of statistical hypotheses.
Philos. Trans. Royal Soc. 231, 289–337.

[45]
Pesarin, F. and L. Salmaso (2010).
Permutation tests for complex data: theory, applications and
John Wiley & Sons.

[46]
Rosset, S., R. Heller, A. Painsky, and E. Aharoni (2020).
Optimal and maximin procedures for multiple testing problems.

[47]
Schweder, T. and E. Spjøtvoll (1982).
Plots of p-values to evaluate many tests simultaneously.
Biometrika 69(3), 493–502.

[48]
Serfling, R. J. (1980).
Approximation theorems of mathematical statistics.
John Wiley & Sons.

[49]
Simes, R. J. (1986).
An improved bonferroni procedure for multiple tests of significance.
Biometrika 73(3), 751–754.

[50]
Steel, R. G. (1959).
A multiple comparison rank sum test: treatments versus control.
Biometrics 15(4), 560–572.

[51]
Storey, J. D. (2002).
A direct approach to false discovery rates.
J. R. Statist. Soc. B 64(3), 479–498.

[52]
Storey, J. D., J. E. Taylor, and D. Siegmund (2004).
Strong control, conservative point estimation and simultaneous
conservative consistency of false discovery rates: A unified approach.
J. R. Statist. Soc. B 66(1), 187–205.

[53]
Tian, J., X. Chen, E. Katsevich, J. J. Goeman, and A. Ramdas (2023).
Large-scale simultaneous inference under dependence.
Scand. J. Stat. 50(2), 750–796.

[54]
van Zwet, W. R. (2011).
Erich lehmann on rank methods.
In Selected Works of EL Lehmann, pp.  331–338. Springer.

[55]
Vesely, A., L. Finos, and J. J. Goeman (2023).
Permutation-based true discovery guarantee by sum tests.
J. R. Statist. Soc. B 85(3), 664–683.

[56]
Vovk, V., A. Gammerman, and G. Shafer (2005).
Algorithmic learning in a random world, Volume 29.

[57]
Westfall, P. H. and S. S. Young (1993).
Resampling-based multiple testing: Examples and methods for
p-value adjustment, Volume 279.
John Wiley & Sons.

[58]
Wilcoxon, F. (1945).
Individual comparisons by ranking methods.
Biometrics Bulletin 1(6), 80–83.