coverbirthdays.R

install.packages("expm")
library(expm)

#If you have 5000 people, what are the chances every birthday is represented at least once?
#method taken from http://mathforum.org/library/drmath/view/56679.html , coded for general solution

#number of days over which people are randomly distributed
n = 365

#Set up an n*n transition matrix.
#n+1 possible states, from 0 days filled to all n days filled. 
#Since we know first person will add 1st day, we look at remaining n states
#From state i, i/n chance of repeating a day already filled, 1-(i/n) of covering a new day.
#Fill in i/n going down diagonal, then 1-(i/n) on entry directly to right of it.
transitions=mat.or.vec(n,n)
for (i in 1:(n-1)) {
  transitions[i+(i-1)*n] = i/n
  transitions[n+i+(i-1)*n] = 1-i/n
}
transitions[n*n] = 1

#Multiply matrix by self 4999 times to find effect of adding 4999 new people.
#Read across first row to find probability of ending at each point after starting with 1.
#Row 1, Col 365 gives prob of ending at 365 covered after starting with 1.
transitions%^%4999