I want to calculate the p.adjusted p-value for all genes in my dataset. I have a total of 3294 genes with p.vlue <= 0.05 but when I try to calculate the adjusted p-value all I get is just 1 for all genes. what am I doing wrong?
here me code (shown only for those 3194 genes but it gives me the same results also if I do it on the unfiltered dataset):
#load dataset # not shown
#FILTER pvalue <= 0.05
test_pvalue <- filter(test, test$p_value <= 0.05 ) #3284 total genes with p.value <= 0.05
pvalue <- test_pvalue$p_value #select column of p.value
f1 <- p.adjust(pvalue,method="bonferroni")
and I got this result (just the few 30 rows printed):
What you see is expected. The Bonferroni correction works by multiplying the pvalues with the number of tests, here 3294.
#/ Example (four values, so everything gets multiplied by 4, and 1 stays 1 of course):
p.adjust(c(0.05, 0.0005, 0.1, 1), "bonferroni")
> 0.200 0.002 0.400 1.000
That means that in order to get an adjusted p-value of 0.05 or smaller your uncorrected p-values must be at least 0.05/3294 = 0.00001517911. All larger p-values will turn non-significant if using 0.05 as cutoff. It can well be that your experiment does not have the statistical power to return these kinds of significances. Bonferroni is quite a stringent correction method. You might want to look into the Benjamini-Hochberg method, but if you lack statistical background knowledge it is generally best to either work with a statistician who gives counsel or to use software packages that do the stats work for you internally.
What kind of experiment is that you are analyzing?