[SciPy-user] Using SciPy/NumPy optimization
Alok Singhal
as8ca at virginia.edu
Tue Feb 27 16:46:29 EST 2007
Hi,
Sorry about the long post. I did not ask the original question, but I
tried the solution posted and could not get it to work.
Original question: given a set of data [[x1, y1], ..., [xn, yn]], how
to fit y = a*(1 + b*c*x) ** (-1.0/b) to it, subject to the
constraints:
0 < b <= 5
-1 <= c <= 1
a ~ >= max(x1, ... xn)
On 27/02/07: 14:25, Robert Kern wrote:
> There are a number of constrained optimizers in scipy.optimize .
> scipy.optimize.fmin_tnc seems most appropriate for simple bounds like you have.
> In order to get a function to minimize that depends on data, I usually like to
> use a class:
>
> import numpy as np
> from scipy.optimize import fmin_tnc
>
> class LossFunction(object):
> def __init__(self, x, y):
> self.x = x
> self.y = y
>
> def __call__(self, abc):
> """ A function suitable for passing to the fmin() minimizers.
> """
> a, b, c = abc
> y = a*(1.0 + b*c*self.x) ** (-1.0/b)
> dy = self.y - y
> return dy*dy
>
> x = np.array([...])
> y = np.array([...])
> lf = LossFunction(x, y)
> abc0 = np.array([x.max(), 2.5, 0.0]) # or whatever
> retcode, nfeval, abc_optimal = fmin_tnc(lf, abc0,
> bounds=[(None, None), (0., 5.), (-1., 1.)])
I used optimize.leastsq to do the fitting, and it worked well me (see
http://www.scipy.org/Cookbook/FittingData). But I was trying the
above method, and unfortunately it doesn't work as presented. When I
define the variables x and y by:
from scipy.optimize import fmin_tnc
from scipy import rand
from numpy import mgrid, log
def f(abc, x):
return abc[0]*(1.0+abc[1]*abc[2]*x)**(-1.0/abc[1])
a = 15.0
b = 2.5
c = 0.3
abc = (a, b, c)
num_points = 151
x = mgrid[1:10:num_points*1j]
y = f(abc, x) + rand(num_points) - 0.5
abc0 = [x.max(), 2.5, 0.0]
lf = LossFunction(x, y)
retcode, nfeval, abc_optimal = fmin_tnc(lf, abc0, bounds=[(None, None), (0., 5.), (-1., 1.)])
The call to fmin_tnc gives me an error:
File "/usr/lib/python2.4/site-packages/scipy/optimize/tnc.py", line 191, in func_and_grad
f, g = func(x, *args)
ValueError: too many values to unpack
Looking at the documentation and the source file, the first parameter
to fmin_tnc should be a function that returns the function value and
the gradient vector (with respect to the parameters being estimated).
Also, it seems that I can tell fmin_tnc to estimate the gradient
itself, if I specify approx_grad as True:
retcode, nfeval, abc_optimal = fmin_tnc(lf, abc0, bounds=[(None, None), (0., 5.), (-1., 1.)], approx_grad=True)
But this gives me an error saying:
File "/usr/lib/python2.4/site-packages/scipy/optimize/optimize.py", line 576, in approx_fprime
grad[k] = (apply(f,(xk+ei,)+args) - f0)/epsilon
ValueError: setting an array element with a sequence.
Then I tried redefining the __call__ in LossFunction to return the
gradient as well:
class LossFunction(object):
def __init__(self, x, y):
self.x = x
self.y = y
def __call__(self, abc):
a, b, c = abc
y = a*(1.0 + b*c*self.x) ** (-1.0/b)
dy = self.y - y
d = [0, 0, 0]
d[0] = (1.0 + b*c*self.x) ** (-1.0/b)
d[1] = y * (log(a*(1.0+b*c*self.x))/b**2 - c*self.x/(b*(1.0+b*c*self.x)))
d[2] = -a/b*(1.0 + b*c*self.x) ** (-1.0/b) * b*self.x
return dy*dy, d
When I call this function instead, I get another error:
File "/usr/lib/python2.4/site-packages/scipy/optimize/tnc.py", line 221, in fmin_tnc
fmin, ftol, rescale)
ValueError: Bad return value from minimized function.
At this point, I gave up, but I would like to understand what I am
doing wrong. To me, it seems like fmin_tnc will find the values of
parameters (a, b, c in this case) that minimize a given function that
depends *only* on a, b, c. I don't understand how to make it work so
as to fit existing data. I am sure I am missing something, but I
don't know what that is.
Thanks for all your help,
Alok
--
Alok Singhal * *
Graduate Student, dept. of Astronomy * * *
University of Virginia
http://www.astro.virginia.edu/~as8ca/ * *
More information about the SciPy-User
mailing list