# [SciPy-User] scipy.optimize.leastsq question

Ralph Kube ralphkube@googlemail....
Mon Jun 28 07:36:38 CDT 2010

```Hello people,
I am having a problem using the leastsq routine. My goal is to
determine three parameters r_i, r_s and ppw so that the residuals
to a model function a(r_i, r_s, ppw) to a measurement are minimal.
When I call the leastsq routine with a good guess of starting values, it
iterates 6 times without changing the vales of the initial parameters
and then exits without an error.
The function a is very complicated and expensive to evaluate. Some
evaluation is done by using the subprocess module of python. Can this
pose a problem for the leastsq routine?

This is in the main routine:

import numpy as N

for t_idx, t in enumerate(time_var):

r_i = 300.
r_s = 1.0
ppw=1e-6
sza = 70.
wl = N.arange(300., 3001., 1.)

albedo_true = compute_albedo(r_i, r_s, ppw, sza, wl)
# This emulates the measurement data
albedo_meas = albedo_true + 0.01*N.random.randn(len(wl))

print 'Optimizing albedo'
p0 = [2.*r_i, 1.4*r_s, 4.*ppw]
plsq2 = leastsq(albedo_residual, p0, args=(albedo_meas, sza,
wl))
print '... done: ', plsq2[0][0], plsq2[0][1], plsq2[0][2]
albedo_model = compute_albedo(plsq2[0][0], plsq2[0][1], plsq2[0][2],
sza, wl)

The residual function:
def albedo_residual(p, y, sza, wvl):
r_i, r_s, ppw = p
albedo = compute_albedo(r_i, r_s, ppw, sza, wvl)
err = albedo - y
print 'Albedo for  r_i = %4.0f, r_s = %4.2f, ppw = %3.2e \
Residual squared: %5f' % (r_i, r_s, ppw, N.sum(err**2))

return err

The definition of the function a(r_i, r_s, ppw)
def compute_albedo(radius_ice, radius_soot, ppw, sza, wvl):

The output is:
Optimizing albedo
Albedo for r_i =  600, r_s = 1.40, ppw = 4.00e-06 		Residual squared:
0.973819
Albedo for r_i =  600, r_s = 1.40, ppw = 4.00e-06 		Residual squared:
0.973819
Albedo for r_i =  600, r_s = 1.40, ppw = 4.00e-06 		Residual squared:
0.973819
Albedo for r_i =  600, r_s = 1.40, ppw = 4.00e-06 		Residual squared:
0.973819
Albedo for r_i =  600, r_s = 1.40, ppw = 4.00e-06 		Residual squared:
0.973819
Albedo for r_i =  600, r_s = 1.40, ppw = 4.00e-06 		Residual squared:
0.973819
... done:  600.0 1.4 4e-06

To check for errors, I implemented the example code from
http://www.tau.ac.il/~kineret/amit/scipy_tutorial/ in my code and it
runs successfully.

I would be glad for any suggestion.

Cheers, Ralph
```

More information about the SciPy-User mailing list