Improving the Benjamini-Hochberg Procedure for Discrete Tests

TL;DR

This paper introduces new modifications to the Benjamini-Hochberg procedure that account for data discreteness and signal strength, providing theoretical FDR control guarantees and demonstrating improved power in genomic data analysis.

Contribution

It develops theoretically guaranteed, discrete-aware FDR procedures that adapt to signal proportion, filling a gap in existing methods.

Findings

New FDR procedures control the false discovery rate under independence.

The methods show increased power in simulations and real genomic datasets.

The approach combines discreteness and signal proportion adaptively.

Abstract

To find interesting items in genome-wide association studies or next generation sequencing data, a crucial point is to design powerful false discovery rate (FDR) controlling procedures that suitably combine discrete tests (typically binomial or Fisher tests). In particular, recent research has been striving for appropriate modifications of the classical Benjamini-Hochberg (BH) step-up procedure that accommodate discreteness. However, despite an important number of attempts, these procedures did not come with theoretical guarantees. The present paper contributes to fill the gap: it presents new modifications of the BH procedure that incorporate the discrete structure of the data and provably control the FDR for any fixed number of null hypotheses (under independence). Markedly, our FDR controlling methodology allows to incorporate simultaneously the discreteness and the quantity of signal of the data (corresponding therefore to a so-called $\pi\_0$-adaptive procedure). The power advantage of the new methods is demonstrated in a numerical experiment and for some appropriate real data sets.

Figures

Tables

Table 1. Table 1: Number of rejections (discoveries) for the pharmacovigilance and Arabidopsis methylation data.

Data set	[BH]	[DBH-SU]	[Heyse]	[A-DBH-SU]	[A-DBH-SD]
Pharmacovigilance	24	27	27	27	27
Arabidopsis methylation	2097	2358	2379	2446	2453

Table 2. Table 2: Average power of FDR procedures for N = 25 𝑁 25 N=25 (see Section  5 ).

$m$	$m_3$	$m_1$	$q$	[BH]	[Heyse]	[DBH-SU]	[A-DBH-SU]	[A-DBH-SD]
800	80	144	0.15	0.0000	0.0004	0.0003	0.0003	0.0004
		144	0.25	0.0004	0.0197	0.0177	0.0177	0.0135
		144	0.4	0.0803	0.4425	0.4247	0.4247	0.4130
		360	0.15	0.0000	0.0007	0.0006	0.0006	0.0007
		360	0.25	0.0004	0.0244	0.0209	0.0209	0.0153
		360	0.4	0.0803	0.4529	0.4509	0.4509	0.4487
		576	0.15	0.0000	0.0009	0.0007	0.0007	0.0008
		576	0.25	0.0004	0.0343	0.0259	0.0259	0.0231
		576	0.4	0.0803	0.5367	0.4741	0.4741	0.4999
	240	112	0.15	0.0000	0.0003	0.0003	0.0003	0.0002
		112	0.25	0.0005	0.0276	0.0249	0.0249	0.0157
		112	0.4	0.2148	0.5365	0.5012	0.5012	0.4951
		280	0.15	0.0000	0.0003	0.0003	0.0003	0.0002
		280	0.25	0.0005	0.0315	0.0272	0.0272	0.0175
		280	0.4	0.2147	0.5758	0.5536	0.5536	0.5495
		448	0.15	0.0000	0.0005	0.0003	0.0003	0.0004
		448	0.25	0.0005	0.0372	0.0308	0.0308	0.0207
		448	0.4	0.2145	0.5920	0.5741	0.5741	0.5775
	640	32	0.15	0.0000	0.0002	0.0002	0.0002	0.0002
		32	0.25	0.0010	0.0378	0.0341	0.0341	0.0174
		32	0.4	0.4243	0.6174	0.5955	0.6828	0.6621
		80	0.15	0.0000	0.0002	0.0002	0.0002	0.0002
		80	0.25	0.0010	0.0388	0.0347	0.0347	0.0179
		80	0.4	0.4242	0.6282	0.6128	0.6841	0.6638
		128	0.15	0.0000	0.0002	0.0002	0.0002	0.0002
		128	0.25	0.0010	0.0400	0.0354	0.0354	0.0183
		128	0.4	0.4240	0.6353	0.6265	0.6854	0.6656
2000	200	360	0.15	0.0000	0.0002	0.0002	0.0002	0.0002
		360	0.25	0.0001	0.0156	0.0142	0.0142	0.0100
		360	0.4	0.0730	0.4486	0.4317	0.4317	0.4197
		900	0.15	0.0000	0.0002	0.0002	0.0002	0.0002
		900	0.25	0.0001	0.0192	0.0166	0.0166	0.0125
		900	0.4	0.0730	0.4517	0.4511	0.4511	0.4509
		1440	0.15	0.0000	0.0003	0.0002	0.0002	0.0002
		1440	0.25	0.0001	0.0286	0.0211	0.0211	0.0165
		1440	0.4	0.0730	0.5402	0.4684	0.4684	0.4984
	600	280	0.15	0.0000	0.0002	0.0002	0.0002	0.0002
		280	0.25	0.0001	0.0239	0.0213	0.0213	0.0115
		280	0.4	0.2058	0.5350	0.4988	0.4988	0.4960
		700	0.15	0.0000	0.0002	0.0002	0.0002	0.0002
		700	0.25	0.0001	0.0290	0.0239	0.0239	0.0132
		700	0.4	0.2058	0.5750	0.5590	0.5590	0.5516
		1120	0.15	0.0000	0.0002	0.0002	0.0002	0.0002
		1120	0.25	0.0001	0.0350	0.0283	0.0283	0.0157
		1120	0.4	0.2057	0.5908	0.5739	0.5739	0.5761
	1600	80	0.15	0.0000	0.0001	0.0001	0.0001	0.0000
		80	0.25	0.0003	0.0379	0.0342	0.0342	0.0126
		80	0.4	0.4223	0.6196	0.5928	0.6860	0.6591
		200	0.15	0.0000	0.0001	0.0001	0.0001	0.0000
		200	0.25	0.0003	0.0387	0.0351	0.0351	0.0131
		200	0.4	0.4222	0.6281	0.6152	0.6869	0.6602
		320	0.15	0.0000	0.0001	0.0001	0.0001	0.0000
		320	0.25	0.0003	0.0396	0.0360	0.0360	0.0137
		320	0.4	0.4220	0.6327	0.6278	0.6877	0.6617

Equations

$$\overline{F}(t)=\frac{1}{m}\sum _{i=1}^m{F}_i(t),t\in [0,1],$$

$${\overline{F}}_{\text{SU}}(t)=\frac{1}{m}\sum _{i=1}^m\frac{F_i\left(t\right)}{1-F_i\left({\tau }_m\right)};{\overline{F}}_{\text{SD}}(t)=\frac{1}{m}\sum _{i=1}^m\frac{F_i\left(t\right)}{1-F_i\left(t\right)},t\in [0,1],$$

$$F_i(t)=\underset{P\in \mathcal{P}:i\in {\mathcal{H}}_0(P)}{\sup }{\mathbf{P}}_{X\sim P}(p_i(X)\le t),t\in [0,1],1\le i\le m$$

$$F_i(t)\le t,t\in [0,1],1\le i\le m,$$

$$F_i(t)=t,t\in [0,1],1\le i\le m.$$

$$F_i(t)=t,t\in {\mathcal{A}}_i,1\le i\le m.$$

$$\text{FDR}(R,P)={\mathbf{E}}_{X\sim P}\left[\frac{{\sum }_{i\in {\mathcal{H}}_0(P)}1\{p_i\le \hat{t}(X)\}}{1\vee {\sum }_{i=1}^{m}1\{p_i\le \hat{t}(X)\}}\right],P\in \mathcal{P}.$$

$${\tau }_k=\left((1-\lambda )\frac{\alpha k}{m-k+1}\right)\wedge \lambda ,1\le k\le m;$$

$${\tau }_k=\frac{\alpha k}{m-(1-\alpha )k+1},1\le k\le m;$$

$${\tau }_k=\max \{t\in \mathcal{A}:\overline{F}(t)\le \alpha k/m\},1\le k\le m;$$

$${\tau }_m$$

$${\tau }_k$$

$${\tau }_k=\max \left\{t\in \mathcal{A}:\frac{1}{m}\sum _{i=1}^m\frac{F_i\left(t\right)}{1-F_i\left(t\right)}\le \alpha k/m\right\},1\le k\le m.$$

$${\tau }_k$$

$${\tau }_k=\max \left\{t\in \mathcal{A}:{\left(\frac{F\left(t\right)}{1-F\left(t\right)}\right)}_{(1)}+\cdots +{\left(\frac{F\left(t\right)}{1-F\left(t\right)}\right)}_{(m-k+1)}\le \alpha k\right\},1\le k\le m,$$

$$\text{FDR}(\mathbf{SU}(\tau ),P)$$

$$\text{FDR}(\mathbf{SD}(\tau ),P)$$

$$\text{FDR}(\mathbf{SU}(\tau ),P)$$

$$\text{FDR}(\mathbf{SD}(\tau ),P)$$

$$\text{FDR}(\mathbf{SU}(\tau ),P)$$

$$\text{FDR}(\mathbf{SD}(\tau ),P)$$

$$\underset{1\le k\le m}{\max }\frac{m{\tau }_k}{3k}\left(\frac{1}{1-{\tau }_k}+\frac{2}{1-2{\tau }_k}\right),$$

$$\Psi (\lambda )=\min \left(\lambda ,\underset{1\le k\le m}{\max }\left(\frac{1}{k}\sum _{i=1}^m\frac{F_i\left(\lambda k/m\right)}{1-F_i\left(\lambda \right)}\right)\right).$$

$${\tau }_k$$

$$\text{FDR}(\mathbf{SU}(\tau ))$$

$$\sum _{i\in {\mathcal{H}}_0}\mathbf{E}\left(\frac{F_i\left({\tau }_{{\hat{k}}^{\prime -i}+1}\right)}{{\hat{k}}^{\prime -i}+1}\right)\le \sum _{i=1}^m\mathbf{E}\left(\frac{F_i\left({\tau }_{{\hat{k}}^{\prime -i}+1}\right)}{{\hat{k}}^{\prime -i}+1}\right)\le \sum _{i=1}^m\underset{1\le k\le m}{\max }\frac{F_i({\tau }_k)}{k}.$$

$$\sum _{i\in {\mathcal{H}}_0}\mathbf{E}\left(\frac{F_i\left({\tau }_{{\hat{k}}^{\prime -i}+1}\right)}{{\hat{k}}^{\prime -i}+1}\right)$$

$$=\sum _{i\in {\mathcal{H}}_0}\mathbf{E}\left(\frac{F_i\left({\tau }_{{\hat{k}}^{\prime }+1}\right)}{1-F_i\left({\tau }_m\right)}\frac{1\{p_i>{\tau }_m\}}{{\hat{k}}^{\prime }+1}1\{{\hat{k}}^{\prime }+1\le m\}\right),$$

$$\mathbf{E}\left(\sum _{i=1}^m\frac{F_i\left({\tau }_{{\hat{k}}^{\prime }+1}\right)}{1-F_i\left({\tau }_m\right)}\frac{1\{p_i>{\tau }_{{\hat{k}}^{\prime }+1}\}}{{\hat{k}}^{\prime }+1}1\{{\hat{k}}^{\prime }+1\le m\}\right)\le \underset{0\le k\le m-1}{\max }\underset{\begin{array}{c}A\subset \{1,\dots ,m\}\\ |A|=m-k\end{array}}{\max }\sum _{i\in A}\frac{F_i\left({\tau }_{k+1}\right)}{1-F_i\left({\tau }_m\right)}\frac{1}{k+1},$$

$$\text{FDR}(\mathbf{SD}(\tau ))=\sum _{i\in {\mathcal{H}}_0}\mathbf{E}\left(\frac{1\{p_i\le {\tau }_{\tilde{k}}\}}{\tilde{k}}\right)=\sum _{i\in {\mathcal{H}}_0}\mathbf{E}\left(\frac{1\{p_i\le {\tau }_{{\tilde{k}}^{-i}+1}\}}{{\tilde{k}}^{\prime -i}+1}\right)\le \sum _{i\in {\mathcal{H}}_0}\mathbf{E}\left(\frac{F_i\left({\tau }_{{\tilde{k}}^{-i}+1}\right)}{{\tilde{k}}^{-i}+1}\right),$$

$$\sum _{i\in {\mathcal{H}}_0}\mathbf{E}\left(\frac{F_i\left({\tau }_{{\tilde{k}}^{-i}+1}\right)}{{\tilde{k}}^{-i}+1}\right)$$

$$\le \mathbf{E}\left(\sum _{i=1}^m\frac{F_i\left({\tau }_{\tilde{k}+1}\right)}{1-F_i\left({\tau }_{\tilde{k}+1}\right)}\frac{1\{p_i>{\tau }_{\tilde{k}+1}\}}{\tilde{k}+1}1\{\tilde{k}+1\le m\}\right)$$

$$\le \mathbf{E}\left(\underset{0\le k\le m-1}{\max }\underset{\begin{array}{c}A\subset \{1,\dots ,m\}\\ |A|=m-k\end{array}}{\max }\sum _{i\in A}\frac{F_i\left({\tau }_{k+1}\right)}{1-F_i\left({\tau }_{k+1}\right)}\frac{1}{k+1}\right),$$

$$\{1\le i\le m:p_i\le {\tau }_k\}$$

Full text

false discovery rate,

discrete hypothesis testing,

type I error rate control,

adaptive procedure,

step-up algorithm,

step-down algorithm,

and

and

Multiple testing procedures are now routinely used to find significant items in massive and complex data. An important focus has been given to methods controlling the false discovery rate (FDR) because this scalable type I error rate “survives” to high dimension. Since the original procedure of Benjamini and Hochberg (1995), much effort has been undertaken to design FDR controlling procedures that adapt to various underlying structures of the data, such as the quantity of signal, the signal strength and the dependencies, among others.

In this work, we deal with adaptation to discrete data. This type of data arises in many relevant applications, in particular when data are represented by counts. Examples can be found in clinical studies (see e.g. Westfall and Wolfinger (1997)), genome-wide association studies (GWAS) (see e.g. Dickhaus et al. (2012)) and next generation sequencing data (NGS) (see e.g. Chen and Doerge (2015b)). It is well known (see e.g. Westfall and Wolfinger (1997)) that using discrete test statistics can generate a severe power loss, already at the stage of the single tests. A consequence is that using “blindly” the BH procedure with discrete $p$-values will control the FDR in a too conservative manner. Therefore, more powerful procedures that avoid this conservatism are much sought after in applications, see for instance Karp et al. (2016), van den Broek et al. (2015) and Dickhaus et al. (2012).

In the literature, building multiple testing procedures that take into account the discreteness of the test statistics has a long history that can be traced back to Tukey and Mantel (1980). Some null hypotheses can be a priori excluded from the study because the corresponding tests are unable to produce sufficiently small $p$-values. This results in a multiplicity reduction that should increase the power. While this idea has been exploited in Tarone (1990) and in a more general manner in Westfall and Wolfinger (1997) for family-wise error rate, an attempt has been made for FDR later in Gilbert (2005). More recently, Heyse (2011) has proposed a more powerful solution, relying on the following averaged cumulative distribution function (c.d.f.):

[TABLE]

where each $F_{i}$ corresponds to the c.d.f. of the $i$-th test $p$-value under the null hypothesis. To illustrate the potential benefit of using $\overline{F}$, Figure 1 displays this function for the pharmacovigilance data from Heller and Gur (2011) (see Section 4 for more details). The critical values of the Heyse procedure can be obtained by inverting $\overline{F}$ at the values $\alpha k/m,1\leq k\leq m$. Thus, the smaller the $\overline{F}$-values, the larger the critical values. Here, Heyse critical values improve the BH critical values roughly by a factor 3, thereby yielding a potentially strong rejection enhancement. Furthermore, since the functions $F_{i}$’s are known, so is $\overline{F}$. Hence, the user has a good prior idea of the improvements reachable by this discrete approach. Unfortunately, the Heyse procedure does not rigorously control the FDR in general; counter-examples are provided in Heller and Gur (2011) and Döhler (2016).

Meanwhile, different solutions have been explored by modifying directly the $p$-values, either by randomization (see Habiger (2015) and references therein), or by shrinking them to build so-called midP-values (see Heller and Gur (2011) and references therein). Other approaches incorporate discreteness to obtain less conservative FDR estimates, see, e.g., Pounds and Cheng (2006), or by combining grouping and weighting approaches, see Chen and Doerge (2015b).

Overall, although many new procedures have been proposed in the literature, only few of them have been proved to achieve a rigorous FDR control under standard conditions, especially in the finite sample case. To the best of our knowledge, we can only refer to the discretised version of the procedure of Benjamini and Liu (1999) introduced by Heller and Gur (2011) and to the asymptotic work of Ferreira (2007). Chen and Doerge (2015b) sum up the status quo by noting that ’…how to derive better FDR procedures in the discrete paradigm remains an urgent but still unresolved problem.’ This paper offers a solution by presenting new procedures that achieve both theoretical validity and good practical performance.

The paper is organised as follows: after having precisely defined the setting in Section 1, we introduce in Section 2 new procedures relying on the following modifications of the $\overline{F}$ function

[TABLE]

(with the convention $1/0=+\infty$), where an appropriate choice of $\tau_{m}$ is made. To feel how light these modifications are, Figure 1 displays these functions and shows they are very close to the original $\overline{F}$ for small values of $t$. In addition, we also introduce more powerful “adaptive” versions, meaning that the derived critical values are designed in a way that “implicitly estimate” the overall proportion of true null hypotheses and thus may outperform the original Heyse procedure. Next, in Section 3, we establish rigorous FDR control of the corresponding non-adaptive and adaptive procedures under standard conditions. Our proofs rely on new bounds on FDR that generalise some prominent results of the multiple testing literature. These bounds are the main mathematical contributions of the paper and are interesting in their own right, beyond the discrete setting. Also, to explore in detail the improvement of our procedures, we analyse both real and simulated data in Sections 4 and 5. Finally, while the proofs are given in appendix (together with some additional procedures), complementary results are provided in Appendix C.

1 Preliminaries

1.1 General model

Let us observe a random variable $X$, defined on a probabilistic space and valued in an observation space $(\mathcal{X},\mathfrak{X})$. We consider a set $\mathcal{P}$ of possible distributions for the distribution of $X$ and we denote the true one by $P$. We assume that $m$ null hypotheses $H_{0,i}$, $1\leq i\leq m$, are available for $P$ and we denote the corresponding set of true null hypotheses by \mathcal{H}_{0}(P)=\{1\leq i\leq m\>:\>\mbox{H_{0,i}$is satisfied by$P}\}. We also denote by $m_{0}(P)=|\mathcal{H}_{0}(P)|$ the number of true nulls.

We assume that the user has at hand a set of $p$-values to test each null, that is, a set of random variables $\{p_{i}(X),1\leq i\leq m\}$, valued in $[0,1]$. Throughout the paper, we also make the important (but classical) assumption that the $p$-values $p_{i}(X)$, $1\leq i\leq m$, are mutually independent.

Now, we denote $\mathcal{F}=\{F_{i},1\leq i\leq m\}$, where for each $i\in\{1,\dots,m\}$, the function

[TABLE]

is assumed to be known. Note that we necessarily have $F_{i}(\cdot)$ non decreasing, $F_{i}(t)\in[0,1]$, $F_{i}(1)=1$ and we add the technical condition $F_{i}(0)=0$. Loosely, each $F_{i}$ corresponds to the cumulative distribution of $p_{i}$ under the null. Above, we have taken the supremum to cover the case where the null hypothesis is composite: in that situation, each $F_{i}$ is adjusted according to the least favorable configuration within the null $H_{0,i}$.

Here are some conditions on $\mathcal{F}$ that will be useful to compare some of the studied procedures (these conditions are not assumed in our results unless explicitly mentioned):

[TABLE]

Condition (2) ensures that the $p$-values have marginals stochastically lower-bounded by a uniform variable under the null, called a super-uniform distribution in the sequel. This is the classical setting which is used in most of the work dealing with FDR controlling theory, see, e.g., Benjamini and Hochberg (1995). Condition (3) is more restrictive: if each null hypothesis is a singleton, it is equivalent to the $p$-values having uniform marginals under the null.

1.2 Discrete and continuous modelling

In order to describe the overall support of $p$-value distributions we assume one of the two following situations to be at hand throughout the paper (except in Section 3 which is written in a more general manner):

• •

Continuous case: for all $i\in\{1,\dots,m\}$, $F_{i}$ is continuous. In that case, we let $\mathcal{A}_{i}=[0,1]$, $1\leq i\leq m$ and $\mathcal{A}=\cup_{i=1}^{m}\mathcal{A}_{i}=[0,1]$ is the overall $p$-value support.

• •

Discrete case: each $p$-value $p_{i}$ (both under the null and alternative) takes values in some finite set $\mathcal{A}_{i}=\{a_{i,k},0\leq k\leq K_{i}\}$, where $K_{i}\geq 0$ and $(a_{i,k})_{0\leq k\leq K_{i}}\in[0,1]^{K_{i}+1}$ is an increasing sequence (with $a_{i,0}=0$, $a_{i,K_{i}}=1$). We denote $\mathcal{A}=\cup_{i=1}^{m}\mathcal{A}_{i}$ the overall $p$-value support.

The continuous setting is typically valid in situations where the $p$-values are calibrated from test statistics having a continuous distribution under the null. In this situation, (3) is often satisfied. The discrete setting typically arises in situations where the $p$-values are calibrated from test statistics having a finitely supported distribution under the null. In this situation, we generally have that (3) holds true only on the support of $F_{i}$, that is,

[TABLE]

In the discrete framework, let us underline that while (4) will typically hold, the equality $F_{i}(t)=t$, $t\in\mathcal{A}$ will fail in general because $\mathcal{A}$ contains points of $\mathcal{A}_{j}$ for $j\neq i$. Then $\overline{F}(t)$ defined by (1) will be smaller than $t$ in general (see Figure 1), which is exactly the property that we want to exploit in this paper.

To illustrate the above framework, we provide below two simple examples (for more advanced examples, see for instance Chen and Doerge (2015b)).

Example 1.1* (Gaussian testing).*

Observe $X=(X_{i})_{1\leq i\leq m}$ with independent coordinates and marginals $X_{i}\sim\mathcal{N}(\mu_{i},1)$, where $\mu_{i}\in\mathbb{R}$ is the parameter of interest, $1\leq i\leq m$. In that situation, a possible hypothesis testing problem is to consider the nulls $H_{0,i}:$ “$\mu_{i}\leq 0$” against $H_{1,i}:$ “$\mu_{i}>0$”. Then $p_{i}(X)=1-\Phi(X_{i})$, $1\leq i\leq m$, is a family of $p$-values satisfying (3) (where $\Phi$ denotes the c.d.f. of a standard Gaussian variable).

Example 1.2* (Binomial testing).*

Observe $X=(X_{i})_{1\leq i\leq m}$ with independent coordinates and marginals $X_{i}\sim\mathcal{B}(n_{i},\theta_{i})$, where $n_{i}\geq 1$ is known and $\theta_{i}\in(0,1)$ is the parameter of interest, $1\leq i\leq m$. In that situation, a possible hypothesis testing problem is to consider the nulls $H_{0,i}:$ “$\theta_{i}\leq 1/2$” against $H_{1,i}:$ “$\theta_{i}>1/2$”. Then $p_{i}(X)=T_{i}(X_{i})$, $1\leq i\leq m$, define a family of $p$-values where $T_{i}(x)=2^{-n_{i}}\sum_{j=0}^{x}{n_{i}\choose j}$ is the upper-tail distribution function of a binomial distribution of parameters $(n_{i},1/2)$. The support of the $p$-values under the null and alternative is covered by letting $K_{i}=n_{i}+1$ and $a_{i,k}=2^{-n_{i}}\sum_{j=0}^{k-1}{n_{i}\choose j}$, $1\leq k\leq K_{i}$. We merely check in that case that (3) is violated while (2) and (4) hold.

1.3 Step-wise procedures

First define a critical value sequence as any nondecreasing sequence ${\tau}=(\tau_{k})_{1\leq k\leq m}\in[0,1]^{m}$ (with $\tau_{0}=0$ by convention).

The step-up procedure of critical value sequence $\tau$, denoted by ${\mathbf{SU}}(\tau)$, rejects the $i$-th hypothesis if $p_{i}\leq\tau_{\hat{k}}$, with $\widehat{k}=\max\{k\in\{0,1,...,m\}\>:\>p_{(k)}\leq{\tau}_{k}\},$ where $p_{(1)}\leq p_{(2)}\leq...\leq p_{(m)}$ denote the ordered $p$-values (with the convention $p_{(0)}=0$).

The step-down procedure of critical value sequence $\tau$, denoted by ${\mathbf{SD}}(\tau)$, rejects the $i$-th hypothesis if $p_{i}\leq\tau_{\tilde{k}}$, with $\widetilde{k}=\max\{k\in\{0,1,...,m\}\>:\>\forall k^{\prime}\leq k,\>p_{(k^{\prime})}\leq{\tau}_{k^{\prime}}\}.$ It is straightforward to check that, for the same set of critical values, the step-up version always rejects more hypotheses than the step-down version. More comments and illustrations on step-wise procedures can be found in Blanchard et al. (2014) and Dickhaus (2014), among others.

1.4 False discovery rate

We measure the quantity of false positives of a step-up (resp. step-down) procedure by using the false discovery rate (FDR), introduced and popularised by Benjamini and Hochberg (1995), which is defined as the averaged proportion of errors among the rejected hypotheses. More formally, for some procedure $R$ rejecting the $i$-th hypothesis if $p_{i}\leq\hat{t}(X)$ (for some threshold $\hat{t}(X)$), we let

[TABLE]

The main contribution of this work is to propose procedures that control the FDR at a prescribed level $\alpha$ and that incorporate the knowledge of the $F_{i}$’s in a way that increases the number of discoveries.

2 Procedures

In this section we briefly review some existing methods for FDR control and introduce our new procedures.

2.1 Existing methods

We use the following methods as starting points for constructing new procedures.

[BH]: the seminal procedure proposed in Benjamini and Hochberg (1995), corresponding to the step-up procedure ${\mathbf{SU}}(\tau)$, with critical values $\tau_{k}=\alpha k/m$, $1\leq k\leq m$;

• -

[BR-$\lambda$]: an adaptive version of the BH procedure that was proposed in Blanchard and Roquain (2009), corresponding to the step-up procedure ${\mathbf{SU}}(\tau)$, with critical values

[TABLE]

• -

[GBS]: an adaptive version of the BH procedure that has been proposed in Gavrilov et al. (2009), corresponding to the step-down procedure ${\mathbf{SD}}(\tau)$, with critical values

[TABLE]

• -

[Heyse]: the step-up procedure ${\mathbf{SU}}(\tau)$ using critical values given by

[TABLE]

where $\overline{F}$ is defined by (1). This procedure was proposed in Heyse (2011).

The rationale behind the critical values of [BR-$\lambda$] and [GBS] is that they are intended to mimic the oracle critical values $\tau_{k}=\alpha k/m_{0}(P)$, $1\leq k\leq m$, which are less conservative than those of [BH] when $m_{0}(P)/m$ is not close to $1$, see, e.g., Benjamini et al. (2006); Blanchard and Roquain (2009) for more details on adaptive procedures.

Let us now comment on [Heyse]. First, in the continuous setting where (2) holds, $\overline{F}(t)\leq t$, $t\in[0,1]$, and thus the critical values given by (8) satisfy $\tau_{k}\geq\alpha k/m$, $1\leq k\leq m$, which means that [Heyse] rejects at least as many hypotheses as [BH]. When (3) additionally holds, we have $\overline{F}(t)=t$, $t\in[0,1]$, and the two critical value sequences are the same. Second, in the discrete setting where (2) holds, $\mathcal{A}$ is finite and $\tau_{k}$ is not necessarily greater than $\alpha k/m$ anymore. However, [Heyse] is also less conservative (or equal) than [BH] in the latter case, as stated in the following result (proved in Appendix C for completeness).

Lemma 2.1**.**

Consider the model of Section 1.1 assuming (2), both in the continuous and discrete setting described in Section 1.2. Then the set of nulls rejected by [Heyse] is larger than the one of [BH] (almost surely). Furthermore, under (4), these two rejection sets are equal (almost surely) if $F_{i}=F_{j}$ for all $i\neq j$.

The equality case of Lemma 2.1 was provided in Proposition 2.3 of Heller and Gur (2011), who presented it as a limitation of Heyse procedure in the discrete case. However, we argue that the condition $F_{i}=F_{j}$ for all $i\neq j$ is a somehow extreme configuration which is rarely met in practice (in the discrete case). More typically, the $F_{i}$’s have an heterogeneous structure implying that $\overline{F}(t)$ is smaller than $t$ (see Figure 1). This entails that [Heyse] can substantially improve [BH] (see Figure 2).

While [Heyse] incorporates the knowledge of the $F_{i}$’s in a natural way (see also Remark 2.2 below), it is not correctly calibrated for a rigorous FDR control: as shown in Heller and Gur (2011); Döhler (2016), it fails to control the FDR in general. We propose suitable modifications of [Heyse] in the next sections.

Remark 2.2* (Empirical Bayes point of view on the Heyse procedure).*

We claim that [Heyse] corresponds to a suitable empirical Bayes procedure. To see this, consider the “binomial example” of Section 1.2, but assume now that the counts $n_{1},\dots,n_{m}$ are observed from a sample $N_{1},\dots,N_{m}$ i.i.d. of a priori distribution $\nu$. Unconditionally, the $p$-values $p_{i}$, $i\in\mathcal{H}_{0}$, are thus i.i.d. with c.d.f. $\bar{F}_{0}=\sum_{n\geq 0}\nu(\{n\})F_{0,n},$ where $F_{0,n}$ is the c.d.f. jumping at each $x_{k,n}=2^{-n}\sum_{j=0}^{k-1}{n\choose j}$ with $F_{0,n}(x_{k,n})=x_{k,n}$, $1\leq k\leq n+1$. This suggests to normalise the $p$-values $p_{i}$ as $\bar{F}_{0}(p_{i})$ which leads to the step-up procedure with critical values $\tau_{k}=\max\{t\>:\>\bar{F}_{0}(t)\leq\alpha k/m\}$. Following an empirical Bayes approach, the prior $\nu$ can be estimated by $\hat{\nu}(\{n\})=m^{-1}\sum_{i=1}^{m}\mathbf{1}_{\{N_{i}=n\}}$, which gives rise to the estimator of $\bar{F}_{0}$ defined by $\hat{\bar{F}}_{0}=\sum_{n\geq 0}\hat{\nu}(\{n\})F_{0,n}=m^{-1}\sum_{i=1}^{m}F_{0,N_{i}},$ which is equal to $\overline{F}$ given by (1). Hence, the corresponding (empirical Bayes) step-up procedure reduces to [Heyse].

2.2 Two new methods

We now present two procedures that aim at correcting [Heyse] :

[DBH-SU]: the step-up procedure ${\mathbf{SU}}(\tau)$ using the critical values defined in the following way:

[TABLE]

• -

[DBH-SD]: the step-down procedure ${\mathbf{SD}}(\tau)$ using the critical values defined in the following way :

[TABLE]

[DBH-SU] can be seen as a correction of [Heyse]: the correction term in the critical values (10) lies in the additional denominator $1-F_{i}\left(\tau_{m}\right)$. A consequence is that [DBH-SU] can be more conservative than [BH]. However, the magnitude of this phenomenon is always small, as the next lemma shows (proved in Appendix C for completeness).

Lemma 2.3**.**

Under the conditions of Lemma 2.1, the set of nulls rejected by [DBH-SU] contains the one of [BH] taken at level $\alpha/(1+\alpha)$ (almost surely).

For [DBH-SD], the following result can be established.

Lemma 2.4**.**

Under the conditions of Lemma 2.1, the set of nulls rejected by [DBH-SD] contains the one of the step-down procedure with critical values $(\alpha k/m)/(1+\alpha k/m)$, $1\leq k\leq m$ (almost surely).

From (10) and (11) it is clear that the critical values of [DBH-SD] are always at least as large as those for [DBH-SU]. However, since the step-up direction is more powerful than the step-down direction (see Section 1.3) neither of the two generally dominates the other one.

Remark 2.5*.*

We may ask whether we can construct a uniform improvement of [BH] that incorporate the $F_{i}$’s. There is indeed such a procedure (see procedure [RBH] in Appendix A.1 for more details). However, the improvement brought by the $F_{i}$’s information is less substantial than for [DBH-SU], so we have chosen to omit [RBH] from the main stream of the paper.

2.3 Adaptive versions

In this section, we define adaptive versions of [DBH-SU] and [DBH-SD] in the following way:

[A-DBH-SU]: the step-up procedure ${\mathbf{SU}}(\tau)$ using the critical values defined in the following way: $\tau_{m}$ as in (9) and for $1\leq k\leq m-1$,

[TABLE]

where each $\left(\frac{F\left(t\right)}{1-F\left(\tau_{m}\right)}\right)_{(j)}$ denotes the $j$-th largest elements of the set $\left\{\frac{F_{i}\left(t\right)}{1-F_{i}\left(\tau_{m}\right)},1\leq i\leq m\right\}$.

• -

[A-DBH-SD]: the step-down procedure ${\mathbf{SD}}(\tau)$ using the critical values defined in the following way :

[TABLE]

where each $\left(\frac{F\left(t\right)}{1-F\left(t\right)}\right)_{(j)}$ denotes the $j$-th largest elements of the set $\left\{\frac{F_{i}\left(t\right)}{1-F_{i}\left(t\right)},1\leq i\leq m\right\}$.

Note that the critical values of [A-DBH-SU] and [A-DBH-SD] are clearly larger than or equal to those of their non-adaptive counterparts [DBH-SU] and [DBH-SD], respectively. This means that the adaptive versions are always less conservative.

The following result establishes a connection of the adaptive procedures to the [BR-$\lambda$] and [GBS] procedures (proved in Appendix C for completeness).

Lemma 2.6**.**

Under the conditions of Lemma 2.1, the following holds:

• (i)

the set of nulls rejected by [A-DBH-SU] contains the one of [BR-$\lambda$] (almost surely), where $\lambda$ is taken equal to (9);

• (ii)

the set of nulls rejected by [A-DBH-SD] contains the one of [GBS] (almost surely);

The above lemma ensures that the user can incorporate the knowledge of the $F_{i}$’s in adaptive procedures with a “no loss” guarantee with respect to [BR] and [GBS]. This is a somehow striking fact, coming loosely from a “fortunate marriage” between the proof technics of discreteness theory and adaptation theory.

Remark 2.7*.*

We may ask whether we can build a procedure that is a uniform improvement of [BR-$\lambda$], for any fixed value of $\lambda\in(0,1)$. We propose a solution in Appendix A.2, called [DBR-$\lambda$]. It does not improve uniformly [DBH-SU], but is an interesting variant of [A-DBH-SU].

3 New FDR bounds

In this section, we present new FDR bounds which are the main mathematical contributions of this paper and that are of independent interest. They generalise some classical bounds from super-uniform null distributions to arbitrary heterogeneous (not necessarily discrete) null distributions, and immediately yield FDR control of our new procedures.

3.1 Results

First, remember that the model of Section 1.1 basically only assume independence between the $p$-values (and not super-uniformity of the null distribution). The following result holds.

Theorem 3.1**.**

In the model of Section 1.1, for any critical values $\tau_{k}$, $1\leq k\leq m$ and for all $P\in\mathcal{P}$, we have

[TABLE]

The proof of Theorem 3.1 is deferred to Appendix B. It combines several techniques: the first tool is an expression of the FDR introduced by Ferreira (2007) (step-up case) and Roquain and Villers (2011) (step-down case). A second idea comes from the work Blanchard and Roquain (2009) (step-up case) and Gavrilov et al. (2009) (step-down case), which introduced a new term (here, the denominator $(1-F_{i}(\cdot))$) to make the adaptive argument works fine. Finally, another inspiration is the study of Roquain and van de Wiel (2009) and Döhler (2016) that allowed to deal with heterogeneous FDR thresholding. Let us underline that the obtained proof is especially concise, which means that these different techniques fit together perfectly well, which is perhaps surprising at first glance, see Appendix B.3.

Next, let us note that taking the maximum over the subset $A$ in (14) and (15) allows us to adapt to the unknown number of true null hypotheses: loosely, if $k-1$ is the number of rejections, $A$ corresponds to the acceptation set (hence of cardinality $m-k+1$), which “estimates” $\mathcal{H}_{0}$ and thus the sums in (14) and (15) are indexed by a set “close” to the unknown set $\mathcal{H}_{0}$. Taking the maximum then corresponds to account for the least favorable possible $\mathcal{H}_{0}$.

Finally, let us underline again that the above bounds do not use the super-uniformity of the $F_{i}$’s which makes them quite general and flexible tools. As a case in point, consider mid-$p$-values which were introduced by Lancaster (1961) and are sometimes used for analysing discrete data (see e.g. Karp et al. (2016)). These $p$-values are no longer super-uniform under the null hypotheses, however our theorem can accomodate such distributions in a natural way to still yield valid FDR controlling procedures. In addition, note that our bounds can be useful outside the discrete setting, when the $F_{i}$’s are continuous but with flat parts, see the (toy) Example 3.3 below.

3.2 Rationale and relation to previous work

Let us now give some intuition behind these bounds by showing how it allows to cover previous work in the literature.

First, assuming the super-uniformity $F_{i}(t)\leq t$ for all $t$ and $i$, then these bounds entail

[TABLE]

which immediately recover the fact that [BH], [BR-$\lambda$] (with $\tau_{m}=\lambda$) and [GBS] all control the FDR at level $\alpha$. To this respect, bounds (16), (17) and (18) encompass Theorem 3.1 of Benjamini and Hochberg (1995), Theorem 9 of Blanchard and Roquain (2009) and Theorem 1.1 of Gavrilov et al. (2009), respectively.

Second, by removing the adaptative part of the bounds, that is, by replacing $A$ by $\{1,\dots,m\}$, we obtain the simpler but more conservative bounds

[TABLE]

where $\overline{F}_{\textnormal{SU}}$ and $\overline{F}_{\textnormal{SD}}$ are defined in the introduction, see also Figure 1. These variants illustrate perhaps more intuitively how the Heyse-type procedures take advantage of the heterogeneous structure: if some of the $F_{i}$’s are really small, they will not contribute much into $\overline{F}_{\textnormal{SU}}$ (or $\overline{F}_{\textnormal{SD}}$), offering some additional room for the other $F_{j}$’s.

Finally, these bounds immediately imply that our new procedures enjoy the desired FDR controlling property.

Corollary 3.2**.**

*In the model of Section 1.1, both in the continuous and discrete setting described in Section 1.2, the procedures [DBH-SU]; [DBH-SD]; [A-DBH-SU]; [A-DBH-SD] all control the FDR at level $\alpha$. *

Example 3.3*.*

Assume that the hypotheses are structured in $3$ non-overlapping groups $S_{1}$, $S_{2}$ and $S_{3}$, each of cardinality $m/3$ (assumed to be an integer). Assume that $F_{i}(x)$ is equal to $x$ if $i\in S_{1}$, [math] if $i\in S_{2}$, and $F_{i}(x)=2x$ ($x\in[0,1/4]$); $1/2$ ($x\in[1/4,3/4]$); $2x-1$ ($x\in[3/4,1]$), if $i\in S_{3}$. Then the bound (20) becomes (for $\tau_{k}\leq 1/4$):

[TABLE]

which entails a new step-down FDR controlling procedure by taking $\tau_{k}$ such that the above expression is equal to $\alpha$. For $\alpha$ small, we get $\tau_{k}\approx\alpha k/m$ which yields to a procedure close to a step-down version of [BH]. It thus controls the FDR even though the super-uniformity of the $p$-values is violated. In particular, this illustrates that our methodology exceeds the scope of discrete tests.

4 Empirical data

To illustrate the performance of FDR-controlling procedures for discrete data, we analyse two benchmark data sets which have also been used in previous publications. In what follows, our main goal is to compare the performance of the new procedures [DBH-SU], [A-DBH-SU] and [A-DBH-SD] to the classical [BH] procedure. As a further benchmark we also include [Heyse] in the analysis. All analyses were performed using the R language for statistical computing (R Core Team, 2016).

4.1 Pharmacovigilance data

This data set is derived from a database for reporting, investigating and monitoring adverse drug reactions due to the Medicines and Healthcare products Regulatory Agency in the United Kingdom. It contains the number of reported cases of amnesia as well as the total number of adverse events reported for each of the $m=2446$ drugs in the database. For more details we refer to Heller and Gur (2011) and to the accompanying R-package ’discreteMTP’ (Heller et al., 2012), which also contains the data. Heller and Gur (2011) investigate the association between reports of amnesia and suspected drugs by performing for each drug a Fisher’s exact test (one-sided) for testing association between the drug and amnesia while adjusting for multiplicity by using several (discrete) FDR procedures.

4.2 Next generation sequencing data

We also revisit the next generation sequencing (NGS) count data analysed by Chen and Doerge (2015b), to which we also refer for more details. More specifically, we reanalyse the methylation data set for cytosines of Arabidopsis in Lister et al. (2008) which is part of the R-package ’fdrDiscreteNull’ (Chen and Doerge, 2015a). This data set contains the counts for a biological entity under two different biological conditions or treatments. Following Chen and Doerge (2015b), $m=7421$ genes whose treatment-wise total counts are positive but row-total counts are no greater than 100 are analysed using the exact binomial test, see Chen and Doerge (2015b).

4.3 Results

Table 1 summarises the number of discoveries for the pharmacovigilance and NGS data when using the respective FDR procedures at level $\alpha=0.05$.

Compared to the classical FDR controlling procedures, the new procedures are able to detect three additional candidates linking amnesia and drugs in the pharmacovigilance data. Note also that for this data, they reject the same number of hypotheses as [Heyse], even though [Heyse] is not correctly calibrated for FDR control. For the Arabidopsis data, the new procedures improve considerably on [BH]. Moreover, there is a clear separation between the adaptive and non-adaptive procedures.

Figure 2 illustrates graphically the data and the critical constants of the involved multiple testing procedures.

In particular, the benefit of taking discreteness into account becomes more apparent: for the pharmacovigilance data, the discrete critical values are considerably (by a factor of $2.5-3.5$) larger than their respective classical counterparts. This leads to more powerful procedures. For the NGS data, we can observe quite clearly that the [DBH-SU] critical constants are dominated by the [A-DBH-SU] constants, as explained in Section 2. This leads to roughly 100 additional rejections. Again, the discrete critical values are considerably larger than their respective classical counterparts. In 2.2 we mentioned that the correction factor $1-F_{i}(\tau_{m})$, introduced for guaranteeing FDR control of [DBH-SU], may lead to a procedure which is more conservative than [BH]. However, Figure 2 shows that – at least for the data sets considered here – this risk is by far compensated by the benefit of taking discreteness adequately into account.

5 Simulation study

We now investigate the power of the procedures from the previous section in a simulation study similar to those described in Gilbert (2005), Heller and Gur (2011) and Döhler (2016). Again, we focus on comparing the performance of the new discrete procedures to [BH].

5.1 Simulated Scenarios

We simulate a two-sample problem in which a vector of $m$ independent binary responses (“adverse events”) is observed for each subject in two groups, where each group consists of $N=25$ subjects. Then, the goal is to simultaneously test the $m$ null hypotheses $H_{0i}:$ “$p_{1i}=p_{2i}$”, $i=1,\ldots,m$, where $p_{1i}$ and $p_{2i}$ are the success probabilities for the $i$th binary response in group 1 and 2, respectively. We take $m=800,2000$ where $m=m_{1}+m_{2}+m_{3}$ and data are generated so that the response is $Bernoulli(0.01)$ at $m_{1}$ positions for both groups, $Bernoulli(0.10)$ at $m_{2}$ positions for both groups and $Bernoulli(0.10)$ at $m_{3}$ positions for group 1 and $Bernoulli(q)$ at $m_{3}$ positions for group 2 where $q=0.15,0.25,0.4$ represents weak, moderate and strong effects respectively. The null hypothesis is true for the $m_{1}$ and $m_{2}$ positions while the alternative hypothesis is true for the $m_{3}$ positions. We also take different configurations for the proportion of false null hypotheses, $m_{3}$ is set to be $10\%$, $30\%$ and $80\%$ of the value of $m$, which represents small, intermediate and large proportion of effects (the proportion of true nulls $\pi_{0}$ is $0.9$, $0.7$, $0.2$, respectively). Then, $m_{1}$ is set to be $20\%$, $50\%$ and $80\%$ of the number of true nulls (that is, $m-m_{3}$) and $m_{2}$ is taken accordingly as $m-m_{1}-m_{3}$.

For each of the 54 possible parameter configurations specified by $m,m_{3},m_{1}$ and $q$, $10000$ Monte Carlo trials are performed, that is, $10000$ data sets are generated and for each data set, an unadjusted two-sided $p$-value from Fisher’s exact test is computed for each of the $m$ positions, and the multiple testing procedures mentioned above are applied at level $\alpha=0.05$. The power of each procedure was estimated as the fraction of the $m_{3}$ false null hypotheses that were rejected, averaged over the $10000$ simulations. For random number generation the R-function rbinom was used. The two-sided $p$-values from Fisher’s exact test were computed using the R-function fisher.test.

5.2 Results

We have computed the (average) power of the five procedures under investigation in all the scenarios (see Table 1 in Appendix C for the full display). For weak and moderate effects, i.e. $q=0.15$ and $q=0.25$, none of the procedure possesses relevant power. For strong effects, the results are summarised in Figure 3. (Since the power of the discrete procedures is slightly increasing in $m_{1}$ for fixed $m_{3}$ and $q$, we present – in order to avoid over-optimism – the configuration with smallest $m_{1}$).

The results are consistent with the findings of the previous section: The new discrete procedures are considerably more powerful than their classical counterparts and perform roughly similarly for small and intermediate proportions of alternatives. When the proportion of alternatives is large, the benefit of using adaptive procedures – especially the [A-DBH-SU] procedure – is clearly visible in Figure 3.

6 Conclusion and discussion

In this paper, we provided new bounds for the FDR of step-up and step-down procedures that use discrete test statistics. This allowed to define a new class of multiple testing procedures that provably control the FDR while they incorporated the discreteness of the tests statistics in a convenient way. We have shown that our approach can be seen as correcting and improving the approach of Heyse (2011): while it ensures a theoretical control, it can also make more rejections when the signal amplitude is strong enough.

Let us also mention that our procedures are fully usable in practice. An R-package is in preparation and it will be presented in a companion paper (which will also deal with computational aspects).

Finally, this paper opens several directions for future research, especially by trying to extend our arguments to other frameworks. For instance, it may be worth to relax the independence requirements. To this respect, we believe that our procedures will inheritate the behavior of BH procedure: while the FDR control is likely to be maintained under “realistic” dependence, formally proving such a result is probably a challenging problem.

Acknowledgements

This work has been supported by the CNRS (PEPS FaSciDo) and the French grant ANR-16-CE40-0019.

Appendix A Additional procedures

A.1 A rescaled BH procedure

The procedure [RBH] (rescaled-BH) is defined as the step-up procedure using the critical values $\tau_{k}=\lambda_{\alpha}k/m$, $1\leq k\leq m$, where $\lambda_{\alpha}=\max\{\lambda\in[0,1]\>:\>\Psi(\lambda_{\alpha})\leq\alpha\}$ for

[TABLE]

The following result is straightforward from Theorem 3.1 (SU part).

Corollary A.1**.**

In the model of section 1.1 with the additional assumption (2), we have $\forall P\in\mathcal{P}$, $\textnormal{FDR}(\textnormal{RBH},P)\leq\alpha$.

Moreover, if $\alpha$ is such that the equality $\Psi(\lambda_{\alpha})=\alpha$ holds true, then $\lambda_{\alpha}\geq\Psi(\lambda_{\alpha})=\alpha$ and [RBH] always dominates [BH] in terms of critical values and therefore rejects at least as many hypotheses.

A.2 A discrete BR procedure

The procedure [DBR-$\lambda$] (discrete BR) is defined as the step-up procedure ${\mathbf{SU}}(\tau)$ using the critical values defined in the following way $\tau_{m}=\max\left\{t\in\mathcal{A}\>:\>t\leq\left((1-\lambda)m\alpha\right)\wedge\lambda\right\}$ (that is, the discretised last critical value of [BR-$\lambda$]) and

[TABLE]

where each $\left(F\left(t\right)\right)_{(j)}$ denotes the $j$-th largest elements of the set $\left\{F_{i}\left(t\right),1\leq i\leq m\right\}$. The following result is straightforward from Theorem 3.1 (SU part).

Corollary A.2**.**

In the model of section 1.1 with the additional assumption (2), we have $\forall P\in\mathcal{P}$, $\textnormal{FDR}(\textnormal{DBR},P)\leq\alpha$. Moreover, the set of nulls rejected by [DBR-$\lambda$] is larger than the one of [BR-$\lambda$] (almost surely), with equality (almost surely) under (4) and $F_{i}=F_{j}$ for all $i\neq j$.

Appendix B Proof of Theorem 3.1

B.1 Lemmas for step-down and step-up procedures

Let us introduce the following modifications of ${\mathbf{SU}}(\tau)$ :

• •

${\mathbf{SU}}^{\prime}(\tau)={\mathbf{SU}}(\tau^{\prime})$ the step-up with $m$ critical values defined by $(\tau^{\prime}_{1},\dots,\tau^{\prime}_{m})=(\tau_{2},\dots,\tau_{m},\tau_{m})$;

• •

for some given index $i\in\{1,\dots,m\}$, ${\mathbf{SU}}^{\prime-i}(\tau)={\mathbf{SU}}(\tau^{\prime-i})$ the step-up with $m-1$ critical values defined by $(\tau^{\prime-i}_{1},\dots,\tau^{\prime-i}_{m-1})=(\tau_{2},\dots,\tau_{m})$ and restricted to the $p$-values of the set $\{p_{j},j\neq i\}$.

The following lemma holds (classical from Ferreira and Zwinderman (2006) and proved in Appendix C):

Lemma B.1**.**

For all $i\in\{1,\dots,m\}$, the following assertions are equivalent: (i) $p_{i}\leq{\tau}_{\hat{k}}$; (ii) $p_{i}\leq{\tau}_{\hat{k}^{\prime-i}+1}$; (iii) $\widehat{k}^{\prime-i}+1=\widehat{k}$, where $\widehat{k}^{\prime-i}$ denotes the number of rejected hypotheses of the procedure ${\mathbf{SU}}^{\prime-i}(\tau)$. Moreover, we have $\{p_{i}>{\tau}_{m}\}\subset\{\widehat{k}^{\prime}=\widehat{k}^{\prime-i}\}$, where $\widehat{k}^{\prime}$ denotes the number of rejected hypotheses of the procedure ${\mathbf{SU}}^{\prime}(\tau)$.

Let us introduce the following modifications of ${\mathbf{SD}}(\tau)$ :

• •

for some given index $i\in\{1,\dots,m\}$, ${\mathbf{SD}}^{-i}(\tau)={\mathbf{SD}}(\tau^{-i})$ the step-down procedure with $m-1$ critical values defined by $(\tau^{-i}_{1},\dots,\tau^{-i}_{m-1})=(\tau_{1},\dots,\tau_{m-1})$ and restricted to the $p$-values of the set $\{p_{j},j\neq i\}$.

• •

for some given index $i\in\{1,\dots,m\}$, ${\mathbf{SD}}^{\prime-i}(\tau)={\mathbf{SD}}(\tau^{\prime-i})$ the step-down procedure with $m-1$ critical values defined by $(\tau^{\prime-i}_{1},\dots,\tau^{\prime-i}_{m-1})=(\tau_{2},\dots,\tau_{m})$ and restricted to the $p$-values of the set $\{p_{j},j\neq i\}$.

The following lemma holds (classical from Ferreira and Zwinderman (2006) and proved in Appendix C):

Lemma B.2**.**

For all $i\in\{1,\dots,m\}$, the following assertions are equivalent: (i) $p_{i}\leq{\tau}_{\tilde{k}}$; (ii) $p_{i}\leq{\tau}_{\tilde{k}+1}$; (iii) $p_{i}\leq{\tau}_{\tilde{k}^{-i}+1}$; (iv) $\tilde{k}^{\prime-i}+1=\tilde{k}$, where $\widetilde{k}^{-i}$ is the number of rejections of ${\mathbf{SD}}^{-i}(\tau)$ and $\widetilde{k}^{\prime-i}$ is the number of rejections of ${\mathbf{SD}}^{\prime-i}(\tau)$. Moreover, we have $\{p_{i}>{\tau}_{\tilde{k}^{-i}+1}\}\subset\{\widetilde{k}=\widetilde{k}^{-i}\}$.

B.2 Proof of Theorem 3.1, step-up part

By using Lemma B.1 $(ii)$ and $(iii)$ and independence, we easily obtain

[TABLE]

where the last expectation is taken only with respect to $(p_{j},j\neq i)$. Now, on the one hand,

[TABLE]

Next, on the other hand, by using again the independence,

[TABLE]

where the latter equality comes from the last assertion of Lemma B.1. Now, since ${\tau}_{\hat{k}^{\prime}+1}\leq\tau_{m}$, we have that the last display is smaller or equal to

[TABLE]

by taking the maximum over all the possible realizations of the set $A=\{1\leq i\leq m\>:\>p_{i}>{\tau}_{\hat{k}^{\prime}+1}\}=\{1\leq i\leq m\>:\>p_{i}>{\tau^{\prime}}_{\hat{k}^{\prime}}\}$ which is the index set corresponding to the non-rejected null hypotheses of $\mathbf{SU}(\tau^{\prime})$ (the latter being by definition of cardinality $m-\hat{k}^{\prime}$). This concludes the proof.

B.3 Proof of Theorem 3.1, step-down part

It is similar to the step-up case, but relies now on Lemma B.2 $(iii)$ and $(iv)$ (and still independence): we have

[TABLE]

which gives the first part of the bound. Next, by using independence and the last assertion of Lemma B.2, we obtain

[TABLE]

because $\{1\leq i\leq m\>:\>p_{i}>{\tau}_{\tilde{k}+1}\}$ is equal to $\{1\leq i\leq m\>:\>p_{i}>{\tau}_{\tilde{k}}\}$, which is the set of non-rejected hypotheses of $\mathbf{SD}(\tau)$. Since $\mathbf{SD}(\tau)$ rejects exactly $\tilde{k}$ hypotheses, the proof is completed.

Appendix C Supplement

C.1 Proofs for lemmas comparing procedures

The lemmas presented here rely on the fact that, there is almost surely no $p$-value in $[0,1]\backslash\mathcal{A}$ (both in the continuous and discrete cases). All symbols “$=$” or “$\subset$” are intended to be valid almost surely in this section.

A result which will be extensively used in the proofs of this section is the following one : for $p$-values valued in the set $\mathcal{A}$, then the step-up procedure with critical values $\tau_{k}$, $1\leq k\leq m$, has the same rejection set as the step-up procedure with critical values $\xi_{k}=\max\left\{t\in\mathcal{A}\>:\>t\leq\tau_{k}\right\}$, $1\leq k\leq m$. This fact comes from the simple following observation : for all $k$,

[TABLE]

The $\xi_{k}$’s are called the “effective” critical values of $\mathbf{SD}(\tau)$ or $\mathbf{SU}(\tau)$ in the sequel.

C.1.1 Proof of Lemma 2.1

The effective critical values of the BH procedure are given by $\xi_{k}=\max\left\{t\in\mathcal{A}\>:\>t\leq\alpha k/m\right\}$, $1\leq k\leq m$. If (2) holds, then $\overline{F}(t)\leq t$ and each $\xi_{k}$ is clearly smaller than the $k$-th critical values of [Heyse]. This implies that the rejection set of [Heyse] is larger than the one of [BH]. Conversely, under (4) and if $F_{i}=F_{j}=\overline{F}$ for all $i\neq j$, we always have $\overline{F}(t)=F_{i}(t)=t$ for $t\in\mathcal{A}$. This implies that the $\xi_{k}$’s are the critical values of [Heyse] and shows the reversed inclusion.

C.1.2 Proof of Lemmas 2.3 and 2.4

Let $\tau_{k}$, $1\leq k\leq m$, be the critical values of [DBH-SU]. Let $\xi_{k}=\max\left\{t\in\mathcal{A}\>:\>t\leq\frac{\alpha}{1+\alpha}\frac{k}{m}\right\}$ be the effective critical values of the [BH] procedure at level $\alpha/(1+\alpha)$. Now, for all $t\in[0,1]$, we have by (2),

[TABLE]

where the last inequality follows from the definition of $\tau_{m}$. Thus we have $\overline{F}_{\textnormal{SU}}(\xi_{m})\leq\alpha$, which in turn implies $\xi_{m}\leq\tau_{m}$. Additionally, the bound (21) yields for $1\leq k<m$

[TABLE]

where we used that $\xi_{m}\leq\tau_{m}$. This proves Lemma 2. The proof of Lemma 3 is analogue and is left to the reader.

C.1.3 Proof of Lemma 2.6

Let us first focus on the case $(i)$ and denote by $\tau_{k}$, $1\leq k\leq m$, the critical values of [A-DBH-SU]. From (2), we have for $1\leq k\leq m-1$,

[TABLE]

which correspond to the effective critical values of [BR-$\lambda$] with $\lambda=\tau_{m}$. Now consider the case $(ii)$ and denote again by $\tau_{k}$, $1\leq k\leq m$, the critical values of [A-DBH-SD]. From (2), we have for $1\leq k\leq m$,

[TABLE]

which correspond to the effective critical values of [GBS]. This implies the result.

C.2 Proofs of technical lemmas for step-down and step-up procedures

C.2.1 Proof of Lemma B.1

First note that for any step-up procedure

[TABLE]

which is sometimes more handy, because this definition avoids to rely explicitly on the order statistics of the $p$-values.

Now, it is not difficult to check that $\widehat{k}^{\prime-i}\geq\widehat{k}-1$ always holds: this comes from the inequality

[TABLE]

because $\tau^{\prime-i}_{\ell-1}=\tau_{\ell}$ for $\ell\in\{2,\dots,m\}$ (note that we can assume without loss of generality $\widehat{k}\geq 1$ here). This means that $(i)$ implies $(ii)$. Now, when $p_{i}\leq{\tau}_{\hat{k}^{\prime-i}+1}$, we have

[TABLE]

which implies $\hat{k}^{\prime-i}+1\leq\sum_{j=1}^{m}{\mathbf{1}\{p_{j}\leq\tau_{\hat{k}^{\prime-i}+1}\}}$ and thus $\widehat{k}^{\prime-i}+1\leq\widehat{k}$. Since, again, $\widehat{k}^{\prime-i}\geq\widehat{k}-1$ always holds, we have $\widehat{k}^{\prime-i}+1=\widehat{k}$. Hence, $(ii)$ implies $(iii)$. Now, if $\widehat{k}^{\prime-i}+1=\widehat{k}$, we have

[TABLE]

by definition of $\tau^{\prime-i}$, which gives that $(iii)$ implies $(i)$. Now, to prove the last statement, we first note that $\widehat{k}^{\prime}\geq\widehat{k}^{\prime-i}$ always holds. Furthermore, if $p_{i}>{\tau}_{m}$ let us prove $\widehat{k}^{\prime}\leq\widehat{k}^{\prime-i}$. First, $\widehat{k}^{\prime}=m$ is impossible because $p_{i}$ is above ${\tau}_{m}$ and thus $p_{i}$ cannot be rejected by $\mathbf{SU}^{\prime}(\tau)$. Hence, $\widehat{k}^{\prime}\leq m-1$ and thus $\tau^{\prime-i}_{\hat{k}^{\prime}}$ is well defined. Now, since $p_{i}>{\tau}_{m}$, we obtain

[TABLE]

which implies $\widehat{k}^{\prime}\leq\widehat{k}^{\prime-i}$ by definition of ${\mathbf{SU}}^{\prime-i}(\tau)$.

C.2.2 Proof of Lemma B.2

First note that for any step-down procedure

[TABLE]

Now, we check that $\widetilde{k}^{\prime-i}+1\geq\widetilde{k}$ always holds. Since $\sum_{j\neq i}{\mathbf{1}\{p_{j}\leq\tau^{\prime-i}_{\tilde{k}^{\prime-i}+1}\}}<\tilde{k}^{\prime-i}+1$, we have

[TABLE]

which gives $\widetilde{k}<\widetilde{k}^{\prime-i}+2$ by definition of $\widetilde{k}$ and thus $\widetilde{k}\leq\widetilde{k}^{\prime-i}+1$. Next, if $p_{i}\leq{\tau}_{\tilde{k}}$, we have

[TABLE]

so that $\widetilde{k}>\widetilde{k}^{\prime-i}$ and thus $\widetilde{k}\geq\widetilde{k}^{\prime-i}+1$. This proves that $(i)$ implies $(iv)$. Next, if $p_{i}>{\tau}_{\tilde{k}^{-i}+1}$, then

[TABLE]

which entails $\widetilde{k}<\widetilde{k}^{-i}+1$ and thus $\widetilde{k}\leq\widetilde{k}^{-i}$. This proves $\widetilde{k}\neq\widetilde{k}^{\prime-i}+1$. Hence, $(iv)$ implies $(iii)$. The fact that $(iii)$ implies $(ii)$ is obvious because $\widetilde{k}\geq\widetilde{k}^{-i}$ always holds. Finally, we merely check that $\widetilde{k}$ is such that

[TABLE]

which means that the set of $p$-values rejected at threshold $\tau_{\tilde{k}}$ is the same as the set of $p$-values rejected at threshold $\tau_{\tilde{k}+1}$. This gives that $(ii)$ implies $(i)$. For the last assertion, it has been proved in the above reasoning while showing that $(iv)$ implies $(iii)$.

C.3 Table for the simulations