Why Compare the Two Adjustments?
In a genomic analysis, the same family of tests can produce two different adjusted values: the Benjamini–Hochberg (BH) adjusted p-values and Storey q-values. Comparing them is useful because it separates two ingredients in an FDR analysis: the ordered observed p-values and the estimate of how many tests are truly null. BH uses the full number of tests as its scaling factor; Storey’s method can use an estimate of the null proportion below one.
That distinction can matter for a discovery list. For a fixed set of p-values and a fixed Storey estimate, q-values may be smaller than the corresponding BH-adjusted p-values. But the difference is not evidence by itself that one analysis is correct and the other is not. It depends on the estimate, the assumptions behind the methods, and the chosen decision threshold.
What Is Being Estimated?
The hypothesis family must be the same for both calculations: it contains the same \(m\) tests, with the same p-values and the same ordering. In the BH calculation, the multiplier is \(m\). Storey’s approach estimates the proportion \(\pi_0=m_0/m\) of tests whose null hypotheses are true, where \(m_0\) is the unknown number of true nulls. A common estimate at a chosen \(\lambda\in[0,1)\) is
Under a true null hypothesis, a valid continuous p-value is uniform on \([0,1]\), so its probability of exceeding \(\lambda\) is \(1-\lambda\). The denominator therefore represents the expected number of null p-values in that upper tail if all tests were null; the observed tail count is used to estimate the null proportion. The estimate is random and can be sensitive to \(\lambda\). In particular, estimating \(\pi_0\) below one does not establish that the estimate is accurate.
The standard BH FDR guarantee applies under its stated conditions, including independent valid p-values. Storey’s estimated q-values also rely on valid p-values and a suitable null-proportion estimate; using a smaller estimate does not, by itself, establish a guarantee under every dependence pattern or finite-sample setting. Here the comparison is an arithmetic comparison of two adjustments on a common list, not a claim that their guarantees are interchangeable.
The Exact Relationship for a Fixed Estimate
When both calculations use the same ordered p-values and Storey’s method uses one fixed \(\hat{\pi}_0\le 1\), its q-values are a rescaling of BH-adjusted p-values. To see this, write \(A_i=\min_{j\ge i}mp_{(j)}/j\). Then \(p^{BH}_{(i)}=\min(1,A_i)\), while \(q_{(i)}=\min(1,\hat{\pi}_0 A_i)\). Since \(\hat{\pi}_0\le1\),
This identity also explains why taking the running minimum does not break the scaling relationship: multiplying every candidate \(mp_{(j)}/j\) by the same nonnegative constant multiplies their minimum by that constant. If \(\hat{\pi}_0=1\), the two adjusted values coincide. If \(0<\hat{\pi}_0<1\), every Storey q-value is that fraction of its BH-adjusted counterpart. The ranks and ordering remain unchanged; the adjusted values and, at some thresholds, the rejection set can differ.
Worked Example: Calculate Both Adjustments on One List
Suppose ten hypothetical tests have the following ordered p-values: \(0.001, 0.004, 0.009, 0.015, 0.025, 0.040, 0.31, 0.48, 0.62, 0.81\). At \(\lambda=0.50\), exactly two p-values exceed \(\lambda\), so \(\hat{\pi}_0=\min\{1,2/[10(1-0.50)]\}=0.40\). All p-values are treated as valid for this illustration, and the same ten-test family is used for both methods.
For each rank, first calculate \(mp_{(i)}/i\), then take the minimum of that value and all later-rank values, capped at 1. In this list those ratios increase with rank, so the running minimum from rank \(i\) onward is the ratio at rank \(i\). The calculations are:
| Rank \(i\) | \(p_{(i)}\) | BH adjustment \(10p_{(i)}/i\) | Storey q-value \(0.40(10p_{(i)}/i)\) |
|---|---|---|---|
| 1 | 0.001 | 0.0100 | 0.00400 |
| 2 | 0.004 | 0.0200 | 0.00800 |
| 3 | 0.009 | 0.0300 | 0.0120 |
| 4 | 0.015 | 0.0375 | 0.0150 |
| 5 | 0.025 | 0.0500 | 0.0200 |
| 6 | 0.040 | 0.0667 | 0.0267 |
| 7 | 0.31 | 0.4429 | 0.1771 |
| 8 | 0.48 | 0.6000 | 0.2400 |
| 9 | 0.62 | 0.6889 | 0.2756 |
| 10 | 0.81 | 0.8100 | 0.3240 |
For example, at rank 7 the unadjusted BH ratio is \(10(0.31)/7=0.442857142\ldots\). The Storey calculation using the unrounded ratio is \(0.40\times[10(0.31)/7]=0.177142857\ldots\), which rounds to 0.1771. Rounding the ratio before multiplying can give a slightly different final digit.
At a target level \(\alpha=0.05\), BH rejects ranks 1 through 5 because their BH-adjusted values are at most 0.05. Storey’s q-values are at most 0.05 through rank 6, so the Storey calculation includes six tests. The difference comes from the chosen \(\hat{\pi}_0=0.40\), not from any change in the input p-values or their ranks.
How the Estimate Changes the Comparison
The estimate of \(\pi_0\) depends on \(\lambda\). A larger \(\lambda\) uses a narrower upper tail, so fewer p-values contribute to its count and the estimate may be more variable. The next calculation uses the same hypothetical p-value list to show how changing \(\lambda\) changes the Storey q-values without changing BH.
Worked Example: A Second Lambda Gives a Different Scaling Factor
Keep the ten p-values above but use \(\lambda=0.60\). Exactly two p-values exceed 0.60, so \(\hat{\pi}_0=\min\{1,2/[10(1-0.60)]\}=0.50\). At rank 7, the Storey q-value is calculated from the unrounded ratio:
At \(\lambda=0.50\), the corresponding value was \(0.40\times[10(0.31)/7]=0.177142857\ldots\approx0.1771\). The BH-adjusted value at that rank remains \(10(0.31)/7=0.442857142\ldots\approx0.4429\) in both comparisons. The two Storey estimates produce different q-values because they use different estimated null proportions, not because the observed p-value at rank 7 changed.
At \(\alpha=0.05\), both estimates include rank 6 and exclude rank 7 in this list: their rank-6 q-values are \(0.50(10(0.040)/6)=0.0333\) and \(0.40(10(0.040)/6)=0.0267\), while both rank-7 values exceed 0.05. Thus a change in q-values need not change the rejection set at every threshold.
Worked Example: The Decision Can Depend on the Target Level
Use the original estimate \(\hat{\pi}_0=0.40\) and the same list, but set the target level to \(\alpha=0.20\). At rank 7, the Storey q-value is approximately 0.1771, which is no larger than 0.20. At rank 8, the q-value is \(0.40[10(0.48)/8]=0.24\), which is larger than 0.20. Therefore Storey’s q-value rule includes ranks 1 through 7.
For BH, the rank-6 adjusted value is \(10(0.040)/6=0.0667\), while the rank-7 value is \(10(0.31)/7=0.4429\). The BH-adjusted values at ranks 7 and above exceed 0.20, so BH includes ranks 1 through 6. The extra rank in the Storey set follows directly from the smaller multiplier. It does not mean that the probability rank 7 is a false discovery is 0.1771.
Interpretation, Diagnostics, and Limitations
A BH-adjusted p-value is the smallest BH target level at which the corresponding test would be rejected by the step-up rule. A Storey q-value is the minimum estimated FDR over thresholds that include that test, under the chosen estimation approach. Both are useful for ranking discoveries and applying a specified threshold, but a q-value is not the posterior probability that its particular null hypothesis is true, nor the probability that this particular result is false.
A practical comparison should make the following checks:
- Confirm the family. Both methods must use the same \(m\) tests and the same p-values. Filtering tests before one calculation but not the other invalidates a direct comparison.
- Confirm the ordering and monotonicity. Sort the p-values once. Adjusted values should be nondecreasing with rank after the reverse cumulative minimum is applied.
- Confirm the null-proportion estimate. Record the chosen \(\lambda\), the number of p-values exceeding it, and any cap at 1. A value below one should be reported as an estimate, not as a known fact.
- Check for numerical consistency. With one fixed estimate, each Storey q-value should equal \(\hat{\pi}_0\) times the corresponding BH-adjusted value, up to rounding.
- Consider sensitivity. If reasonable choices of \(\lambda\) produce substantially different estimates or discovery sets, that instability should qualify the interpretation.
A common mistake is to say, “Storey found fewer nulls, so each reported discovery has only a 4% chance of being false.” This is incorrect: the estimated proportion \(\hat{\pi}_0\) is a family-level ingredient in an estimated FDR calculation, not a per-test posterior probability. A second mistake is to compare Storey q-values from one family with BH adjustments from another and attribute the difference to \(\hat{\pi}_0\). The correction is to hold the test family fixed, state the Storey estimate and its tuning choice, and interpret any thresholded set as a set selected by that procedure.
Storey’s estimate can be noisy, especially when the upper tail contains few p-values; its value can also depend on model misspecification or non-uniform null p-values. In a genomic application, the calibration and dependence structure of the tests matter to error-rate guarantees. The numerical scaling identity remains valid for a fixed list and estimate, but it does not settle whether the estimate is appropriate or whether a particular procedure’s theoretical conditions hold.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions.
- For \(m=8\), suppose \(\hat{\pi}_0=0.75\) and the BH-adjusted value at one rank is 0.24. What is the Storey q-value at that rank?
- In a family of 20 tests, four p-values exceed \(\lambda=0.80\). Calculate the Storey estimate using the stated formula.
- For the hypothetical list in this tutorial, why do BH and Storey give the same rank ordering even when their adjusted values differ?
- At \(\alpha=0.05\), the rank-7 Storey q-value in the main example is 0.1771. Is rank 7 selected, and does this mean its null hypothesis has a 17.71% probability of being true?
Solutions. (1) The fixed-estimate scaling relation gives \(0.75(0.24)=0.18\). (2) \(\hat{\pi}_0=\min\{1,4/[20(1-0.80)]\}=\min\{1,4/4\}=1\). (3) Both adjustments are nondecreasing functions of the same ordered p-values, and multiplying all BH-adjusted values by the same positive constant preserves their order; it changes their scale, not their rank ordering. (4) Rank 7 is not selected at \(\alpha=0.05\), since \(0.1771>0.05\). The q-value is not the probability that this individual null hypothesis is true; it is a minimum estimated family-level FDR for thresholds that include the test.