Curve Fitting using MATLAB
Aim:
- To generate a code in order to produce a curve which is a best fit for the original data set plot.
- To determine the Goodness of curve fit or how well the equation is representing the data points.
Introduction:
Curve fitting is a process in which we generate a curve or mathematical function which tries to match with the original data set or which possess the closest proximity to the real series of data.
Curve fitting is used to -
- To carry out system simulation and optimisation.
- To determine values of properties of intermediate points.
- To do trend analysis.
- For Hypothesis testing.
What is Best fit and Perfect fit?
Perfect fit is the curve which passes through all the data points and represents the real data.Basically it is a curve fit without any errors.
A best fit is used to represent the scattered data points ,and show what trend it follows.
Here our concern is Best fit-
Let us generate a curve which has linear and cubic order of polynomial.
Application of functions -
load('data') = it loads the 'data' file on which we have to work.
polyfit = polyfit(x,y,n)= polyfit finds the co-efficients of a polynomial p(x) of degree n that fits the data p[x(i)] to y(i) in a least square sense.
polyval = polyval(x,p) evaluates a polynomial for the given set of x values.
The link for data is given as follows-
https://drive.google.com/file/d/1vWi9Vj3AdQf0gqUK1u9jwS8zg92BHGR0/view
Code for Linear Curve Fit is as follows -
The code functions are explained in the code comments.
% Linear Curve Fit
clear all
close all
clc
cp_data = load('data'); %this commands loads the data already stored in directory.
temperature = cp_data(:,1); %loads the temparature data from column 1
cp = cp_data(:,2); %loads the specific heat data from column 2
%code for curve fit
% going by equation - a*T + b, where a and b are coefficients and T is a variable.
%co_effs = polyfit(temperature,cp,1);
%predicted_cp = polyval(co_effs,temperature);
%polyfit = this command generates the co-efficients of polynomial
%which is used to produce a curve.
%polyval = this command evaluates a polynomial for a given set of x values.
% linear curve fit
co_effs_linear = polyfit(temperature,cp,1);
predicted_cp_linear = polyval(co_effs_linear,temperature);
plot(temperature,predicted_cp_linear,'linewidth',1,'color','r');
hold on
plot(temperature, cp, 'linewidth',1,'color','g')
xlabel('Temperature[K]')
ylabel('Specific Heat[KJ/Kmol-K]')
legend('Linear curve fit','original data set')
grid on
title('Linear Curve fit')
n = 3200
k = 1
%for linear fit
error_linear = cp - predicted_cp_linear;
SSE_LINEAR = sum(error_linear.^2);
SSR_LINEAR = sum((predicted_cp_linear - mean(cp)).^2);
SST_LINEAR = SSR_LINEAR + SSE_LINEAR;
RSQUARE_LINEAR = SSR_LINEAR/SST_LINEAR;
RMSE_LINEAR = sqrt(SSE_LINEAR/n);
ADJ_RSQUARE = 1-[(1-RSQUARE_LINEAR)*(n-1)/(n-k-1)];
Code for Cubic Curve Fit is as follows-
% Cubic Curve Fit
clear all
close all
clc
cp_data = load('data');
temperature = cp_data(:,1);
cp = cp_data(:,2);
plot(temperature, cp, 'linewidth',1,'color','g')
hold on
%code for curve fit
% going by equation- a*T^2 + b*T + c
%co_effs = polyfit(temperature,cp,3);
%predicted_cp = polyval(co_effs,temperature);
% cubic curve fit
co_effs_cubic = polyfit(temperature,cp,3);
predicted_cp_cubic = polyval(co_effs_cubic,temperature);
plot(temperature,predicted_cp_cubic,'linewidth',1,'color','r');
xlabel('Temperature[K]')
ylabel('Specific Heat[KJ/Kmol-K]')
legend('Cubic curve fit','original data set')
grid on
title('Cubic Curve fit')
n = 3200
k = 1
%for cubic fit
error_cubic = cp - predicted_cp_cubic;
SSE_CUBIC = sum(error_cubic.^2);
SSR_CUBIC = sum((predicted_cp_cubic - mean(cp)).^2);
SST_CUBIC = SSR_CUBIC + SSE_CUBIC;
RSQUARE_CUBIC = SSR_CUBIC/SST_CUBIC;
RMSE_CUBIC = sqrt(SSE_CUBIC/n);
ADJ_RSQUARE = 1-[(1-RSQUARE_CUBIC)*(n-1)/(n-k-1)];As we can see that our curve fit tries to match with the original data set plot.Still it is not the best curve (RSQUARE < 1). Also if we keep on increasing the order of polynomial,it will throw an error in command window. As we increase the order of polynomial, the curve gets very oscillatory.
Error thrown in command window-
In order to eliminate it , we will use Splitwise Polynomial method to generate a plot and divide the data set in 4 folds,keeping the order of polynomial as 3.
% Splitwise Polynomial Method (Data in 4 parts)
clear all
close all
clc
cp_data = load('data');
temperature = cp_data(:,1);
cp = cp_data(:,2);
%code for curve fit
% going by equation- a*T^2 + b*T + c
%co_effs = polyfit(temperature,cp,2);
%predicted_cp = polyval(co_effs,temperature);
%for linear curve fit
co_effs_0 = polyfit(temperature,cp,1);
predicted_cp_0 = polyval(co_effs_0,temperature);
plot(temperature,predicted_cp_0,'linewidth',1,'color','g');
hold on
%
%co_effs = polyfit(temperature,cp,3);
%predicted_cp_c = polyval(co_effs,temperature);
%applying splitwise polynomial method to cubic polynomial
temperature_1 = cp_data(1:800,1);
cp_1 = cp_data(1:800,2);
co_effs_1 = polyfit(temperature_1,cp_1,3);
predicted_cp_1 = polyval(co_effs_1,temperature_1);
temperature_2 = cp_data(801:1600,1);
cp_2 = cp_data(801:1600,2);
co_effs_2 = polyfit(temperature_2,cp_2,3);
predicted_cp_2 = polyval(co_effs_2,temperature_2);
temperature_3 = cp_data(1601:2400,1);
cp_3 = cp_data(1601:2400,2);
co_effs_3 = polyfit(temperature_3,cp_3,3);
predicted_cp_3 = polyval(co_effs_3,temperature_3);
temperature_4 = cp_data(2401:3200,1);
cp_4 = cp_data(2401:3200,2);
co_effs_4 = polyfit(temperature_4,cp_4,3);
predicted_cp_4 = polyval(co_effs_4,temperature_4);
predicted_cp = [predicted_cp_1;predicted_cp_2;predicted_cp_3;predicted_cp_4];
plot(temperature,predicted_cp,'linewidth',1,'color','b');
hold on
plot(temperature, cp, 'linewidth',1,'color','r')
xlabel('Temperature[K]')
ylabel('Specific Heat[KJ/Kmol-K]')
legend('Linear Curve fit','Splitwise 4','Original Data Set')
grid on
title('Splitwise polynomial method by dividing data in 4 parts')
% to check the goodness of fit
n = 3200
k = 1
%for linear fit
error_linear = cp - predicted_cp_0;
SSE_LINEAR = sum(error_linear.^2);
SSR_LINEAR = sum((predicted_cp_0 - mean(cp)).^2);
SST_LINEAR = SSR_LINEAR + SSE_LINEAR;
RSQUARE_LINEAR = SSR_LINEAR/SST_LINEAR;
RMSE_LINEAR = sqrt(SSE_LINEAR/n);
%for cubic fit
error_cubic = cp - predicted_cp;
SSE_CUBIC = sum(error_cubic.^2);
SSR_CUBIC = sum((predicted_cp - mean(cp)).^2);
SST_CUBIC = SSR_CUBIC + SSE_CUBIC;
RSQUARE_CUBIC = SSR_CUBIC/SST_CUBIC;
RMSE_CUBIC = sqrt(SSE_CUBIC/n);
ADJ_RSQUARE = 1-[(1-RSQUARE_CUBIC)*(n-1)/(n-k-1)];Let us also check the plot by dividing the plot in 8 parts.
% Splitwise Polynomial (Data in 8 parts)
cp_data = load('data');
temperature = cp_data(:,1);
cp = cp_data(:,2);
% linear curve fit
co_effs_0 = polyfit(temperature,cp,1);
predicted_cp_0 = polyval(co_effs_0,temperature);
plot(temperature,predicted_cp_0,'linewidth',2,'color','g');
hold on
%
co_effs = polyfit(temperature,cp,3);
predicted_cp_c = polyval(co_effs,temperature);
%code for curve fit
% going by polynomial - a*T^3 + b*T^2 + c*T + d
%co_effs = polyfit(temperature,cp,3);
%predicted_cp = polyval(co_effs,temperature);
% applying splitwise polynomial method dividing the data in 8 parts.
% using cubic polynomial
temperature_1 = cp_data(1:400,1);
cp_1 = cp_data(1:400,2);
co_effs_1 = polyfit(temperature_1,cp_1,3);
predicted_cp_1 = polyval(co_effs_1,temperature_1);
temperature_2 = cp_data(401:800,1);
cp_2 = cp_data(401:800,2);
co_effs_2 = polyfit(temperature_2,cp_2,3);
predicted_cp_2 = polyval(co_effs_2,temperature_2);
temperature_3 = cp_data(801:1200,1);
cp_3 = cp_data(801:1200,2);
co_effs_3 = polyfit(temperature_3,cp_3,3);
predicted_cp_3 = polyval(co_effs_3,temperature_3);
temperature_4 = cp_data(1201:1600,1);
cp_4 = cp_data(1201:1600,2);
co_effs_4 = polyfit(temperature_4,cp_4,3);
predicted_cp_4 = polyval(co_effs_4,temperature_4);
temperature_5 = cp_data(1601:2000,1);
cp_5 = cp_data(1601:2000,2);
co_effs_5 = polyfit(temperature_5,cp_5,3);
predicted_cp_5 = polyval(co_effs_5,temperature_5);
temperature_6 = cp_data(2001:2400,1);
cp_6 = cp_data(2001:2400,2);
co_effs_6 = polyfit(temperature_6,cp_6,3);
predicted_cp_6 = polyval(co_effs_6,temperature_6);
temperature_7 = cp_data(2401:2800,1);
cp_7 = cp_data(2401:2800,2);
co_effs_7 = polyfit(temperature_7,cp_7,3);
predicted_cp_7 = polyval(co_effs_7,temperature_7);
temperature_8 = cp_data(2801:3200,1);
cp_8 = cp_data(2801:3200,2);
co_effs_8 = polyfit(temperature_8,cp_8,3);
predicted_cp_8 = polyval(co_effs_8,temperature_8);
predicted_cp = [predicted_cp_1;predicted_cp_2;predicted_cp_3;
predicted_cp_4;predicted_cp_5;predicted_cp_6;
predicted_cp_7;predicted_cp_8];
plot(temperature,predicted_cp,'linewidth',3,'color','b');
hold on
plot(temperature, cp, 'linewidth',2,'color','r')
xlabel('Temperature[K]')
ylabel('Specific Heat[KJ/Kmol-K]')
legend('linear curve fit','cubic curve fit','original')
grid on
title('Splitwise_8 Cubic Curve fit')
% to check the goodness of fit
n = 3200;
k = 1;
%for linear fit
error_linear = cp - predicted_cp_0;
SSE_LINEAR = sum(error_linear.^2);
SSR_LINEAR = sum((predicted_cp_0 - mean(cp)).^2);
SST_LINEAR = SSR_LINEAR + SSE_LINEAR;
RSQUARE_LINEAR = SSR_LINEAR/SST_LINEAR;
RMSE_LINEAR = sqrt(SSE_LINEAR/n);
%for cubic fit
error_cubic = cp - predicted_cp;
SSE_CUBIC = sum(error_cubic.^2);
SSR_CUBIC = sum((predicted_cp - mean(cp)).^2);
SST_CUBIC = SSR_CUBIC + SSE_CUBIC;
RSQUARE_CUBIC = SSR_CUBIC/SST_CUBIC;
RMSE_CUBIC = sqrt(SSE_CUBIC/n);
ADJ_RSQUARE = 1-((1-RSQUARE_CUBIC)*(n-1)/(n-k-1));Now we can see that we have got the best possible fit,visually.
So there are certain parameters on which we can judge the Goodness of fit,mathematically.
Let us assume a curve fit of function `y = f(x)` and let `(x1,y1),(x2,y2),....(x_n,y_n)` be the data points.
and let `(x_i,y_i)` be the instantaneous data point from the original data set. Mean is the Arithmetic mean of data provided.n is the total no of data points.
the following is given as -
ERROR = `abs[y(i) - "value of f(x) at x(i)"] = abs[y(i) - f(x_i)]`
sum of square due to errors= SSE
SSE = `sum_(i = 1)^n(ER ROR) = [y(i) - f(x_i)]^2`
SSR = `sum_(i=1)^n[f(x_i) - (Mean)]^2`
SST = SSR + SSE
`R^2 = (SSR)/(SST)`
Goodness of the fit is determined by following parameters=
- The sum of squares due to error (SSE)
- R SQUARE (R^2)
- Adjusted R SQUARE
- Root Mean Squared Error (RMSE).
SSE is the sum of squared deviations of actual values from predicted values.
SSR is the sum of squared deviations of predicted values (predicted using regression) from the mean value.
R SQUARE is the square of the correlation between the response values and the predicted response values.
It is also called the square of the multiple correlation coefficient and the coefficient of multiple determination.
RMSE is known as the fit standard error and the standard error of the regression.
It is an estimate of the standard deviation of the random component in the data.
Adjusted R-square = Adjusted R-square disregards infinitesimally small increments to R-square,
and only considers the variables which actually contribute to the output of the system in a drastic way.
How to improve the cubic fit?
One method is to increase the order of polynomial until the curve fit matches with the original data plot.
The other is to split the data in parts by splitwise method so that the resultant curve matches with the given original plot.
For a best fit,the R SQUARE value has to be 1 and RMSE value has to be minimum.
Using these four quantities Matlab/Octave effectively deduces the goodness of fit, typically if R-square is greater than 0.95 the fit is assumed to be good.
Let us have a look at the workspace data of various graphs mentioned and deduce which data represented for a best fit.
Linear Curve Fit Data-
Cubic Curve Fit Data-
Splitwise polynomial method curve fit by dividing data in 4 parts-
Splitwise polynomial method curve fit by dividing data in 8 parts-
Since the Splitwise polynomial method used to divide the data in 8 parts represents that
R SQUARE value is 1 and RMSE value is minimum, in comparison to other data.
Hence the Cubic Curve Fit using Splitwise polynomial in 8 parts is the best curve fit.
REFERENCES -
https://projects.skill-lync.com/projects/Curve-Fitting-and-criterion-to-choose-best-fit-41860
NPTEL :: Mechanical Engineering - Design and Optimization of Energy systems