问题
我正在寻找尾部(1e-10 and 1 - 1e-10)中正态分布的高精度值,因为我正在使用的R包将超出此范围的任何数字设置为这些值,然后调用qnormqt函数。
我注意到的是,在看尾部时,R中的qnorm实现不是对称的。这对我来说是非常令人惊讶的,因为众所周知这种分布是对称的,而且我已经看到其他对称语言的实现。我已经检查了qt函数,它的尾部也不对称。
以下是qnorm函数的结果:

x       qnorm(x)                qnorm(1-x)              qnorm(1-x) + qnorm(x)
1e-2    -2.3263478740408408     2.3263478740408408      0.0 (i.e < machine epsilon)
1e-3    -3.0902323061678132     3.0902323061678132      0.0 (i.e < machine epsilon)
1e-4    -3.71901648545568       3.7190164854557084      2.8421709430404007e-14
1e-5    -4.2648907939228256     4.2648907939238399      1.014299755297543e-12
1e-10   -6.3613409024040557     6.3613408896974208      -1.2706634855419452e-08
很明显,在x的值接近于0或1时,此函数会崩溃。是的,在“正常”使用中这不是问题,但是我正在研究边缘情况,并将很小的概率乘以非常大的值,在这种情况下,错误(1e-08)会变成很大的值。
注意:我已经尝试使用1-x并输入实际数字0.000010.99999来进行此操作,但准确性问题仍然存在。
问题
首先,这是qnormqt实现的已知问题吗?我在文档中找不到任何内容,对于Algorithm AS 241文件中所述的10^-314,p值的算法应该是准确的16位数字。
来自R doc的引用:

如果R代码实现7位数字版本,为什么要声明16位数字呢?还是“准确”但原始算法不对称且错误?
如果R同时实现了Algorithm AS 241的两个版本,我可以打开16位数字的版本吗?
或者,R中是否有更准确的qnorm版本?
或者,解决我的问题的另一种解决方案,我需要在尾部使用高精度以实现分位数功能。
R版
>version
platform       x86_64-w64-mingw32
arch           x86_64
os             mingw32
system         x86_64, mingw32
status
major          3
minor          3.2
year           2016
month          10
day            31
svn rev        71607
language       R
version.string R version 3.3.2 (2016-10-31)
nickname       Sincere Pumpkin Patch

最佳答案

事实证明(正如Spencer Graves在his response中针对R-devel列表服务器上的相同问题所指出的那样)qnorm() 确实按照广告中的指示执行。仅仅是为了在发行版的最上层获得高度准确的结果,您需要利用函数的lower.tail参数。

这样做的方法如下:

options(digits=22)

## For values of p in [0, 0.5], specify lower tail probabilities
qnorm(p = 1e-10)                      ## x: P(X <= x) == 1e-10
# [1] -6.3613409024040557

## For values of p in (0.5, 1], specify upper tail probabilities
qnorm(p = 1e-10, lower.tail=FALSE)    ## x: P(X > x)  == 1e-10     (correct approach)
# [1] 6.3613409024040557
qnorm(p = 1 - 1e-10)                  ## x: P(X <= x) == 1-(1e-1)  (incorrect approach)
# [1] 6.3613408896974208

问题是1-1e-10(例如)易受浮点舍入错误的影响,因此它与1(间隔的上端)的距离与1e-100(间隔的下端)的距离实际上并不相同。当以更熟悉的伪装出现时,潜在的问题(它是R-FAQ 7.31!)变得很明显:
1 - (1 - 1e-10) == 1e-10
## [1] FALSE

最后,这是一个快速的确认,qnorm()可为其帮助文件中声明的值提供准确的(或至少对称的)结果:
qnorm(1e-314)
## [1] -37.906647423565666
qnorm(1e-314, lower.tail=FALSE)
## [1] 37.906647423565666

## With this failing in just the way (and for just the reason) you'd now expect
qnorm(1-1e-314)
# [1] Inf
1 == (1-1e-314)
# [1] TRUE

关于r - 从尾部的qnorm获取高精度值,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/43362644/

10-13 09:14