Linear and Cubic curve fitting
CURVE FITTING:
The objective of curve fitting is to theoretically describe experimental data with a model (function or equation) and to find the parameters associated with this model.
LINEAR LEAST SQUARES METHOD:
Curve fitting finds the values of the coefficients (parameters) which make a function match the data as closely as possible. The best values of the coefficients are the ones that minimize the value of the expression given below:
`sum_(i=0)^n(y(i) - f(x(i))^2`
The above expression is called sum of squared errors (SSE).
For the linear least squares method, the curve fit equation can be only a polynomial of nth order such as the following:
`y= a*x^n + b*x^(n-1) + c*x^(n-2) + ............. + d`
In the above expression, `a,b,c,......,d` are the coefficients/parameters which need to be calculated. `y` denotes the dependent variable and `x` is the independent variable. For the linear least squares method, the dependent variable `y` is always a linear function of al the coefficients. If the curve has `y` as a non linear function of any of the coefficients (for eg. `y = a*e^(b*x)`), then we need to use non-linear least squares method.
The above plot represents a linear curve fit for a set of data points. In the above plot, the value of error is given by:
Error = `abs(y(i) - f(x(i)))`
Since we want the coefficients to be such, that the value of SSE (and thus the error) to be minimum, we can find out the coefficients in the following way:
`[[(delSSE)/(dela) = 0],[(delSSE)/(delb) = 0],[(delSSE)/(delc) = 0],[.],[.],[.],[(delSSE)/(deld) = 0]]`
BEST FIT AND PERFECT FIT:
A curve fit is the best fit for a particular set of data if the error between the predicted and the actual data values i.e. the error is minimised. That being said, if a curve is a best fit, it doesn't necessarily mean that the curve passes through all the original points. If a curve fit passses through all the original data points, then the error in that case is 0 and the curve fit is called perfect fit. In essence, every perfect fit is the best fit, but every best fit is not the perfect fit.
There are many ways to make/determine if a curve fit is the best fit or not based on some parameters which will be discussed shortly after.
GOODNESS OF FIT PARAMETERS:
1). SSE
2). R-Squared(R): It is a measure of how successful the fit is in explaining the variation of the data. The value of R lies between 0 and 1 and is directly proportional to the goodness of the fit i.e. greater the value of R, better the fit is. Usually if a curve gives a R of greater than 0.95, it is considered agood fit.
`R = (SSR)/(SST)`
where, `SSR` is the sum of squares of regression and is given by:
`SSR = sum_(i=0)^n(f(x(i)-M)^2` , where `M` denotes the mean of all the independent variable inputs.
`SST = SSR + SSE`
3). Adjusted R-Squared: This is the normalised version of R-square which neglects the variation of R-Square due to insignificant independent variables
4). Root mean squared error (RMSE): It measures how spread out the original data is around the best fit curve.
`RMSE = sqrt((SSE)/n)` where n is the number of original datapoints.
For this project we will be considering a data of 3200 values of temperature (independent variables) and their corresponding oxygen specific heat values (dependent variables). We will try to attempt to fit two types of curves: Cubic and Linear
LINEAR:
Initially, we will be preparing the data by loading it into MATLAB and will be assigining the Cp and T data to arrays.
%preparing data
ip_data = load('data');%loading the data values into MATLAB
temperature = ip_data(:,1);%temperature values are the first column of the data vector
c_p = ip_data(:,2); %specific heat values are the second column of the data vector
In this, we will be using the polyfit function not to find out just the coefficients of the fit, but also to find out the estimated errors and the mean and standard deviation of the temperature values.
%preparing the linear curve fit
%calculating the coefficients of the curve fit
%creating a structure which gives the error estimates b/w original and fit
%creating a 2 element vector which has the mean and standard deviation of x
%centering temperature at 0 and scaling it to have unit standard deviation
[coeffs s mu] = polyfit(temperature,c_p,1); %curve fit is linear (of order 1)
[c_p_predicted deltaE] = polyval(coeffs,temperature,s,mu);
The array deltaE carries all the errors for each data point. The two element vector mu carries the mean and the standard deviation of the temperature values. Using these values, the polyfit centers the temperature values at 0 and scales it to have unit standard deviation.
To compare the curve fit and the original data points, we can plot the curves on the same graph.
hold on
%plotting the curve fit
plot(temperature,c_p_predicted,'linewidth',3)
%plotting the orginal data curve
plot(temperature,c_p,'linewidth',3)
xlabel('TEMPERATURE')
ylabel('SPECIFIC HEAT')
legend('Curve-Fit','Original Data Curve')
title('LINEAR CURVE FIT')The plot obtained:
To find the goodness of the fit, we can find out the values of the various parameters discussed (SSE, R, RMSE).
%preparing the sum of squared errors (SSE) array
for i = 1:length(c_p_predicted)
a = c_p(i) - c_p_predicted(i);
sse(i) = a^2; %defining the value of each term of array ssr
end
%summation of all elements of sse
SSE = sum(sse)
%preparing the sum of squares of regression (SSR)
%mean of the predicted c_p values for each input value of temperature
mean_fx = sum(c_p_predicted)/length(c_p_predicted);
mean_x = mu(1); %mean of all the temperature value inputs
for i = 1:length(c_p_predicted)
b = c_p_predicted(i) - mean_fx;
ssr(i) = b^2; %defining the value of each element of array ssr
end
%summation of all elements of array ssr
SSR = sum(ssr)
%total sum of squares=sum of squares of regression+sum of squares of error
SST = SSR+SSE
%R-square
r_squared = SSR/SST
%root mean squared error(how much dispersed the original data is about fit)
rmse = (SSE/length(c_p_predicted))^0.5
%average error estimate
error = sum(deltaE)/length(deltaE)The calculated values are:
Since, R<0.95, thus linear fit is not that good.
We can also use the MATLAB inbuilt curve fitting toolbox to judge the goodness of the fit.
Residual Plot:
Residual plot plots the residuals for each data point. Residuals are the quantity which are left when the fit data is substracted from the raw data.
The complete linear curve fit code is given below:
clear all
close all
clc
%preparing data
ip_data = load('data');%loading the data values into MATLAB
temperature = ip_data(:,1);%temperature values are the first column of the data vector
c_p = ip_data(:,2); %specific heat values are the second column of the data vector
%preparing the linear curve fit
%calculating the coefficients of the curve fit
%creating a structure which gives the error estimates b/w original and fit
%creating a 2 element vector which has the mean and standard deviation of x
%centering temperature at 0 and scaling it to have unit standard deviation
[coeffs s mu] = polyfit(temperature,c_p,1); %curve fit is linear (of order 1)
[c_p_predicted deltaE] = polyval(coeffs,temperature,s,mu);
hold on
%plotting the curve fit
plot(temperature,c_p_predicted,'linewidth',3)
%plotting the orginal data curve
plot(temperature,c_p,'linewidth',3)
xlabel('TEMPERATURE')
ylabel('SPECIFIC HEAT')
legend('Curve-Fit','Original Data Curve')
title('LINEAR CURVE FIT')
%preparing the sum of squared errors (SSE) array
for i = 1:length(c_p_predicted)
a = c_p(i) - c_p_predicted(i);
sse(i) = a^2; %defining the value of each term of array ssr
end
%summation of all elements of sse
SSE = sum(sse)
%preparing the sum of squares of regression (SSR)
%mean of the predicted c_p values for each input value of temperature
mean_fx = sum(c_p_predicted)/length(c_p_predicted);
mean_x = mu(1); %mean of all the temperature value inputs
for i = 1:length(c_p_predicted)
b = c_p_predicted(i) - mean_fx;
ssr(i) = b^2; %defining the value of each element of array ssr
end
%summation of all elements of array ssr
SSR = sum(ssr)
%total sum of squares=sum of squares of regression+sum of squares of error
SST = SSR+SSE
%R-square
r_squared = SSR/SST
%root mean squared error(how much dispersed the original data is about fit)
rmse = (SSE/length(c_p_predicted))^0.5
%average error estimate
error = sum(deltaE)/length(deltaE)CUBIC:
The code and parameters will stay the same forthecubic curve, except now we will have to calculate a total of 4 coefficients for a cubic curve.
The code for the cubic fit is given below:
clear all
close all
clc
%preparing data
ip_data = load('data');
temperature = ip_data(:,1);
c_p = ip_data(:,2);
%preparing the linear curve fit
%calculating the coefficients of the curve fit
%creating a 2 element vector which has the mean and standard deviation of x
%centering temperature at 0 and scaling it to have unit standard deviation
[coeffs s mu] = polyfit(temperature,c_p,3); %curve fit is cubic (of order 3)
[c_p_predicted deltaE] = polyval(coeffs,temperature,s,mu);
hold on
%plotting the curve fit
plot(temperature,c_p_predicted,'linewidth',3)
%plotting the orginal data curve
plot(temperature,c_p,'linewidth',3)
xlabel('TEMPERATURE')
ylabel('SPECIFIC HEAT')
legend('Curve-Fit','Original Data Curve')
title('CUBIC CURVE FIT')
%preparing the sum of squared errors (SSE) array
for i = 1:length(c_p_predicted)
a = c_p(i) - c_p_predicted(i);
sse(i) = a^2; %defining the value of each term of array ssr
end
%summation of all elements of sse
SSE = sum(sse)
%preparing the sum of squares of regression (SSR)
%mean of the predicted c_p values for each input value of temperature
mean_fx = sum(c_p_predicted)/length(c_p_predicted);
mean_x = mu(1);
for i = 1:length(c_p_predicted)
b = c_p_predicted(i) - mean_fx;
ssr(i) = b^2; %defining the value of each element of array ssr
end
%summation of all elements of array ssr
SSR = sum(ssr)
%total sum of squares=sum of squares of regression+sum of squares of error
SST = SSR+SSE
%R-square
r_squared = SSR/SST
%root mean squared error(how much dispersed the original data is about fit)
rmse = (SSE/length(c_p_predicted))^0.5
%average error estimate
error = sum(deltaE)/length(deltaE)
The plot obtained is:
The calculated parameters are:
We can validate the results with the help of curve fitting tool:
Since, R>0.95, the cubic fit can be considered a good fit.
Residual Plot:
OBTAINING THE BEST FIT:
The first method to increase the goodness of a fit is by increasing the order of the curve fit polynomial.
We can observe the graphs for the increasing order polynomials in order to justify this:
4th order polynomial curve fit:
Goodness of fit:
5th order polynomial curve fit:
Goodness of fit:
6th order polynomial curve fit:
Goodness of fit:
7th order polynomial curve fit:
Goodness of fit:
8th order polynomial fit:
Goodness of fit:
9th order polynomial fit:
Goodness of fit:
As seen above, the plots for the 8th and 9th order curve fits depict perfect curves i.e. there is no difference between the predicted and the actual specific heat values (error = 0). This shows us that, by increasing the order of the polynomial, we can increase the goodness of a fit and obtain the best fit.
Another good way of increasing the goodness of the fit is by making the use of split curve fits. We can use polynomials with multiple orders for specific different ranges of temperature in order to make the fit more alligned to the original data.
In this case, for the first 550 values of temperature (300K to 850K), we can keep the curve fit as a linear polynomial, and for rest of the temperature ranges, we can keep the curve fit as a cubic polynomial.
The code for the split curve fit is given below:
clear all
close all
clc
%preparing data
cp_data = load('data'); %reading data from the data input file
temperature = cp_data(:,1);
temperature1 = cp_data(1:550,1);
temperature2 = cp_data(551:end,1);
cp = cp_data(:,2);
%curve fitting
%assume cp is dependant on tempearture in form of a nth order ploynomial
%calculating the coeffs of the polynomial from the data
coeffs1 = polyfit(temperature1,cp(1:550),1); %second order polynomial
coeffs2 = polyfit(temperature2,cp(551:end),3); %fifth order polynomial
predicted_cp1 = polyval(coeffs1,temperature1); %calculating predicted cp (Y) from the value of the poynomial and temperature (X)
predicted_cp2 = polyval(coeffs2,temperature2);
%comparing with original data and cp value
%plotting
plot(temperature,cp,'linewidth',6,'color','r')
hold on
plot(temperature1,predicted_cp1,'linewidth',3,'color','b')
plot(temperature2,predicted_cp2,'linewidth',3,'color','b')
xlabel('TEMPERATURE')
ylabel('SPECIFIC HEAT')
title('SPLIT CURVES')
legend('CURVE-FIT','ORIGINAL DATA SET')The plot obtained by split curve fitting is:
As seen in the above figure, the curve fit (blue) allignes to the original data curve (red) much better as compared to the linear or cubic plots.
REFERENCES:
http://collum.chem.cornell.edu/documents/Intro_Curve_Fitting.pdf
https://www.statisticshowto.datasciencecentral.com/rmse/
https://terpconnect.umd.edu/~toh/models/CurveFitting.html
https://ipython-books.github.io/93-fitting-a-function-to-data-with-nonlinear-least-squares/
https://projects.skill-lync.com/projects/Curve-Fitting-and-criterion-to-choose-best-fit-41860