Hi, I tried to calculate the formula for the birthday problem (the probability that at least two people out of a group of n people share the same birthday) But the factorial-function allows me only to calculate factorials up to 170. So is there a way to push that limit? to solve this formula: (factorial(365) / factorial((365-23))) / (365^23) (n=23)
J?rg Gro? wrote:> Hi, > > I tried to calculate the formula for the birthday problem > (the probability that at least two people out of a group of n people > share the same birthday) > > But the factorial-function allows me only to calculate factorials up to > 170. > > > So is there a way to push that limit? > > to solve this formula: > > (factorial(365) / factorial((365-23))) / (365^23)Obviously you can easily rewrite this formula to: prod(343:365) / (365^23) or factorial(23) * choose(365, 23) / (365^23) Uwe Ligges> (n=23) > > ______________________________________________ > R-help at r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code.
On 28-Sep-08 17:51:55, Uwe Ligges wrote:> J?rg Gro? wrote: >> Hi, >> I tried to calculate the formula for the birthday problem >> (the probability that at least two people out of a group of >> n people share the same birthday) >> >> But the factorial-function allows me only to calculate >> factorials up to 170. >> >> So is there a way to push that limit? >> >> to solve this formula: >> >> (factorial(365) / factorial((365-23))) / (365^23) > > Obviously you can easily rewrite this formula to: > > prod(343:365) / (365^23) > > or > > factorial(23) * choose(365, 23) / (365^23) > > Uwe Ligges > >> (n=23)I would put it in an even "safer" form: n <- 23 prod( ((365-(n-1)):365)/rep(365,n) ) In other word: It evaluates (343/365)*(344/365)* ... *(365/365) 365^N --> "Inf" if N > 120, whereas n<-150 prod( ((365-(n-1)):365)/rep(365,n) ) # [1] 2.451222e-16 Best wishes, Ted. -------------------------------------------------------------------- E-Mail: (Ted Harding) <Ted.Harding at manchester.ac.uk> Fax-to-email: +44 (0)870 094 0861 Date: 28-Sep-08 Time: 19:09:40 ------------------------------ XFMail ------------------------------
On Sun, 28 Sep 2008, J?rg Gro? wrote:> Hi, > > I tried to calculate the formula for the birthday problem > (the probability that at least two people out of a group of n people share > the same birthday) > > But the factorial-function allows me only to calculate factorials up to 170. >You might want to look at the pbirthday() and qbirthday() functions and the examples on their help page. The function implements Diaconis & Mosteller's approximation, which works for different values of two as well as of n and 365. One of the examples compares to the exact calculation, which it computes as x2<- 1-sapply(10:100, function(n)prod((365:(365-n+1))/rep(365,n))) This is obviously slow compared to the formula with factorials, but it does only take two milliseconds on my laptop, so it's not *that* slow. The phrase "premature optimization" springs to mind. -thomas
Em Dom, 2008-09-28 ?s 19:43 +0200, J?rg Gro? escreveu:> Hi, > > I tried to calculate the formula for the birthday problem > (the probability that at least two people out of a group of n people > share the same birthday) > > But the factorial-function allows me only to calculate factorials up > to 170. > > > So is there a way to push that limit? > > to solve this formula: > > (factorial(365) / factorial((365-23))) / (365^23) > > (n=23)Log experession n<-23 exp(sum(log(1:365))-sum(log(1:(365-n)))-n*log(365)) [1] 0.4927028 -- Bernardo Rangel Tura, M.D,MPH,Ph.D National Institute of Cardiology Brazil