>>>>> "ken" == ken knoblauch <[EMAIL PROTECTED]>
>>>>> on Mon, 23 Jan 2006 09:43:28 +0100 writes:
ken> Actually, since NaN's are also detected in na.action
ken> operations, a simpler fix might just be to use the
ken> na.rm = TRUE option of min
ken> upper <- min(n^k/(c^(k - 1)), 1, na.rm = TRUE)
Well, I liked your first fix better -- thank you for it! --
since it's always good practice to formulate such as to avoid
overflow when possible.
All things considered, I think I'd go for
upper <- min( exp(k * log(n) - (k-1) * log(c)), 1, na.rm = TRUE)
Martin
Ken> Recent news articles concerning an article from The
Ken> Lancet with fabricated data indicate that in the sample
Ken> containing some 900 or so patients, more than 200 had the
Ken> same birthday. I was curious and tried out the p and q
Ken> birthday functions but pbirthday could not handle 250
Ken> coincidences with n = 1000. The calculation of upper
Ken> prior to using uniroot produces NaN,
Ken> upper<-min(n^k/(c^(k-1)),1)
Ken> I was able to get it to work by using logs, however, as
Ken> in the following version
>> function(n, classes = 365, coincident = 2){
>> k <- coincident
>> c <- classes
>> if (coincident < 2) return(1)
>> if (coincident > n) return(0)
>> if (n > classes * (coincident - 1)) return(1)
>> eps <- 1e-14
>> if (qbirthday(1 - eps, classes, coincident) <= n)
>> return(1 - eps)
>> f <- function(p) qbirthday(p, c, k) - n
>> lower <- 0
>> upper <- min( exp( k * log(n) - (k-1) * log(c) ), 1 )
>> nmin <- uniroot(f, c(lower, upper), tol = eps)
>> nmin$root
>> }
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
