Use mutual information to find a cutoff separating two distributions

mutinfo.R

## Use mutual information to define a value separating two
## distributions.

entropy <- function(x,y){
  ## shannon entropy of x or joint entropy of x and y
  if(missing(y)){
    freqs <- table(x)/length(x)
  }else{
    stopifnot(length(x) == length(y))
    freqs <- table(paste(x,y))/length(x)
  }
  -sum(freqs*log(freqs))
}

mutinfo <- function(x,y){
  entropy(x) + entropy(y) - entropy(x,y)
}

## some pretend data
data <- data.frame(
                   V = c(rnorm(100)+5, rnorm(100)+7),
                   W = c(rep('within', 100), rep('between', 100))
                   )

## mutual information at various values of d
cuts <- sort(unique(round(data$V, 1)))
vals <- with(data, sapply(cuts, function(d) mutinfo(W, V < d)))

## D is the point of maximum mutual information - it divides "withins"
## from "betweens"
D <- cuts[which.max(vals)]

par(mfrow=c(1,2))
plot(cuts, vals, ylab="mutual information")
abline(v=D, col='red', new=FALSE)
boxplot(V ~ W, data=data)
abline(h=D, col='red', new=FALSE)