curve fit a custom polynomial

6 ビュー (過去 30 日間)
Randy Chen
Randy Chen 2021 年 10 月 26 日
コメント済み: Star Strider 2021 年 10 月 27 日
I have the following 2nd order polynomial in the r-z coordinates:
Right now I have four sets of coordinats (r,z), how should I do a curve fit such that I can get an expression for Z in terms of r?
z = [-6.41 -12.4 2.143 102];
r = [13.58 15.7636 12.96 46.6];

採用された回答

Star Strider
Star Strider 2021 年 10 月 27 日
One approach —
syms A B C D r z
sf = A*r^2 + B*z + C*r + D == 0
sf = 
sfiso = isolate(sf, z)
sfiso = 
zfcn = matlabFunction(rhs(sfiso), 'Vars',{[A B C D],r})
zfcn = function_handle with value:
@(in1,r)-(in1(:,4)+in1(:,3).*r+in1(:,1).*r.^2)./in1(:,2)
z = [-6.41 -12.4 2.143 102];
r = [13.58 15.7636 12.96 46.6];
B0 = rand(1,4);
nlm = fitnlm(r, z, zfcn, B0)
Warning: The model is overparameterized, and model parameters are not identifiable. You will not be able to compute confidence or prediction intervals, and you should use caution in making predictions.
nlm =
Nonlinear regression model: y ~ F(in1,r) Estimated Coefficients: Estimate SE tStat pValue ________ ________ _______ __________ b1 -0.61998 0.070975 -8.7353 0.072563 b2 2.4979 0.87073 2.8688 0.21353 b3 29.339 1.4571 20.135 0.031591 b4 -275.65 0.16252 -1696.1 0.00037534 Number of observations: 4, Error degrees of freedom: 1 Root Mean Squared Error: 3.89 R-Squared: 0.998, Adjusted R-Squared 0.995 F-statistic vs. constant model: 290, p-value = 0.0415
Bv = nlm.Coefficients.Estimate;
pv = nlm.Coefficients.pValue;
Out = table({'A';'B';'C';'D'},Bv,pv, 'VariableNames',{'Parameter','Value','p-Value'})
Out = 4×3 table
Parameter Value p-Value _________ ________ __________ {'A'} -0.61998 0.072563 {'B'} 2.4979 0.21353 {'C'} 29.339 0.031591 {'D'} -275.65 0.00037534
rv = linspace(min(r), max(r));
zv = predict(nlm, rv(:));
figure
plot(r, z, 'pg')
hold on
plot(rv, zv, '-r')
hold off
grid
xlabel('r')
ylabel('z')
legend('Data','Model Fit', 'Location','best')
The Warning was thrown because the number of parameters are not less than the number of data pairs.
Experiment to get different results.
.
  2 件のコメント
Randy Chen
Randy Chen 2021 年 10 月 27 日
Thank you! I was trying out the code myself, but this error occured:
syms A B C D r z ;
sf = A*r^2 + B*z + C*r + D == 0;
sfiso = isolate(sf,z);
zfcn = matlabFunction(rhs(sfiso), 'Vars',{[A B C D],r});
z = [12.96 46.6 13.5 46.5188]
r = [2.143 102 2.41814 101.5];
B0 = rand(1,4);
nlm = fitnlm(r, z, zfcn, B0);
Bv = nlm.Coefficients.Estimate;
pv = nlm.Coefficients.pValue;
Out = table({'A';'B';'C';'D'},Bv,pv, 'VariableNames',{'Parameter','Value','p-Value'});
rv = linspace(min(r), max(r));
zv = predict(nlm, rv(:));
figure
plot(r, z, 'pg')
hold on
plot(rv, zv, '-r')
hold off
grid
xlabel('r')
ylabel('z')
legend('Data','Model Fit', 'Location','best')
Unrecognized function or variable 'charcmd'.
Error in sym/isolate (line 108)
throw(CaughtMException);
Error in hw10 (line 4)
sfiso = isolate(sf,z);
Star Strider
Star Strider 2021 年 10 月 27 日
As always, my pleasure!
The isolate function was introduced in R2017a.
Use solve instead —
syms A B C D r z
sf = A*r^2 + B*z + C*r + D == 0
sf = 
sfiso = solve(sf, z)
sfiso = 
zfcn = matlabFunction(sfiso, 'Vars',{[A B C D],r})
zfcn = function_handle with value:
@(in1,r)-(in1(:,4)+in1(:,3).*r+in1(:,1).*r.^2)./in1(:,2)
z = [-6.41 -12.4 2.143 102];
r = [13.58 15.7636 12.96 46.6];
B0 = rand(1,4);
nlm = fitnlm(r, z, zfcn, B0)
Warning: The model is overparameterized, and model parameters are not identifiable. You will not be able to compute confidence or prediction intervals, and you should use caution in making predictions.
nlm =
Nonlinear regression model: y ~ F(in1,r) Estimated Coefficients: Estimate SE tStat pValue ________ ________ _______ __________ b1 -0.725 0.082997 -8.7353 0.072563 b2 2.9211 1.0182 2.8688 0.21353 b3 34.308 1.7039 20.135 0.031591 b4 -322.35 0.19005 -1696.1 0.00037534 Number of observations: 4, Error degrees of freedom: 1 Root Mean Squared Error: 3.89 R-Squared: 0.998, Adjusted R-Squared 0.995 F-statistic vs. constant model: 290, p-value = 0.0415
Bv = nlm.Coefficients.Estimate;
pv = nlm.Coefficients.pValue;
Out = table({'A';'B';'C';'D'},Bv,pv, 'VariableNames',{'Parameter','Value','p-Value'})
Out = 4×3 table
Parameter Value p-Value _________ _______ __________ {'A'} -0.725 0.072563 {'B'} 2.9211 0.21353 {'C'} 34.308 0.031591 {'D'} -322.35 0.00037534
rv = linspace(min(r), max(r));
zv = predict(nlm, rv(:));
figure
plot(r, z, 'pg')
hold on
plot(rv, zv, '-r')
hold off
grid
xlabel('r')
ylabel('z')
legend('Data','Model Fit', 'Location','best')
I like isolate because of the output format, and some of its other characteristics.
.

サインインしてコメントする。

その他の回答 (1 件)

Rik
Rik 2021 年 10 月 27 日
You have two options: rewrite your equation to be a pure quadratic and use polyfit, or use a function like fit or fminsearch on this shape.
If you have trouble implementing either of these two, feel free to comment with what you tried.
  2 件のコメント
Randy Chen
Randy Chen 2021 年 10 月 27 日
I'm not sure how to rewrite that expression to be pure quadratic. Is it Bz = -Ar^2-Cr-D? But why do I have to rewrite it before using polyfit?
Rik
Rik 2021 年 10 月 27 日
Polyfit will fit a pure polynomial of the form
f(x)=p(1)*x^n +p(2)*x^(n-1) ... +p(n)*x +p(n+1)
That means you can determine the values of -A/B, -C/B, and -D/B with polyfit.
As you may conclude from this: there is no unique solution for your setup, unless you have other restrictions to the values you haven't told yet.
z = [-6.41 -12.4 2.143 102];
r = [13.58 15.7636 12.96 46.6];
p=polyfit(r,z,2)
p = 1×3
0.2482 -11.7452 110.3524

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeLinear and Nonlinear Regression についてさらに検索

タグ

Community Treasure Hunt

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

Start Hunting!

Translated by