Extrapolation: How do I extrapolate a curve with minimum error.

조회 수: 9 (최근 30일)
Anand Ra
Anand Ra 2021년 7월 7일
답변: Scott MacKenzie 2021년 7월 7일
I would like to extrapolate the following data set in the interval of 300 seconds, assuming it reaches equilibrium ( =1) after 3600seconds from the data with minimum error. Kindly requesting help. Thanks!
  • I currentlly have data collected upto 18600 seconds in the interval of 300 seconds).
  • From most the answers that I have browsed I understand the first step to use any cruve fitting tool is to use polyfit and use interp
  • I am not sure how to proceed to extrapolate using interp.
Below is my code generated from fitting using the curvefitting tool app - to monitor the residual and other fitting stats
% The time data
t1 = [0:300:18600]';
% processed data
y_obs = [
0
1.6216e-01
2.9813e-01
4.0805e-01
4.9338e-01
5.5928e-01
6.0894e-01
6.5506e-01
6.8876e-01
7.1773e-01
7.4284e-01
7.6433e-01
7.8448e-01
7.9912e-01
8.1484e-01
8.2950e-01
8.3760e-01
8.4707e-01
8.5767e-01
8.6299e-01
8.6862e-01
8.7433e-01
8.7879e-01
8.8736e-01
8.9152e-01
8.9903e-01
9.0343e-01
9.0984e-01
9.1344e-01
9.1615e-01
9.2046e-01
9.2313e-01
9.2640e-01
9.2992e-01
9.3134e-01
9.3242e-01
9.3616e-01
9.3986e-01
9.4201e-01
9.4145e-01
9.4434e-01
9.4252e-01
9.4131e-01
9.4249e-01
9.4283e-01
9.4395e-01
9.4355e-01
9.4690e-01
9.4919e-01
9.5255e-01
9.5626e-01
9.6282e-01
9.6283e-01
9.6672e-01
9.6439e-01
9.6887e-01
9.7099e-01
9.7529e-01
9.8034e-01
9.8302e-01
9.8494e-01
9.8828e-01
9.8769e-01
];
%% Fit: 'untitled fit 1'.
[xData, yData] = prepareCurveData( t1, y_obs );
% Set up fittype and options.
ft = fittype( 'exp2' );
opts = fitoptions( 'Method', 'NonlinearLeastSquares' );
opts.Display = 'Off';
opts.StartPoint = [0.854210120329305 7.68324724853675e-06 -0.838821705600288 -0.000648039892367098];
% Fit model to data.
[fitresult, gof] = fit( xData, yData, ft, opts );
% Create a figure for the plots.
figure( 'Name', 'untitled fit 1' );
% Plot fit with data.
subplot( 2, 1, 1 );
h = plot( fitresult, xData, yData );
legend( h, 'y_obs vs. t1', 'untitled fit 1', 'Location', 'NorthEast', 'Interpreter', 'none' );
% Label axes
xlabel( 't1', 'Interpreter', 'none' );
ylabel( 'y_obs', 'Interpreter', 'none' );
grid on
% Plot residuals.
subplot( 2, 1, 2 );
h = plot( fitresult, xData, yData, 'residuals' );
legend( h, 'untitled fit 1 - residuals', 'Zero Line', 'Location', 'NorthEast', 'Interpreter', 'none' );
% Label axes
xlabel( 't1', 'Interpreter', 'none' );
ylabel( 'y_obs', 'Interpreter', 'none' );
grid on
  댓글 수: 4
Anand Ra
Anand Ra 2021년 7월 7일
That was funny!
Btw I see "interp" is also used to extrapolate. Is there any other function or solver that I can use to extrapolate my data?
Star Strider
Star Strider 2021년 7월 7일
I like the Galsworthy quote, too!
The interp1 function allows extrapolation. To do that it is necessary to specify an extrapolation method.
You will need to determine the correct method that gives the result you want.
(I generally do not suggest extrapolating, since the behaviour of the data beyond the region-of-fit is entirely unknown. If the model suggests an asymptote, that could be appropriate, however in general, extrapolation is not advisable.)
.

댓글을 달려면 로그인하십시오.

채택된 답변

Mathieu NOE
Mathieu NOE 2021년 7월 7일
hello
my suggestion for exponential curve fitting and extrapolation (adjust the new time vector t2 according to your needs)
% The time data
t1 = [0:300:18600]';
% processed data
y_obs = [
0
1.6216e-01
2.9813e-01
4.0805e-01
4.9338e-01
5.5928e-01
6.0894e-01
6.5506e-01
6.8876e-01
7.1773e-01
7.4284e-01
7.6433e-01
7.8448e-01
7.9912e-01
8.1484e-01
8.2950e-01
8.3760e-01
8.4707e-01
8.5767e-01
8.6299e-01
8.6862e-01
8.7433e-01
8.7879e-01
8.8736e-01
8.9152e-01
8.9903e-01
9.0343e-01
9.0984e-01
9.1344e-01
9.1615e-01
9.2046e-01
9.2313e-01
9.2640e-01
9.2992e-01
9.3134e-01
9.3242e-01
9.3616e-01
9.3986e-01
9.4201e-01
9.4145e-01
9.4434e-01
9.4252e-01
9.4131e-01
9.4249e-01
9.4283e-01
9.4395e-01
9.4355e-01
9.4690e-01
9.4919e-01
9.5255e-01
9.5626e-01
9.6282e-01
9.6283e-01
9.6672e-01
9.6439e-01
9.6887e-01
9.7099e-01
9.7529e-01
9.8034e-01
9.8302e-01
9.8494e-01
9.8828e-01
9.8769e-01
];
f = @(a,b,c,x) a.*(1-exp(-(x.*b).^c));
obj_fun = @(params) norm(f(params(1), params(2), params(3),t1)-y_obs);
sol = fminsearch(obj_fun, [y_obs(end),0,1]);
a_sol = sol(1)
b_sol = sol(2)
c_sol = sol(3)
% y_fit = f(a_sol, b_sol,c_sol, t1);
% extrapolation
t2 = [0:300:18600+30*300]';
y_extra = f(a_sol, b_sol,c_sol, t2);
figure;
plot(t1, y_obs, '+', 'MarkerSize', 10, 'LineWidth', 2)
hold on
plot(t2, y_extra);grid on
xlabel('time t');
ylabel('y');

추가 답변 (1개)

Scott MacKenzie
Scott MacKenzie 2021년 7월 7일
Here's what I put together. The extrapolation uses the 4-coefficient model returned by fit in your code:
xx = 1:1000:40000;
yy = fitresult.a * exp(fitresult.b*xx) + fitresult.c*exp(fitresult.d*xx);
% Create a figure for the plots.
figure( 'Name', 'untitled fit 1' );
% Plot fit with data.
h = plot( fitresult, xData, yData );
hold on;
plot(xx, yy, 'linewidth', 1); % extrapolation from model
set(gca,'ylim', [0 1.2]);
legend( h, 'y_obs vs. t1', 'untitled fit 1', 'Location', 'SouthEast', 'Interpreter', 'none' );
xlabel( 't1', 'Interpreter', 'none' );
ylabel( 'y_obs', 'Interpreter', 'none' );
grid on

카테고리

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

Community Treasure Hunt

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

Start Hunting!

Translated by