r - 小数点——R语言中0的概率值

Question

如何处理R中的p值？

我期望非常低的 p 值，例如:

1.00E-80

我需要-log10

-log10(1.00E-80)

-log10(0) 是 Inf，但 Inf 也有四舍五入的感觉。

但似乎在 1.00E-308 之后，R 产生 0。

1/10^308  
[1] 1e-308

 1/10^309 
[1] 0

lm 函数显示的 p 值的精度是否与截止点 1e-308 相同，或者它只是设计为我们需要一个截止点，我需要考虑一个不同的截止点 - 例如 1e-100(例如)将 0 替换为 <1e-100。

Answer 1

有多种可能的答案——哪一个最有用取决于上下文:

• 在一般情况下，R 确实无法存储比 .Machine$double.xmin 更接近零的浮点值。 ，它因平台而异，但通常(正如您发现的那样)大约为 1e-308 .如果你真的需要处理这么小的数字并且找不到直接在对数刻度上工作的方法，你需要搜索 Stack Overflow 或 R wiki 来寻找处理任意/扩展精度值的方法(但你可能应该尝试在对数尺度上工作——这会更容易)
• 在许多情况下，R 实际上在内部计算(自然)对数刻度上的 p 值，并且可以根据要求返回对数值，而不是在给出答案之前对它们取幂。例如，dnorm(-100,log=TRUE)给出-5000.919。您可以通过除以 log10 直接转换为 log10 比例(无需取幂然后使用 log(10) ) : dnorm(-100,log=TRUE)/log(10) =-2171，这太小而无法用 float 表示。对于 p*** (累积分布函数)函数，使用 log.p=TRUE而不是 log=TRUE . (这一点在很大程度上取决于您的特定上下文。即使您没有使用内置的 R 函数，您也可以找到一种方法来提取对数刻度的结果。)
• 在某些情况下，R 将 p 值结果呈现为 <2.2e-16即使已知更精确的值:(t1 <- t.test(rnorm(10,100),rnorm(10,80)))

打印

....
t = 56.2902, df = 17.904, p-value < 2.2e-16

但您仍然可以从结果中提取精确的 p 值

> t1$p.value
[1] 1.856174e-18

(在许多情况下，此行为由 format.pval() 函数控制)

这一切如何与 lm 一起工作的例证:

d <- data.frame(x=rep(1:5,each=10))
set.seed(101)
d$y <- rnorm(50,mean=d$x,sd=0.0001)
lm1 <- lm(y~x,data=d)

summary(lm1)将斜率的 p 值打印为 <2.2e-16 , 但如果我们使用 coef(summary(lm1)) (不使用 p 值格式)，我们可以看到该值为 9.690173e-203。

更极端的情况:

set.seed(101); d$y <- rnorm(50,mean=d$x,sd=1e-7)
lm2 <- lm(y~x,data=d)
coef(summary(lm2))

显示 p 值实际上已下溢为零。然而，我们仍然可以在对数尺度上得到答案:

tval <- coef(summary(lm2))["x","t value"]
2*pt(abs(tval),df=48,lower.tail=FALSE,log.p=TRUE)/log(10)

给出 -692.62(您可以在前面的示例中检查此方法，其中 p 值没有溢出并看到您得到与摘要中打印的相同的答案)。