0

I am trying to get a simple fit to my data of an exponential decay of the form a*(x-x0)**b, where I know a and b must be negative. So if you plot it on a log-log plot, I should see a linear trend for the obtained data.

As such, I'm giving scipy.optimize initial guesses where a and b are negative, but it keeps ignoring them and giving me the error,

OptimizeWarning: Covariance of the parameters could not be estimated

.. and giving me values for a and b that are positive. It then also does not give me an exponential decay, but a parabola that bottoms out and begins to increase.

I have tried many different guesses as to the initial parameters over a large range of values (one such is in the code below), but none worked without giving me the nonsensical return and the error. This has made me start to wonder if my code is wrong, or if there's just some obvious way to get good initial guesses into the code that won't be rejected.

import math
import numpy as np
import sys
import matplotlib.pyplot as plt
import scipy as sp
import scipy.optimize
from scipy.optimize import curve_fit
import numpy.polynomial.polynomial as poly

x= [1987, 1993.85, 2003, 2010.45, 2009.3, 2019.4]

t= [31, 8.6, 4.84, 1.96, 3.9, 1.875]

def model_func(x, a, b, x0):
    return (a*(x-x0)**b)

# curve fit
p0 = (-.0005,-.0005,100)
opt, pcov = curve_fit(model_func, x, t,p0)
a, b, x0 = opt

# test result
x2 = np.linspace(1980, 2020, 100)
y2 = model_func(x2, a, b,x0)

coefs, cov = poly.polyfit(x, t, 2,full=True)
ffit = poly.polyval(x2, coefs)

plt.loglog(x,t,'.')
plt.loglog(x2, ffit,'--', color="#1f77b4")

print('S = (',coefs[0],'*(t-',coefs[2],')^',coefs[1])
kb3hts
  • 165
  • 1
  • 4
  • 9
  • 1
    What is `poly` in your code? Can you add all the required imports so that it runs without being edited? – Warren Weckesser Jul 23 '19 at 22:41
  • I sometimes do it using `scipy`'s `differential_evolution` and the later fine-tune. You can see an example in [this post (last example)](https://stackoverflow.com/a/52484882/1534017). – Cleb Jul 24 '19 at 09:59
  • Sorry, I forgot the imports! I'm adding them into the original post now, at the top. – kb3hts Jul 24 '19 at 15:18

0 Answers0