# Method of p value correction for a set of dependent p values

I just got a list of 100 p-values, some of which are associated. For details, these 100 p-values could be separated into 7 groups, each group was calculated from a expression dataset. However, these datasets had 10%~60% overlap among each other, that is the pvalues inside a dataset were independent, but were not across the dataset, which bring the whole 100 p-values not independent.

So is there any solutions you know to directly and simply solve this? Some ref. pointed out to use the BH adjustment inside the R package, while choosing a loose criteria, e.g. 0.2, to define the significance level.

tmp$V7 = p.adjust(tmp$V6,method="BY")
tmpV6[1:5]
[1] 0.040461040 0.001250268 0.037407030 0.009361665 0.006866937
tmpV7[1:5]
[1] 0.31423224 0.08659423 0.29570166 0.14456543 0.13611251
tmp$V7 = p.adjust(tmp$V6,method="BH")
tmp$V6[1:5] [1] 0.040461040 0.001250268 0.037407030 0.009361665 0.006866937 tmp$V7[1:5]
[1] 0.05252837 0.01447545 0.04943072 0.02416616 0.02275313


Please take a look at this, tmp$V6 is the original p-value, while tmp$V7 is after correction.

However, it seems that after BY correction, I get even worse result than BH (or fdr) method.

-

You should definitely correct for multiple comparisons! There are a huge number of strategies for this, many of which are suitable for dependent data. The Benjamini–Hochberg–Yekutieli procedure is computationally simple approach to setting the false discovery rate and works under arbitary dependence assumptions.

In brief, you first perform all $m$ of your tests, to get a list of p-values $p_1 \ldots p_m.$ Next, you find the largest value of $k$ such that: $$p_k \le \alpha \cdot \frac{k}{m \cdot c(m)}$$ where $c(m)=1$ if the tests are independent or positively correlated, or for arbitrary dependence $c(m) = \sum_i^m \frac{1}{i} \approx \ln(m) + \gamma$, where $\gamma$ is the Euler–Mascheroni constant (0.5772...).

You then accept as significant all the $p$-values less than or equal to $p_k$, and reject the rest as non-significant.

Note that this is somewhat more liberal (i.e., has better power, but more Type I errors) than something which controls the Familywise Error Rate. However, people typically use the same $\alpha$ value (0.05, 0.01, etc). Setting $\alpha=0.2$ seems like it would be asking for a lot of spurious calls!

tmp$V7 = p.adjust(tmp$V6,method="BY") tmpV6[1:5] [1] 0.040461040 0.001250268 0.037407030 0.009361665 0.006866937 tmpV7[1:5] [1] 0.31423224 0.08659423 0.29570166 0.14456543 0.13611251 tmp$V7 = p.adjust(tmp$V6,method="BH") tmp$V6[1:5] [1] 0.040461040 0.001250268 0.037407030 0.009361665 0.006866937 tmp$V7[1:5] [1] 0.05252837 0.01447545 0.04943072 0.02416616 0.02275313 plz take a look at this, tmp$V6 is the original pvalue, while tmp$V7 is after correction. However, it seems that after BY correction, I get even worse result than BH(or fdr) method – JF.JIANG Feb 21 '13 at 18:18