scipy curve_fit strange result matplotlib

Viewed 36

inbefore my code:

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

xFit= np.arange(0.0, 20, 0.01)
dev_Fluenz1= np.array([68.9, 21.81, 9.38, 3.73])
dev_Fluenz2= np.array([137.68 , 42.34, 18.75, 7.47 ])
dev_Fluenz3= np.array([80.34, 23.82 , 10.06, 3.76 ])
dev_Fluenz4= np.array([203.7, 61.67 , 10.06, 10.33 ])
dev_Fluenz5= np.array([135.74, 46.23 , 19.42 , 11.21 ])


dev_Fluenz6= np.array([382.83, 112.95, 50.02, 14.95])
dev_Fluenz7= np.array([382.45 , 117.62 , 50.01, 14.95 ])
dev_Fluenz8= np.array([147.32 , 43.67 , 17.88 , 5.01 ])
dev_Fluenz9= np.array([282.91, 85.64 , 35.63 , 13.77 ])
dev_Fluenz10= np.array([150.82, 51.37 , 20.5 , 11.21 ])


dev_x2= np.array([2,5,10,20])

plt.plot(dev_x2,dev_Fluenz1, 'bo')
popt, pcov = curve_fit(func, dev_x2, dev_Fluenz1)
plt.plot(xFit, func(xFit,*popt),color='b', linestyle='--',label=f'Reales DOE Sa <= 0,3 Fluenz 1 J/cm**2 \u03bcm F(x) = {round(popt[0])} * e^({round(popt[1])}*x)')

results in : enter image description here

changing xFit= np.arange(0.0, 20, 0.01) to xFit= np.arange(0.0, 200, 0.01)

and

`dev_x2= np.array([2,5,10,20])` to `dev_x2= np.array([20,50,100,200])`

results in enter image description here

Why?

Thank you for your help!

1 Answers

curve_fit takes an initial guess as a starting point of the fit, see https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html.

If you do not provide it, as you did not, the default is used:

p0array_like, optional Initial guess for the parameters (length N). If None, then the initial values will all be 1 (if the number of parameters for the function can be determined using introspection, otherwise a ValueError is raised).

For your second test case that blows up, it matters because with the default initial value of the exponent b in func, for the given range of dev_x2`, you get very large values and the optimizer cannot find a good solution

So you need to provide a sensible starting value. For example, note the p0 = [... bit below

xFit= np.arange(0.0, 200, 0.01)
dev_x2= np.array([20,50,100,200])

plt.plot(dev_x2,dev_Fluenz1, 'bo')
popt, pcov = curve_fit(func, dev_x2, dev_Fluenz1, p0 = [100,  -0.01])
plt.plot(xFit, func(xFit,*popt),color='b', linestyle='--',label=f'Reales DOE Sa <= 0,3 Fluenz 1 J/cm**2 \u03bcm F(x) = {round(popt[0])} * e^({round(popt[1])}*x)')

you get a sensible result: enter image description here

Related