Curve fitting tool with custom equation, odd power polynomial

I am trying to fit an odd power polynomial of the form (f(x) = p1*x^3 + p3*x) to a data set of with 38 points, which clearly shows a smooth cubic polynomial. I am using the custom equation type in the curve fitting tool to do this. I understand the custom fit is sensitive to the starting parameters, so I first run the standard 3rd order polynomial fit to my data. The fit includes even power terms, but I extract the odd order coefficients and use those as a starting point in my custom fit. It does not work well at all.

As a test, I plugged in the full 3rd order polynomial equation into the custom equation (i.e. f(x) = p1*x^3 + p2*x^2 + p3*x + p4), and set my starting values to the same order of magnitude solved for in the standard polynomial fit, and the custom equation fit still did not fit well at all. How do I fix this?

Alternatively, how can I best fit an odd power polynomial?

Answer 1

[Matt J]

It seems like overkill to use the curve fit toolbox just to fit a 1D polynomial. And it should not depend on any starting parameter guess. For a polynomial fit, the coefficients should be the unique solution to a set of linear equations. Below is what I use. For you, the usage would be,

>> polyfitc(x,y,[3 1])

function p=polyfitc(x,y,nvec)
%Simple 1D polynomial fitting with particular coefficients constrained to
%zero
%
%    p=polyfitc(x,y,nvec)
%
%in: 
%
%  nvec: A vector of integer exponents present in the polynomial.
%    x: x data
%    y: y data
%
%out:
%
%   p: vector of polynomial coefficients. Order of polynomial is max(nvec).

           A=bsxfun(@power,x(:),nvec(:).');
           [QQ,RR]=qr(A,0);
           coeffs = RR\( QQ'*y(:) );

            p=zeros(1,max(nvec)+1);
           p(nvec+1)=coeffs;
           p=p(end:-1:1);

[John]

Indeed, it's overkill, which is why I don't understand how badly the curve fitting toolbox fails.

Your code seems to work quite well, however. Thanks for that!

[dpb]

I've not used it enough to know but it seems you're not using a LSQ routine that way to even have to have any estimates. Surely there's a way to fit an OLS problem w/ custom terms? If not, for this problem it simply isn't the right toolset.

But, that doesn't help fix up your data if the model doesn't fit it well.

Answer 2

[dpb]

That to me indicates the odd-only terms aren't a very good model for your data, first off. Did you look at the interval estimates around the even coefficients? Are the coefficients estimates small in magnitude?

Anyway, to fit a model, just use the left-divide operation and set up the model explicitly.

X=[x.^3 zeros(size(x)) x ones(size(x))];

where x represents the vector x of your observations. Then

p=X\y;

The above includes an intercept; for no intercept use zeros() for the last column, too, instead of ones or simply omit it.

If your x data have sizable range you'll get better numerical stability if you center and standardize first; if so remember the coefficients are in the transformed space.

[John]

To give some context, let me show a couple screen shots. First is the polyfit

As we can see, it fits well. Now let's try the custom equation with just odd power terms and use order of magnitude of the coefficients solved here as a starting point.

As you can see, pretty awful. Clearly, I am using the tool the wrong way, or something about this data makes it hard to fit using the custom equation.

In any case, the explicit method you described, which was also given by the previous commenter Matt J, works perfectly in my case.

Thanks for the feedback!

[Matt J]

If you want errors lower than 1e-6, then you probably need to set TolFun substantially lower than 1e-6. This might drive the algorithm to use more iterations, so MaxIter and MaxFunEvals might need to be raised as well.

Finally, not sure what the Robust='Off' flag is doing, but maybe it's a good idea to turn it on.

[dpb]

I can't read it well enough to tell what the actual fitted coefficients/errors are but the problem I see is that you're using a nonlinear least squares solution instead of OLS. If you fix that it'll estimate the same coefficients.

Just out of curiosity, what if you use the OLS estimates as the NLLSQ estimates?

ADDENDUM:

I see there's no way to specify a custom equation for an OLS problem--seems a significant oversight/limitation in the tool...you might also try switching to Marquardt-Levenberg just to see.

I used a set of data that were created from an odd-poly and the tool worked fine without me ever specifying an initial guess. That might be another thing to try??? Again, this would seem all to be just for the exercise...

[Marc]

Ok the fit is "bad" but is it wrong? Your fit to the full polynomial suggests a significantly sized p2 term. So maybe this is as good as your going to get?

Do you have the stats toolbox? Try regress().