Add C++ log-likelihood function to speed up optimization in dev_beta_binom()

[nehill197]

I wrote a non-vectorized C++ function to calculate the beta binomial log-likelihood, and used this in the call to optimize() within dev_beta_binom().
Microbenchmarking shows that optimizing using the C++ function is much faster than the previous version that optimized using the extras::log_lik_beta_binom() function.

Requires importing Rcpp.

[joethorley]

@krlmlr - I've added you as a reviewer for discussion while we are working in Vancouver.

[krlmlr]

Thanks, Nicole. Before adding native code to a package, can we explore how fast we can get with pure R? Something along the following lines might give a decent speed-up already:

log_lik_beta_binom <- function(x, size = 1, prob = 0.5, theta = 0) {
  alpha <- prob * 2 * (1 / theta)
  beta <- (1 - prob) * 2 * (1 / theta)
  lbeta_binom <- lgamma(size + 1) - lgamma(x + 1) - lgamma(size - x + 1) +
    lgamma(x + alpha) + lgamma(size - x + beta) - lgamma(size + alpha + beta) +
    lgamma(alpha + beta) - lgamma(alpha) - lgamma(beta)
  
  # FIXME: Rename bol to a more evocative name
  bol <- !is.na(x + size + prob + theta)
  if (length(bol) == 1) {
    # FIXME: Add fast path for scalar inputs
    ...
  } else if (length(bol) == 0) {
    numeric(0)
  } else {
    lbeta_binom[bol & ((x == 0 & prob == 0) | (x == size & prob == 1))] <- 0
    lbeta_binom[bol & x != 0 & prob == 0] <- -Inf
    lbeta_binom[bol & x != size & prob == 1] <- -Inf
    lbeta_binom[bol & x > size] <- -Inf
    bol_theta <- !is.na(theta)
    lbeta_binom[bol_theta & theta < 0] <- NaN
    use_binom <- bol_theta & theta == 0
    if (any(use_binom)) {
      lbinom <- log_lik_binom(x, size = size, prob = prob)
      lbeta_binom[use_binom] <- lbinom[use_binom]
    }
    lbeta_binom
  }
}

[nehill197]

Closing this as we decided to go with speeding up the log-likelihood function in R, among other things, see #69