Python program to perform curve fitting for a Temperature v/s Cp data
AIM : Python program to perform Curve Fitting using given input data.
The program should achieve the following objectives:
- Coding to fit a linear and cubic polynomial for the given Temperature v/s Cp data.
- Plotting the linear and cubic fit curves along with the raw data points.
BACKGROUND : Curve fitting is the way we model or represent a data spread by assigning a ‘best fit‘ function (curve) along the entire range. Ideally, it will capture the trend in the data and allow us to make predictions of how the data series will behave in the future.
It is process of manipulating the free variables to try and make the equations reproduce the known curve as closely as possible.
If you know the underlying equations are correct, then this may be quite a reasonable way of identifying the value of the free variables. If the underlying equations are not correct, it will produce a solution which apparently fits the data quite well, but is completely bogus. For this reason, it is only valid in situations where you know the equations are correct. Often, but not always, this is the case.
Because "curve fitting" can make totally incorrect sets of equations superficially look correct, it is often used as a pejorative term indicating the underlying equations are wrong. Incorrect use of curve fitting can be detected by checking the validity of solutions which are outside the range of data used to curve fit - the solutions don't extrapolate properly.
Now, with respect to Curve Fitting, the following questions are generic with respect to the code written below and need to be answered while performing the curve fitting analysis.
- How to make a curve fit perfectly?
Upon finding the curve which satisfies the data points there are 4 quantities which help us measure the goodness of fit criteria or how well the equations are representing the data points. These are –
- Sum of Squares due to Error(SSE) : This statistic measures the total deviation of the response values from the fit to the response values. It is also called the summed square of residuals and is usually labelled as SSE.
`SSE= sum_(i=1)^n{w_i(y_i=f_i)^2}`
Here yi is the observed data value and fi is the predicted value from the fit. wi is the weighting applied to each data point, usually wi = 1. A value closer to 0 indicates that the model has a smaller random error component, and that the fit will be more useful for prediction.
- R-Square : This statistic measures how successful the fit is in explaining the variation of the data. Put another way, 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.
R-square is defined as
`R- Square= 1-{[sum_(i=1)^n{w_i(y_i-f_i)^2}]}/{[sum_(i=1)^n{w_i(y_i-y_a)^2}}`
`R-Square= 1-frac{SSE}{SST}`
Here fi is the predicted value from the fit, ya is the mean of the observed data yi is the observed data value. wi is the weighting applied to each data point, usually wi=1. SSE is the sum of squares due to error and SST is the total sum of squares. R-square can take on any value between 0 and 1, with a value closer to 1 indicating that a greater proportion of variance is accounted for by the model. For example, an R-square value of 0.8234 means that the fit explains 82.34% of the total variation in the data about the average.
If you increase the number of fitted coefficients in your model, R-square will increase although the fit may not improve in a practical sense. To avoid this situation, you should use the degrees of freedom adjusted R-square statistic described below.
Note that it is possible to get a negative R-square for equations that do not contain a constant term. Because R-square is defined as the proportion of variance explained by the fit, if the fit is actually worse than just fitting a horizontal line then R-square is negative. In this case, R-square cannot be interpreted as the square of a correlation. Such situations indicate that a constant term should be added to the model.
- Degree of Freedom Adjusted R-square : This statistic uses the R-square statistic defined above, and adjusts it based on the residual degrees of freedom. The residual degrees of freedom is defined as the number of response values n minus the number of fit coefficients m estimated from the response values.
v = n-m
v indicates the number of independent pieces of information involving the n data points that are required to calculate the sum of squares. Note that if parameters are bounded and one or more of the estimates are at their bounds, then those estimates are regarded as fixed. The degrees of freedom is increased by the number of such parameters.
The adjusted R-square statistic is generally the best indicator of the fit quality when you compare two models that are nested – that is, a series of models each of which adds additional coefficients to the previous model.
`adjusted R-Square=1-frac{SSE(n-1)}{SST(v)}`
The adjusted R-square statistic can take on any value less than or equal to 1, with a value closer to 1 indicating a better fit. Negative values can occur when the model contains terms that do not help to predict the response.
- Root Mean Squared Error : This statistic is also 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, and is defined as
`RMSE=s=(MSE)^frac{1}{2}`
where MSE is the mean square error or the residual mean square
`MSE=frac{SSE}{v}`
Just as with SSE, an MSE value closer to 0 indicates a fit that is more useful for prediction.
2. What does popt and pcov mean?
(a) popt is the array that refers to the optimal values of the coefficients so that the sum of the squared residuals is minimised.
(b) pcov is the estimated covariance(measure of how much two random variables vary together) of popt or the output coefficients.
3. What does np.array(temperature) do?
np.array(temperature) will pass the temperature as numpy array to the functions where the values of the output coefficients is to be calculated.
3. What does the * in *popt mean?
The * in *popt helps to refer the values of the coefficients. Instead of passing the coefficients in the function, we are passing popt which already contains the optimised values of all the coefficients.
SOLUTION ALGORITHM
- We will begin by reading the data for Temperature v/s Specific Heat from an external text file.
- The temperature and the specific heat will be stored in separate variables and a curve will be plotted with them which will be actual curve.
- Now, first we will assume a linear function to represent the data and with the help of curve_fit command we will generate the coefficients of the linear function.
- Then we will also check the fitness characteristics of the curve by using fit_cp.
- After this, we will compare the curve generated with the predicted values with the actual curve to determine whether the curve fits perfectly.
- Steps from 2 to 5 will then be repeated for an assumed cubic function.
- Finally we will compare the results.
- The code for the entire program will be written using Sublime Text Editor.
PYTHON CODE :
"""
Python program to perform curve fitting for a given data set
Script written by Abhishek Arora
"""
import math
import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit
#defining the function to read data from the file
def read_data():
temperature =[]
cp= []
#reading the data line after line
for line in open('data','r'):
#using the split command to break the line where a comma is encountered.
#Python will take the line,split it where comma occurs and it will give the array containing 2 entries
values= line.split(',')
temperature.append(float(values[0]))
cp.append(float(values[1]))
return temperature,cp
def func_lin(t,a,b):
return a*t + b
def func_cubic(t,a,b,c,d):
return a*pow(t,3) + b*pow(t,2) + c*t + d
#main program
#reading the temperature and cp from the external file
temperature,cp=read_data()
#performing curve fitting
#for a linear function
popt,pcov= curve_fit(func_lin,temperature,cp)
#checking the fitness of the curve
fit_cp= func_lin(np.array(temperature),*popt)
#plotting the results
plt.figure(1)
plt.plot(temperature,cp,color="blue")
plt.plot(temperature,fit_cp,color="red")
plt.legend(['actual data', 'curve fit linear'])
plt.xlabel('temperature(K)')
plt.ylabel('Cp')
plt.show()
#performing curve fitting
#for a cubic function
popt,pcov= curve_fit(func_cubic,temperature,cp)
#checking the fitness of the curve
fit_cp= func_cubic(np.array(temperature),*popt)
#plotting the results
plt.figure(2)
plt.plot(temperature,cp,color="blue")
plt.plot(temperature,fit_cp,color="red")
plt.legend(['actual data', 'curve fit cubic'])
plt.xlabel('temperature(K)')
plt.ylabel('Cp')
plt.show()OUTPUT
The following output is obtained after running the code :
*visually it is evident that the curve for a cubic polynomial fits more perfectly than a linear one.
CONCLUSION
- From the output, it is clear that the curve for cubic function fits more perfectly as compared to the linear function.
- It hence can be generalized for the given data set that as we increase the order of the function representing the given data set, the accuracy is enhanced. However, it does not necessarily imply for all the other data sets.
- The accuracy of the cubic function can be further enhanced by using several other techniques such as performing split curve fitting or performing regresssion analysis. Regression analysis is generally performed to find the best fit for a given set of data points. Regression is the process of analyzing the relationship between different variables and building a mathematical model that can be used to predict the value of a variable(label) based on the values of remaining variables(features).
REFERENCES