Search code examples
pythoncurve-fittingscipy-optimize

scipy.curve_fit() does not work for an exponential


I have four points: sigma2 = [15, 8, 5.5, 3.21], obtained for points [100, 300, 500, 1000]. I want to fit an exponential curve to this, and then measure the value for 800.

However, when I try to do so, I simply get a horizontal line.

Here is my code:


sigma2 = [15, 8, 5.5, 3.21]

def exp_fit(x, a, b, c):
    return a * np.exp(-b * x) + c

x_data = np.asarray([100, 300, 500, 1000], dtype=np.float64)
y_data = np.asarray(sigma2, dtype=np.float64)

param_exp, _ = curve_fit(exp_fit, x_data, y_data)

x_test = np.arange(100, 1100, 100)
print(exp_fit(x_test, *param_exp))
plt.plot(x_test, exp_fit(x_test, *param_exp), 'r-',
         label='fit: a=%5.3f, b=%5.3f, c=%5.3f' % tuple(param_exp))

Any help would be greatly appreciated. I realize it is a simple thing, but I can't quite put my finger on the error I am making.

I did try to negate the initial value of x as suggested in other posts, but it still does not work. I tried to bound the exponential as follows:

    param_exp, _ = curve_fit(exp_fit, x_data, y_data, bounds=(0, [0.001, 0.001, 0.001]))

and now I see a straight line with a negative slope.


Solution

  • SciPy's curve fitting methods rely on gradient information - in other words, if adjusting a parameter doesn't seem to do anything, they won't use it. For that reason, it's important to provide good initial guesses for the optimizers.

    The default guess provided for b here is 1, which results in the exponential term being near-zero.

    >>> b = 1
    >>> np.exp(-b * x_data)
    array([3.72007598e-044, 5.14820022e-131, 7.12457641e-218, 0.00000000e+000])
    

    The following initial guess seems to work well:

    c0 = a0 = 1
    b0 = np.log(np.mean(y_data)) / np.mean(x_data)
    p0 = a0, b0, c0
    param_exp, _ = curve_fit(exp_fit, x_data, y_data, p0=p0)