On Fri, 14 Sep 2007, [EMAIL PROTECTED] wrote:
I believe that this may be more appropriate here in r-devel than in r-help.
The normal hazard function, or reciprocal Mill's Ratio, may be obtained
in R as dnorm(z)/(1 - pnorm(z)) or, better, as dnorm(z)/pnorm(-z) for
small values of z. The latter formula breaks dowm numerically for me
(running R 2.4.1 under Windows XP 5.1 SP 2) for values of z near 37.4
or greater.
OK,
mr <- function(z )dnorm( z )/ ( pnorm(z,lower=FALSE) )
mr.2 <- function(z) exp(dnorm( z, log=TRUE ) - (pnorm(z, lower=FALSE, log=TRUE
)))
mr(seq(10,100,by=10)) # breaks before 40
[1] 10.09809 20.04975 30.03326 NaN NaN NaN NaN NaN
NaN NaN
mr.2(seq(10,100,by=10)) #seems robust
[1] 10.09809 20.04975 30.03326 40.02497 50.01998 60.01666 70.01428
80.01250 90.01111 100.01000
mr.2(1e4)
[1] 10000
mr.2(1e5)
[1] 99999.97
mr.2(1e6)
[1] 999980.2
mr.2(1e7)
[1] 9990923
mr.2(1e8) # Oh well!
[1] 65659969
I guess you get large than 1e4, you should just use the asymptotic result.
HTH,
Chuck
Looking at the pnorm documentation I see that it is based on Cody (1993)
and thence, going one step further back, on Cody (1969). Consulting
Cody (1969) I see that the algorithm for pnorm(z) [or actually erf(z)]
is actually based on rational function approximations for the
reciprocal Mill's ratio itself, as I rather expected.
I wonder if anyone has dug out a function for the reciprocal Mill's
ratio out of the pnorm() code? Anticipating the obvious response I don't
believe that this would be one of the things I might be good at!
Murray Jorgensen
References
Cody, W. D. (1993)
Algorithm 715: SPECFUN A portable FORTRAN package of special function
routines and test drivers.
ACM Transactions on Mathematical Software 19, 2232.
Cody, W. D. (1969)
Rational Chebyshev Approximations for the Error Function
Mathematics of Computation, Vol. 23, No. 107. (Jul., 1969), pp. 631-637.
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Charles C. Berry (858) 534-2098
Dept of Family/Preventive Medicine
E mailto:[EMAIL PROTECTED] UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
