Thread Subject: weighted linear LS fit to a 3rd or 4th order polynomial

I want to perform a weighted linear LS fit using an N-th order polynomial
fitting model.

In particular I want to fit to 3rd and 4th order polynomials, and I want to the
LS procedure to accept a user-defined weighting function
    w(k) = W( d(k)^2 ),
where d(k)^2 is the squared deviation from the fitted line for datum k.

Although weighting in this sense is mentioned in the LLS section of the Help
for the Curve Fitting toolbox, I couldn't actually find a Matlab function that
will do this. I may have missed it because my grasp of linear algebra is weak.

If there is a suitable function, I'd be grateful for a rudimentary example of
how to use it with N = 4 and user-specified W.

Thanks

ross

ross wrote:
>
>
> I want to perform a weighted linear LS fit using an N-th order
> polynomial
> fitting model.
>
> In particular I want to fit to 3rd and 4th order polynomials, and I
> want to the
> LS procedure to accept a user-defined weighting function
> w(k) = W( d(k)^2 ),
> where d(k)^2 is the squared deviation from the fitted line for
> datum k.
>
> Although weighting in this sense is mentioned in the LLS section of
> the Help
> for the Curve Fitting toolbox, I couldn't actually find a Matlab
> function that
> will do this. I may have missed it because my grasp of linear
> algebra is weak.
>
> If there is a suitable function, I'd be grateful for a rudimentary
> example of
> how to use it with N = 4 and user-specified W.

I assume this is for a robust fitting
scheme. If so, then look at robustfit
from the stats toolbox.

If you don't have that TB, then you
will want to use a solver like
fminsearch, fminunc, fmincon, etc.,
to do the work. For example,

Consider the simple problem

n = 100;
x = 2*rand(n,1) - 1;
y0 = exp(x);

% make some noise with a
% non-gaussian distribution
noise = randn(n,1)/2;
noise = sign(noise).*abs(noise).^4;

y = y0 + noise;

plot(x,y0,'ro',x,y,'b+')

% Fit this curve with a
% fourth order polynomial model
%
% y = a1+a2*x+a3*x^2+a4*x^3+a5*x^4
%
% using a robust fitting scheme.
% First, use polyfit with no weights.
% Note how poorly polyfit does on
% the noisy data.
polyfit(x,y0,4)
ans =
  0.043738 0.17625 0.49911 0.99789 1

polyfit(x,y,4)
ans =
  -1.0048 0.47634 1.8633 0.63957 0.69082

% Can we do better?
Wfun = @(R) erf(R);
obj = @(c) sum(Wfun(y-polyval(c,x)).^2);

coef = fminsearch(obj,[1 1 1 1 1])
coef =
 -0.054147 0.11587 0.65318 1.0641 0.98598

These estiamtes are considerably
better.

HTH,
John D'Errico

"John D'Errico" <woodchips@rochester.rr.com> wrote:

>I assume this is for a robust fitting
>scheme. If so, then look at robustfit
>from the stats toolbox.

but that only fits linear models, n'est-ce pas?
the bisquare function used in robustfit is for weighting.

>If you don't have that TB, then you
>will want to use a solver like
>fminsearch, fminunc, fmincon, etc.,
>to do the work. For example,

From the Help descriptions, it seems to follow that:
-- fminunc seems better than fminsearch because the model is polynomial and the
derivs are easy to calculate..
-- fminunc seems better than fmincon because there are no constraints.
I will have a play around to see if this works out in practice.

These nonlinear optimising functions are very liberating for people like me who
are mathematically challenged. It's no longer necessary to stick to functions
with simple analytical form for the fit model or the weighting.

Thanks for your help, John.

ross

Hey All,

I have a similar issue. I was tempted to write my own little
algorithm, but that will take to much time. I remember there
was a function called 'fit' and you could use a linear fit
and give it weights. In addition to being able to do a
weighted linear fit you could also access the confidence
interval, and get an estimate for the uncertainty in the fit
parameters. My question is the following: "Where is that
function? and What can I use in place of it?"

Much obliged,

Nate

ross <ross@nospam.com> wrote in message
<45701b99@news.eftel.com>...
> "John D'Errico" <woodchips@rochester.rr.com> wrote:
>
> >I assume this...

"Nathan Orloff" <orloffREMOVETHIS@nist.gov> wrote in message
news:ffb4v3$8fq$1@fred.mathworks.com...
> Hey All,
>
> I have a similar issue. I was tempted to write my own little
> algorithm, but that will take to much time. I remember there
> was a function called 'fit' and you could use a linear fit
> and give it weights. In addition to being able to do a
> weighted linear fit you could also access the confidence
> interval, and get an estimate for the uncertainty in the fit
> parameters. My question is the following: "Where is that
> function? and What can I use in place of it?"
>
> Much obliged,
>
> Nate
>
> ross <ross@nospam.com> wrote in message
> <45701b99@news.eftel.com>...
>> "John D'Errico" <woodchips@rochester.rr.com> wrote:
>>
>> >I assume this...
>

There are number of related functions in the statistics toolbox.

http://www.mathworks.com/access/helpdesk/help/toolbox/stats/
click on "Category" then click on "Regression Analysis" you will see a
number of functions releated to least squares fit. For example "robustfit"
and others.

Nasser

I actually ended up using lscov(A,b,W). It did the job.
Thanks everyone.