Re: GSoC 2015: Optimization Package: Non-linear and constrained least squares lsqcurvefit, lsqlin, lsqnonlin

[Olaf Till]

Hi AsmaA,

On Mon, Feb 23, 2015 at 02:10:29PM -0800, AsmaA wrote:
> Hi,
> 
> I am interested in the GSoC project for implementing Nonlinear and
> constrained least squares optimization problem in Octave [1]. 

I wasn't aware of the cited ([1]) paragraph on the wiki. I reproduce
it here for convenience:

-- start quote

  Nonlinear and constrained least squares

  The Optimization package is missing the functions lsqcurvefit,
  lsqlin, lsqnonlin to conveniently solve least-squares problems that
  are nonlinear and/or constrained. There are free implementations of
  the needed algorithms in other languages, such as minpack in Fortran
  and levmar in C. This project would link to or port these
  implementations and develop Matlab-compatible Octave wrappers.

  Mentor: Nir Krakauer

-- end quote

Although it is true that we havn't functions called 'lsqcurvefit' or
'lsqnonlin' in the optim package, the rest of the paragraph (when was
it written?) does not adequately describe the situation. The wrapper
functions 'nonlin_residmin' and 'nonlin_curvefit' already provide
access to backends for least squares nonlinear constrained
optimization. The wrappers provide a convenient interface similarly to
'lsqnonlin' and 'lsqcurvefit', but are differently named because there
are some differences (not by chance, but by decision). Additionally,
they provide some functionality which 'lsqnonlin' and 'lsqcurvefit'
don't provide.

> I am studying the Levenberg-Marquardt algorithm from [2].

This seems to be a general introduction to LM. An actual algorithm
must also, among others, take measures to be numerically stable. It is
usually best to start from one of the several already existing
algorithms (I know this is your intention). In the optim package we
have, among others, an SVD-based algorithm.

> I have taken a
> look at the levmar.c library [3].

On this web-page no non-linear constraints seem to be mentioned, which
would make the new algorithm inferior to the existing algorithm with
respect to this. But it may still be worth having the new one.

> It provides wrapper mex files to call it
> from MATLAB. The FAQ indicates the author heard that the wrappers are
> already Octave compatible. But he isn't sure.
> 
> I have quite a bit of experience with MATLAB but not mex/Octave.
> 
> My next step is to learn mex, call levmar.c from MATLAB and run a few
> examples between lsqcurvefit/lsqnonlin and levmar.c from MATLAB and
> comparing data. Then call mex from Octave and do similar.
> 
> At this point, I'd like a little feedback if I'm proceeding in the right
> direction. And any pointers would be appreciated :).

A new, alternative, optimization algorithm can be useful. There might
be problems better solvable by one algorithm than by the other. An
example problem for such a case in favour of the new algorithm would
be good, however.

Such a new, alternative, algorithm should be integrated as a backend
of 'nonlin_residmin' and 'nonlin_curvefit' (the latter is actually
only a convenience-wrapper to the former).

'nonlin_residmin' provides configurable access to functions computing
the Jacobian, so any new Jacobian-function should be treated
separately from the new algorithm.

There is a somewhat outdated description of the backend interface of
'nonlin_residmin' in the optimization package. More can be learned
from studying an existing backend. It is probably also advisable to
co-ordinate the integration of a new backend with me.

If anyone wants to exactly mimic the interface of 'lsqnonlin' and
'lsqcurvefit' (and personally I'd think this would not be the right
way ... sorry) this can either be done by writing them as seperate
wrappers for the existing backends, or by writing them as wrappers to
'nonlin_residmin' or 'nonlin_curvefit'. It would probably be wrong to
write a new algorithm only for them.

I'd think a project 'officially' supported by the Octave community
should not be based on mex-files (use m-code and/or oct-files).

I'm not sure, but it could be worth rewriting an algorithm of a
C-library in m-code (or as an oct-file) so that Octaves native
interfaces to e.g. decompostions and matrix operations are used, and
modifications are more easily possible. And: to my experience, the
objective- or model-function of the optimization problem is often so
expensive that writing least squares optimizing backends completely in
C would barely give a speed advantage. For expensive operations like
decompositions, m-code calls into Octave anyway, which interfaces with
compiled code to perform them.

Olaf

> References:
> 
> [1]
> http://wiki.octave.org/Summer_of_Code_Project_Ideas#Nonlinear_and_constrained_least_squares
> 
> [2] http://www.imm.dtu.dk/pubdb/views/edoc_download.php/3215/pdf/imm3215.pdf
> 
> [3] http://users.ics.forth.gr/~lourakis/levmar/index.html

-- 
public key id EAFE0591, e.g. on x-hkp://pool.sks-keyservers.net

signature.asc
Description: Digital signature