[R] how to use mle with a defined function

Ben Bolker bolker at ufl.edu
Tue Jul 3 17:10:26 CEST 2007


Lin Pan <linpan1975 <at> yahoo.com> writes:

> 
> 
> Hi all,
> 
> I am trying to use mle() to find a self-defined function. Here is my
> function:
> 
> test <- function(a=0.1, b=0.1, c=0.001, e=0.2){
> 
> # omega is the known covariance matrix, Y is the response vector, X is the
> explanatory matrix
> 
> odet = unlist(determinant(omega))[1]
> 
> # do cholesky decomposition
> 
> C = chol(omega)
> 
> # transform data
> 
> U = t(C)%*%Y
> WW=t(C)%*%X
> 
> beta = lm(U~W)$coef
> 
> Z=Y-X%*%beta
> V=solve(t(C), Z)
> 
> 0.5*odet + 0.5*(t(V)%*%V)
> 
> }
> 
> and I am trying to call mle() to calculate the maximum likelihood estimates
> for function (0.5*odet+0.5*(t(V)%*%V)) by
> 
> result = mle(test, method="Nelder-Mead")
> 
> But I get the following error message:
> 
> Error in optim(start, f, method = method, hessian = TRUE, ...) : 
>         (list) object cannot be coerced to 'double'


  Can you give a self-contained example (e.g., make up a small
matrix?)  There are a few places where I don't understand what
you're doing:

  - is the WW above a typo for W?

  - if omega is fixed, why calculate its (log) determinant
every time the function is called?  (this shouldn't change
the answer, just slow things down)
 
  - is X a single vector or a design matrix?
it seems like you might want lm(U~W-1) instead of lm(U~W)

  - it doesn't look like your parameters are used in the
function at all -- are they really parameters for calculating
omega?  [because of this, I can get as far as computing
the Hessian, and then I get "system is exactly singular"]

   [ off topic: somewhat disturbingly, the CAPTCHA
word for posting from Gmane was mother f***ers ... ]



More information about the R-help mailing list