# How to avoid NaN in optim?

Suppose, I have a function and want to optimize it. But if I use optim() which gives warnings(). How can I avoid these warnings of NaN?

    myfun<-function(par, x){

f<- sum(x)*length(x)+sum(log(gamma(par))*x)+1
return(-f)
}
optim(0.1, myfun, x=c(1,5,4,7,8,5,6,5,45,8))

$par [1] 4.203895e-46$value
[1] -10762.39
$counts function gradient 502 NA$convergence
[1] 1
$message NULL There were 50 or more warnings (use warnings() to see the first 50) warnings() 2: In log(gamma(par)) : NaNs produced 3: In log(gamma(par)) : NaNs produced 4: In log(gamma(par)) : NaNs produced  - ## migrated from stats.stackexchange.comJan 16 '13 at 9:36 This question came from our site for people interested in statistics, machine learning, data analysis, data mining, and data visualization. would this not be a better question to ask stackoverflow? – ACD Jan 16 '13 at 4:43 ## 3 Answers As AdamO said, you have an issue with non-positive values. Purely for optimization purposes, exponentiating to make sure things are non-negative and adding some salt to avoid zero generally does the trick. So... something like: myfun<-function(par, x){ par <- exp(par) + 10^-10 f<- sum(x)*length(x)+sum(log(gamma(par))*x)+1 return(-f) }  Then none of the two optimizers will give you any trouble. optim( log(0.1), myfun, x=c(1,5,4,7,8,5,6,5,45,8), method="BFGS") optim( log(0.1), myfun, x=c(1,5,4,7,8,5,6,5,45,8), method="CG") # I logged because I exponentiate in the function.  Basically you have a constrained optimization problem and you want to express it as an unconstrained one. Therefore you exponentiate your solution space to make sure it is non-negative one. Just a word of caution: Exponentiation can occasionally lead to problems if you test for very large values, eg. par=113; would lead to evaluating "gamma(exp(113))" which is actually quite a large number. :) - The problem is that your optimal solution is a boundary value (at 0, and negative values of par aren't allowed). You can verify this by plotting myfun with the given values of$x\$. I don't think general optimizers are well configured to handle such issues.

-

Use lgamma(par) instead of log(gamma(par)).

-
Generally speaking on stack exchange, longer answers are preferred. To improve your answer you could add things such as more detail as to what the underlying problem with log(gamma(par)) is, why lgamma(par) should fix (or work around) the problem, and how to avoid this sort of problem in the future. –  Mark Booth Jan 21 '13 at 16:20