[Regression] Piecewise Regression with numpy

[Regression] Piecewise Regression with numpy

Piecewise Regression이란 하나의 식으로 전체 구간을 설명하기 어려울 경우 구간에 따라 모델을 다르게 쓰는 기법이다.

우선 데이터를 생성한다.

import numpy as np
import matplotlib.pyplot as plt
from scipy import optimize

np.random.seed(0)
x = np.random.normal(0, 1, 1000) * 10
# 조건에 맞는 값을 특정 다른 값으로 변환하기
# np.where(조건, 조건에 맞을 때 값, 조건과 다른 때 값) if 문과 비슷함. 그냥 where은 조건에 맞는 인덱스 반환
y = np.where(x < -15, -2 * x + 3, np.where(x < 10, x + 48, -4 * x + 98)) + np.random.normal(0, 3, 1000)
plt.scatter(x, y, s = 5, color = u'b', marker = '.', label = 'scatter plt')

원래는 random array가 매번 바뀌지만 np.random.seed는 괄호 안의 값에 따라 동일한 랜덤 array를 생성하도록 한다. 

a = np.array([1, 2, 3, 4, 5]) 의 배열이 있을때, np.where(a < 3)을 하면 조건에 맞게 1과 2의 인덱스 값(0, 1)을 반환한다.

3개의 parameter를 이용하는 방식도 있다.

np.where(a < 3, 10, a) 를 하면 [10, 10, 3, 4, 5]를 반환한다. 즉, 조건에 맞으면 그 배열 값을 10으로 바꾸고 아니면 원래 값으로 둔다는 식이다. if 문과 비슷하다.

x와 y를 scatter plot으로 보면 아래와 같다. 

위의 데이터에 np.piecewise를 적용하기 전에 쉬운 예로 이해를 돕고자 한다.

x = np.linspace(-2.5, 2.5, 6)
print(x)  # [-2.5 -1.5 -0.5  0.5  1.5  2.5]

condlist = [x < 0, x >= 0]
print(condlist) ## [array([True, True, True, False, False, False]), array([False, False, False, True, True, True,])]

t = np.piecewise(x, [x < 0, x >= 0], [lambda x: -3*x, lambda x: 3*x])
plt.plot(range(len(x)), t)

나누고자 하는 구간이 2개이면 condlist의 조건도 2개이고 2개의 array를 반환한다. 그리고 그 condlist가 첫번째 lambda와 두번째 lambda의 input으로 입력되고, True면 계산되고 False면 0을 반환하여 두 lambda 함수의 합을 최종적으로 np.piecewise의 값으로 반환한다.

하지만 우리는 위에서 0으로 두었던 분기점을 정하지 않았기 때문에 데이터에 맞춰서 최적의 분기점을 찾도록 할 것이다.

그러면 scipy.optimize의 목적함수로 쓰일 piecewise_linear 함수를 정의해야 한다. 

condlist는 3구간으로 나올수 있도록 하고(분기점의 수만큼 조건식이 필요) 그에 대응될 lambda 함수의 모양을 짠다. lambda 함수가 이어지도록 되어있는데 k1*x + b + k2*(x-x0)에서 앞의 k1*x + b 는 분기점에서의 y값(새로운 상수항)이 되므로 이렇게 쓰지 않으면 에러가 난다.

def piecewise_linear(x, x0, x1, b, k1, k2, k3):
    condlist = [x < x0, (x >= x0) & (x < x1), x >= x1]
    funclist = [lambda x: k1*x + b, lambda x: k1*x + b + k2*(x-x0), lambda x: k1*x + b + k2*(x-x0) + k3*(x - x1)]
    return np.piecewise(x, condlist, funclist)

p , e = optimize.curve_fit(piecewise_linear, x, y)

이 과정을 거쳐 데이터에 가장 적합한(에러가 최소화되는) x0, x1, b, k1, k2, k3이 정해지고 아래와 같이 선이 그려진다.

xd = np.linspace(-30, 30, 1000)
plt.plot(x, y, "o")
plt.plot(xd, piecewise_linear(xd, *p))