const-ae
commented
7 years ago Hi sidderb,
thanks for the pull request, this is really addressing an important issue.
Unfortunately your suggested way is not ideal, because it can lead to unexpected behavior when compared to what the user would get if he ran the analysis manually.
The problem is that you call p.adjust on just a subset of the p-values with the fixed n. But this can lead to more conservative estimates of the FDR than necessary. I have created a small example that should make the point clear:
p_values <- c(0.001, 0.002, 0.1, 0.3, 0.9)
p.adjust(p_values, method = "fdr")
#> [1] 0.005 0.005 0.1666667 0.375 0.9
sapply(p_values, function(p) p.adjust(p, method = "fdr", n = length(p_values)))
#> [1] 0.005 0.010 0.500 1.000 1.000As you can see calling the p.adjust method on the list of p-values is not the same as calling it in a loop on each one individually and setting n=5.
To get around this problem one would need to refactor the calculation of the p-values and first collect all of them across groups and comparisons, and only then call p.adjust. Alternatively one could use the bonferroni correction, for which I know that the correction of a single p-value does not depend on the other p-values.
Correcting for multiple testing is an important issue in many fields where multiple tests are often ran simultaneously. I often find I want plots generated with ggsignif to reflect the multiple testing corrected p-values I generate during my analysis. Hence, I have added a new parameter, FDR, which if TRUE (default) causes the p_value generated by 'test' to be corrected for multiple testing using p.adjust(p_value, method=fdr, n=length(comparisons)).