3

I have two lists that I am trying to do an exponential fit of form y=a*e^(bx) between. I am using an approach similar to the second answer from here but the results are not matching what I know to be true from testing with excel. Here is my code:

import numpy as np
from scipy.optimize import curve_fit
exp_constants = [62.5, 87.5, 112.5, 137.5, 162.5, 187.5, 212.5, 237.5, 262.5, 287.5]
means = [211.94, 139.30, 80.09, 48.29, 26.94, 12.12, 3.99, 1.02, 0.09, 0.02]
def func(x1, a, b):
  return a * np.exp(b * x1)
popt, pcov = curve_fit(func, exp_constants, means)

When returning popt[0] and popt[1] I get 3.222e-127 and 1.0 respectively. However, when checking with excel the correct exponential equation should be y=7231.3e^(-0.04x). I am not very familiar with the curve_fit approach, is there something that I am missing in my code or a better approach to getting the correct exponential fit?

Edit: Here is the plot that is made with the following code:

plt.figure()
plt.plot(exp_constants, means, 'ko', label="Data")
plt.plot(exp_constants, func(exp_constants, *popt), 'r-', label="Fitted Curve")
plt.legend()
plt.show

plot

Community
  • 1
  • 1
mbreezy
  • 99
  • 1
  • 5

1 Answers1

6

I guess the problem is that you do not provide an initial guess for the parameters, so as per the manual, curve_fit uses [1, 1] as a guess. The optimization might then get stuck at some local minimum. One other thing you should do is to change your xdata and ydata lists to numpy arrays, as shown by this answer:

import numpy as np
from scipy.optimize import curve_fit
exp_constants = np.array([62.5, 87.5, 112.5, 137.5, 162.5, 187.5, 212.5,
                          237.5, 262.5, 287.5])
means = np.array([211.94, 139.30, 80.09, 48.29, 26.94, 12.12, 3.99,
                  1.02, 0.09, 0.02])
def func(x1, a, b):
  return a * np.exp(b * x1)
guess = [100, -0.1]
popt, pcov = curve_fit(func, exp_constants, means, p0 = guess)

The exact value of the guess is not important, but you should probably have at least the order of magnitude and the signs right, so that the optimization can converge to the optimal value. I just used some random numbers close to the 'correct answer' you mentioned. When you don't know what to guess, you can do a polyfit(xdata, log(ydata), 1) and some basic math to get an initial value, as shown by this answer to the question you linked.

Quick plot:

x = np.linspace(exp_constants[0], exp_constants[-1], 1000)
plt.plot(exp_constants, means, 'ko', x, popt[0]*np.exp(popt[1]*x), 'r')
plt.show()

Result:

enter image description here

Community
  • 1
  • 1
Bas Swinckels
  • 18,095
  • 3
  • 45
  • 62