0

I have been learning optimization methods for a few days now. The following code that I wrote returned a RuntimeWarning.

import numpy as np
from scipy.optimize import minimize

def func(a, x):
    return 1 + (x - 0.5) * a

def log_like(a, x):
    sum1 = 0 
    for i in range(len(x)):
        sum1 += np.log(func(a, x[i]))
    return sum1
    
def log_like_prime(a, x):
    sum1 = 0
    for i in range(len(x)):
        sum1 += (x[i] - 0.5) / (1 + (x[i] - 0.5) * a)
    return sum1

def log_like_prime2(a, x):
    sum1 = 0
    for i in range(len(x)):
        sum1 += -(x[i] - 0.5) ** 2.0 / (1 + (x[i] - 0.5) * a) ** 2.0
    return sum1

x = [0.89, 0.03, 0.50, 0.36, 0.49]
a = -1
a_opt = minimize(
    log_like, a, args=(x,), method="Newton-CG",
    jac=log_like_prime, hess=log_like_prime2
)
print(a_opt)

Returns the following error:

fun: array([0.03194467])
     jac: array([0.18690836])
 message: 'Warning: Desired error not necessarily achieved due to precision loss.'
    nfev: 21
    nhev: 1
     nit: 0
    njev: 21
  status: 2
 success: False
       x: array([-1.])
py:17: RuntimeWarning: invalid value encountered in log
  sum1 += np.log(func(a, x[i]))
py:17: RuntimeWarning: invalid value encountered in log
  sum1 += np.log(func(a, x[i]))

It should not be returning an invalid value as for the given values of x = [0.89, 0.03, 0.50, 0.36, 0.49], no negative value must be returned by the function inside the logarithm part. I cannot understand why such a problem is there.

mkrieger1
  • 19,194
  • 5
  • 54
  • 65
gautam bhuyan
  • 21
  • 1
  • 4
  • Are you sure you want to minimize the function instead of maximizing it? There's no local minimum for your choice of `x`. – joni Oct 01 '21 at 13:53

1 Answers1

0

Your initial guess for a, -1, gives strictly positive values func(a, x[i]) for all i.

However scipy's minimization algorithm changes the value of a (that's how it finds a better value). Here, it goes towards larger negative values, eventually resulting in negative values for func(a, x[i]) and thus the log error.

To debug (I agree that scipy is not friendly here) try writing np.seterr(all='raise') at the beginning of your script, then the warning raised in the log will transform into an error, and typing %debug% in the next ipython cell, you'll be able to check the value of a`:

~/miniconda3/lib/python3.8/site-packages/scipy/optimize/_differentiable_functions.py in fun_wrapped(x)
     68         def fun_wrapped(x):
     69             self.nfev += 1
---> 70             return fun(x, *args)
     71 
     72         def update_fun():

~/tmp.py in log_like(a)
     12     sum1 = 0
     13     for i in range(len(x)):
---> 14         sum1 += np.log(func(a, x[i]))
     15     return sum1
     16 

FloatingPointError: invalid value encountered in log

In [21]: %debug
> /home/mathurin/tmp.py(14)log_like()
     12     sum1 = 0
     13     for i in range(len(x)):
---> 14         sum1 += np.log(func(a, x[i]))
     15     return sum1
     16 

ipdb> a
a = array([-2.7761329])
P. Camilleri
  • 12,664
  • 7
  • 41
  • 76