Dear S+ers,
In order to find how accurately S+ calculates binomial probabilities I
wrote the following program (for S-Plus for Windows, Version 3.1
Release 3):
> acc.bin
function()
{
dens <- numeric(18)
half <- numeric(18)
for(i in 1:18) {
dens[i] <- dbinom(10^i - 1, 2 * 10^i - 1, 0.5)
half[i] <- pbinom(10^i - 1, 2 * 10^i - 1, 0.5)
}
cbind(dens, half)
}
> acc.bin()
dens half
[1,] 1.761971e-001 5.000000e-001
[2,] 5.634848e-002 5.000000e-001
[3,] 1.783901e-002 5.000000e-001
[4,] 5.641825e-003 5.000000e-001
[5,] 1.784122e-003 5.000000e-001
[6,] 5.641895e-004 5.000000e-001
[7,] 1.784124e-004 5.000000e-001
[8,] 5.641894e-005 5.000000e-001
[9,] 1.784122e-005 4.999994e-001
[10,] 5.642081e-006 5.000312e-001
[11,] 1.784552e-006 5.001274e-001
[12,] 5.680232e-007 5.014349e-001
[13,] 1.699281e-007 4.762227e-001
[14,] 7.795084e-008 6.854453e-001
[15,] 3.160881e-010 1.431820e-004
[16,] 2.515439e+030 -1.257719e+046
[17,] 5.537519e+305 -Inf
[18,] 0.000000e+000 1.000000e+000
Up to size=10^12 everything is OK, but for bigger values of size the
accumulated errors make the result unreliable.
1) How one can find the Fortran codes for the functions dbinom and
pbinom?
2) Do these functions use the identities of the kind
P(Sn<=m)=1-Ip(m+1,n-m) where Ix(a1, a2) is the cdf at point x of the
beta(a1,a2) random variable?
Remigijus
-----------------------------------------------------------------------
This message was distributed by s-news@wubios.wustl.edu. To unsubscribe
send e-mail to s-news-request@wubios.wustl.edu with the BODY of the
message: unsubscribe s-news
|