Curve fitting
Aim:- To fit a linear and cubic polynomial for the given Cp data and measure goodness of fit.
Theory:-
The process of constructing a curve using the available data and finding its equation is known as curve fitting.Its very helpful in predicting the future values of any data set.
There are 4 quantities which help us to measure the goodness of fit
1. The sum of squares due to error (SSE)
`SSE ("sum of squared errors") = sum_(i=1)^n"Error(i)"^2`
`SSE = sum_(i=1)^n(Y(i)-f(x(i)))^2`
2. R-square
`SSR = sum_(i=1)^n(f(x(i))-Mean)^2`
`SST=SSR+SSE`
`R^2=(SSR)/(SST)`
3. Adjusted R-square
`R_(adj)^2=1-(1-"R"^2)[(n-1)/(n-k-1)]`
4. Root mean squared error (RMSE)
`RMSE = sqrt((sum_(i=1)^n(Y(i)-f(x(i)))^2)/n) = sqrt((SSE)/n)`
Code:-
#Importing libraries
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
#1st order fit function
def or1(t,a,b):
return a*t + b
#3rd order fit function
def or3(t,a,b,c,d):
return a*pow(t,3) + b*pow(t,2) + c*t + d
#5th order fit function
def or5(t,a,b,c,d,e,f):
return a*pow(t,5) + b*pow(t,4) + c*pow(t,3) + d*pow(t,2) + e*t + f
#7th order fit function
def or7(t,a,b,c,d,e,f,g,h):
return a*pow(t,7) + b*pow(t,6) + c*pow(t,5) + d*pow(t,4) + e*pow(t,3) + f*pow(t,2) + g*t + h
#8th order fit function
def or8(t,a,b,c,d,e,f,g,h,i):
return a*pow(t,8) + b*pow(t,7) + c*pow(t,6) + d*pow(t,5) + e*pow(t,4) + f*pow(t,3) + g*pow(t,2) + h*t + i
#reading thermodynamic data file
def read_file():
temperature=[]
cp=[]
for line in open('data','r'):
values=line.split(',')
temperature.append(float(values[0]))
cp.append(float(values[1]))
return [temperature,cp]
#goodness evaluation function
def gof(fit_cp,popt):
#R_square
SSR=np.sum(pow((np.array(np.array(fit_cp)-np.mean(cp))),2))
SSE=np.sum(pow(np.array(np.array(cp)-np.array(fit_cp)),2))
SST=SSR+SSE
Rsq=SSR/SST
#Adjusted R_square
n=len(cp)
k=0
arsq=1-((1-Rsq)*(n-1)/(n-k-1))
#Root mean squared error(RMSE)
rmse=np.sqrt(SSE/n)
print("{:2d}\t{:21.6f}\t{:12.4f}\t{:12.4f}\t{:19.4f}".format(len(popt)-1,SSE,Rsq,arsq,rmse))
#main program
temperature,cp=read_file()
#getting coefficients and covarience of coefficients
popt1, pcov1 = curve_fit(or1,temperature,cp)
popt3, pcov3 = curve_fit(or3,temperature,cp)
popt5, pcov5 = curve_fit(or5,temperature,cp)
popt7, pcov7 = curve_fit(or7,temperature,cp)
popt8, pcov8 = curve_fit(or8,temperature,cp)
#fitting curve
fit_cp1=or1(np.array(temperature), *popt1)
fit_cp3=or3(np.array(temperature), *popt3)
fit_cp5=or5(np.array(temperature), *popt5)
fit_cp7=or7(np.array(temperature), *popt7)
fit_cp8=or8(np.array(temperature), *popt8)
#plotting the curves
plt.figure()
plt.plot(temperature,cp,color="r",lw=2)
plt.plot(temperature,fit_cp1,color="g",lw=2)
plt.plot(temperature,fit_cp3,color="b",lw=2)
plt.legend(['Actual data','Linear fit','Cubic fit'])
plt.xlabel('Temperature[K]',fontsize=12)
plt.ylabel('Specific Heat',fontsize=12)
plt.title('Curve Fitting',fontsize=15)
#plotting best fit
plt.figure()
plt.plot(temperature,cp,color="r",lw=2)
plt.plot(temperature,fit_cp8,color="b",lw=2)
plt.legend(['Actual data','8th order fit'])
plt.xlabel('Temperature[K]',fontsize=12)
plt.ylabel('Specific Heat',fontsize=12)
plt.title('Perfect fit',fontsize=15)
# Tabulation of results
print('\t\t\t\t\tFitness characteristics / Goodness of fit')
print('--------------------------------------------------------------------------------------')
print('Order\tSum of Squared Error\tR-square\tAdjuste R-square\tRoot Mean Square Error')
print('--------------------------------------------------------------------------------------')
#goodness evaluation
gof(fit_cp1,popt1)
gof(fit_cp3,popt3)
gof(fit_cp5,popt5)
gof(fit_cp7,popt7)
gof(fit_cp8,popt8)
plt.show()
Explanation:-
- First the required libraries are imported. And the functions of different degree of polynomials are defined.
- Defining the function 'read_file' which is further called for reading the thermodynamic data file line by line using the 'for' loop with 'line' as a counter variable and storing the data i.e temperature and specific heat into two arrays 'temperature' and 'cp' respectively after splitting each line using the command 'line.split(',')' where line is the variable having string format containing the data from each line and comma(,).
- Then we define the function for goodness evaluation of fit. In which we calculate all 4 quantities and print them. And the file is read.
- Now the coefficients and covarience of coefficients are calculated using curve_fit function and curves are fitted for different degree of polynomials.
- Then the plots are plotted with corresponding labels, title and legend.
- At last the results are printed and the evaluation of fits is done.
popt:-
First output argument i.e 'popt' is a 2D array which consists the optimal (best fit) values of coefficients of a polynomial, used for curve fitting process, for which the sum of squared residuals for a given data is minimized. 'popt' can be used to extract the coefficients of a polynomial.
pcov :-
Second output argument i.e 'pcov' is also a 2D array which gives an estimated values of covariance of coefficients (popt) of a polynomial. The diagonal elements of an array/matrix gives the value of covariance of a coefficient with itself which is also called as variance of that coefficent. Other non-diagonal elements of the matrix gives the covariance of a coefficient with other coefficient. It is a square matrix (no. of rows = no. of columns). This array/matrix indicates the uncertainities and correlation between the coefficients of a polynomial and it is mostly useful when the data has uncertainities. 'pcov' can be used to find the standard deviation (square root of the value of diagonal elements) of a given dataset.
np.array(temperature):-
The command 'np.array(temperature)' is used for accessing the individual elements of the temperature array with the help of numpy module. Numpy module/library is imported as 'np' and used to perform mathematical operations like addition, subtraction, multiplication, division, etc. on each element of an array or between the arrays.
* in *popt mean:-
The star i.e '*' in '*popt' unpacks the 'popt' array so that the optimized coefficients (parameters) values can be accessed from the array 'popt'. In other words, '*' (* means all) aids in referring to all the values stored in the array 'popt'.
What needs to be done in order to make the curve fit perfect ?
To make a curve fit perfect the error estimates should be considered. The perfect fit is a fit where it should match the original curve and shouldn't have any scaling and centering errors in that particular degree of the polynomial. The R-square value must be equal to 1. This can be done by increasing the degree of polynomial.
Output:-
Figure 1
We can see that the cubic fit has less error compared to the linear fit.
Figure 2
Result
In the above figure the results of the goodness of fit for different degree of polynomials are tabulated. As the degree of polynomial increases the R-square value increses and reaches to value of 1. Hence, the fit is perfect fit.
Hence, we can conclude that the curve which is of 8th degree polynomial is a perfect fit as it has R-square value of 1 and we can also see in figure 2 that it exactly coincides with the actual data.