About implementation of lsqnonlin

Hi

I have time-course data variable 'prodconcnumbers', size 18000X22 meaning 18000 data points for each of the 22 molecules in my model.

I have 30 parameters to be fit.

I invoke fitting using lsqnonlin as follows:

karbitraryvalues (My initial value parameters)=row vector of some 30 starting values.

lb and ub were again assigned bounds.

lsqnonlinoptions=optimset('Algorithm','Levenberg-Marquardt');

[kestimates,errorperspecies]= lsqnonlin(@lsqnonlin_fun, karbitraryvalues,lb,ub,lsqnonlinoptions);

So, first I built my objfunction as follows:

for eachspecies from 1 to 22

errorpermolecule(eachspecies)=prodconcnumbers(:,eachspecies)-estimatedconcnumbers(:,eachspecies);

end

So, in conclusion, my obj.function I want to minimize errorperspecies is of length 22 where as I have 30 parameters to b estimated using this lsqnonlin.

I am getting the following error:

    ??? Error using ==> lsqncommon at 102
    Levenberg-Marquardt and Gauss-Newton algorithms do not handle bound constraints and trust-region-reflective algorithm requires at least as
    many equations as variables; aborting.

What mistake am I doing? Appreciate your time and help.

[Matt J]

lsqcurvefit might be simpler here than lsqnonlin. Your data comparison is simple subtraction.

[Matt J]

Abdelwahab Afifi's comment moved here:

How you define your objective function interms of the variables' vector of length 30 ?

Answer 1

[Matt J]

Since you have bound constraints, there was no point choosing the Levenberg-Marquardt algorithm, since it will not handle those. When it saw that you were applying bounds, lsqnonlin defaulted back to the trust-region-reflective algorithm.

The other problem was that you are trying to fit 30 parameters with only 22 data points. The trust-region-reflective algorithm cannot process under-determined problems, but even if it could, you would get poor results because your problem is undetermined.

Find more data points with which to do the fit.

[John D'Errico]

+1. I'll add that in general, you need at least as many pieces of data as you have parameters to fit, but often you will be far happier if you have quite a bit more.

With fewer data points then parameters, the problem will typically have infinitely many solutions.

[SreeHarish Muppirisetty]

Matt and John D'Errico

Thanks a lot for inputs. I am slightly confused about the term data points. In my problem, I have concentrations of 22 different molecules over time measured from 0s to 18000s for every 1s. This data structure 'prodconcnumbers' (of size 18000X22) is my data. Do you mean to say that I should have 30 molecules instead of the current 22 to fit these 30 parameters? Thanks for your time again.

[Matt J]

The code you've shown suggests that your objective function is only generating one error term per molecule for a total of 22 error measurements. If you really have 18000*22 measurements, couldn't you be returning the error for each of those?

[SreeHarish Muppirisetty]

Thank you, that is what I already started doing with implementation of 'lsqcurve' which you suggested, i.e., calculating error at each of 18000 time points for all the 22 species. So, the order of my error (obj.function value) will be also 18000X22 like that of data. Will see how it works. Will close this question later.