## #######################
## Rat survival example ##
##########################
x <- c(52,10,40,104,51,27,146,30,46)
n <- length(x)

## Empirical likelihood calculations
nn <- 99
mu <- seq(min(x)+1e-2, max(x)-1e-2, length=nn)
l <- numeric(nn)
for (i in 1:nn) {
  f <- function(lambda, theta) {
    sum((x-theta)/(1+lambda*(x-theta)))
  }
  lambda <- uniroot(f, interval=(1-n)/(n*(range(x) - mu[i])), theta=mu[i])$root
  w <- 1/(n*(1+lambda*(x-mu[i])))
  l[i] <- sum(log(n*w))
}

## Slide 9
plot(mu,exp(l),type="l",xlab=expression(mu),ylab=expression(R(mu)))
plot(mu,2*l,type="l",xlab=expression(mu),ylab=expression("2log"~R(mu)),ylim=c(-20,0))
abline(h=-qchisq(.95,1))

## Slide 10
ind <- which(mu > mean(x))
upper <- approx(l[ind], mu[ind], -qchisq(.95, 1)/2)$y
ind <- which(mu < mean(x))
lower <- approx(l[ind], mu[ind], -qchisq(.95, 1)/2)$y
emplik <- c(lower,upper)