Week 5 - Curve fitting

Aim

The aim of this program is to fit a linear and cubic polynomial for the given data between heat capacity and temperature.

Program

The linear polynomial and cubic polynomial are defined in 'func_1' and 'func_3'. The data file is read in a function called 'readfile' which returns the values of temp and heat capacity as arrays. The file is read using a for loop and the values are obtained using 'split()' command with comma as the delimiter.

The 'curvefit' command is imported from 'scipy.optimize' and is used to obtain the coefficients for the curve fit. New values of heat capacity are calculated using the existing temperature array and the coefficients for the curve fit. A function to calculate the R-squared value is defined. An R-squared value above 0.95 deems the fit to be good.The original data, linear curve fit and cubic curve fit are plotted together in one figure.

1. 'popt' stores the coefficents of the respective model function and 'pcov' is 2D array storing teh estimated covariance of 'popt'
2. np.array(temperature) coverts the list 'temperature' into a numpy array which takes much less space and its performance is better. Also it can be operated on directly.
3. The '*' in '*popt' implies popt is an iterable and elements from 'popt' are considered as additional arguments.
4. The linear is fit is not good as its R-squared value is 0.92 and from the plot it is visible that it doesnt explain the original data. The cubic fit is good with an R-squared value of 0.99 and it explains much of teh original data.
5. To make the fit perfect higher order polynomials can be considered or various other curves are available for curve fitting instead of just polynomials.

Code

import math
import numpy as np
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
#linear function
def func_1(t,a,b):
	return a*t+b
#cubic function
def func_3(t,a,b,c,d):
	return a*pow(t,3)+b*pow(t,2)+c*t+d
#function to read file
def readfile():
	temp=[]
	cp=[]
	for line in open("data","r"):
		values=line.split(",")
		temp.append(float(values[0]))
		cp.append(float(values[1]))
	return [temp,cp]

[temp,cp]=readfile()
[popt1,pcov1]=curve_fit(func_1,temp,cp) #curve fit linear
[popt2,pcov2]=curve_fit(func_3,temp,cp) #curve fit cubic
new_cp1=func_1(np.array(temp),*popt1)
new_cp2=func_3(np.array(temp),*popt2)
#function to caluclate r square value 
def r_squared(cp,new_cp):
	residuals=cp-new_cp
	ss_res=np.sum(residuals**2)
	ss_tot=np.sum((cp-np.mean(cp))**2)
	r_square=1-(ss_res/ss_tot)
	return r_square
r1=r_squared(cp,new_cp1)
r2=r_squared(cp,new_cp2)
print(f"r squared value for linear fit={r1}")
print(f"r squared value for cubic fit={r2}")
#plotting
plt.plot(temp,cp,color="red")
plt.plot(temp,new_cp1,color="blue")
plt.plot(temp,new_cp2,color="green")
plt.xlabel("Temperature")
plt.ylabel("Cp")
plt.legend(["original data","curve fit linear","curve fit cubic"])
plt.show()

Results