Non-linear Regression with Python

비선형모델을 시계열 데이터에 fitting해달라는 요구를 받고 프로그램을 만들었다.
그런데 막상 다 만들고나서야 scipy에서 커브피팅 모듈을 제공한다는걸 알게되었다.
심지어 내꺼보다 더 잘돼

어쨌든 눈물을 머금고 기록을 남겨본다.

1. Scipy를 이용한 비선형회귀

정말 뭐 별거 없다.

  • 함수를 정의하고
  • 피팅 대상이 될 x값, y값을 넣어주면
  • 그 결과가 공분산과 함께 잘 나온다.
  • 잘 안된다 싶으면 초기값을 잘 설정해주거나 iteration을 높여준다.

코드로 보면 정말 간단하다.

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

weight, cov = curve_fit(function, x, y, [0, 0], maxfev=100000)

[0, 0]은 초기값, maxfev는 최대 iteration이다.
가끔 복잡한 모델을 사용하면 함수 콜이 너무 많아져서 계산이 불가능하다는 에러가 나온다.

RuntimeError: Optimal parameters not found: Number of calls to function has reached maxfev = 1000

이럴땐 초기값을 잘 주거나, maxfev를 높여주면 해결된다.
왜 이런 에러가 뜨는지는 모듈을 뜯어봐야 알 것 같다.
느낌상으론 maximum iteration을 잡아둔것같은데, 내 구현과 얼마나 다른지 확인하는 겸 뜯어봐야겠다.

아 어쨋든 잘된다.

물론 피팅이 덜 된것은 있는데, 이건 알고리즘 문제가 아니라 전처리의 문제이니 RANSAC이나 평활화로 해결하면 될 일이다.

2. 경사하강법

회귀문제는 1. 수학적(해석적) 접근이나 2. 수치해석적 접근의 두 가지 방법으로 해결한다.
수학적 접근이라 함은 최대우도추정(Maximum Likelihood Estimation, MLE)이나 최소제곱법(Ordinary Least Squares, OLS) 등의 방법으로 최적해를 계산하는 방법을 말하는 것이고
수치해석적 접근이라 함은 경사하강법(Gradient Descent, GD)과 같은 방법으로 근사해를 계산하는 방법을 말하는 것이다.

수학적 접근방법의 가장 큰 장점은 빠르고, 수학적으로 optimal하다는 것이 증명되어있다는 점인데
전제조건이나 연산속도때문에 MLE나 OLS를 사용하지 못하는 경우가 발생할 수 있다.
예를 들어 Loss가 Non-convex하거나 수십만차원을 가진 행렬의 역행렬을 구해야한다거나 하는 경우에는 수치해석적으로 근사해를 구하는 것이 좋다.

간단하게 OLS를 살펴보자

  • Loss가 Convex하다는 가정하에, Loss를 MSE로 설정한다.
  • Loss의 최솟값을 구하기 위해 각 weight의 gradient가 0이 되는 점을 구한다.

가장 좋은 예시가 최소제곱법(Ordinary Least Squares, OLS)을 이용하는 선형회귀 문제이다.
예시는 위키만 봐도 충분하다.

반면 경사하강법은 Loss가 낮아지는 곳으로 계속 나아가는 방법이다.
이는 다시말해, 현 지점에서 구한 Gradient의 역방향으로 계속 나아가는 것이다.
이 때 중요한 것은

  • 얼만큼씩 나아갈것인가 (step size)
  • 나아가는 방향에 일관성을 둘 것인가 (momentum)

이다.
물론 이런것들에 대한 고려를 하지 않아도 근사해는 구할 수 있다.
다만 얼마나 빠르게 근사할 수 있는가만 다를 뿐이다.
최적화에 대한 연구는 이미 충분히 진행됐고, 아직도 진행중이다.
보통 Adam을 많이 사용한다.

일단 주어진 모델을 바탕으로 Loss의 gradient를 계산해야한다.
이건 뭐.. 계산이니 직접 하자.
예를 들어 다음과 같은 지저분한 과정을 거쳐야 한다.
(아래 예시는 Loss가 아닌, model에 대한 gradient이다)

이렇게 구해낸 gradient를 바탕으로 iteration을 반복하며 Loss를 계속 줄여나가면 된다.
최적화는 이 gradient를 어떻게 update하는지에 따라 달려있다.
필자도 Fully GD를 쓰다가 RMSProp으로 바꿔봤더니 수렴하는 속도가 엄청 빨라졌다.

다만 이런 수치해석적 방법은 하이퍼파라미터가 늘어난다는 문제가 있다.
예를 들어

  • gradient를 얼마나 update시킬 것인지 (learning rate)
  • gradient update에 있어 최근 gradient를 얼마나 반영할 것인지 (decay)
  • gradient에 가중평균을 사용한다면 그 가중치는 어떻게 할 것인지 (gamma)
  • 초기값은 어떻게 선정할 것인지
    (Non-convex에서 초기값 선정의 중요성은 이전 포스팅에서도 언급한 바 있다.)

등의 추가적으로 고려해야 할 사항들이 있다.
뭐 여튼 필자가 구현한것도 안되지는 않았는데, 시간도 오래걸리고 generalization 또한 약해서 매번 값이 다르게 나왔다.
scipy는 어떻게 구현했는지 꼭 확인해봐야겠다.

*추가
확인해보니 scipy의 least_squares 함수는 Levenberg-Marquardt 알고리즘을 통해 함수식을 피팅한다.
워낙 유명한 알고리즘인데 단순하게만 말하면 뉴턴-랩슨법과 경사도 하강법을 조합하여 지역 최적해를 찾아내는 방법이다.
지역 최적해를 찾는다는 말은 즉 초기값에 매우 민감하다는 것을 의미한다.
iteration만 변화를 줄 것이 아니라, 초기값에도 변화를 주어서 Random Walk시키면 더 좋은 결과가 나올 것이다.

댓글 남기기