Curve fitting (Python)

1. AIM : To write codes in Python to perform curve fitting.

2. OBJECTIVE :

1. To wite codes to fit a linear and cubic polynomial for the Cp data.
2. To plot the linear and cubic fit curves along with the raw data points.
3. To measure the fitness characteristics for both the curves.

The link for the Cp vs Temperature data file is given below-

https://drive.google.com/file/d/1vWi9Vj3AdQf0gqUK1u9jwS8zg92BHGR0/view?usp=sharing

3. THEORY :

Curve fitting: - Curve fitting is the process of constructing a curve, or mathematical function, that has the best fit to a series of data points, possibly subject to constraints. Curve fitting can involve either interpolation, where an exact fit to the data is required, or smoothing, in which a "smooth" function is constructed that approximately fits the data. A related topic is regression analysis, which focuses more on questions of statistical interference such as how much uncertainty is present in a curve that is fit to data observed with random errors. Fitted curves can be used as an aid for data visualization, to infer values of a function where no data are available, and to summarize the relationships among two or more variables. Extrapolation refers to the use of a fitted curve beyond the range of the observed data, and is subject to a degree of uncertainty since it may reflect the method used to construct the curve as much as it reflects the observed data.

Curve fitting is the way we model or represent a data spread by assigning the best fit function (curve) along with 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.

Types of curve fitting include:

Interpolation, where you discover a function that is an exact fit to the data points. Since this assumes no measurement error, it has limited applicability to real-life scenarios.

Smoothing is when we find a function that is an approximate fit to the data points, but we give room for error and we allow our actual points to be near, but not necessarily on the line; given the error is minimized overall.

Linear and Polynomial Curve fitting :

(i) Linear curve fitting, or linear regression, is when the data is fit to a straight line. Although there might be some curve to your data, a straight line provides a reasonable enough fit to make predictions.

Since the equation of a generic straight line is always given by f(x)= a x + b, the question becomes: what a and b will give us the best fit line for our data?

Considering the vertical distance from each point to a prospective line as an error, and summing them up over our range, gives us a concrete number that expresses how far from ‘best’ the prospective line is.

A line that provides a minimum error can be considered the best straight line.

Since it’s the distance from our points to the line we’re interested in—whether it is positive or negative distance is not relevant—we square the distance in our error calculations. This also allows us to weight greater errors more heavily. So this method is called the least square approach.

(ii) Polynomial curve fitting is when we fit our data to the graph of a polynomial function. The same least-squares method can be used to find the polynomial, of a given degree, that has a minimum total error.

To choose the best fit for the curve the following four parameters help us to measure the goodness of fit criteria or how well the equations are representing the given datapoints:

1. The sum of squares due to error (SSE),
2. R-square (R2),
3. Adjusted R-squrae,
4. Root mean squared error (RMSE)

1.) The sum of square due to error (SSE): This statistic measures the total deviation of the responses values from the fit to the response values. It is also called the summed square of residuals and is usually labelled as SSE.

Here, TSS = Total sum of squares = `∑(y_i-y_(mean) )^2`

RSS = Residual sum of squares = `∑(y_i-y_(fit) )^2`

ESS = Explained sum of squares = `∑(y_(fit)-y_(mean) )^2`

2.) R-square:

R-square = `"RSS" /"TSS"` or `1-("ESS" /"TSS" )`

3.) Adjusted R-square: Adjusted R-square is a normalizing variation of R-square which neglects the variation of R-square due to insignificant independent variables.

`("R"_(adjusted))^2 =1- ((1-"R" ^2 )("N" -1))/("N" -"p" -"1" )`

Where,

R² = sample R-square

P = Number of predictors

N = Total sample size.

Where,

n = number of data points

k = number of independent variables

4.) Root mean squared error (RMSE): This is also a statistical way of measuring how well the curve fit is. As the name suggests it has something to do with the Squaring of the errors, then summing them and then taking the average and finally take the root.

It is a very similar measure as to the SSE (Sum of the Square of the Error) which was discussed earlier in the report. Mathematically,

`RSME = √((∑_(i=1)^N[y(i) -f(x(i))]^2 )/N)`

Therefore we can see that:

`RSME = √((SSE)/N)` . As it is evident that the smaller the RSME, the better the curve fit.

4.BODY OF CONTENT :

Solution steps :-

1. At first, necessary modules such as math,matplotlib,NumPy,scipy are imported to perform their respective functions.
2. Two functions are defined which return linear and cubic equations.
3. The following data is loaded into Python by placing it in the right directory and a function read_file is created. For loop is used to run each line by using 'open' command with 'r' to read it.
4. From the data obtained, two individual column matrices of temperature and specific heat are derived.
5. For curve fitting, 'popt' and 'pcov' are used as output for storing parameters obtained by 'curve_fit' command which takes input as the required function(linear or cubic) , temperature, and Cp. Now, another variable named 'fit_cp' stores the new values of Cp calculated by using temperature array.
6. Plots of both raw dataset and predicted dataset are made and compared to find fitness characteristics for both the curves for both linear and cubic polynomial.
7. To measure the fitness characteristics for both the curves, various parameters such as SSE,SSR,SST,R2,RMSE are calculated. The plot which has highest R2 or closest to 1 is considered to be a good fit.

Code :

# A program to measure fitness characteristics for the linear and cubic polynomial for the Cp data:
import math
import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit

# Curve fit function
# Linear function
def function_0(t,a,b):
	return a*t + b
# Cubic function
def function_2(t,a,b,c,d):
	return a*pow(t,3) + b*pow(t,2) + c*t +d


# 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]

# Main Program
temperature , cp = read_file()

popt, pcov = curve_fit(function_0,temperature,cp)
fit_cp_l = function_0(np.array(temperature), *popt)
plt.figure(1)
plt.plot(temperature,cp,'k--')
plt.plot(temperature,fit_cp_l,color='red',linewidth = 1)
plt.legend(['Actual data','Curve fit'])
plt.xlabel('Temperature (K)')
plt.ylabel('Cp')
plt.title('Linear Curve fitting')
plt.show()


popt, pcov = curve_fit(function_2,temperature,cp)
fit_cp_c = function_2(np.array(temperature), *popt)
plt.figure(3)
plt.plot(temperature,cp,'k--')
plt.plot(temperature,fit_cp_c,color='red',linewidth = 1)
plt.legend(['Actual data','Curve fit'])
plt.xlabel('Temperature (K)')
plt.ylabel('Cp')
plt.title('Cubic Curve fitting')
plt.show()


"""
Measuring the fitness characetristics of the curves

	SSE (sum of error squared)
	SSR (sum of squares of the regression)
	SST = SSE + SSR
	R^2 = SSR/SST
	RMSE = (SSE/l)^0.5

	Finding mean of Cp data
	mean = ( sum of all elements)/(total number of elements)
"""

s = np.sum(cp)
l = np.size(cp)
m = s/l
print('Mean of all cp :',m)
print('')

# For Linear curve fit
linear_sse = 0
linear_ssr = 0
for i in range(l):
	linear_error = abs((np.sum((cp[i] - fit_cp_l[i]))))
	linear_sse = linear_sse+pow(linear_error,2)
	linear_ssr = linear_ssr+ np.sum(pow((fit_cp_l[i] - m),2))


linear_sst = linear_sse + linear_ssr
print('linear_sst :',linear_sst)
linear_R2 = linear_ssr/linear_sstprint('linear_R2 :',linear_R2)
linear_rmse = pow((linear_sse/l),0.5)
print('linear_RMSE :',linear_rmse)
print('')


# For Cubic curve fit
cubic_sse = 0
cubic_ssr = 0
for j in range(l):
	cubic_error = abs(np.sum((cp[i] - fit_cp_c[i])))
	cubic_sse = cubic_sse+(pow(cubic_error,2)) 
	cubic_ssr = cubic_ssr+np.sum(pow((fit_cp_c[i] - m),2))


cubic_sst = cubic_sse + cubic_ssr
print('cubic_sst :',cubic_sst)
cubic_R2= cubic_ssr/cubic_sst
print('cubic_R2 :',cubic_R2)
cubic_rmse = pow((cubic_sse/l),0.5)
print('cubic_RMSE :',cubic_rmse)

5. RESULTS :

Outputs of the above program are as follows -

For linear Curve fitting, the plot obtained as,

For Cubic Curve fitting, the plot obtained as,

6. Definitions and Functionalities:

Modules used : matplot.lib, scipy, numpy, imageio.

Functions : -

6.1 curve_fir()

The curve_fit is a function in the scipy.optimize (Optimisation and Root finding) library of scipy module. It is essentially a non-linear least square fit tool.

pcov and pcov: The curve_fit function returns two arrays. The first one, popt is a one-dimensional array that contains the coefficients of the determined function from curve_fit. The second array, pcov, is a 2-dimensional array, or a matrix, that contains the estimated covariance of the popt.

Covariance is the variance between the two variables. Thus, popt 2D array contains the covariances of each value (coefficient) in the popt list with the other values in the list that the curve_fit algorithm calculates. The scipy.optimize calculates these with the least squares. Now, the diagonal of this matrix gives the covariance of a coefficient with itself, which is simply the variance of the coefficient. And the (i,j) and (j,i) values are the same and give the covariance between the two (i,j)th or (j,i)th co-efficients.

6.2 np.array()

This is a method or function from the NumPy library which creates an array with the required arguments. In this assignment, we need to pass the temperature list as an array.

In the above code, if 'func' is the function to which we want the curve to be fitted, then 'fit' is assigned to calculate the value of the function at a particular value of temp with the required coefficients (in this case p1). However, we want the values at all 'temp' values from the function. So, we cannot simply pass poly_function(temp,*p1).

6.3 ployfit()

The polyfit() function from the NumPy module is another curve fitting tool which is essentially a least squares polynomial fit. Unlike the curve_fit() function, the polyfit() function doesn't require the definition of the function of the curve we wish to fit, but, it simply asks for the degree of the polynomial to which we want to fit the data to. Thus, curve_fit() can be applied when we have a definition of the function at hand (which may be exponential, logarithmic, linear or even a polynomial) and polyfit() is exclusively used to fit some polynomial to the data.