TAGS :Viewed: 16 - Published at: a few seconds ago

[ 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:

Data with erf fitting plot

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.

enter image description here


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')