Re: [Help-gsl] Help with Geodesic Acceleration

[Patrick Alken]

Hello,

  In the 1D case, if you have n residuals f_i, where:

f_i = y_i - A*exp(-(x - B)^2 / (2 C^2))

then the fvv vector is the second directional derivative of f in the
direction of v. The first directional derivative is given by:

D_v f = (v . del) f

The second directional derivative is:

fvv = D_v^2 f = (v . del)^2 f = sum_{ij} v_i v_j d_i d_j f

So in order to compute fvv you need to evaluate the sum above. The first
step is to evaluate all second derivatives of your residuals with
respect to all combinations of model parameters, ie:

d^2 f_i / dA^2
d^2 f_i / dA / dB
d^2 f_i / dA / dC
...

The math gets a little messy so I've omitted it here. But there are 9
second derivatives you need (although some are symmetric so you can save
some work). Once you have all 9 second derivatives, then you just need
the sum:

(fvv)_k = sum_{ij} v_i v_j d_i d_j f_k

where the vector v is provided by the library to the fvv function. k
runs from 1 to n. i and j runs from 1 to 3. In the 1D case, v is a
3-vector (3 model parameters). For your 2D case, v will have 5 elements,
and you'll also need 25 second derivatives computed.

I have attached an example program I am using for fitting 1D gaussians.
Take a look at the fvv function in the attached file and see if you can
work out the corresponding function for the 2D case. It will be very
messy just due to the nature of taking derivatives of gaussians,
unfortunately.

Although in my experience it is better to have an analytic function for
the geodesic acceleration term, rather than use finite differences. I've
found that the finite difference approximation is often not accurate
enough, unfortunately.

You might also try experimenting with the dogleg methods if LM is not
performing well on your problem. Sometimes dogleg can converge quickly
while LM does not, and vice versa - so it may not be necessary to do all
the effort of computing these second derivatives!

Good luck!
Patrick

On 10/28/2016 04:50 PM, Hauser Fabian wrote:
> Dear whom it may concern,
>
>
> I want to use the Geodesic Acceleration (GSL V2.2.1) in combination with the 
> Levenberg-Marquardt? Algorithm to minimize a 2D Gaussian function. The 
> Problem is I have no idea who to develop the formula for f_vv. The source 
> code 'test_gaussian.c' and the documentation couldn't help me.
>
>
> Can someone explain me (mathematically) how develop any formal for Geodesic 
> Acceleration for this approach?
>
>
> Thanks in advance,
>
> Fabian Hauser

gaussfit.c
Description: Text Data