1

I want to program the maximum likelihood of a gamma distribution in R; until now I have done the following:

library(stats4)
x<-scan("http://www.cmc.edu/pages/faculty/MONeill/Math152/Handouts/gamma-arrivals.txt")
loglike2<-function(LL){
alpha<-LL$a
beta<-LL$b
(alpha-1)*sum(log(x))-n*alpha*log(beta)-n*lgamma(alpha)}

mle(loglike2,start=list(a=0.5,b=0.5))

but when I want to run it, the following message appear:

Error in mle(loglike2, start = list(a = 0.5, b = 0.5)) : 
  some named arguments in 'start' are not arguments to the supplied log-likelihood

What am I doing wrong?

halfer
  • 19,824
  • 17
  • 99
  • 186
Little
  • 3,363
  • 10
  • 45
  • 74

2 Answers2

0

From the error message it sounds like mle needs to be able to see the variable names listed in start= in the function call itself.

loglike2<-function(a, b){
    alpha<-a
    beta<-b
    (alpha-1)*sum(log(x))-n*alpha*log(beta)-n*lgamma(alpha)
}

mle(loglike2,start=list(a=0.5,b=0.5))

If that doesn't work you should post a reproducible example with all variables defined and also explicitly indicate which package the mle function is coming from.

Community
  • 1
  • 1
MrFlick
  • 195,160
  • 17
  • 277
  • 295
  • thanks @MrFlick, I have updated the question with data to work on – Little Feb 22 '15 at 14:21
  • You still haven't indicated what package you are getting the `mle` function from. Plus the `n` variable is still undefined. – MrFlick Feb 22 '15 at 14:23
  • sorry my bad, now it is included – Little Feb 22 '15 at 14:25
  • I've updated my answer. Looks like `mle` didn't want a call, it actually wanted a function object. You need to make sure the parameter of your function match the names you pass in `start=`. You can't pass a list because `mle` wont look "inside" your function to see how that list is being used, – MrFlick Feb 22 '15 at 14:29
  • thanks, but still got an error: Error in optim(start, f, method = method, hessian = TRUE, ...) : non-finite finite-difference value [2] – Little Feb 22 '15 at 14:31
  • Can you explain what is n ? – Metrics Feb 22 '15 at 16:00
0

The error message is unfortunately criptic because it indicates mising values owing to the fact that alpha and gamma have to be positive and mle optimizes over the real numbers. Hence, you need to transfomt the vector over which the function is being optimized, like so:

library(stats4)
x<-scan("http://www.cmc.edu/pages/faculty/MONeill/Math152/Handouts/gamma-arrivals.txt")

loglike<-function(alpha,beta){
    (alpha-1)*sum(log(x))-n*alpha*log(beta)-n*lgamma(alpha)
}

fit  <-  mle(function(alpha,beta)
                 # transfrom the parameters so they are positive
                 loglike(exp(alpha),exp(beta)),
             start=list(alpha=log(.5),beta=log(.5)))

# of course you would have to exponentiate the estimates too.
exp(coef(fit1))

note that the error now is that you are using n in loglike() which you have not defined. If you define n, then you get an error stating Lapack routine dgesv: system is exactly singular: U[1,1] = 0. which is caused either by a not very good guess for the start value of alpha and beta or (more likely) that loglike() does not have a minima (I think your deleted post from last night had a slightly different formula which I was able to get working, but not able to respond to b/c the post was deleted...)

FYI, if you want to inspect the alpha and beta parameters that cause the errors, you can use scoping assignment to post the most recently called parameters to the environment in which loglike() is defined as in:

loglike<-function(alpha,beta){
    g <<- c(alpha,beta)
    (alpha-1)*sum(log(x))-n*alpha*log(beta)-n*lgamma(alpha)
}
Jthorpe
  • 9,756
  • 2
  • 49
  • 64