p_orig =
Problem with polyfit command (R2015a)
이전 댓글 표시
Below is a small script for using the polyfit command but surprisingly the last command gives me a completely wrong polynomial p. I don't understand why. Thanks in advance for your help.
Below the script is the response from my version of MATLAB (R2015a).
%--------------------------------------------------------------------------------------------------------------
%butta_sto_test
%
clear all
clc
x=[1:9]
y=[5,6,10,20,28,33,34,36,42]
p=polyfit(x,y,1)
[p,S]=polyfit(x,y,1)
[p,S,mu]=polyfit(x,y,1)
%----------------------------------------------------------------------------------------------------------------------
x =
1 2 3 4 5 6 7 8 9
y =
5 6 10 20 28 33 34 36 42
p =
4.9833 -1.1389
p =
4.9833 -1.1389
S =
R: [2x2 double]
df: 7
normr: 8.4581
p =
13.6474 23.7778
S =
R: [2x2 double]
df: 7
normr: 8.4581
mu =
5.0000
2.7386
댓글 수: 1
Dyuman Joshi
2024년 5월 30일
You are, most likely, not using the appropriate syntax of polyval() to evaluate the polyfit output obtained, as Matt has shown below.
채택된 답변
The result from polyfit is correct. I suspect you are simply not using polyval properly to evaluate the fit:
x=[1:9] ;
y=[5,6,10,20,28,33,34,36,42];
[p,S,mu]=polyfit(x,y,1);
xx=linspace(1,9);
plot(x,y,'o',xx,polyval(p,xx,S,mu)); legend Data Line-Fit
댓글 수: 8
Dyuman Joshi
2024년 5월 30일
OP was probably not using the syntax you used to evaluate the polynomial obtained with that particular syntax.
T Hafid
2024년 5월 30일
When you call polyfit with the syntax,
[p,S,mu]=polyfit(x,y,1);
the fit is done after pre-normalizing x, which you can read more about in the polyfit documentation. The coefficients that polyfit gives you will be different to reflect this pre-normalization.
After applying polyfit, you commonly want to evaluate the polynomial at new locations xnew. For this, it is common to use polyval. If you use polyval with the input syntax I showed you, it will evaluate the polynomial in a way that accounts for the normalization. That way, you do not have to apply the normalization to xnew manually.
You understand the normalization of the polynomial obtained from your second call to polyfit ?
x=[1:9] ;
y=[5,6,10,20,28,33,34,36,42];
p=polyfit(x,y,1);
[p,S]=polyfit(x,y,1)
[p,S,mu]=polyfit(x,y,1);
syms z
p_orig = p(1)*((z-mean(x))/std(x))+p(2)
fliplr(double(coeffs(p_orig,z)))
T Hafid
2024년 5월 30일
The "two" polynomials p1 and p3 returned from polyfit define a unique polynomial p - only their polynomial coefficients a_i and b_i differ because of the different basis functions used in their development.
p1 is written in the standard basis {1,X,X^2,...,X^n} as
p1(X) = sum_{i=0}^{i=n} a_i * X^i
while p3 is written in the basis {1,((X-mean(xdata))/std(xdata)),((X-mean(xdata))/std(xdata))^2,...,((X-mean(xdata))/std(xdata))^n} as
p3(X) = sum_{i=0}^{i=n} b_i * ((X-mean(xdata))/std(xdata))^i
It would be less confusing if the third form of the polyfit command also provided the same polynomial p as the first two forms.
The math operations required to undo the normalization would result in floating point errors in the unnormalized coefficients. The floating point errors could then negatively impact numerical accuracy when the polynomial is evaluated later.
추가 답변 (1개)
From the polyfit documentation page: "[p,S,mu] = polyfit(x,y,n) performs centering and scaling to improve the numerical properties of both the polynomial and the fitting algorithm. This syntax additionally returns mu, which is a two-element vector with centering and scaling values. mu(1) is mean(x), and mu(2) is std(x). Using these values, polyfit centers x at zero and scales it to have unit standard deviation,"
If you call polyfit with three outputs, p is not a polynomial in x. It is a polynomial in the centered and scaled 
.
.xdata=[1:9];
y=[5,6,10,20,28,33,34,36,42];
p = polyfit(xdata, y, 1);
Let's look at p symbolically.
psym = poly2sym(p);
polynomialInX = vpa(psym, 5)
Now let's look at the polynomial in the centered and scaled 
.
.[p, ~, mu] = polyfit(xdata, y, 1);
syms xhat
polynomialInXhat = vpa(poly2sym(p, xhat), 5)
These look different. But what happens if we substitute the expression for into polynomialInXhat?
into polynomialInXhat?syms x
vpa(subs(polynomialInXhat, xhat, (x-mu(1))/mu(2)), 5)
That looks the same as polynomialInX. What if we evaluate both polynomials, polynomialInX at the unscaled X data and polynomialInXhat at the scaled X data?
valueUnscaled = vpa(subs(polynomialInX, x, xdata), 5)
valueScaled = vpa(subs(polynomialInXhat, xhat, (xdata-mu(1))./mu(2)), 5)
The difference doesn't really matter that much for the 1st degree polynomial and the small magnitude x data you're using. But suppose you were doing something that required you to take the fourth power of a year, like if you were trying to fit the census data to the population?
load census
format longg
pUnscaled = polyfit(cdate, pop, 4)
That leading coefficient is tiny becuase you're working with large numbers when you raise 2020 to the fourth power. That's why you receive a warning.
2020^4
But if you'd centered and scaled the years from 1900 to 2020:
[pScaled, ~, mu] = polyfit(cdate, pop, 4)
Now you're taking powers of numbers on the order of:
normalizedYears = normalize(cdate, 'center', mu(1), 'scale', mu(2))
and those numbers aren't nearly as large.
normalizedYears(end)^4
I'd much rather work with 6.7 than 16649664160000 and a leading coefficient near 0.7 rather than 4e-8.
카테고리
도움말 센터 및 File Exchange에서 Polynomials에 대해 자세히 알아보기
Translated by 
