CODEWITHSUNDEEP
python
Avik Mukherjee
Avik Mukherjee

Reputation: 1

Python Data Fitting

The fit that I am getting

I am getting a horrible fit when I am trying to fit a parabola to this data. I am initially making a histogram of the data which is the position of an object and then plotting the negative log values of the histogram bin counts to the position using a parabola fit. the code I am using is this:

time,pos=postime()
plt.plot(time, pos)
poslen=len(pos)
plt.xlabel('Time')
plt.ylabel('Positions')
plt.show()


n,bins,patches = plt.hist(pos,bins=100)
n=n.tolist()
plt.show()
l=len(bins)   
s=len(n)
posx=[]
i=0
j=0
pbin=[]
sig=[]
while j < (l-1): 
    pbin.append((bins[j]+bins[j+1])/2)
    j=j+1

while i < s:

    if n[i]==0:
        pbin[i]=0
    else:
        sig.append(np.power(1/n[i],2))
        n[i]=n[i]/poslen

        n[i]=np.log(n[i])
        n[i]=n[i]*(-1)
    i=i+1


n[:]=[y for y in n if y != 0] 
pbin[:]=[y for y in pbin if y != 0]       

from scipy.optimize import curve_fit
def parabola(x, a , b):
    return a * (np.power(x,2)) + b

popt, pcov = curve_fit(parabola, pbin, n)
print popt

plt.plot(pbin,n)
plt.plot(pbin, parabola(pbin, *popt), 'r-')

Upvotes: 0

Views: 122

Answers (2)

Da Qi
Da Qi

Reputation: 635

You can also use matrix inversion.

import numpy as np
import matplotlib.pyplot as plt 
x = np.linspace(-5,5,100)
Y = (np.power(x,2) + np.random.normal(0,1,x.shape)).reshape(-1,1)

X = np.c_[np.ones(x.shape), x, np.power(x,2)]
A = np.linalg.inv(X.transpose().dot(X)).dot(X.transpose().dot(Y))
Yp = X.dot(A)

fig = plt.figure()
ax = fig.add_subplot()
plt.plot(x,Y,'o',alpha=0.5)
plt.plot(x,Yp)
plt.show()

The matrix form is

X*A=Y
A=(Xt*X)-1*Xt*Y

You can have a better idea here if needed. It does not always work out and you may want to apply some form of regularization.

Upvotes: 0

Gerges
Gerges

Reputation: 6509

I am not sure why you are computing the histogram... But here is a working example which does not require histogram computation.

import numpy as np
from scipy.optimize import curve_fit
from matplotlib import pyplot


time_ = np.arange(-5, 5, 0.1)
pos = time_**2 + np.random.rand(len(time_))*5

def parabola(x, a, b):
    return a * (np.power(x, 2)) + b

popt, pcov = curve_fit(parabola, time_, pos)
yfit = parabola(time_, *popt)

pyplot.plot(time_, pos, 'o')
pyplot.plot(time_, yfit)

enter image description here

Also, if your time_ vector is not uniformly sampled, and you want it to be uniformly sampled for the fit, you can do: fittime_ = np.linsapce(np.min(time_), np.max(time_)) and then yfit = parabola(fittime_, *popt).

Upvotes: 1

Related Questions