p.adjust for bonferroni p value correction does not appear to be working correctly:> p <- runif(50) > > p[1] 0.08280254 0.08955706 0.19754389 0.52812033 0.68907772 0.21849696 0.02774784 0.23923562 0.03482480 0.76437481 0.87236155 0.76438604 [13] 0.37432032 0.89630318 0.01626565 0.08152060 0.55715478 0.47736921 0.77968275 0.17388127 0.37212900 0.18363170 0.51655538 0.14526733 [25] 0.60870820 0.13752392 0.92412799 0.73045115 0.89887453 0.33744577 0.84966571 0.97797283 0.20571554 0.29115022 0.75928867 0.12929511 [37] 0.64923057 0.68168196 0.19311014 0.83818106 0.85592243 0.56276287 0.80822911 0.53377044 0.44691466 0.59198417 0.36114259 0.35431768 [49] 0.10381694 0.57738429> p.adjust(p, method = "bonferroni", n=50)[1] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 [15] 0.8132824 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 [29] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 [43] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 for the one corrected pvalue that is less than 1, it seems to have been divided by 0.02 and not 50. I can do this one manually but would be nice if it worked:) [[alternative HTML version deleted]]
> p.adjust for bonferroni p value correction does not appear to be working > correctly:You should re-check what a Bonferroni correction does, or at least reboot your intuition (mine needs rebooting all the time). All the p-values should _increase_ by a factor of n, with a ceiling of 1.0. Dividing by n would imply incorrectly that individual events have become less probable as the number increases. The result you have obtained is what is supposed to happen. S Ellison> > p <- runif(50) > > > > p > [1] 0.08280254 0.08955706 0.19754389 0.52812033 0.68907772 0.21849696 > 0.02774784 0.23923562 0.03482480 0.76437481 0.87236155 0.76438604 [13] > 0.37432032 0.89630318 0.01626565 0.08152060 0.55715478 0.47736921 > 0.77968275 0.17388127 0.37212900 0.18363170 0.51655538 0.14526733 [25] > 0.60870820 0.13752392 0.92412799 0.73045115 0.89887453 0.33744577 > 0.84966571 0.97797283 0.20571554 0.29115022 0.75928867 0.12929511 [37] > 0.64923057 0.68168196 0.19311014 0.83818106 0.85592243 0.56276287 > 0.80822911 0.53377044 0.44691466 0.59198417 0.36114259 0.35431768 [49] > 0.10381694 0.57738429 > > > p.adjust(p, method = "bonferroni", n=50) > [1] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > 1.0000000 > 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > [15] 0.8132824 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > 1.0000000 > 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > [29] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > 1.0000000 > 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > [43] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > 1.0000000 > 1.0000000 > > for the one corrected pvalue that is less than 1, it seems to have been divided > by 0.02 and not 50. > > I can do this one manually but would be nice if it worked:) > > [[alternative HTML version deleted]] > > ______________________________________________ > R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code.******************************************************************* This email and any attachments are confidential. Any use...{{dropped:8}}
Part of the confusion is that Bonferroni is often described as an adjustment to the significance level (sig-level) not the p-value. For example, to evaluate p-values when there are 8 tests, we compare the p-value to the sig-level/8 so the sig-level decreases. The p.adjust() function adjusts the p-value instead of the sig-level so we multiply the p-value by the number of tests (p-value * 8). As a result the p-values increase. ------------------------------------- David L Carlson Department of Anthropology Texas A&M University College Station, TX 77840-4352 -----Original Message----- From: R-help [mailto:r-help-bounces at r-project.org] On Behalf Of S Ellison Sent: Wednesday, August 3, 2016 6:05 AM To: Alicia Ellis; r-help at r-project.org Subject: Re: [R] p.adjust not working correctly?> p.adjust for bonferroni p value correction does not appear to be working > correctly:You should re-check what a Bonferroni correction does, or at least reboot your intuition (mine needs rebooting all the time). All the p-values should _increase_ by a factor of n, with a ceiling of 1.0. Dividing by n would imply incorrectly that individual events have become less probable as the number increases. The result you have obtained is what is supposed to happen. S Ellison> > p <- runif(50) > > > > p > [1] 0.08280254 0.08955706 0.19754389 0.52812033 0.68907772 0.21849696 > 0.02774784 0.23923562 0.03482480 0.76437481 0.87236155 0.76438604 [13] > 0.37432032 0.89630318 0.01626565 0.08152060 0.55715478 0.47736921 > 0.77968275 0.17388127 0.37212900 0.18363170 0.51655538 0.14526733 [25] > 0.60870820 0.13752392 0.92412799 0.73045115 0.89887453 0.33744577 > 0.84966571 0.97797283 0.20571554 0.29115022 0.75928867 0.12929511 [37] > 0.64923057 0.68168196 0.19311014 0.83818106 0.85592243 0.56276287 > 0.80822911 0.53377044 0.44691466 0.59198417 0.36114259 0.35431768 [49] > 0.10381694 0.57738429 > > > p.adjust(p, method = "bonferroni", n=50) > [1] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > 1.0000000 > 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > [15] 0.8132824 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > 1.0000000 > 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > [29] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > 1.0000000 > 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > [43] 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 1.0000000 > 1.0000000 > 1.0000000 > > for the one corrected pvalue that is less than 1, it seems to have been divided > by 0.02 and not 50. > > I can do this one manually but would be nice if it worked:) > > [[alternative HTML version deleted]] > > ______________________________________________ > R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code.******************************************************************* This email and any attachments are confidential. Any use...{{dropped:9}}