CODEWITHSUNDEEP
roptimizationtypescoerce
user3907725
user3907725

Reputation: 11

Error while using optim function on likelihood equation on R

I'm a beginner to R and ran into an error while trying to use the optim function.

I have a likelihood equation that I'd like to maximize, so I've implemented the following code:

>datafile=read.delim("clipboard")

> log.lik=function(theta, x, mu){
+ b=theta[1]
+ beta=theta[2]
+ tau=theta[3]
+ result=function(b,beta, tau){(-sum(x)/b)-sum(log(-beta*(x-tau)))-sum(log(integrate(exp(-x/b)/(1+exp(-beta(x-tau)))), lower=1500, upper=Inf))}
+ return(result)
+ }
> estimate=c(1,1,1)
> model=optim(par=estimate, fn=log.lik, control=list(fnscale=-1), x=datafile, mu=1500)

Everything works until the optim function at which I get the following error message: Error in optim(par = estimate, fn = log.lik, control = list(fnscale = -1), : cannot coerce type 'closure' to vector of type 'double'

Would anyone have an idea of what might be the issue here? Any help would be appreciated!

The data file is just one column of simulated financial losses in a csv format. When I output the datafile variable here's a sample of what I get:

       X1946774
1      34949037
2     734018898
3     393502463
4     388573133
5      93213300
6      74982868
7      55322550
8      10828207
9       4530577
10      3786748
11      2041762
12    342745985
13    292313639
14    259569928
15    143871771
16     53691635
17     24489644
18     20506718
19     14281945

Edited code incorporating changes from comments:

 > log.lik=function(theta,x,mu){
    + b=theta[1]
    + beta=theta[2]
    + tau=theta[3]
    + integrand<-function(x,b,beta,tau){exp(-x/b)/(1+exp(-beta*(x-tau)))}
    + result<-(-sum(x)/b)-sum(log(-beta*(x-tau)))-sum(log(integrate(integrand, lower=mu, upper=Inf)))
    + return(result)
    + }
    > model=optim(par=estimate, fn=log.lik, control=list(fnscale=-1), x=datafile, mu=1500)

Upvotes: 0

Views: 833

Answers (1)

jlhoward
jlhoward

Reputation: 59355

Too long for a comment.

First, you're still using integrate(...) incorrectly. This function returns a list (read the documentation!!). The $value element of this list is the integral. So to find the integral of f from a to b, use:

integrate(f,a,b,...)$value

Sadly, this is the least of your problems. Fundamentally, you're taking too many shortcuts. You can't just throw together some code - you need to pay attention to the math.

For instance, did you plot the values of your integrand(...) function for initial values of theta, in the range (mu,Inf)?? If you had, you would have seen that the integrand is 0 in this range, because mu=1500 and b=1, and exp(-1500/1) is numerically 0; hence the integral is 0, and the log of the integral is undefined. Also, your objective function includes the term log(-beta*(x-tau)), but for beta=tau=1, -beta*(x-tau) < 0 for all of the x in your dataset, and again the log is undefined.

For the record, I didn't downvote your question (because I find the practice offensive...), but you really do need to work on understanding if your log-likelihood function is correct, and when you've done that, take a close look at your initial estimates.

Upvotes: 1

Related Questions