问题
我正在寻找尾部正态分布的高精度值(1e-10 and 1 - 1e-10)
,因为我使用的 R 包将任何超出此范围的数字设置为这些值,然后调用qnorm
and qt
功能。
我注意到的是qnorm
从尾部来看,R 中的实现并不对称。这对我来说非常令人惊讶,因为众所周知,这个分布是对称的,并且我已经看到其他语言中的对称实现。我已经检查过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.00001
and 0.99999
并且准确性问题仍然存在。
问题
首先,这是一个已知问题吗?qnorm
and qt
实施?我在文档中找不到任何内容,该算法对于 p 值应该是准确的 16 位数字10^-314
如中所述算法 AS 241 https://www.jstor.org/stable/2347330?seq=1#page_scan_tab_contents paper.
引用 R 文档:
Wichura, M. J. (1988) 算法 AS 241:正态分布的百分点。应用统计学,37, 477–484。
它提供高达约 16 位数字的精确结果。
如果R代码实现了7位数字版本,为什么它声称是16位数字?或者是“准确”但原始算法不对称且错误?
如果 R 确实实现了两个版本算法 AS 241 https://www.jstor.org/stable/2347330?seq=1#page_scan_tab_contents我可以打开 16 位版本吗?
或者有没有更准确的版本qnorm
在 R 中?
或者,我的问题的另一种解决方案,我需要分位数函数尾部的高精度。
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