# [SciPy-user] Using SciPy/NumPy optimization

Robert Kern robert.kern@gmail....
Tue Feb 27 14:25:59 CST 2007

```Brandon Nuttall wrote:
> Folks,
>
> The usual confessions: I am relatively new to Python programming and SciPy.
> I have a problem I'm looking for some help in solving.
>
> I have a list of data pairs [[x1,y1], [x2,y2], ..., [xn,yn]] and am trying
> to find the best fit of that data to an equation:
>
>       y = a*(1+b*c*y)^(-1/b)

Presumably, you mean

y = a*(1 + b*c*x) ** (-1.0/b)

to correct a typo and use Python notation.

> The parameters, b and c, are constrained:
>       1) 0<b<=5
>       2) -1<=c<=1
>
> Parameter a is only weakly constrained in that x is usually >= max(x1, x2,
> ... , xn).
>
> My "best fit" goal is either to minimize the root mean square deviation (or
> consequently maximize the r-square value).

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.)])

--
Robert Kern

"I have come to believe that the whole world is an enigma, a harmless enigma