With the growing availability of large-scale datasets, decision-making in healthcare, information technology, and government agencies is increasingly driven by data analyses. This data-driven paradigm, however, comes with great risk if the databases contain sensitive information of individuals such as health records or financial data. Without appropriate adjustments, statistical analysis applied to these databases can lead to privacy violation. For example, Homer et al. demonstrate that, under certain conditions, it is possible to determine whether an individual with a known genotype is in a genome-wide association study (GWAS) even when only minor allele frequencies are revealed [38]. Such privacy issues have serious implications: at best, individuals and agencies are discouraged from sharing their data for research purposes due to the concern of privacy leakage, impeding scientific progress [41]; at worst, potential adversaries could make use of sensitive information to jeopardize the social foundations of liberal democracy [64].
Being able to conduct data analysis in a way that preserves privacy, therefore, is key to removing barriers to scientific research while preventing breaches of personal data. First introduced by Dwork et al. [25], differential privacy (Definition 2.1) has put private data analysis on a rigorous foundation. A differentially private algorithm is required to hide the presence or absence of any individual or small group of individuals, the intuition being that an adversary unable to tell whether or not a given individual is even a member of the database surely cannot glean information specific to this individual. In computer science, considerable efforts have been made to develop private data release mechanisms [25,47,2] and private machine learning algorithms under differential privacy constraints, for example, boosting [35], empirical risk minimization [15], private PAC learning [4], and deep learning [1,12]. On the statistical front, differential privacy has been added to and incorporated into many statistical methods in areas of robust statistics [24], nonparametric density estimation [63], hypothesis testing [61,31], finite-sample confidence intervals [40], functional data analysis [34], network data analysis [39], and linear regression [43,62].
In this paper, we provide the first differentially private multiple testing procedure. The problem of multiple testing arises in many privacy-sensitive applications such as a GWAS, where a large number of single-nucleotide polymorphisms (SNPs) are tested simultaneously for an association with a disease and the hope is to control some error rate for the significant SNPs. Perhaps the most popular error rate is the false discovery rate (FDR), which, roughly speaking, is the expected fraction of erroneously rejected hypotheses among all rejected hypotheses. This notion of type I error rate was introduced in the seminal work of Benjamini and Hochberg [5], along with the Benjamini-Hochberg procedure (BHq) that controls the FDR under certain conditions. This procedure is detailed in Algorithm 1.
Our interest in privacy-preserving FDR control arose as a group of researchers showed how to use one-way marginals, specifically, allele frequency statistics, together with the DNA of a target individual and allele frequency statistics for the general population, to determine the target's presence or absence in the study [38]. In response, the US National Institutes of Health and the Wellcome Trust changed the access policy to statistics of this type in the studies they fund. Although differential privacy has been shown to permit nontrivial estimates of very large numbers of statistical queries [10,36], the errors introduced in these techniques are -and must be [14,29] -too large for the (typical) setting, where the number of alleles exceeds the square of the number of data subjects.
Our procedure, which is referred to as PrivateBHq henceforth (Algorithm 4), is derived by recognizing the iterative nature of the BHq procedure and making each iteration differentially private. PrivateBHq provides unconditional end-to-end privacy. For now, regarding p-values as functions of a dataset, our proof of the privacy guarantees of PrivateBHq relies on a new definition of sensitivity that is tailor-made for p-values (Definition 2.5). Loosely speaking, this definition evaluates how insensitive p-values are to perturbations of any individual record in the database. All computations satisfy the definition for some choice of the privacy parameters, but not all choices of these parameters yield useful results when we enforce privacy. Popular examples of p-values are described in terms of these privacy parameters in Section 2.
Another contribution of this paper lies in our proof of the FDR control of PrivateBHq and beyond. In short, all existing proof strategies for FDR control are invalid for PrivateBHq. Thus, a new technique for proving FDR control is needed. To this end, we (1) Develop a novel proof of FDR control for a class of multiple testing procedures, including the original (non-private) BHq and many of its variants -a proof requiring different assumptions than those found in the vast literature on this topic (see Section 3) -and (2) Relate the FDR control and power properties of PrivateBHq to the corresponding properties of the non-private version.
The outline of the remainder of the paper is as follows. The next two subsections elucidate the two contributions, namely developing PrivateBHq and proving FDR control for a class of procedures, and the following subsection consolidates privacy and inferential properties together for PrivateBHq. To make this paper self-contained, in Section 2 we give a brief introduction to differential privacy, followed by the complete development of the PrivateBHq procedure. Section 3 is devoted to establishing FDR control of a broad class of multiple testing procedures and, as an application, Section 4 proves FDR control of PrivateBHq and argues its power as well. The paper is concluded by a discussion in Section 5.
Input: nominal level 0 < q < 1 and p-values p 1 , . . . , p m Output: a set of rejected hypotheses 1: sort the p-values in increasing order:
2: for j = m to 1 do reject p (1) , . . . , p (j) and halt
end if 8: end for
In words, the BHq (step-up) procedure finds the largest j such that p (j ) ≤ qj /m and rejects all p-values below qj /m.
1.1. Making BHq private. The original BHq is our starting point in developing the PrivateBHq procedure. The original procedure is non-private because the data of a single individual can affect the p-values of all hypotheses simultaneously, possibly changing the outcome of the BHq procedure dramatically.
To make the BHq private, for now we need two facts about differential privacy: (1) differential privacy is closed under composition, permitting us to bound the cumulative privacy loss over multiple differentially private computations. This allows us to build complex differentially private algorithms from simple differentially private primitives, and (2) we will make use of the well-known Report Noisy Max (respectively, Report Noisy Min) primitive [26], in which appropriately distributed fresh random noise is added to the result of each computation, and the index of the computation yielding the maximum (respectively, minimum) noisy value is returned. By returning only one index the procedure allows us to pay an accuracy price for a single computation rather than all computations.
A natural approach to obtaining a private version of BHq is by repeated use of Report Noisy Max: Starting with j = m and decreasing: use Report Noisy Max to find the (approximately) largest p-value; estimate that p-value and, if the estimate is above a certain more conservative critical value than qj/m, accept the corresponding null hypothesis, remove it from consideration, and repeat. Once a hypothesis is found with its p-value below the threshold, reject all the remaining hypotheses. The principal difficulty with this approach is that every iteration of the algorithm incurs a privacy loss, which can be mitigated only by increasing the magnitude of the noise used by Report Noisy Max. Since each iteration corresponds to the failure of rejecting a null hypothesis, this step-up procedure is paying in privacy precisely for all null hypotheses accepted, which are by definition not the "interesting" ones. Moreover, recognizing that most null hypotheses in a typical GWAS would be accepted, it is fundamentally difficult to preserve information content while protecting individual privacy by emulating the step-up procedure.
Instead of starting with the largest p-value and considering the values in decreasing order, another approach is to start with the smallest p-value and consider the values in increasing order, rejecting hypotheses one by one until we find a p-value above some threshold. This widely studied variant is called the BHq step-down procedure, which, in contrast to the aforementioned BHq step-up procedure, finds the largest j such that p (i) ≤ qi/m for all i ≤ j and then rejects p (1) , . . . , p (j) . Their definitions reveal that the step-down procedure shall be more conservative than its step-up counterpart. This variant, however, can assume less stringent critical values than the BHq critical values while still offering FDR control, often allowing more discoveries than the step-up counterpart [32].
If we make the natural modifications to the step-down procedure using Report Noisy Min, also known as the Private Min (Algorithm 2), instead of Report Noisy Max, then we pay a privacy cost only for nulls rejected in favor of the corresponding alternative hypotheses, which by definition are the "interesting" ones. Since the driving application of BHq is to select promising directions for future investigation that have a decent chance of panning out, we can view its outcome as advice for allocating resources. Thus, a procedure that finds a relatively small number of high-quality hypotheses, still achieving FDR control, may be as useful as a procedure that finds a much larger set.
A new technique for proving FDR control. While various techniques have been developed in the literature for proving FDR control, they are not applicable to privacypreserving procedures. Any privacy-preserving procedure is necessarily randomized. Consequently, the jth most significant noisy p-value may not necessarily correspond to the jth most significant true p-value. Even worse, PrivateBHq may compare a noisy p-value to a critical value with a different rank and, as an inevitable result, a larger p-value may be rejected while a smaller p-value is accepted. This is in stark contrast to the (non-private) BHq and most of its variants, which reject p-values that are contiguous in sorted order.
These facts about the PrivateBHq procedure destroy some crucial properties for proving FDR control in existing approaches. For example, it is not clear how to adapt the elegant martingale technique for FDR control, proposed by Storey, Taylor, and Siegmund [57]. In essence, this approach is to construct an empirical process indexed by a threshold under which a p-value is rejected. In the case of PrivateBHq, unfortunately, no such threshold exists for singling out p-values for declaring significance. Another technique that appears frequently in the FDR control literature (see, for example, [7,53,30,50,11,37]) is based on a crucial property of BHq: provided that a p-value is rejected, the effective threshold for declaring significance is completely determined by the remaining p-values. Unfortunately, this property is not satisfied by PrivateBHq either.
To pursue a new strategy for PrivateBHq, we observe that, although PrivateBHq might skip some of the minimum p-values, nevertheless it preserves a key property with high probability: if R rejections are made, the largest rejected p-value is roughly upper bounded by qR/m. This motivates us to give the following definition. Definition 1.1. Given any cutoffs 0 < q 1 ≤ q 2 ≤ • • • ≤ q m , a multiple testing procedure is said to be compliant with {q j } m j=1 , if all rejected p-values are always bounded above by q R , where R is the number of rejections.
In the case of no rejections (R = 0), as a convention, the (non-existent) rejected p-value is considered to be bounded above by q R . Compliance is an instance of a more general condition termed self-consistency [8,9], which, roughly speaking, requires that any rejected p-value be upper bounded by a general function of the total number of rejections. Interestingly, the compliance condition as a special instance of self-consistency has not been considered in the literature. Here, we prefer to use the compliance condition as self-consistency further allows a procedure to incorporate prior information about each hypothesis into the cutoffs, which is beyond the scope of this paper. Using the BHq critical values {qj/m} as the cutoffs (referred to as BHq-compliance henceforth), however, our condition is sufficiently general to cover many classical multiple testing procedures, including both the step-down and step-up procedures, the generalized step-up-step-down procedures [60,51] and particularly the PrivateBHq procedure (Proposition 4.1). The compliance condition is solely determined by the number of rejected p-values and the size of the largest one, without requiring that each rejected p-value be below its associated critical value. As a consequence, this condition permits skipping the smallest p-values and this is well-suited for differentially private procedures.
As revealed by this work, FDR control, roughly speaking, is a consequence of BHqcompliance together with the independence with the null condition (Definition 1.2). As such, our finding offers more than expected, applying to far more examples than PrivateBHq. In detail, we consider a generalized FDR [52,54] defined as
where V denotes the number of true null hypotheses that are falsely rejected (false discoveries).
The present paper primarily focuses on the case of k ≥ 2 and, whenever clear from the context, the term FDR control in this paper stands for FDR k control. Note that FDR k reduces to the usual FDR if the positive integer k is set to 1. This slightly relaxed FDR permits no more than k -1 false discoveries without any penalty, trading off for more power improvement while still maintaining a meaningful interpretation of the rejected hypotheses.
The difference between the original FDR and FDR k becomes negligible if the number of discoveries R is large. As an aside, we remark that the compliance condition is not satisfied by the FDR k -controlling procedures developed in [52,54]. Now we introduce the independence with the null condition, which is concerned with the distribution of the p-values. This condition is satisfied by the three examples in Appendix B. Definition 1.2. A set of m test statistics are said to satisfy a condition referred to as independence within the null, or IWN for short, if the true null test statistics are jointly independent.
More elaboration on this new condition is carried out following Theorem 2 below. With the two preparatory definitions in place, we offer the following theorem. Let m 0 denote the number of true null hypotheses and π 0 := m 0 /m be the true null proportion. Theorem 1 . If the test statistics obey the IWN condition, then any procedure that is compliant with the BHq critical values {qj/m} m j=1 must satisfy
for every k ≥ 2, where C k is a universal constant.
We immediately obtain the following corollary. This bound involves an additional factor C k , compared with the usual bound π 0 q in the FDR literature. Explicitly, letting {ξ j } ∞ j=1 be i. i. d. exponential random variables with mean 1, the constant is given as
For example, C 2 ≈ 2.41, C 3 ≈ 1.85, C 10 ≈ 1.32, and C k tends to 1 as k → ∞. In particular, C 1 defined in (1.2) is infinite, and this is exactly why Theorem 1 does not apply to the usual FDR. Theorem 1 is optimal for all k ≥ 2 as we show next.
Theorem 2 . Given any C < C k , if q is sufficiently small and m is sufficiently large, then there exists a BHq-compliant procedure applied to a set of IWN p-values such that FDR k > Cq.
In the literature, existing FDR-controlling procedures often assume independence between the true null and false null test statistics (see [5,6]) or certain sophisticated correlation structures between these two sets of test statistics, such as the positive regression dependent on subset (PRDS) property [7,42,55] (see also [51,9]). Roughly speaking, the PRDS property holds if the test statistics exhibit certain positive dependence on each true null test statistic. In particular, the dependence between true and false nulls cannot be arbitrary. For the sake of completeness, we emphasize that the literature has considered a few cases for FDR control with an arbitrary correlation between the two sets of test statistics [7,8], but, unfortunately, the associated procedures are often extremely conservative. As a well-known example, Benjamini and Yekutieli showed in Theorem 1.3 of [7] that the BHq procedure gives FDR control using critical values at level q/(1
log m + 0.577) in place of q. In fact, BHq with this log-factor correction could be even more conservative than the Bonferroni method [45].
In contrast, Theorem 1 makes no assumptions regarding the dependence between the true nulls and false nulls while still controlling the FDR up to a small multiplicative factor, as the IWN condition is concerned only with the true nulls. As such, Theorem 1 is a contribution of independent interest to the vast FDR literature. Notably, the dependence can even be "adversarial" in the sense that the false null p-values can even be constructed as arbitrary functions of the true null p-values. This provides positive evidence toward understanding the robustness of the BHq procedure observed in a wide range of theoretical and empirical studies [56,33,16].
In this section, we first introduce the differential privacy machinery at a minimal level and then focus on developing the PrivateBHq procedure.
2.1. Preliminaries on differential privacy. A database D = (d 1 , d 2 , . . . , d n ) ∈ X n consists of n data items (for example, health records of n individuals), where X is a sample universe. Data items need not be independent (for example, health records of siblings). Two databases D, D = (d 1 , d 2 , . . . , d n ) are said to be neighbors, or adjacent, if they differ only in one data item. That is, there is exactly one j such that d j = d j . A (randomized) mechanism M is an algorithm that takes a database as input and releases some (randomized) response of interest. We denote by range(M) the collection of all possible outputs of the mechanism M. In the context of genome-wide association studies, a database D records genotypes of individuals, and M, for example, is a mechanism that releases the minor allele frequency of a SNP plus some random noise.
Differential privacy, now sometimes called pure differential privacy, was defined and first constructed in [25]. The relaxation defined next is sometimes referred to as approximate differential privacy. Definition 2.1 (Differential Privacy [25,23]). A (randomized) mechanism M is (ε, δ)differentially private for some nonnegative ε, δ if for all adjacent databases D, D and for any measurable event S ⊂ range(M),
Pure differential privacy is the special case where δ = 0. In the definition above, both databases D, D are fixed and the probabilities are taken over the randomness of the mechanism M. The parameters ε and δ measure the desired privacy protection. With small ε and δ, this definition states that the likelihood of the released response is indifferent to changing a single individual in the database, thus leaking little indication of whether a particular individual is in the database even if all the other individuals are known. This provides strong privacy protection for each individual in and outside the database.
To report a statistic f = f (D) in a differentially private manner, it is necessary to randomize the mechanism. As its name suggests, the Laplace mechanism preserves privacy by perturbing f with noise generated from the Laplace distribution Lap(λ), whose probability density is exp(-|x|/λ)/(2λ). The scale λ > 0 should be calibrated to the sensitivity of the statistic f , defined as follows.
Definition 2.2. Let f be a real or vector valued function that takes as input a database. The sensitivity of f , denoted as ∆f , is the supremum of f (D) -f (D ) 1 over all adjacent D, D , where • 1 denotes the 1 norm.
Formally, for any function f that maps databases to R r for some positive integer r, we have the following result.
Lemma 2.3 Laplace Mechanism [25]. The Laplace mechanism M L that outputs
preserves (ε, 0)-differential privacy, where Z j are i. i. d. draws from Lap(∆f /ε).
Intuitively, sensitivity quantifies the effect of any individual in the dataset on the outcome of the analysis. In this mechanism, Laplace noise with magnitude proportional to the sensitivity has the effect of masking the characteristics of any individual, thereby preserving privacy.
A simple algorithm that integrates the Laplace mechanism is the Private Min, which is better known as the Report Noisy Min in the literature [26] and will be the building block of PrivateBHq, introduced in Section 2.3. Consider a collection of scalar functions f 1 , . . . , f m . The Private Min adds Laplace noise to each f j and then reports the smallest noisy count (with fresh noise added) and its index. A formal description of Private Min is given in Algorithm 2. The following lemma concerns its privacy property.
Lemma 2.4 . The Private Min, as detailed in Algorithm 2, is (ε, 0)-differentially private.
A peek at the proof of this well known lemma, which for completeness appears in the appendix, reveals that reporting each of j and f j (D) + Z is (ε/2, 0)-differentially private, hence leading to a total privacy loss of (ε/2, 0) + (ε/2, 0) = (ε, 0). Here we have used the simple fact that differential privacy loss adds up under the composition of sequential mechanisms, that is, the union of the outputs of a sequence of mechanisms that each preserve (ε j , δ j )-differential privacy is ( ε j , δ j )-differentially private [25]. As an aside, the Advanced Composition Theorem [27] (see Lemma 2.9 in Section 2.4), provides a much tighter bound on this privacy degradation.
Input: database D, functions f 1 , . . . , f m each with sensitivity at most ∆, and privacy parameter ε Output: index j and approximation to f j (D)
1: for j = 1 to m do 2:
set f j = f j (D) + Z j , where Z j is independently sampled from Lap(2∆/ε) 3: end for 4: return j = argmin j f j and f j (D) + Z, where Z is a fresh draw from Lap(2∆/ε) Looking ahead, and omitting some technicalities, PrivateBHq will operate on differentially private approximations to the logarithms of p-values, returned by multiple invocations of Private Min. Since differential privacy is closed under post-processing [26], any subsequent computation on these differentially privately obtained values can never increase privacy loss. Thus, PrivateBHq is indeed differentially private, for all p-value functions. Its statistical properties will depend on the kinds of p-value computations that are performed, which we turn to next.
Multiple testing procedures ubiquitously act on a set of p-values that are computed by functions that operate on databases. A p-value in our context is frequently referred to as the function on databases for computing the p-value instead of its numerical value, in contrast with the vast statistical literature that often does not distinguish between the function that maps a database to a p-value and the result of the mapping.
We now consider making p-value computations private as the first step toward developing a private multiple testing procedure. In many important p-value computations (see Example 2.6), a larger p-value is affected more in magnitude by the change of a single data item than a smaller p-value. As a result, directly adding noise to the p-values may overprotect privacy and completely overwhelm signals in small p-values. This would inevitably lead to significant detection power loss as the smallest p-values are more likely to correspond to promising hypotheses.
Our solution will be to (very carefully) work with the logarithms of the p-values. This strategy is motivated by the observation that, although the (additive) sensitivity of a p-value may vary greatly, oftentimes the relative change (that is, the ratio) of a p-value on two neighboring databases is very stable, regardless of the magnitude of the p-value, unless it is extremely small. In light of this observation, the sensitivity of a p-value, that is, the worst-case change due to the replacement of an individual in the database, is best measured multiplicatively. Below, η and ν are nonnegative. Our PrivateBHq algorithm will make explicit use of both parameters in ensuring privacy. The parameter ν is introduced in recognition of the fact that a very small p-value may jump or fall by a relatively large multiplicative factor between adjacent databases. This parameter is normally much less than the Bonferroni level q/m (see, for example, [19]), resulting in essentially no power loss for truncating p-values at ν. A p-value can satisfy different pairs of (η, ν)-multiplicative sensitivities. In short, the two parameters η and ν exhibit a certain trade-off relationship in the sense that one can increase (resp. decrease) η and decrease (resp. increase) ν in a careful way such that a p-value still satisfies this condition. Every p-value satisfies (η, ν)-sensitivity for some values of the parameters. Moreover, given p-value functions p 1 , p 2 with multiplicative sensitivities (η 1 , ν 1 ) and (η 2 , ν 2 ) respectively, it is immediate that both functions satisfy (max{η 1 , η 2 }, max{ν 1 , ν 2 })-sensitivity, so given a collection of p-values there always exist η, ν so that all of the p-values in the collection are (η, ν)-sensitive.
Given an (η, ν)-sensitivity p-value function p and a database D, we work with the logarithmic mapping θ(D; p, ν) = log max{ν, p(D)} This statistic satisfies θ(D) -η ≤ θ(D ) ≤ θ(D) + η for all neighboring databases D, D . In other words, θ has an additive sensitivity bounded by η. Hence, Lemma 2.2 ensures that adding Laplace noise Lap(η/ε) to θ(D) preserves (ε, 0)-differential privacy.
We will see below via examples that two large and important classes of p-value computations are (η, ν)-sensitive for some small η and ν, with rigorous proofs given in Appendix A; as a consequence of this, preserving privacy for these p-values only requires a small amount of noise, leading to negligible accuracy loss. Recall that m denotes the total number of hypotheses.
Example 2.6 (Binomial Distribution). Suppose the n individuals in D are, respectively, associated with n i. i. d. Bernoulli variables ξ 1 , . . . , ξ n , each of which takes the value 1 with probability α and the value 0 otherwise. Let T denote the sum. A p-value p(D) for testing
where t is the realization of T on the database D. Denote by t the counterpart of t on a neighboring database D . Without loss of generality, assume t = t + 1. The difference between the two p-values, |p(D) -p (D)| = 1 2 n n t , attains its maximum at t = n/2 or (n + 1)/2 ( x denotes the greatest integer that is less than or equal to x) and decays rapidly as t deviates from n/2. This implies that additive sensitivity is not a good measure of the variability of this p-value construction.
Instead, we fix a (very) small ν and denote by η the maximum of log p(D) p(D ) subject to the constraint p(D ) ≥ ν. The p-value by definition is (η, ν)-sensitive. To evaluate η, observe that the log-likelihood ratio To test H 0 : λ = 1 against the alternative H 1 : λ > 1, we consider the p-value p(D) = P λ=1 (T ≥ t), where t is the realization of T (note that the value t differs at most by A between adjacent databases). With the same notations as in Example 2.6, this p-value is (η, ν)-multiplicatively sensitive with ν = m -1-c and η log n n for any small constant c > 0. Similarly, the analysis applies to the case of a Gaussian distribution. In short, consider i. i. d. random variables ξ 1 , . . . , ξ n drawn from the normal distribution N (µ, 1) truncated at -A and A, which has density e
Using the same proof strategy as for the exponential distribution, one can show that this p-value strategy is (η, ν)-multiplicatively sensitive with some ν = m -1-c and η log n n . We remark that (η, ν)-sensitivity is a worst-case guarantee on the sensitivity of a p-value function. Only the interpretation of the p-value requires the i.i.d. assumption. Regarding the above-mentioned two examples, the asymptotic expressions of the privacy parameter η can be easily made precise.
As seen in both examples, the parameter η vanishes roughly at the rate O(n -1/2 ), implying that less noise is required for privacy protection as the sample size becomes larger. This appealing feature is impossible without the restriction p ≥ ν for some appropriate choice of ν. Specifically, in the absence of this constraint, or equivalently by setting ν = 0, we shall have η = n + 1 in the first example and η = ∞ in the second, requiring a vast or even an infinite amount of multiplicative noise for preserving privacy. This would completely dilute any signal of interest. To be complete, we note that not all p-value computations necessarily lead to vanishing η and ν as n → ∞. An example from [61,65] considers a privacy-preserving release of χ 2 -statistics computed from allelic contingency tables. For the sake of simplicity, here we consider 2 × 2 contingency tables with n/2 cases and n/2 controls:
Table 1. Two neighboring allelic contingency tables.
In the case of a fixed a > 5, one can show that the two p-values computed from the two tables neither differ by a negligible factor nor both tend to zero as n → ∞. This fact is elaborated in detail in the appendix.
2.3. Developing PrivateBHq. The PrivateBHq procedure (Algorithm 4) is the sequential composition of Algorithm 3, which we refer to as the peeling mechanism, denoted as peeling. In a little more detail, given (non-private) p-value functions p 1 , . . . , p m and a prescribed number of invocations m ≤ m, PrivateBHq first applies Private Min m times to the logarithms of the p-values, "peeling off" and removing from further consideration the approximately smallest element with each new invocation of Private Min. These m pre-selected hypotheses are thought of as promising hypotheses. In particular, the number m as an upper bound on the total number of discoveries shall be much less than m. This recognizes that, in many application scenarios, much fewer are truly significant in an ocean of mediocre hypotheses.
During the peeling procedure, in order to keep track of indices within the original set, peeling removes a function from further consideration by redefining it to be +∞, ensuring that it will not be returned by future invocations of the Private Min. The Laplace noise scale λ shall be chosen to adjust for the privacy protection target, factoring in the multiplicative sensitivities of p 1 , . . . , p m and the number of invocations m . let (i j , fi j (D)) be returned by Private Min applied to (D, f 1 , . . . , f m ) with Laplace noise scale λ 3:
With m hypotheses yielded by peeling in place, PrivateBHq supplies quantities in logarithmic scale instead of, in the conventional setting, the m raw p-values and critical values to the (step-up) BHq procedure. This difference however does not affect the way BHq proceeds. To be concrete, BHq first orders the noisy values θi 1 , . . . , θi m as θ(i 1 ) ≤ • • • ≤ θ(i m ) , and then rejects any corresponding hypotheses if θi j is below max{γ j : θ(i j ) ≤ γ j }, with the convention that max ∅ = -∞. As we will see in Section 4, the cutoffs γ 1 , . . . , γ m are chosen specifically to ensure FDR control of PrivateBHq; roughly speaking, γ j is slightly below the logarithm of the corresponding BHq critical value qj/m, where the gap between the two accounts for the multiplicative sensitivity of the p-values and the uncertainty brought by the Laplace mechanism.
Input: database D, parameters ε, δ, η, ν, (η, ν)-multiplicatively sensitive p-value functions p 1 , . . . , p m , number of invocations m , Laplace noise scale λ = λ(ε, δ, η, m ) , and cutoffs Lemma 2.8 [23,63]. Let M be an (ε, δ)-differentially private mechanism and g be any (measurable) function. Then g(M) also preserves (ε, δ)-differential privacy.
This lemma implicitly assumes the range of the mechanism M falls into the domain of g. In our context, taking g to be step-up BHq, Lemma 2.8 shows that it suffices to establish the differential privacy property of peeling. By construction, each θ j has sensitivity no more than η. Lemma 2.4 then immediately ensures that the Private Min, which is invoked sequentially m times in PrivateBHq, guarantees on its own (2η/λ, 0)-differential privacy. Making use of the fact that, at worst "(ε, δ)'s add up" (see the discussion right below Lemma 2.4), one can conclude that the peeling mechanism is (2m η/λ, 0)-differentially private. Equivalently, to achieve (ε, 0)-differential privacy for peeling, and therefore also for PrivateBHq, we can set the Laplace noise scale to be λ = 2m η/ε. In this way, the noise level grows linearly with m . Surprisingly, we can trade a little bit of δ for a significant improvement on ε, as shown by the lemma below. Lemma 2.9 Advanced Composition [27]. For all ε, δ ≥ 0 and δ > 0, running l mechanisms sequentially that are each (ε, δ)-differentially private preserves (ε 2l log(1/δ ) + lε(e ε -1), lδ + δ )-differential privacy.
This lemma holds no matter how each mechanism adaptively depends on information released by prior mechanisms. Taking δ = 0 in Lemma 2.9, we easily obtain the main theorem of this section, with its proof deferred to the appendix. This theorem shows adding Laplace noise with scale of order roughly O( √ m ) is sufficient for protecting privacy of PrivateBHq.
Theorem 3 . Let η, ν be chosen so that all the p-value functions input to PrivateBHq are (η, ν)-sensitive. Given ε ≤ 0.5, δ ≤ 0.1 and m ≥ 10, PrivateBHq with Laplace noise scale λ = η 10m log(1/δ)/ε, or larger, is (ε, δ)-differentially private.
We remark that the constraints on ε, δ, and m are used to optimize the constants for practical use.
The main purpose of this section is to prove Theorem 1. The proof strategy contains two novel elements: an upper bound on FDR k involving only true null p-values (Equation (3.1) below) and a backward submartingale that allows us to use a martingale maximal inequality. In addition, this section attempts to obtain the optimal constant C k for Theorem 1 in Section 3.2, where we give some intuition behind Theorem 2, and considers a new variant of the FDR in Section 3.3.
Throughout the section, we focus on an arbitrary BHq-compliant procedure. That is, any p-value rejected by the procedure is not greater than qR/m, where R denotes the total number of rejections.
3.1. Controlling FDR k . In this subsection, we prove Theorem 1. However, the proof presented here does not seek to optimize the constant C k in Theorem 1. We consider
which gives FDR k ≡ E FDP k by taking expectation. The following upper bound on the FDP k for k ≥ 2 of the BHq-compliant procedure serves as the basis for our analysis:
are the order statistics of the m 0 true null p-values. To prove (3.1), denote by V the number of false rejections. If V ≤ k -1, (3.1) holds since FDP k = 0. Otherwise, the largest rejected true null p-value is at least p 0 (V ) and, therefore, one must have p 0 (V ) ≤ qR/m due to the compliance condition. As a consequence, we get
.
(3.
2)
The IWN condition imposed in Theorem 1 ensures the joint independence of the true null pvalues, each of which is, by definition, stochastically larger than or equal to U (0, 1). Thus, the ordered true null p-values can be replaced by the order statistics
of m 0 i. i. d. uniform random variables on (0, 1), while (3.2) remains true in the expectation sense (recall that π 0 = m 0 /m):
.
Therefore, Theorem 1 follows from the lemma below.
denote the order statistics of n i. i. d. uniform variables on (0, 1). There exists an absolute constant c k such that
The proof of this lemma starts by recognizing a well-known representation in law for uniform order statistics:
where T j = ξ 1 + • • • + ξ j and ξ 1 , . . . , ξ n+1 are i. i. d. exponential random variables with mean 1. Writing
Intuitively, the maximum is likely to be attained at some small index j as W j /n is close to 1 for a large value of j, due to the law of large numbers. This intuition can be indeed made rigorous by the fact that W 1 , . . . , W n+1 is a backward submartingale, as shown by the following lemma.
Lemma 3.2 . With respect to the filtration F j := σ(T j , T j+1 , . . . , T n+1 ) for j = 1, . . . , n + 1, the stochastic process W 1 , . . . , W n+1 is a backward submartingale. That is, E(W j |F j+1 ) ≥ W j+1 for j = 1, . . . , n.
The proof of Lemma 3.2 is deferred to the appendix. Next, we apply this lemma to prove (3.4) (hence Lemma 3.1 follows immediately) using the following martingale maximal inequality (for a proof, see pages 71-73 of [48]).
Proof of Lemma 3.1. Since Lemma 3.2 asserts that W j /n is a backward submartingale, Lemma 3.3 concludes
To complete the proof, it suffices to show that for a fixed k the expectation above involving k/(nU (k) ) is uniformly bounded for all n ≥ k. To this end, observe that U (k) is distributed as Beta(k, n + 1 -k), and this allows us to evaluate the expectation as
To obtain an upper bound that is independent of n, it suffices to show that n k B(k, n + 1 -k) has a lower bound depending only on k. Indeed, this is the case:
≥ (k -1)!.
Optimizing the bounds. The constant C k in Theorem 1 matters from a practical perspective. This section is aimed at finding the optimal constants for all k ≥ 2. Compared with what has been performed in Section 3.1, this improvement is based on a delicate property about the expectation in (3.4), as detailed by the following lemma.
for n ≥ k ≥ 2, where U (j) 's are the order statistics of n i. i. d. uniform variables on (0, 1). Then,
The monotonicity in Lemma 3.4 reveals that the optimal C k in (3.4) takes the form (recall that
Note that C k does not seem to admit a closed-form expression. Nevertheless, this optimal constant can be easily computed via simulations. While relegating the full proof of Theorem 2 to Appendix A, here we provide a proof sketch based on the construction of a BHq-compliant procedure and a set of p-values satisfying the IWN condition to show the optimality of C k . Explicitly, let the true null p-values be m 0 i. i. d. uniform variables U 1 , . . . , U m 0 between 0 and 1, and let all the m -m 0 false null pvalues be 0. Denote by j the index k ≤ j ≤ m 0 that maximizes j/U (j) . The BHq-compliant procedure rejects the j smallest true null p-values and any max{ mU (j ) /q -j , 0} of the false null p-values ( x denotes the least integer that is greater than or equal to x), which by construction are all 0. This procedure is compliant (self-consistent) but not nonincreasing (a procedure is called nonincreasing if it never rejects more if some p-value gets larger), so the FDR-controlling results in [9] do not apply to our case. Taking q sufficiently small and assuming that m -m 0 is sufficiently large, we get FDP k ≈ qj /(mU (j ) ) with high probability. Consequently, we get
which tends to C k q by taking m 0 → ∞ and m 0 /m → 1.
For the moment, suppose the limit can be taken under the expectation in (3.5). As such, the optimal constant for FDR k is
where the last equality results from applying the strong law of large numbers to T n /n. Recognizing that the integrable random variable max k≤j<∞ j/T j decreases to 1 almost surely as k increases to infinity, Lebesgue's dominated convergence theorem readily asserts that
, where o k (1) denotes a sequence of numbers tending to 0 as k → ∞. This is formally stated in the proposition below, where we consider a sequence of multiple testing problems indexed by l such that both m l , k l → ∞ as l → ∞.
Proposition 3.5 . Under the assumptions of Theorem 1, as k → ∞, we have
To make the derivation of the optimal C k above rigorous, we must validate (3.6). In fact, the Vitali convergence theorem together with the following lemma ensures that the limit lim n→∞ and expectation E can be interchanged. Lemma 3.6 . For a fixed k ≥ 2, the sequence of random variables max k≤j≤n jT n+1 nT j are uniformly integrable for n ≥ k.
While the proof of Lemma 3.6 is deferred to the appendix, the proof of Lemma 3.4 is given below.
the order statistics of n + 1 i. i. d. uniform random variables on (0, 1). Then,
are distributed the same as the order statistics of n i. i. d. uniform random variables on (0, 1) and, moreover, are independent of U (n+1) . Making use of this fact, we get
Since the density of U (n+1) is (n + 1)x n for 0 < x < 1, we readily see that
This completes the last step in certifying The (tiny) shaded band illustrates the 99%-coverage confidence interval for each k using normal approximation. Now, we turn to numerically evaluate C k using the expression (3.6). Although the distribution of each j/T j admits an analytical expression, it is however not clear how to calculate the distribution of the maximum of j/T j over j. In view of this difficulty, we resort to Monte Carlo simulations, and Figure 1 presents the results that are averaged over 10 4 independent replicates. For instance, C 2 ≈ 2.41, C 3 ≈ 1.85, C 4 ≈ 1.65, C 5 ≈ 1.54, and C 25 ≈ 1.18. In passing, we remark that the estimated values of C k as a function of k are fairly accurate as indicated by the uniformly short widths of the confidence intervals for all k.
Controlling FDR k . To further leverage the martingale-based proof idea, we consider a variant of the FDR defined as
which includes the usual FDR as an example by taking k = 1. This relaxed FDR differs insignificantly from the usual FDR if a large number of discoveries are expected, which is often the case in modern multiple testing applications such as genome-wide association studies. For the moment, we do not intend to advocate the use of this new FDR definition in practice as it is clear that future investigation is needed.
In the following, we aim to prove Theorem 4, a counterpart of Theorem 1 for the FDR k . A similarity between the two theorems lies in that their proofs both make use of martingale arguments. That being said, the bound on the FDR k in Theorem 1 cannot carry over to the FDR k because FDR k ≤ FDR k . Theorem 4 . If the test statistics obey the IWN condition, then any BHq-compliant procedure satisfies
for any k ≥ 1.
A number of remarks are as follows. This theorem allows us to take k = 1, thus giving a bound on the usual FDR: FDR ≤ q + 2 √ q. For example, we can set q = 0.0024 if the FDR is aimed to be controlled at 10%. Such a bound is not available in Theorem 1. For completeness, the bound for k = 1 might not be sharp since Doob's 2 martingale maximal inequality used in the proof of Theorem 4 is generally not sharp. Indeed, this bound can be improved using a careful treatment of (3.1) (see [59]). For k ≥ 2, the bound here is larger than that in Theorem 1, namely 2/ √ qk ≥ C k -1, due to the optimality of C k and the fact FDR k ≥ FDR k . The following proof actually establishes a stronger bound, π 0 q + 2 π 0 q/k, on the FDR k . Recall that π 0 is the true null proportion m 0 /m. Proof of Theorem 4. Due to the compliance condition, the number of false discoveries satisfies
Thus, we get an upper bound on FDP := V R (with the convention 0/0 = 0) that takes the following form:
Consequently, we get
Similar to what has been argued in Section 3.1, the inequality (3.7) still holds if all true null p-values are replaced by m 0 i. i. d. uniform variables U 1 , . . . , U m 0 on (0, 1). This observation shows that it suffices to prove
To show (3.8), denote by V j = #{1 ≤ i ≤ m 0 : U i ≤ qj/m} and Y j = V j /j. Conditional on Y j+1 , for every i ∈ {1 ≤ i ≤ m 0 : U i ≤ q(j + 1)/m} the random variable U i is uniformly distributed on [0, q(j + 1)/m]. Hence, the conditional expectation of V j given Y j+1 is
In words, Y j is a backward martingale and, as a consequence, (Y j -qm 0 /m) + is a backward submartingale. This fact allows us to apply Doob's 2 martingale maximal inequality to (Y j -qm 0 /m) + , yielding
Using Jensen's inequality, the left-hand side of (3.8) satisfies
As an application of Theorem 1, this section considers FDR control and power of PrivateBHq. Throughout this process, we take the assumptions of Theorem 3 as given. That is, we assume that each p i is (η, ν)-sensitive and the parameters satisfy ε ≤ 0.5, δ ≤ 0.1, and m ≥ 10. From Theorem 3, PrivateBHq preserves (ε, δ)-differential privacy, and for brevity this fact will not be reiterated in this section.
The proposition below demonstrates that the PrivateBHq is indeed compliant by making the cutoffs {γ j } in Algorithm 4 slightly more stringent than the logarithms of the BHq critical values. Proposition 4.1 . For any 0 < q < 1, use the cutoffs
for j = 1, . . . , m in PrivateBHq. Under the assumptions of Theorem 3, this procedure is compliant with the BHq critical values qj/m with probability at least 1 -0.1q.
As a remark, the first term log qj m in (4.1) corresponds to the non-private cutoff and the second term -η √ 10m log(1/δ) log(6m /q) ε is used to handle the added noise. Notably, the constant 0.1 above can be replaced by any positive constant provided that the second term is appropriately scaled. The proof of Proposition 4.1 is given later after Theorem 5.
The compliance condition shown in Proposition 4.1 together with Theorem 1 implies FDR control of PrivateBHq. More precisely, letting C denote the event that the rejected p-values are compliant, we have
for every k ≥ 2. As such, to control the FDR at level, say 10% (a common level used in practice), we can set q = 0.1/(C k + 0.1) in PrivateBHq. This proves the following theorem.
Theorem 5 . Under the same assumptions as in Proposition 4.1 and if the test statistics satisfy the IWN condition, the PrivateBHq procedure gives
To prove Proposition 4.1, we first present a simple lemma that gives a concentration bound on Laplace random variables, and its proof can be found in the appendix. Lemma 4.2 . Let Z 1 , . . . , Z n be i. i. d. Lap(λ) random variables. For any 0 < α < 1, the following two statements are true:
(1) With probability at least 1 -α, all Z j are larger than -λ log n 2α .
(2) With probability at least 1 -α, all |Z j | are smaller than λ log n α . Proof of Proposition 4.1. Let θi j = log max{ν, p i j }+Z i j be yielded by peeling in Algorithm 4, where Z i j follows Lap(λ) for j = 1, . . . , m . The parameter λ = η 10m log(1/δ)/ε is as in Theorem 3. Taking α = 0.1q, Lemma 4.2 shows that
uniformly for j = 1, . . . , m with probability at least 1 -0.1q. Next, we show that on the event (4.2), PrivateBHq is compliant. Denote by R Pt the number of rejections made by this procedure. If θi j is rejected, it must satisfy log max{ν,
Plugging (4.2) into this display gives log max{ν, p i j } ≤ log qR Pt m .
Thus, p i j ≤ qR Pt /m for all rejected p i j on the event (4.2), which happens with probability at least 1 -0.1q. This completes the proof.
Next, Theorem 6 shows that the PrivateBHq procedure with a slightly inflated nominal level is at least as powerful as the BHq step-down procedure. The proofs of this theorem and its corollary are deferred to the appendix. To state the theorem, let R SD denote the number of rejections made by the (non-private) step-down procedure. Theorem 6 . Fix q and assume ν ≤ q/m. Under the assumptions of Theorem 5, run the PrivateBHq procedure at level
and the BHq step-down procedure at level q. Then, the numbers of rejections satisfy
with probability tending to one as m → ∞.
When R SD ≥ m and the event (4.3) happens, PrivateBHq must reject all p-values passing through peeling. In the case where non-null p-values are significant enough to pass through peeling, this fact suggests that PrivateBHq achieves high power. This high-power property, however, is appealing if q is only slightly larger than q or, put more simply, the number 24η m log(1/δ) log m/ε is small. With regard to Examples 2.6 and 2.7, this is equivalent to having a sufficiently large sample size n. The following corollary formalizes this point.
Corollary 4.3 . In Examples 2.6 and 2.7, fix ε, δ and assume m ≤ min{n 1-c , m} for constant c > 0. Under the assumptions of Theorem 6, the claims of both Theorems 5 and 6 hold as m, n → ∞ if PrivateBHq is performed at level (1 + c )q for a sufficiently small constant c > 0.
4.1. Empirical evaluation. In this subsection, we evaluate the price paid for privacy in terms of FDR control and power in the PrivateBHq procedure. The aim is to provide a better picture of how much detection power would be compromised due to privacy guarantees for FDR control. For completeness, this comparison includes a private version of Bonferroni's method, which is referred to as PrivateBonf in this paper. PrivateBonf is perhaps the simplest baseline for private multiple hypothesis testing. This procedure is detailed as follows. As in Algorithm 4, let p 1 , . . . , p m be (η, ν)-sensitive p-values and set θ j = log max{p j , ν} for all j. The PrivateBonf procedure adds independent Lap( λ) noise to all θ j , where λ = η 10m log(1/δ)/(2ε), and rejects those with noisy counts below log q m -η 10m log(1/δ) log(5m/q) 2ε .
The following result is concerned with privacy and family-wise error rate (FWER) control of PrivateBonf. Note that the FWER denotes the probability that at least one false positive is made. The proof is deferred to the appendix. As an aside, the privacy guarantee in this result might be improved by using the sparse vector technique [36] and this is left for future investigation. Proposition 4.4 . Under the assumptions of Theorem 3, the following two statements are true:
(1) PrivateBonf is (ε, δ)-differentially private;
(2) PrivateBonf satisfies FWER ≤ 1.1q. q q q q q q q q q q 0.00
PrivateBHq PrivateBonf q q q q q q 0.00
PrivateBHq PrivateBonf q q q q q q 0.00 µ FDR q BHq PrivateBHq PrivateBonf q q q q q q q q 0.00 Figures 2 and3 present, respectively, the FDR and the power of PrivateBHq, the (non-private) BHq step-up procedure, and PrivateBonf by simulations. Unless specified, we set m = 10 5 , m = 100, q = 0.1, η = 10 -4 , ν = 0.5q/m, ε = 0.5, and δ = 0.001. To construct the p-values, we let p i = Φ(ξ i -µ) for i = 1, . . . , m 1 and p m 1 +1 , . . . , p m be i.i.d. uniform variables on (0, 1), where Φ is the CDF of N (0, 1), ξ 1 , . . . , ξ m 1 are i.i.d. copies of N (0, 1), and the default values of µ and m 1 (the number of true effects) are set to 4 and 100, respectively. In summary, the FDR of PrivateBHq is empirically controlled at q in almost all scenarios, though Theorem 5 is only concerned with FDR k for k ≥ 2. Moreover, PrivateBonf is uniformly the least powerful among the three procedures. This is not surprising given that PrivateBonf is inherently developed for FWER control. Looking closely, the performance of PrivateBHq is comparable to that of BHq when ε is not too small and η is not too large. Notably, the power of PrivateBHq deteriorates when the number of true effects exceeds 150, which is due to the truncation of the PrivateBHq procedure at m = 100.
For completeness, we refer interested readers to a set of numerical comparisons of an FDR k -controlling procedure [52,54], PrivateBHq, and PrivateBonf in the appendix.
This paper has developed a privacy-preserving multiple testing procedure termed PrivateBHq for FDR control. On the privacy side, we propose a new notion of sensitivity tailored to p-values and recognize the sequential nature of the BHq (step-down) procedure so as to keep q q q q q q q q q q 0.00
PrivateBHq PrivateBonf q q q q q q 0.00 µ Power q BHq PrivateBHq PrivateBonf q q q q q q q q 0.00 PrivateBHq efficient under the differential privacy constraint. Differential privacy of this whole pipeline follows from the composition nature of differential privacy. On the statistical side, as a major contribution of the paper, it is proved that a large class of multiple testing procedures, including the step-up, step-down, and PrivateBHq procedures, control the FDR k only provided the joint independence of the true null test statistics. A novel aspect of this result lies in the absence of any assumption on the dependence between the true nulls and false nulls. Notably, some recent progress has been made along this direction using the the FDR-linking technique [59]. Looking forward, our work raises a number of open questions. First, it would be interesting to take into account prior knowledge, such as the importance of hypotheses and beliefs about which are true nulls, into the design of a differentially private procedure. Second, it would be of interest to develop private procedures for control of other popular error rates such as the q-value [58]. Moreover, it is natural to wonder if the bound in Theorem 1 can improve by imposing some structure on the dependence between the true null and false null test statistics. Third, recognizing the vital importance of peeling in our PrivateBHq, an interesting direction is to investigate alternatives to peeling, such as the oneshot approach to the problem of private top-k selection [49]. Last, it would be interesting to consider other notions of privacy such as concentrated differential privacy and Gaussian differential privacy [28,13,18].
Finally, we wish to make a connection to a remarkable property of differential privacy: it protects against false discoveries due to adaptive data analysis, where an analysis is informed by prior interactions with the same database [22,21,3]. Adaptivity is ubiquitous in practice as the analyst is often not clear a priori what are the right questions to ask about a database. In the multiple testing context, this issue arises when hypotheses are adaptively selected based on prior discoveries. A question of great interest is to develop a multiple testing procedure that continues to preserve privacy in the presence of adaptivity.
C. Dwork is Editor-in-Chief of this Journal. As per the Journal's policies, she did not participate in the editorial process, and all aspects of the editorial workflow were blinded (not accessible) to her. into (A.2), we get
Therefore, we can assume t ≤ 7n/8. Provided that n/2 ≤ t ≤ 7n/8, we can apply Littlewood's theorem [44,46].
, where Φ(x) and φ(x) are the cumulative distribution function and density function of N (0, 1) respectively, this theorem gives
where
Next, we consider 1 2 n n t -1 By Stirling's formula and using the fact that t = (0.5 + o(1))n, we get
Thus, we get
Thus, with ν = m -1-c , we can choose η In particular,
Through exponential tilting, we get
Using saddle point approximation, we get
2θκ (θ a )nθ a Thus, we have
Next, we evaluate κ (θ a ) and θ a . Denote by µ and σ 2 the mean and variance of ζ, respectively. We get θ a = o(1) and κ (θ a ) = Var θa (ζ) = σ 2 + o(1). In particular, from (A.4) we get
.
which together with (A.5) yields
To evaluate the ratio P (na -A ≤ T < na) P (T ≥ na) , it remains to approximate P (na -A ≤ T < na). We use the local central limit theorem to do this. Explicitly, using the local central limit theorem, we get
≤ e -n(aθa-κ(θa)) E θa e θaA 1(na -A ≤ T < na)
= e θaA e -n(aθa-κ(θa)) P θa (na
where σ a is the standard deviation of ζ tilted at θ a . That is,
Finally, combing (A.6) and (A.7) gives
As such, we can choose
An example of p-value computation from [61,65]. Now we show that the contingency table example in Section 2.2 does not give p-values that are (η, ν)-sensitive with some η, ν → 0 even if n → ∞. In particular, we focus on two adjacent tables as shown in Table 1. The χ 2 -statistic of the left table is χ 2 L = 0 because a × (n/2 -a) -a × (n/2 -a) = 0 and, as a consequence, the corresponding p-value is p L ≈ P(χ 2 1 ≥ χ 2 L ) = 1. Next, for the right table the statistic equals
Now, assuming 5 ≤ a n, we get
Thus, in this example both p R and p L are bounded below away from 0 and the ratio 2Φ -1 √ 2a+1 does not tend to 1 as n → ∞. As a consequence of this, it is impossible to have both vanishing η and ν for this p-value computation.
Proof of Theorem 3. The PrivateBHq procedure acts on the intermediate results
provided by the peeling. Hence, Lemma 2.8 implies that it suffices to establish the (ε, δ)differential privacy for peeling as a part of PrivateBHq. By Lemma 2.4, each Private Min in peeling is ( Proof of Lemma 3.2. The proof is similar to that of Example 5.6.1 in [20]. By scaling, assume that ξ i are exponential random variables with parameter 1, i.e, E ξ i = 1. Note that W j is measurable with respect to F j . In the proof, we first consider the conditional expectation E(W -1 j |F j+1 ), then return to E(W j |F j+1 ) by applying Jensen's inequality. Specifically, we have
because T m+1 is measurable in F j+1 . Next, observe that by symmetry we get
To complete the proof, note that Jensen's inequality asserts that
as desired.
Proof of Theorem 2. In addition to the proof sketch in Section 3.2, it remains to show that 1 q E j j + max{ mU (j ) /q -j , 0}
→ C k as q → 0, m → ∞, m -m 0 → ∞ and m 0 /m → 1. Since j = O P (1) and mU (j ) is bounded below away from 0 as m 0 → ∞ and m 0 /m → 1, one can show that 1 max{ mU (j ) /q -j , 0} = mU (j ) /q -j + O P (1).
1 One needs to ensure that m -m0 is larger than max{ mU (j ) /q -j , 0} with high probability. Thus, q should tend to 0 slowly as m -m0 → ∞.
Thus, we get 1 q • j j + max{ mU (j ) /q -j , 0} = j mU (j )
+ o P (1)
= max k≤j≤m 0 j mU (j)
+ o P (1).
By assumption, we have
To complete the proof, the last step is to show that 1 q • j j + max{ mU (j ) /q -j , 0} is bounded by an integrable random variable. To this end, we note that if mU (j ) /q ≥ j , then 1 q • j j + max{ mU (j ) /q -j , 0} ≤ j mU (j ) ,
and otherwise q ≥ mU (j ) /j , yielding
.
Note that j mU (j )
is bounded by an integrable random variable by resorting the representation using T j .
Proof of Lemma 3.6. We first prove the case where k ≥ 3. The uniform integrability follows if we show
(A.9)
As proved by Lemma 3.2, jT n+1 nT j is a backward submartingale for j = k, k + 1, . . . , n. Thus, by Doob's maximal inequality, we get
So, the proof would be completed if we verify
To this end, note that nT j are both uniformly integrable, it is sufficiently to show the uniform integrability of 2(T n+1 -T 2 ) nT 2
for n ≥ 2. To this end, note that
which is finite. The proof is complete.
Proof of Lemma 4.2. We first consider part one. Note that
Hence, taking a union bound, we get
The proof of part two is the follows the same reasoning except using
Proof of Theorem 6. In the proof below, we replace the assumption on the nominal level with the relaxed assumption that q ≥ 6m -1.5 . Let 0 < α, α < 1 be specified later. Denote by R SD = min{R SD , m } and let p j 1 , . Then, it must hold that
where on the event (A.12) Z i l ≥ -λ log m 2 α . Therefore, we get To this end, it suffices to set α = m -0.014 and α = m -0.029 /2. Since both α, α → 0 as m → ∞, we see that (A.15) holds with probability tending to one. Recognizing (A.15), to reject all of these R SD hypotheses using PrivateBHq, it is sufficient to have
Since q ≥ 6m -1.5 , we get √ 10 log 6m q ≤ √ 10 log 6m q ≤ √ 10 log 6m 6m -1.5 < 8 log m. Hence, (A.16) is implied by
which is in fact an equality by assumption. Thus, the proof is complete.
Proof of Corollary 4.3. A careful look at the proof of Theorem 5 reveals that the event in Proposition 4.1 holds with probability at least 1 -q/12. As such, the bound on the FDR k in Theorem 5 can be strengthened to (C k + 1/12)q.
In light of the above, we set c such that
Then, PrivateBHq at level (1 + c )q controls the FDR k at level (C k + 0.1)q as ensured by Theorem 5.
It remains to prove that the claim of Theorem 6 also holds. To this end, we only need to show that q that is given in the statement of Theorem 6 is less than (1 + c )q for sufficiently large m, n. That is,
In both examples, η = n -0.5+o (1) and log m = log poly(n) = n o (1) . Thus, we have e 24η √ m log(1/δ) log m ε → 0 due to m ≤ n 1-c . Therefore, in words, PrivateBHq at level (1 + c )q should be at least as powerful as the truncated BHq step-down procedure with probability tending to one.
Proof of Proposition 4.4. We first prove that PrivateBonf is (ε, δ)-differentially private. To this end, we start by observing that each θ j corrupted by Lap( λ) noise is ε -differentially private, where
Using the Advanced Composition Theorem, therefore, it suffices to show that
Under the assumptions of Theorem 3, we have e ε -1 ≤ 1.034ε
Thus, the proof would be completed once we show
This inequality can be easily verified. Next, we turn to show the second statement. As with the proof of Theorem 5, we only need to show that, with probability at least 1 -0.1q, all noisy counts θ j + Lap( λ) with p j > q/m are above log q m -η 10m log(1/δ) log(5m/q) 2ε .
This statement is implied if m i. i. d. Lap( λ) noise terms are all above -η √ 10m log(1/δ) log(5m/q) 2ε with probability at least 1 -0.1q. This claim is true by invoking Lemma 4.2.
B.1. Comparions of FDR k -controlling procedures with others. In this section, we follow the setting of Figures 2 and3, and compare the FDR k and power of a step-up procedure [52,54], PrivateBHq, and PrivateBonf. Specifically, we consider the case of k = 3, and use the critical values
for j = 1, . . . , m, where k = 3 (see Eqn. (5.3) in [54]). For simplicity, we refer to this procedure as BHq k . Figures 4 and5 display the results. q q q q q q q q q q 0.00
PrivateBHq PrivateBonf q q q q q q 0.00
PrivateBHq PrivateBonf q q q q q q 0.00
PrivateBHq PrivateBonf q q q q q q q q 0.00 Example B.1 (Multivariate Normal). Consider observing X ∼ N (µ, Σ). The covariance Σ is constructed as follows: Σ ii = 1 for all 1 ≤ i ≤ m, Σ ij = 0 if i = j and both i, j ∈ N 0 or i, j ∈ N 1 , and Σ ij = -1/ √ m 0 m 1 if one of i, j belongs to N 0 and the other belongs to N 1 (if this value is set to be smaller than -1/ √ m 0 m 1 , the covariance Σ is not positive semidefinite). The distribution of X satisfies the IWN condition and, therefore, Theorem 1 guarantees FDR control of the BHq procedure used to test µ i = 0 against the one-sided q q q q q q q q q q 0.00
PrivateBHq PrivateBonf q q q q q q 0.00
PrivateBHq PrivateBonf q q q q q q q q 0.00 alternative µ i > 0. In contrast, the results of [7] are not applicable because the PRDS property does not hold due to -1/ √ m 0 m 1 < 0. Furthermore, Theorem 1 is still valid for testing against the two-sided alternatives µ i = 0. In general, the PRDS property is not satisfied for two-sided tests (see discussion in Section 3.1 of [7]). Figure 6 presents the empirical FDR, FDR 2 , and FDR 5 of the BHq procedure for both one-sided and two-sided alternatives in this example and, in addition, the bound C k θ 0 q in Theorem 1 for k = 2, 5 in dashed lines. As predicted by Theorem 1, the empirical FDR 2 and FDR 5 are indeed below their corresponding dashed lines. In fact, the empirical values are much below the bounds in Theorem 1, which are derived by assuming a least favorable dependence structure between the nulls and non-nulls. This pattern is also observed in the following two plots. Moreover, these empirical error rates decrease eventually as the number of true effects m 1 increases, which reflects the presence of the true null proportion θ 0 in the bound C k θ 0 q. Notably, this bound can be smaller than the nominal level q if m 1 is sufficiently large. q q q q q q q q q q 0.00 (a) one-sided tests q q q q q q q q q q 0.00 Example B.2 (Multivariate t-Distribution with Different Denominators). Consider observing i. i. d. vectors X (1) , . . . , X (n) from N (µ, Σ), where both µ and Σ are the same as the previous example. To test µ i = 0 against µ i > 0 or µ i = 0, we use the t-test statistics
, where Xi = (X (1)
/n for i = 1, . . . , m. As earlier, Theorem 1 applies to this example, as opposed to the existing FDR literature, which fails to ensure FDR control of the BHq procedure in this example.
Numerical results for Example B.2 are displayed in Figure 7. The setup follows Example B.1, with n being set to 10. While the behavior of the BHq procedure in Figure 6 basically remains the same in the present plot, we wish to point out that the effect of the true null proportion θ 0 is more pronounced in the present simulation study and the three error rates coincide exactly once m 1 exceeds 100 as the BHq in this setting always rejects a substantial number of true nulls. The latter shows the difference between the FDR and FDR k is inconsequential in this example.
Example B.3 (Multivariate Normal with Block-Diagonal Covariance). Consider bivariate normal variables X i , X i with means µ i = 0 and µ i = 0, respectively, for i = 1, . . . , m. Let q q q q q q q q q q 0.00 (a) one-sided tests q q q q q q q q q q 0.00 Var X i = Var X i = 1 for all i and the m pairs (X i , X i ) be jointly independent. Thus, the 2m normal variables exhibit a diagonal-block covariance matrix that is formed by m 2 × 2 blocks on the diagonal. The correlation corr(X i , X i ) within every block varies from -1 to -0.1.
Note that there are m true nulls among the 2m hypotheses and, therefore, θ 0 = 0.5. The IWN condition is satisfied because all true nulls are located in different blocks. Consequently, the BHq procedure maintains FDR k control in this example by applying Theorem 1, as opposed to existing results in the literature, which to our knowledge are not capable of confirming the FDR control for this example. Moreover, the usual FDR control follows from Theorem 1 as a corollary, whose proof can be found in the appendix. As an appealing feature of this result, the dependence within each block can be arbitrary and even be different across blocks.
Corollary B.4 . Fix 0 < q < 1. Assume that {1 ≤ i ≤ m : µ i ≥ c 1 }/m ≥ c 2 for positive constants c 1 , c 2 in Example B.3. For both one-sided and two-sided alternatives, the BHq procedure controls the usual FDR in an asymptotic sense. That is, as m → ∞, we get FDR ≤ (1 + o m (1))q.
The numerical results are summarized in Figure 8. Note that the three FDR variants coincide through the range of within-block correlations. Interestingly, the bound corresponding to k = 5 (the lower dashed line) is below the nominal level q = 0.1.
Proof of Corollary B.4. The proof idea is to apply Theorem 1 or Corollary 1.3 and recognize the number of rejections R in Example B.3 tends to infinity. We only consider the case of one-sided alternatives and the proof of the two-sided alternatives case is very similar. q q q q q q q q q q 0.00 (a) one-sided tests q q q q q q q q q q 0.00 Here, m is set to 5000, µ i is set to 1.5 for all i. Note that the true null proportion θ 0 = 0.5. All points represent the average of 100 independent runs. The correlation between X i and X i is set to be the same across all i, varying from -1 to -0.1.
Assume for the moment that R → ∞ with probability tending to one as m → ∞. Then, we have
which is implied by
Now, we aim to show (B.2) holds with probability tending to one. Denote by A = {i : µ i ≥ c 1 }, which, by assumption, satisfies #A ≥ c 2 m. Consider W i := X i -µ i ∼ N (0, 1) for i ∈ A and let W (1) ≥ • • • ≥ W (#A) be the order statistics. Note that (B.2) simply follows from
which is equivalent to
To prove (B.3), we make two observations: