[ 과제 1 ] Gradient Descent

써먹을 일 없다고 내려두었던 내용들이 누군가에게 도움이 된다면, 

Linear regression 에 대해서는 closed-form 인 pseudo inverse 존재하기 때문에 대부분은 gradient descent 방식을 시도 조차 하지 않습니다. 하지만, nonlinear 로 넘어가면 꽤나 간단하면서 직관적인 최적화 방법이기 때문에 많이 사용하고 있습니다.

제가 이번에 다룰 내용은 polynomial 함수 (y = ax^2 + bx + c), 2차 함수에 대해서 gradient descent 방법을 통해 최적화 하는 방법을 다루도록 하겠습니다. 시작에 앞서 기초적인 내용과 수식을 설명 드리겠습니다.

수식 전개에 앞서서 저랑 약속 한개만 할께요, 모든 벡터는 열벡터 이고 위와 같이 열벡터를 볼드체로 표시하도록 하겠습니다. 혼동하지 마세요~

회사 선배가 갑자기 와서 저에게 데이터를 던져주고  입력으로 말도 안되는 데이터를 주고

"입력이 들어가면 출력이 이런데, 함수를 찾아서 알려줘, 그런데 예전에  이게 2차 함수였던 것 같은데,,"

 모든 힌트가 다 나온 것 같습니다. 주어진 데이터를 관찰 데이터 라고 가정하면 입력값과 출력값 pair 가 n 개 주었졌을 것이고 가정으로는 2차 polynomial 함수 라는 것도 알 수 있었네요. 이를 통해서 오차함수 를 구하면 위와 같이 matrix form 까지 전개가 가능합니다.

이번 주제가 바로 gradient descent 방식을 통한 함수 최적화 라고 말씀드렸죠??

최적화 방법으로 gradient descent 방식을 사용하려면,  Gradient 와 Jacobian 이란 이상한 용어와 마주해야 됩니다.

둘다 정말 간단히 설명 드리면, 미분입니다. Gradient 경우는 multi-variable scalar 함수 에 대해서 각 variable 로 미분을 한 것이고 Jacobian 의 경우는 이를 한단계 더 확장한 multi-variable vector 함수에 대한 미분이라고 생각하시면 됩니다.

다 잊어도 상관없지만, 정말 중요한 것은 둘다 미분이란 것입니다.

아래와 같이 Gradient Descent 방법을 통해서 파라미터를 업데이 해보도록 하겠습니다. 

사용코드는 python 을 통해서 전개 하도록 하겠습니다. 

예전에 , 후배에게 이 내용을 설명해줄 때 넌지시 Jacobian 은 내가 코딩해줄께 라고 생색을 냈는데, 실제 linear regression 에 대한 Jacobian 은 나중에 아시면 너무 쉽습니다. 아니. 큰 의미가 없다고 봐도,,,

 학습 데이터 에 대해서 Gradient Descent 방식을 통해서 예측 값이 점점 다가가는 것을 확인 할 수 있습니다.

아래 python 코드를 완성해 보시면 충분히 어렵지 않다는 것을 알 수 있습니다. Gradient descent 을 사용해서 $f(x) = ax^2 + bx +c$ 를 학습해 보는 코드입니다. 

#!/home/roadking/anaconda3/bin/python
import numpy as np
import numdifftools as nd
import matplotlib.pyplot as plt


def optimizer_gradient(xdata,ydata,learning_rate,max_iter,epsilon=0.0):

    aidx = 0;
    f, axarr = plt.subplots(2,2)
    axarr[0,0].plot(xdata,ydata,'ro')
    axarr[0,0].set_title('original')

    #model 2th-polynomial: ax^2 + bx +c
    #pt = np.random.rand(3,1)
    pt = np.zeros((3,1))

    r_fun = lambda p:(ydata -(p[0,0]*xdata**2+p[1,0]*xdata+p[2,0]))
    Jrp = nd.Jacobian(r_fun)

    for i in range(1,max_iter):
        r = ydata -(pt[0,0]*xdata**2+pt[1,0]*xdata+pt[2,0])
        pt = pt - 2*learning_rate*np.transpose(Jrp(pt)).dot(r)
        error = np.sum(np.absolute(r))
        if(i%500 == 0 or i == max_iter -1):
            print('iter=[{0}], error={1}'.format(i,round(error,4)))
            aidx +=1
            axarr[int(aidx/2),aidx%2].plot(xdata,ydata,'ro')
            axarr[int(aidx/2),aidx%2].plot(xdata,pt[0,0]*xdata**2+pt[1,0]*xdata+pt[2,0])
            axarr[int(aidx/2),aidx%2].set_title('iter[{0}]'.format(i))
    print(pt)
    plt.setp([a.get_xticklabels() for a in axarr[0, :]], visible=False)
    plt.setp([a.get_yticklabels() for a in axarr[:, 1]], visible=False)
    plt.show()

if __name__ == '__main__':
    xdata = np.reshape(np.arange(0,1,0.1),(-1,1))
    ydata = 2.1*xdata**2 -1.5*xdata + 0.5
    optimizer_gradient(xdata,ydata,0.05,1500)