[R] Contour plot of bivariate log-likelihood function

[Dale Steele]

Thanks to prior help from the list, I've made progress in adapting a 
univariate log-likelihood function to a bivariate one for Box-Cox 
transformations.  However, I'm having trouble plotting values of the 
log-likelihood function (z) for given range of lambda(1) and lambda(2).

The current result of my function (below) is z.  How can I adapt it to 
output the values lambda(1), lambda(2) (for each row of all_lambda) and 
z and then plot using contour()?

I'm also stuck on how to allow the variable xL to == log(x) when 
lambda(1)==0 and/or lambda(2)==0.

Thanks.

--Dale


x<-read.table("http://www.stat.lsu.edu/faculty/moser/exst7037/jwdata/T4-1.txt",col.names="x1")
x<-read.table("http://www.stat.lsu.edu/faculty/moser/exst7037/jwdata/T4-5.txt",col.names="x2")
X <- as.matrix(cbind(x1,x2))

  bvboxcox <- function(X, min, max, step)
   {
     seq <- seq(min,max,step)
     all_lambda <- expand.grid(seq,seq)
     all_lambda <- as.matrix(all_lambda)
     c <- nrow(all_lambda)
     res <- numeric(c)

     for (i in seq(1:c)){
       lambda <- all_lambda[i,]
       xL <- (X^lambda - 1)/lambda
       n <- nrow(X)
       t1 <- (-n/2)*log(det(cov(xL)))
       gsum <- apply(X,2,function(x) sum(log(x)))
       t2 <- sum((lambda - 1) * gsum)
       res[i] <- t1+t2
     }
   res
   }

## For example ...

bvboxcox(X, 0.10, 0.17, 0.01)