python SciPy curve_fit with np.exp returns with pcov = inf

Question

I'm trying to optimize a exponential fitting with scipy.optimize.curve_fit. But the result is no good . My code is :

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

# xdata and data is obtain from another dataframe and their type is nparray

xdata =[36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62,63,64,65,66,67,68,69,70 ,71,72]
ydata = [4,4,4,6,6,13,22,22,26,28,38,48,55,65,65,92,112,134,171,210,267,307,353,436,669,669,818,1029,1219,1405,1617,1791,2032,2032,2182,2298,2389]

popt, pcov = curve_fit(func, xdata, ydata)
plt.plot(xdata, func(xdata, *popt), 'r-', label='fit: a=%5.3f, b=%5.3f, c=%5.3f' % tuple(popt))

plt.scatter(xdata, ydata, s=1)
plt.show()

Then I got the result like this:

enter image description here

the result showed that :

pcov = [[inf inf inf] [inf inf inf] [inf inf inf]]
popt = [1  1  611.83784]

I don't know how to make my curve fit well. Can you helo me? Thank you!

Answer 1

Fitting against exponential functions is exceedingly tough because tiny variations in the exponent can make large differences in the result. The optimizer is optimizing across many orders of magnitude, and errors near the origin are not equally weighted compared to errors higher up the curve.

The simplest way to handle this is to convert your exponential data to a line using a transformation:

y' = np.log(y)

Then instead of needing to use the fancier (and slower) curve_fit, you can simply use numpy's polyfit function and fit a line. If you wish, you can transform the data back into linear space for analysis. Here, I've edited your code to do the fit with np.polyfit, and you can see the fit is sensible.

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

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

# xdata and data is obtain from another dataframe and their type is nparray

xdata = np.array([36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62,63,64,65,66,67,68,69,70 ,71,72])
ydata = np.array([4,4,4,6,6,13,22,22,26,28,38,48,55,65,65,92,112,134,171,210,267,307,353,436,669,669,818,1029,1219,1405,1617,1791,2032,2032,2182,2298,2389])

# popt, pcov = curve_fit(func, xdata, ydata)
# plt.plot(xdata, func(xdata, *popt), 'r-', label='fit: a=%5.3f, b=%5.3f, c=%5.3f' % tuple(popt))

# Fit a line (deg=1)
P, pcov = np.polyfit(xdata, np.log(ydata), deg=1, cov=True)
print(pcov)
plt.scatter(xdata, ydata, s=1)
plt.plot(xdata, np.exp(P[0]*xdata + P[1]), 'r-')

plt.legend()
plt.show()

Answer 2

The method is not finding the optimal point. One thing to try is changing the initial guess so that b starts negative, because it looks from your data that b must be negative so that the func fits it decently. Also, from the docs of curve_fit, the initial guess is 1 by default if not specified. A good initial guess is:

popt, pcov = curve_fit(func, xdata, ydata, p0=[1, -0.05, 1])

which gives

popt                                                                                                                                                                                                      
array([ 1.90782987e+00, -1.01639857e-01, -1.73633728e+02])

pcov                                                                                                                                                                                                           
array([[ 1.08960274e+00,  7.93580944e-03, -5.24526701e+01],
       [ 7.93580944e-03,  5.79450721e-05, -3.74693994e-01],
       [-5.24526701e+01, -3.74693994e-01,  3.34388178e+03]])

And the plot