Using Scipy curve_fit with a piece function

Question

I get an optimization warning:

OptimizeWarning: Covariance of the parameters could not be estimated category=OptimizeWarning) 

when trying to match my piecewise function with my data using scipy.optimize.curve_fit . This means that no fit occurs. I can easily substitute a parabola for my data, and I am supplying curve_fit with what, in my opinion, are good initial parameters. Full sample code below. Does anyone know why curve_fit might not mix with np.piecewise ? Or am I making another mistake?

 import numpy as np from scipy.optimize import curve_fit import matplotlib.pyplot as plt def piecewise_linear(x, x0, y0, k1, k2): y = np.piecewise(x, [x < x0, x >= x0], [lambda x:k1*x + y0-k1*x0, lambda x:k2*x + y0-k2*x0]) return y def parabola(x, a, b): y = a * x**2 + b return y x = np.array([-3, -2, -1, 0, 1, 2, 3]) y = np.array([9.15, 5.68, 2.32, 0.00, 2.05, 5.29, 8.62]) popt_piecewise, pcov = curve_fit(piecewise_linear, x, y, p0=[0.1, 0.1, -5, 5]) popt_parabola, pcov = curve_fit(parabola, x, y, p0=[1, 1]) new_x = np.linspace(x.min(), x.max(), 61) fig, ax = plt.subplots() ax.plot(x, y, 'o', ls='') ax.plot(new_x, piecewise_linear(new_x, *popt_piecewise)) ax.plot(new_x, parabola(new_x, *popt_parabola)) ax.set_xlim(-4, 4) ax.set_ylim(-2, 16) 

Answer 1

This is a type problem, you need to change the following line, so x is set as a float:

 x = np.array([-3, -2, -1, 0, 1, 2, 3]).astype(np.float) 

otherwise, piecewise_linear may result in type casting.

To be safe, you could also create starting points:

 popt_piecewise, pcov = curve_fit(piecewise_linear, x, y, p0=[0.1, 0.1, -5., 5.]) 

Answer 2

For completeness, I will point out that fitting a piecewise linear function does not require np.piecewise : any such function can be constructed from absolute values ​​using a multiple of np.abs(x-x0) for each bend. The following gives a good fit to the data:

 def pl(x, x0, a, b, c): y = a*np.abs(x-x0) + b*x + c return y popt_pl, pcov = curve_fit(pl, x, y, p0=[0, 0, 0, 0]) print(pl(x, *popt_pl)) 

The output is close to the original y values:

 [ 8.90899998 5.828 2.74700002 -0.33399996 2.03499998 5.32 8.60500002]