Polynomial fit with centering & scaling

조회 수: 45 (최근 30일)
Sohel Rana
Sohel Rana 2022년 3월 16일
편집: Torsten 2022년 3월 17일
Hello,
When I plotted the following data and applied polyfit with second order, it showed the warning like: Warning: Polynomial is badly conditioned. Add points with distinct X values, reduce the degree of the polynomial, or try centering and scaling as described in HELP POLYFIT. This fitting can also be done using [p, S, mu] = polyfit(x,y,2). However, the fitting coefficients from [p, S, mu] = polyfit(x, y, 2) are totally different from p=polyfit(x,y,2). How can I use the [p, S, mu] = polyfit(x,y,2) and obtain coefficients that correspond to my unscaled system p=polyfit(x,y,2)?
x y
xy = [
1508.936 19.2189
1509.159 54.27226
1509.472 102.8022
1509.852 154.8707
1510.224 199.6335
1510.809333 258.6648
1511.443429 316.4253
1512.112 370.7319
1512.4112 391.8055
1513.118545 441.4118
1514.002667 497.8672
1514.777778 547.8429
1515.6144 598.6195
1516.467 650.4304
1517.3032 701.3778 ]

답변 (3개)

Torsten
Torsten 2022년 3월 17일
편집: Torsten 2022년 3월 17일
Here is an example that shows how to proceed:
x = [1,2,3];
y = [4,9,25];
order = 1;
p = polyfit(x,y,order);
[p1,S,mu] = polyfit(x,y,order);
% Transform coefficients to those of the polynomial p centered at 0
p1 = flip(p1); % Flip order to [p0, ..., pn]
p2 = zeros(1,order+1);
for i = 0:order
for k = i:order
p2(i+1) = p2(i+1) + nchoosek(k, k-i) * p1(k+1)/mu(2)^k * (-mu(1))^(k-i);
end
end
p2 = flip(p2); % Back to original order [pn, ..., p0]
% Compare
p
p2
and here is the link:

Walter Roberson
Walter Roberson 2022년 3월 17일
편집: Walter Roberson 2022년 3월 17일
xy = [
1508.936 19.2189
1509.159 54.27226
1509.472 102.8022
1509.852 154.8707
1510.224 199.6335
1510.809333 258.6648
1511.443429 316.4253
1512.112 370.7319
1512.4112 391.8055
1513.118545 441.4118
1514.002667 497.8672
1514.777778 547.8429
1515.6144 598.6195
1516.467 650.4304
1517.3032 701.3778 ];
x = xy(:,1); y = xy(:,2);
format long g
[p, S, mu] = polyfit(x,y,2)
p = 1×3
-36.6930527200109 231.045137444742 387.978493205332
S = struct with fields:
R: [3×3 double] df: 12 normr: 46.6308023944645
mu = 2×1
1.0e+00 * 1512.38017013333 2.74485447071878
xscaled = (x-mu(1))/mu(2);
pscaled = polyfit(xscaled, y, 2)
pscaled = 1×3
-36.6930527200109 231.045137444742 387.978493205332
syms X t
Xs = (X-mu(1))/mu(2)
Xs = 
punscaled = collect(subs(poly2sym(p, t), t, Xs), X)
punscaled = 
vpa(punscaled, 16)
ans = 
polyfit(x, y, 2)
ans = 1×3
1.0e+00 * -4.8701820734741 14815.3074972764 -11266452.1354997

Torsten
Torsten 2022년 3월 17일
편집: Torsten 2022년 3월 17일
Another option:
xy = [
1508.936 19.2189
1509.159 54.27226
1509.472 102.8022
1509.852 154.8707
1510.224 199.6335
1510.809333 258.6648
1511.443429 316.4253
1512.112 370.7319
1512.4112 391.8055
1513.118545 441.4118
1514.002667 497.8672
1514.777778 547.8429
1515.6144 598.6195
1516.467 650.4304
1517.3032 701.3778 ];
x = xy(:,1); y = xy(:,2);
order = 2;
p = polyfit(x,y,order);
[p1,S,mu] = polyfit(x,y,order);
syms t
p1s = flip(p1)*((t-mu(1))/mu(2)).^(0:order).'; % create polynomial p1
cp1s = coeffs(taylor(p1s,t,0)); % compute Taylor coefficients in t=0
cp1 = double(cp1s); % convert to double
pp = flip(cp1); % flip to original order
p
pp

카테고리

Help CenterFile Exchange에서 Polynomials에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by