1

Can anybody give a formula for the "FDR values"? (I write "FDR values" because strictly speaking, you don't have an FDR value for each test, but an FDR value for the whole study, but nevertheles, you calculate these Q values for each test.) What I am looking is a rigorous formula, not a verbal procedure.

I have the following:

pval=c(0.07,0.08,0.09,0.10)

p.adjust(pval,method="fdr")

0.1 0.1 0.1 0.1

What formulae produce these values?

I thought that the "FDR value" is calculated as \begin{equation} Q=\frac{P_in_{tests}}{i} \end{equation} where $i$ ranks $P_i$ from smallest to largest, but this can't clearly be the case.

stats_boy
  • 46
  • 3
  • The raw, uncorrected $p_{FDR}$s in my example are: 0.28, 0.16, 0.12, 0.10. So, if you correct them from right to left, you get 0.10 for every test. I guess you mean $p_{FDR,i}^{'}=min(p_{FDR,i},p_{FDR,i+1}$, but otherwise your explanation seems to be correct. Thanks! – stats_boy Sep 05 '18 at 11:30

1 Answers1

1

There are many different procedures for correcting $p$-values for the false discovery rate. The one you used above is the Benjamini-Hochberg procedure and works as follows:

  1. Rank all $p$-values from lowest to highest, such that the most significant $p$-value has rank $i=1$, the second most significant rank $i=2$, etc.;
  2. Correct each $p$-value as follows, using the number of tests ($k$) and each $p$-value's rank ($i$):

$$p_{\text{FDR}, i} = \frac{p_{\text{original}, i} \cdot k}{i};$$

  1. Preserve the original ranks of the $p$-values as follows: $$p_{\text{FDR}, i} = \min(p_{\text{FDR}, i}, \, p_{\text{FDR}, i+1})$$

In the help page (?p.adjust), you can see that when you pick method = fdr, this procedure will be used:

Less conservative corrections are also included by [...] Benjamini & Hochberg (1995) ("BH" or its alias "fdr")

Frans Rodenburg
  • 10,376
  • 2
  • 25
  • 58
  • Based on the above discussion, I thought that this was the case. However, I tried to write an R function that does literally that: fdr.own = function(pval){ ord = order(pval) # your step 1 fdr = pval*length(pval) / ord # your step 2 for( i in (length(ord)-1):1 ){ prev = fdr[ord[i+1]] this = fdr[ord[i]] fdr[ord[i]] = min(c(prev, this)) # your step 3 } } – stats_boy Sep 20 '18 at 09:24
  • You're almost there! But to preserve order, you have to do step 3 until all $p_{\text{FDR},i}$ are less than or equal to $p_{\text{FDR},i+1}$. You can do this at once using the cumulative minimum. Check the github page from line #70: https://github.com/SurajGupta/r-source/blob/master/src/library/stats/R/p.adjust.R – Frans Rodenburg Sep 20 '18 at 09:34