[SciPy-user] Another leastsq Jacobian bug
Lin Shao
shao@msg.ucsf....
Sun Jun 24 23:43:02 CDT 2007
Try this:
import numpy as N
import scipy.optimize as O
## Create my data to fit:
x=N.arange(-10,10,dtype=N.float32)
y=(x-1.234)**2+3.456
## I want to fit y=p0*(x-p1)^2+p2
## Now I define my objective
def obj_func(params, xx, yy, mode='col'):
return params[0] * (xx-params[1])**2 + params[2] - yy
## 'mode' is dummy here, because I want to use it for my Jocobian
## Now define my Jacobian
def Jacobian(params,xx,yy,mode='col'):
J = N.empty((len(params),xx.size))
J[0] = (xx-params[1])**2
J[1] = -2*params[0]*(xx-params[1])
J[2] = 1
if mode=='col':
return J
elif mode=='row:
return J.transpose()
else:
raise ValueError, "Unkown mode %s in Jacobian()" % mode
## Hopefully I did my calculus correctly
## First, try without Jocobian given
b=O.leastsq(obj_func, (1.,1.,2.), args=(x, y, 'col'))
print b[0]
## the result is [ 1. 1. 2.]
## obviously a failure, but no warning is returned
## Second, try Jacobian with col_deriv=1
b=O.leastsq(obj_func, (1.,1.,2.), args=(x, y, 'col'), Dfun=Jacobian,
full_output=1, col_deriv=1)
## full_output has to be 1 because of a seg fault bug I reported earlier
print b[0]
## the result is [ 1.00000004 1.23400008 3.45599962]
## Good Job!
## Last, try Jacobian with col_deriv=0
b=O.leastsq(obj_func, (1.,1.,2.), args=(x, y, 'row'), Dfun=Jacobian,
full_output=1, col_deriv=0)
print b[0]
## the result is [ 1.02507777 1.222815 2.2651552 ]
## completely different result and wrong
Sorry it wasn't true what I said in my first report. The leastsq()
result when col_deriv is not 1 is not an exactly the same as the
initial guess. But one does have to agree there's a bug somewhere.
--lin
On 6/23/07, Stefan van der Walt <stefan@sun.ac.za> wrote:
> Hi Lin
>
> On Fri, Jun 22, 2007 at 07:35:05PM -0700, Lin Shao wrote:
> > It seems like when calling leastsq(), if the Jacobian matrix is
> > organized as rows of partial derivatives (w.r.t. each variable), then
> > no optimization is done at all -- the return value is the same as the
> > initial guess. It only works when the matrix is columns of derivatives
> > and set the parameter col_deriv to 1.
>
> It would be helpful if you could provide two short snippets of code to
> illustrate the problems you mention.
>
> Regards
> Stéfan
> _______________________________________________
> SciPy-user mailing list
> SciPy-user@scipy.org
> http://projects.scipy.org/mailman/listinfo/scipy-user
>
More information about the SciPy-user
mailing list