37

我在 R 中发现了 t-tests 和 chi-squared 这个问题，但我认为这个问题通常适用于其他测试。如果我做：

a <- 1:10
b <- 100:110
t.test(a,b)

我得到：t = -64.6472, df = 18.998, p-value < 2.2e-16。我从评论中知道-2.2e-16的值是.Machine$double.eps- 的最小浮点数1 + x != 1，但当然 R 可以表示比这小得多的数字。我还从 R FAQ 知道 R 必须将浮点数舍入到 53 位二进制精度：R FAQ。

几个问题：（1）我将其读取为 53 位精度的二进制数字是否正确，或者 R 中的值< .Machine$double.eps是否计算不准确？(2) 为什么在进行此类计算时，R 没有提供一种方法来显示 p 值的较小值，即使精度有所损失？(3) 有没有办法显示更小的 p 值，即使我失去了一些精度？对于单次测试，2 个十进制有效数字就可以了，对于我要 Bonferroni 正确的值，我需要更多。当我说“失去一些精度”时，我认为 < 53 个二进制数字，但是 (4) 我完全弄错了，任何 p 值< .Machine$double.eps都非常不准确吗？(5) R 只是诚实而其他统计数据包不是吗？

在我的领域，非常小的 p 值是常态，一些例子： http : //www.ncbi.nlm.nih.gov/pubmed/20154341，http: //www.plosgenetics.org/article/info%3Adoi%2F10 .1371%2Fjournal.pgen.1002215这就是为什么我要表示如此小的 p 值。

感谢您的帮助，对于如此曲折的问题感到抱歉。

4

6 回答 6

22

在这里交换答案和评论时，我对几件事感到困惑。

首先，当我尝试 OP 的原始示例时，我得到的p值没有这里讨论的那些小（几个不同的 2.13.x 版本和 R-devel）：

a <- 1:10
b <- 10:20
t.test(a,b)
## data:  a and b 
## t = -6.862, df = 18.998, p-value = 1.513e-06

其次，当我使组之间的差异更大时，我确实得到了@eWizardII 建议的结果：

a <- 1:10
b <- 110:120
(t1 <- t.test(a,b))
# data:  a and b 
# t = -79.0935, df = 18.998, p-value < 2.2e-16
#
> t1$p.value
[1] 2.138461e-25

打印输出的行为t.test是由它的调用驱动的stats:::print.htest（它也被其他统计测试函数调用chisq.test，如 OP 所指出的），它又调用format.pval，它显示的p值小于它的值eps（即.Machine$double.eps默认情况下）为< eps. 我很惊讶地发现自己不同意这些普遍精明的评论者......

最后，虽然担心一个非常小的p值的精确值似乎很愚蠢，但 OP 是正确的，这些值经常被用作生物信息学文献中证据强度的指标——例如，一个人可能会测试 100,000 个候选基因并查看结果p值的分布（搜索“火山图”以获取此类程序的一个示例）。

于 2011-08-14T14:31:31.467 回答
13

两个问题：

1) 1e-16 和 1e-32 的 p 值在统计含义上可能存在什么差异？如果你真的可以证明它是正确的，那么使用记录的值就是要走的路。

2）当您对 R 的数值准确性感兴趣时，为什么要使用维基百科？

R-FAQ 说“其他 [意味着非整数] 数字必须四舍五入到（通常）53 位二进制数字的准确性。” 16 位数字大约是极限。这是在控制台时如何获得准确性的限制：

> .Machine$double.eps
[1] 2.220446e-16

当在 [0,1] 范围内解释时，该数字实际上为零

于 2011-08-07T04:43:29.900 回答
11

您链接到的维基百科页面是 R 不使用的 Decimal64 类型——它使用标准问题双打。

首先，来自.Machine帮助页面的一些定义。

double.eps：最小的正浮点数 'x' 使得 '1 + x != 1'。...通常是“2.220446e-16”。

double.xmin：最小的非零归一化浮点数......通常是'2.225074e-308'。

因此，您可以表示小于 2.2e-16 的数字，但它们的准确性会降低，并且会导致计算问题。尝试一些数字接近最小可表示值的示例。

2e-350 - 1e-350
sqrt(1e-350)

您在评论中提到您想要进行 bonferroni 更正。我建议您改用自己的代码，而不是为此滚动您自己的代码p.adjust(your_p_value, method = "bonferroni")。 pairwise.t.test使用这个。

于 2011-08-07T12:57:58.773 回答
9

尝试这样的事情，t.test(a,b)$p.value看看是否可以为您提供所需的准确性。我相信它与结果的打印有关，而不是与实际存储的计算机值有关，而实际存储的计算机值应该具有必要的精度。

于 2011-08-07T04:23:15.360 回答
5

一些 R 包解决了这个问题。最好的方法是通过包 pspearman。

source("http://www.bioconductor.org/biocLite.R")
biocLite("pspearman")
library("pspearman")
a=c(1:110,110)
b=1:111
out <- spearman.test(a, b, alternative = "greater", approximation="t-distribution")
out$p.value

[1] 3.819961e-294

于 2012-09-27T08:27:34.033 回答
2

最近有同样的问题。统计学家建议：

A <- cor.test(…)
p <- 2* pt(A$statistic,  df = A$parameter, lower.tail=FALSE)
于 2013-12-10T17:05:28.500 回答