## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

## -----------------------------------------------------------------------------
library(multChernoff)

## -----------------------------------------------------------------------------
tailProbBound(x = 20, k = 7, n = 50)

## -----------------------------------------------------------------------------
criticalValue(k = 10, n = 40, p = 0.05)

## -----------------------------------------------------------------------------
library(plyr)

# Corbet butterfly data
# https://en.wikipedia.org/wiki/Unseen_species_problem
corbet.butterfly <- data.frame(j=1:15, n_j=c(118,74,44,24,29,22,20,19,20,15,12,14,6,12,6))
n.butterfly <- Reduce(c, alply(corbet.butterfly, 1, function(.df) rep(.df$j, .df$n_j)))

alpha <- 0.05
n.observed <- c(n.butterfly, 0)    # the last one is the unseen
n <- sum(n.observed)
k <- length(n.observed)
p.observed <- n.observed / n

# critical value
t.alpha <- criticalValue(k, n, p=alpha, verbose = TRUE)
cat(sprintf("critical value = %f\n", t.alpha))

## ----eval = FALSE-------------------------------------------------------------
# # convex program
# p <- Variable(k)
# obj <- p[k]
# constr <- list(p>=0,
#                sum(p) == 1,
#                2 * n * sum(p.observed * (log(p.observed) - log(p))) <= t.alpha)
# prob <- Problem(Maximize(obj), constr)
# result <- solve(prob, verbose=FALSE) # you should try other solvers
# 
# # result
# print(result$status)
# unseen <- result$value
# 
# cat(sprintf("unseen <= %f\n", unseen))