Quadratic to Chaining or Connected Multiple Piecewise in Python

Question

I've been digging in stackoverflow for a while and can't find any example for multiple piecewise curve fitting. I want to convert a quadratic function into multiple chaining (I don't know the exact name of it, but i need every tail connected to the head of the next piecewise, simply "connected") of piecewise function. This is my code so far using scipy.optimize to convert quadratic into 2 pieces of piecewise linear function.

    import scipy.optimize as opt
    import numpy as np
    import copy

    def func_2piecewise(x, m_0, x_1, y_1, m_1):
        y = np.piecewise(x, [x <= x_1, x > x_1],
                            [lambda x:m_0*(x-x_1) + y_1, lambda x:m_1*(x-x_1) + y_1])
        return y

    xmin=0
    xmax=100
    a=0.1
    a0=1
    a00=10
    piece_number=2
    sigma=np.ones(numberOfStep)
    if piece_number==2:
        lower_bounds=[-np.inf,xmin,-np.inf,-np.inf]
        upper_bounds=[np.inf,xmax,np.inf,np.inf]

        w, _ = opt.curve_fit(func_2piecewise, x_sample, y_sample,bounds=(lower_bounds,upper_bounds),sigma=sigma)
        x_0=copy.deepcopy(xmin)
        y_0=func_2piecewise(x_0, *w).tolist()
        [m_0, x_1, y_1, m_1]=w
        result=[x_0,y_0,m_0,x_1,y_1,m_1]

The problem is, I can't implement the same approach for three piecewise (i don't know how to make x_2 > x_1):

    def func_gradients(x_list,y_list):
        len_x_list=len(x_list)
        if len_x_list==1:
            m_list=y_list/x_list
            return m_list
        m_list=[]
        for idx in range(len_x_list-1):
            m_list.append((y_list[idx+1]-y_list[idx])/(x_list[idx+1]-x_list[idx]))
        return m_list
    def func_3piecewise(x, m_0, x_1, y_1, x_2, y_2, m_2):
        y = np.piecewise(x, [x <= x_1, (x > x_1) & (x <= x_2), x > x_2],
                        [lambda x:m_0*(x-x_1) + y_1, lambda x:y_1+(y_2-y_1)*(x-x_1)/(x_2-x_1), lambda x:m_2*(x-x_2) + y_2])
        return y
    if piece_number==3:
        lower_bounds=[-np.inf,xmin,-np.inf,xmin,-np.inf,-np.inf]
        upper_bounds=[np.inf,xmax,np.inf,xmax,np.inf,np.inf]

        w, _ = opt.curve_fit(func_3piecewise, x_sample, y_sample,bounds=(lower_bounds,upper_bounds),sigma=sigma)
        x_0=copy.deepcopy(xmin)
        y_0=func_3piecewise(x_0, *w).tolist()
        [m_0, x_1, y_1, x_2, y_2, m_2]=w
        m_1=func_gradients(x_2-x_1,y_2-y_1)
        result=[x_0,y_0,m_0,x_1,y_1,m_1, x_2, y_2, m_2]

The full code can be seen in pastebin

So, the question is: How to make a chaining (every tail of the piecewise function connected to the head of the next piece, or simply "connected") picewise function in python for general n-pieces? Other algorithm or solver is acceptable.

Edit: I add my result so far for 2 piecewise.

Update: I found that my code (for three pieces) is not working because of a small typo (sorry about this, just tell me if I should delete this question). Now it's working and I update the paste bin. But, if you have a general (flexible, no need to write function for each number variant) function that can generate n number of pieces,I'll gladly accept the answer.

Answer 1

You can parametrize on the distance x2-x1 instead of parametrizing on x2. Because you can give the optimizer bounds, you can set the distance to be greater than 0.

For example, to make a general piecewise-linear function with 4 intervals, define the following:

The points which separate the intervals and x0, x1 and x2. The slopes in the 4 intervals are m0, m1, m2 and m3. The value of the function at x0 is y0.

Define d1 = x1 - x0, d2 = x2 - x1. From here:

x1 = x0 + d1
x2 = x0 + d1 + d2

Then, you have 8 optimization parameters: x0, y0, d1, d2, m0, m1, m2 and m3. By nature of your optimization problem, all except x0 and y0 are non-negative.

Equation for the first interval:

y = m0 * (x - x0) + y0

Equation for the second interval:

y = m1 * (x - x0) + y0

Now you can get the rest of the equations in a recursive way, by applying the previous equation at the rightmost point of its interval. For the x1 point, the value of the function is:

y1 = m1 * d1 + y0

So the third equation is

y =
    m2 * (x - x1) + y1 =
    m2 * (x - x0 - d1) + m1 * d1 + y0

For the x2 point, this gives

y2 = m2 * d2 + y1

So the fourth equation is

y =
    m3 * (x - x2) + y2 =
    m3 * (x - x0 - d1 - d2) + m2 * d2 + m1 * d1 + y0