Cumpleaños

Mi solución al otro problema del cumpleaños

Pues eso, que me piqué —y parte de la culpa la tiene este sujeto— con el otro problema del cumpleaños y he aquí el código —exacto salvo redondeos, no mediante simulaciones— que he usado para resolverlo:

f <- function(n, k = 365, v = NULL){

  if(is.null(v))
    v <- c(1, rep(NA, k))

  res <- 1

  for(j in (k-1):1){
    v[k-j] <- ifelse( is.na(v[k-j]), f(n, k-j, v), v[k-j])
    res    <- res - choose(k,j) * ((k-j)/k)^n * v[k-j]
  }

  res
}

f(2287)
#0.5003708
f(2286)
#0.4994142