1

I've just started learning r language and I am rewriting some Python programs in r. I'm having a problem where r returns me very strange values for the same function when I use r's optim function:

optim(c(10,0.2,0.05),inputs)

and r returns like this:

  Nelder-Mead direct search function minimizer
function value for initial parameters = 590.694928
  Scaled convergence tolerance is 8.80204e-06
Stepsize computed as 1.000000
BUILD              4 99999999999999996863366107917975552.000000 -951.452736
Error in optim(c(10, 0.2, 0.05), inputs, control = list(trace = 3, maxit = 20,  : 
  objective function in optim evaluates to length 0 not 1
In addition: There were 50 or more warnings (use warnings() to see the first 50)

but if I write same function and use scipy.optimize.minimize function in python, python will returns me :

 final_simplex: (array([[1.19418799e+02, 3.70856099e-01, 8.34769615e-03],
       [1.19418799e+02, 3.70856099e-01, 8.34769615e-03],
       [1.19418799e+02, 3.70856099e-01, 8.34769615e-03],
       [1.19418799e+02, 3.70856099e-01, 8.34769615e-03]]), array([-946.69021452, -946.69021452, -946.69021452, -946.69021452]))
           fun: -946.6902145168723
       message: 'Optimization terminated successfully.'
          nfev: 438
           nit: 208
        status: 0
       success: True
             x: array([1.19418799e+02, 3.70856099e-01, 8.34769615e-03])

I have tested that the function in r and python is exactly the same, but the same optimal method return totally different result, anyone know what this is all about?(btw, I also try Matlab, give me same result as python.) I also try different function, like $ x^2-2x+1 $, this kind of function r works well. Actually, I am working with log-likelihood function, so there exists log function, this might be the cause of the error?

Here is my function:

tail_body_likelihood <- function(r,params,pen){
    l = params[1]
    alpha = params[2]
    xmin = params[3]
    C = 2 - exp(-l * xmin)
    L = 0
    r = sort(r)
    f = sum(r < xmin)
    for (i in 1:f) {
        L = L - log(C) + log(l) - l*r[i]
    }
    for (i in (f+1):length(r)){
        L = L - log(C) + log(alpha/xmin) - (alpha+1)*log(r[i]/xmin)
    }
    Fbody = l*exp(-l * xmin) / C
    Ftail = alpha / (C * xmin)
    L = L - pen*(Fbody-Ftail)^2
    return(-L)
}
inputs <- function(x){
    a = tail_body_likelihood(b,x,10)
    return (a)
}

and here is the data.

-----------update-------------------------

I found if use fminsearch with condition:

fminsearch(inputs,c(9.5,2,0.05),lower=c(9,0.005,0.00001),upper = c(200,10,0.1),method='Hooke-Jeeves')

I got the similar result as Python, but it only works with Hooke-Jeeves method.

Yiyun Ding
  • 33
  • 4
  • 2
    It's easier to help you if you include a simple [reproducible example](https://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example) with sample input and desired output that can be used to test and verify possible solutions. Show the R function. – MrFlick Mar 06 '21 at 20:58
  • 1
    There are messages you should carefully analyze: "objective function in optim evaluates to length 0 not 1" and "There were 50 or more warnings (use warnings() to see the first 50)". – Erwin Kalvelagen Mar 06 '21 at 21:03
  • Hi Erwin, I think that is because there is log operation in the function, if I have a negative value in log(), seem r will return a nan, is this correct? – Yiyun Ding Mar 06 '21 at 21:07
  • Hi MrFlick, I have upload the r function and data. – Yiyun Ding Mar 06 '21 at 21:27

0 Answers0