[ Fit erf function to data ]
So I have hysteresis loop. I want to use the erf function to fit it with my data.
A portion of my loop is shown in black on the lower graph.
I am trying to use the scipy.optimize.curve_fit
and scipy.special.erf
function to fit the data with the following code:
import scipy.special
import scipy.optimize
def erfunc(x,a,b):
return mFL*scipy.special.erf((x-a)/(b*np.sqrt(2)))
params,extras = scipy.optimize.curve_fit(erfunc,x_data,y_data)
x_erf = list(range(-3000,3000,1))
y_erf = erfunc(x_erf,params[0],params[1])
mFL
is a constant, a
controls the position of the erf curve and b
the slope of the curve. (To my knowledge)
However, when I plot the obtained x_erf and y_erf data (in blue). I get the following fitting, which is not ideal to say the least:
Is there a way i can get a proper fit?
Edit: Link to data file: https://www.dropbox.com/s/o0uoieg3jkliun7/xydata.csv?dl=0 Params[0] = 1.83289895, Params1 = 0.27837306
Answer 1
I suspect two things are needed for a good fit here. First, I believe you need to add mFL
to your erfunc
function, and second, as suggested by Glostas, you need to specify some initial guesses for your fitting parameters. I created some artificial data in an attempt to replicate your data. The plot on the left is before giving curve_fit
some initial parameters and the plot on the right is after.
Here is the code to reproduce the above plots
import numpy as np
from scipy.special import erf
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
def erfunc(x, mFL, a, b):
return mFL*erf((x-a)/(b*np.sqrt(2)))
x_data = np.linspace(-3000, 3000, 100)
mFL, a, b = 0.0003, 500, 100
y_data = erfunc(x_data, mFL, a, b)
y_noise = np.random.rand(y_data.size) / 1e4
y_noisy_data = y_data + y_noise
params, extras = curve_fit(erfunc, x_data, y_noisy_data)
# supply initial guesses to curve_fit through p0 arg
superior_params, extras = curve_fit(erfunc, x_data, y_noisy_data,
p0=[0.001, 100, 100])
fig = plt.figure()
ax1 = fig.add_subplot(121)
ax2 = fig.add_subplot(122)
ax1.plot(x_data, erfunc(x_data, *params))
ax1.plot(x_data, y_noisy_data, 'k')
ax1.set_title('Before Guesses')
ax2.plot(x_data, erfunc(x_data, *superior_params))
ax2.plot(x_data, y_noisy_data, 'k')
ax2.set_title('After Guesses')