[Haskell-cafe] ANNOUNCE: A Levenberg-Marquardt implementation

[Bas van Dijk]

Dear all,

We like to announce the release of a Haskell binding to Manolis
Lourakis's C levmar library at:

http://www.ics.forth.gr/~lourakis/levmar/

This library implements the Levenberg-Marquardt algorithm which is an
iterative technique that finds a local minimum of a function that is
expressed as the sum of squares of nonlinear functions. It has become
a standard technique for nonlinear least-squares problems and can be
thought of as a combination of steepest descent and the Gauss-Newton
method. When the current solution is far from the correct one, the
algorithm behaves like a steepest descent method: slow, but guaranteed
to converge. When the current solution is close to the correct
solution, it becomes a Gauss-Newton method.

Our binding consists of three packages:

* http://hackage.haskell.org/package/bindings-levmar-0.1

A low-level wrapper around the C library. Note that the C library is
lightly patched so that the functions can be safely called inside
unsafePerformIO. The patched C library is bundled with this package.

* http://hackage.haskell.org/package/levmar-0.1

A high-level wrapper around bindings-levmar. It provides a more
familiar Haskell interface. For example, instead of passing a 'Ptr r'
to the levmar function you can pass a [r]. levmar also provides some
higher-level modules that use some type-level programming to add more
type safety.

* http://hackage.haskell.org/package/levmar-chart-0.1

Finally levmar-chart is a small package for quickly viewing the output
of levmar in a chart.

Unfortunately the documentation of these libraries is not available
from hackage because bindings-levmar won't configure because of a
missing dependency (lapack) on the hackage server. Instead I put the
documentation at the following places:

http://code.haskell.org/bindings-levmar/bindings-levmar-0.1/doc/html/bindings-levmar/index.html
http://code.haskell.org/levmar/levmar-0.1/doc/html/levmar/index.html

Here follows a quick example:

Suppose we have the following model functions:

constant  :: Num r => Model N1 r r
linear    :: Num r => Model N2 r r
quadratic :: Num r => Model N3 r r
cubic     :: Num r => Model N4 r r

constant        a _ = a
linear        b a x = b * x     + constant      a x
quadratic   c b a x = c * x*x   + linear      b a x
cubic     d c b a x = d * x*x*x + quadratic c b a x

And the jacobians:

constantJacob  :: Num r => Jacobian N1 r r
linearJacob    :: Num r => Jacobian N2 r r
quadraticJacob :: Num r => Jacobian N3 r r
cubicJacob     :: Num r => Jacobian N4 r r

constantJacob        _ _ = 1     ::: Nil
linearJacob        _ a x = x     ::: constantJacob      a x
quadraticJacob   _ b a x = x*x   ::: linearJacob      b a x
cubicJacob     _ c b a x = x*x*x ::: quadraticJacob c b a x

Now assume we have some sample data. If you call levmar like this:

levmar cubic
       (Just cubicJacob)
       (-0.05 ::: 0.5 ::: -12 ::: 10 ::: Nil)
       samples
       1000
       defaultOpts
       Nothing
       Nothing
       noLinearConstraints
       Nothing

You get the following fit (using levmar-chart):

http://code.haskell.org/~roelvandijk/code/levmar-chart/cubicFit.png

Note that levmar contains a demo with a lot more examples:

http://code.haskell.org/levmar/Demo.hs

Happy fitting,

Roel and Bas van Dijk