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

Bas van Dijk v.dijk.bas at gmail.com
Thu Sep 10 09:21:40 EDT 2009

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:

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.

A high-level wrapper around bindings-levmar. It provides a more
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.

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:

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):