scipy is not optimizing and returns "Desired error not necessarily achieved due to precision loss"
I copied your example and tried a little bit. Looks like if you stick with BFGS solver, after a few iteration the mu+ alpha * r
will have some negative numbers, and that's how you get the RuntimeWarning.
The easiest fix I can think of is to switch to Nelder Mead solver.
res = minimize(loglikelihood, (0.01, 0.1,0.1), method = 'Nelder-Mead',args = (atimes,))
And it will give you this result:
28.3136498357
status: 0
nfev: 159
success: True
fun: 27.982451280648817
x: array([ 0.01410906, 0.68346023, 0.90837568])
message: 'Optimization terminated successfully.'
nit: 92
Another solution (that worked for me) is to scale your function (and gradients) to values closer to 0. For example, my problem came up when I had to evaluate a log-likelihood of 60k points. This meant that my log-likelihood was a very large number. Conceptually, the log-likelihood was a very very spikey function.
The gradients started off large (to climb this spikey mountain), and then became moderately small, but never less than the default gtol
parameter in the BGFS routine (which is the threshold that all gradients must be below for termination). Also, at this time I had essentially arrived at the correct values (I was using generated data so I knew the true values).
What was happening was that my gradients were approx. 60k * average individual gradient value
, and even if the average individual gradient value
was small, say less than 1e-8, 60k * 1e-8 > gtol
. So I was never satisfying the threshold even though I had arrived at the solution.
Conceptually, because of this very spikey mountain, the algorithm was making small steps, but stepping over the true minimum and never achieved average individual gradient << 1e-8
which implies my gradients never went under gtol
.
Two solutions:
1) Scale your log-likelihood and gradients by a factor, like 1/n
where n
is the number of samples.
2) Scale your gtol
: for example "gtol": 1e-7 * n
Facing the same warning, I solved it by rewriting the log-likelihood function to get log(params)
and log(data)
as arguments, instead of params and data.
Thus, I avoid using np.log()
in the likelihood function or Jacobian, if possible.