MATLAB: Optimizing a (rather large) sum-of-Gaussians – fmincon fails

fmincongaussianMATLABsum

Dear All!
I am attempting to fit a model to data which is made up by a collection of Gaussian-like functions.
each sum term takes a set of 5 parameters and in total there are about 20 sets, so inside my function I just use a loop k=1:20 (the independent variable is T and just passed as a global variable)
function [ SSE, func ] = simpa( params )
for k=1:size(params,2)
f(:,k)= params(k,2).*q.*params(k,4).*Area.*E.*params(k,5).*3e21*0.15*params(k,3).*T.^2.*exp(-params(k,1)./(kB.*T)-kB.*3e21*0.15*params(k,3)./(beta.*params(k,1)).*T.^4.*exp(-params(k,1)./(kB.*T)).*(1-4.*kB.*T./(params(k,1)+20.*(kB.*T./params(k,1)).^2)))
end
and I sum over all k to get the spectrum and then calculate the error
SSE=sum((sum(f,2)-ydata).^2)
and feed this to fmincon.
However, this doesn't work at all. Even if I fine-tune the starting parameters manually pretty close to the data (which actually is not that hard), fmincon is just stuck somewhere in between. It does manage to reduce the sse, but not by much.
Is there any way to (a) optimize my function since it has lots of recurring terms and (b) help fmincon?
Any help is very much appreciated!
Yours, Chris

Best Answer

Because you are minimizing a sum of squares, I suggest that you use lsqnonlin or lsqcurvefit to try to find a minimum. If you take this suggestion, be sure to NOT sum the squares in your objective function.
Usually, using global variables is not a good idea. There are almost always better ways, such as passing T to each summand separately, not as a global.
I suspect that you can scale your problem better. In particular, I would try to get rid of the term 3e21 by appropriate centering and scaling.
It is possible that, along with rescaling, you need to take different finite differences, such as setting the FinDiffRelStep option to 1e-4 or setting the FinDiffType option to 'central'.
Good luck,
Alan Weiss
MATLAB mathematical toolbox documentation