[ML]6. 파이썬으로 Gradient, minibatch Gradient 구현해보기

6 분 소요

0. Review

지난 포스트에서는 Linear Least Square 문제에서의 Gradient와 minibatch gradient, stochastic gradient 에 대해서 알아보았다. 이번포스트는 파이썬을 통해 Gradient , minibatch gradient, stocastic gradient 를 구현해볼 예정이다.

1. 데이터 받기 및 기본적인 작업

첫 포스트 에 언급했던 수업에서 준 Homework의 Data를 이용해 파이썬으로 간단하게 Gradient descent를 이용하였습니다. 데이터 및 기본 틀은 여기에서 다운받을 수 있습니다. 또한, 작성한 코드를 제 깃허브 에서 확인 할 수 있으니(보기 쉽게 실행은 Jupiter Notebook) 참고 바랍니다.(코딩이 아직 익숙하지 않아서, 좀 난잡해 보일 수 있습니다..ㅎ)

Gradient를 실행하기에 앞서서 해야할 작업들 몇가지를 먼저 살펴 봅시다.

  • Feature Scaling :

    Feature이 여러가지가 있고, 각각의 feature이 각각 다른 단위를 가진다고 생각해 봅시다. (ex:몸무게(kg) , 키(cm)) 이렇게 다른 단위를 가지게 된다면, feature마다 다른 범위를 가지게 될 것이고, 결과값이 특정 feature에만 영향을 받게되는 상황이 자주 벌어지게 됩니다. 이를 방지하기 위해 feature의 범위를 제한시키는 것을 feature scaling이라고 합니다. 대표적으로 Min-Max Scaling 이 있고, 특정 값의 위치를 최솟값과 최대값을 기준으로 어디에 있는지를 [0,1]사이에 나타낸 것입니다. 이외에도 다양한 방법이 있는데 데이터마다 적용되는 방식이 조금씩 다르다고 합니다. 더 자세한건 여기를 참고하면 좋을 것 같습니다.

  • Train set and Test Set:

    주어지는 데이터수는 한정되있으므로, 주어진 데이터를 특정 hyperparameter 별로 Gradient 를 시행해 weight을 구하는 train set과, 실제 성능을 구할 때 쓰는 data를 test data으로 나누고, test set을 다시 최종적으로 성과를 확인하는 test set과 hypermeter별로 성능을 확인하는 validation set으로 나눕니다. 통상적으로 약 7:1.5:1.5 ~ 8:1:1 의 비율로 나눈다고 합니다. 원래는 전부다 하는게 맞겠지만, 다 하기는 조금 버거워서, 구현은 그냥 Train set과 test set으로만 나누었습니다.

  • Code

    import numpy as np
    import random
    from sklearn.preprocessing import MinMaxScaler
    from sklearn.model_selection import train_test_split
    
    def featurescale(X_train,X_test):
        ## before doing the gradient, the data must be feature scaled.., 
        ## feature이 여러가지일때, 특정 feature의 범위는 gradient 를 시행할때, 큰 영향을 줄수있다. i.e 키, 몸무게가 feature 인 경우.
        ## 따라서, feature의 크기를 어느정도 조절해주는 도구가 필요한데, 이게 feature Scaling이다.
        ## feature scaling을 함에 있어서, gradient descent의 수렴속도가 훨씬 빨라질 수 잇다. 
        ## feature scaling은 min-max scaling과 standard normal scaling이 있는데, 두가지 모두 자주 사용된다. 
        ## 두가지 방법 은 쓰임세가 살짝 다르긴한데, 대부분의 상황에서 통용될수 있는 min-max scaling을 사용하도록 하겠다.
        scalar=MinMaxScaler()
        scalar.fit(X_train)
        X_train=scalar.transform(X_train)
        X_test=scalar.transform(X_test)
        ##주의!!!: testset의 Scaling은 train set 의 스케일링과 동일하게 진행.
        return X_train,X_test
    
    def splitData(X,y):
        ## 검증을 위해 데이터 X를 train 데이터와 test data로 분리시킴, 통상적으로 약 7:3 언저리 의 비율로 가름
        X_train,X_test,y_train,y_test=train_test_split(X,y,test_size=0.3,random_state=42)
        return X_train,X_test,y_train,y_test
          
    

    직접 구현할 수도 있겠지만, sklearn에 어느정도 익숙해지고 싶어서, sklearn의 툴을 사용하였다. 사용하기 위해서는 위의 import 부분을 실행하고 sklearn에서 적용예시 및 방법들을 보고 코드를 실행하면 된다. 주의할 부분은 feature scaling을 진행할 때, testset에 적용하는 scaling은 training set에 적용된 scaling과 동일하게 진행하는 것이다.

2. Basic Gradient descent 구현

Gradient의 기본적은 내용은 이전포스트들을 참고하면 좋을 것 같습니다. Loss function에는 Square-loss function을 이용하였습니다. Gradient를 구현하는것은 생각보다 간단합니다. 특정 점에서의 gradient를 구해주는 Computegrad function 과, 현재 점에서의 loss를 구해주는 Loss function 두가지만 있으면 됩니다. 자세한 것은 코드를 통해 살펴보겠습니다.

def SquareLossfunction(X,y,theta, l2_reg=0.00):
    ## Square Loss 
    ## X: n* d matrix  : n=number of data, d=number of feature
    ## Y: n*1 matrix: n= number of data
    ## theta: 1* d matrix, d: number of feature
    ## loss =avg(X*thata.T-Y+ l2_reg*(L2norm(theta)))
    ## getting average for every data set? might be slow...->> where SGD Occurs
    ##returns: loss(1*1)
    m=X.shape[0]
    loss_term=np.mean(np.square((np.dot(X,theta)-y)))
    reg_term=np.linalg.norm(theta)
    loss=loss_term+reg_term
    return loss
def computegrad(X,y,theta,l2_reg=0.00):
    ## X: n* d matrix  : n=number of data, d=number of feature
    ## Y: n*1 matrix: n= number of data
    ## returns: d*1 matrix (grad for each feature.)
    m=X.shape[0]
    temp=np.dot(X,theta)-y
    grad_term=(2.0/m)*(np.dot(X.T,temp))
    reg_term=2*l2_reg*theta
    return grad_term  
def gradDescent(X,y,alpha=0.05,num_iter=1000): 
    ## the very basic gradient descent..
    ##things to consider.. iteration, stepsize..
    feat_num=X.shape[1]
    data_num=X.shape[0]
    theta_hist=np.zeros((num_iter+1,feat_num)) ## storing the historical data of theta
    loss_hist = np.zeros(num_iter+1) ## storing loss value to see if gradient is doing well..
    theta_init=np.random.rand(feat_num) ## 통상적으로 0,1 사이에 랜덤하게 생성
    theta_hist[0,:]=theta_init
    loss_hist[0]=SquareLossfunction(X,y,theta_hist[0,:]) 
    ##iteration별로 computegrad, squarelossfunction 계산..
    for i in range(0,num_iter):
      cur_theta=theta_hist[i,:]
      theta_hist[i+1,:]=cur_theta-alpha*computegrad(X,y,cur_theta)
      loss_hist[i+1]=SquareLossfunction(X,y,theta_hist[i+1,:])
    return theta_hist,loss_hist

numpy의 함수중 하나인 np.dot 사용법만 익히면, 이런식으로 매우 간단하게, gradient Descent함수를 구현할 수 있습니다. Overfitting을 방지하기 위해 사용되는 L2-Regularization 포함하여 코딩을 했습니다만, 편의를 위해서 lambda를 0(reg효과 없이)으로 두고 하였습니다. Overfitting 및 regularization 문제는 추후에 더 자세하게 다루도록 하겠습니다.

그럼 실제로 실행해 봅시다. 메인함수나 JupyterNotebook을 열고, 데이터를 불러온 후 실행시켜 봅시다. 실행 코드 및 결과들은 깃허브에도 올려놓았습니다.

df=pandas.read_csv("data.csv") ## Jupiter Notebook이라 간단하게 표현가능, 실제로는 경로를 더 자세하게 표시해줘야함.
y=df['y']
X=df.loc[:,df.columns!='y']
## 데이터 불러오기
X_train,X_test,y_train,y_test=splitData(X,y)
## 성능 평가를 위해 test set 과 train set을 분리.
X_train,X_test=featurescale(X_train,X_test)
## (Min-Max feature scale)
theta_init=np.random.rand(X_train.shape[1]) 
loss=SquareLossfunction(X_train,y_train,theta_init)
num_iter=1000
alpha_samples=[0.1,0.05,0.01,0.001]
plt.figure(figsize=(5, 5))
index=np.arange(0,num_iter+1,1)
for alpha in alpha_samples:
    g_theta_hist,g_loss_hist, _ =gradDescent(X_train,y_train,alpha=alpha,num_iter=num_iter,backtracking=False)
    plt.plot(index,g_loss_hist,label=r"$\alpha$={}".format(alpha))
plt.title("Loss per Iteration")
plt.xlabel("iteration")
plt.ylabel("Loss")
plt.ylim(0,500)
plt.legend()

image

대부분의 $\alpha$에 대해서는 200번정도 안쪽의 iteration에 수렴하는 것을 확인 할 수 있습니다. 하지만, 만약 알파를 0.1이상으로 설정할 경우, Loss가 수렴하지못하고 발산하는 것을 확인 할 수 있습니다…(파란색).

기본적인 gradient descent는 설정되는 alpha값에 매우 민감합니다. 따라서, full-gradient를 적용할때는 backtracking 혹은 다른 Line search 방법들이 필수적으로 요구된다고 합니다.

추가로보면 좋은자료: 여기서 gradient descent 부분을 찾아 alpha를 조정해보면, alpha에 따라 매우 민감하게 변하는 결과값을 확인 할 수 있습니다.

3. Backtracking 구현

앞서 언급했듯이, 일반적인 Gradient descent는 초기값 $\alpha$ 에 매우 민감합니다. Backtracking은 위 같은 상황을 해소시켜줍니다.

코드를 통해 직접 구현해보겠습니다.

  • 코드

def gradDescent(X,y,alpha=0.1,num_iter=1000,backtracking=True): 
    ## the very basic gradient descent..
    ##things to consider.. initialization of theta, iteration..
    feat_num=X.shape[1]
    data_num=X.shape[0]

    theta_hist=np.zeros((num_iter+1,feat_num)) ## storing the historical data of theta
    loss_hist = np.zeros(num_iter+1) ## storing loss value to see if gradient is doing well..
    theta_init=np.random.rand(feat_num) ## 통상적으로 0,1 사이에 랜덤하게 생성
    theta_hist[0,:]=theta_init
    loss_hist[0]=SquareLossfunction(X,y,theta_hist[0,:]) ##이니셜
    bactracknum=0
    ## initialize theta with random number between 0 and 1
    if backtracking:
        alp, beta=0.3,0.9
        for i in range(0,num_iter):
            cur_theta=theta_hist[i,:]
            dx=-computegrad(X,y,cur_theta)
            t=1
            while True:
                n_loss=SquareLossfunction(X,y,cur_theta+t*dx)
                o_loss=SquareLossfunction(X,y,cur_theta)-alp*np.dot(dx,dx)*t
                if n_loss>o_loss:
                    t=0.9*t
                    bactracknum+=1
                else:
                    break
            theta_hist[i+1]=cur_theta-t*computegrad(X,y,cur_theta)
            loss_hist[i+1]=SquareLossfunction(X,y,theta_hist[i+1,:])
    else:
        for i in range(0,num_iter):
            cur_theta=theta_hist[i,:]
            theta_hist[i+1,:]=cur_theta-alpha*computegrad(X,y,cur_theta)
            loss_hist[i+1]=SquareLossfunction(X,y,theta_hist[i+1,:])
    return theta_hist, loss_hist,bactracknum

구현은 간단합니다. while 문에서 적당한 stepsize가 될 때 까지 stepsize에 $\beta$를 곱해서 줄여주고, 원하는 stepsize가 된다면, x를 이동시킨 후 다시 stepsize를 초기화시켜같은 행위를 반복합니다. 설정한 범위 내에서 가장 큰 stepsize를 사용할 수 있기 때문에, alpha를 설정해가면서 쓰는 일반적인 Gradient descent보다 많이 쓰인다고 합니다. +) 보통 $\beta$의 초기값은 0.9언저리에, $\alpha$ 의 초기값은 0.5 언저리로 시작한다고 합니다.

  • 결과 및 비교
    index=np.arange(0,num_iter+1,1)
    alpha_samples=[0.05,0.025,0.01]
    for alpha in alpha_samples:
      if alpha==0.001:backtrack=True
      g_theta_hist,g_loss_hist, _ =gradDescent(X_train,y_train,alpha=alpha,num_iter=num_iter,backtracking=False)
      plt.plot(index,g_loss_hist,label=r"$\alpha$={}".format(alpha))
    b_theta_hist,b_loss_hist,backtracknum=gradDescent(X_train,y_train,num_iter=num_iter,backtracking=True)
    plt.plot(index,b_loss_hist,label=r"with backtracking" )
    plt.title("Loss per Iteration")
    plt.xlabel("iteration")
    plt.ylabel("Loss")
    plt.ylim(0,10)
    plt.legend()
    

image

Backtracking을 시행했을 때, Stepsize를 크게 설정한 경우(아마 0.05언저리가 basic gradient를 했을 때, 수렴할 수 있는 가장 큰 stepsize가 아닐 까 싶습니다)와 비슷한 경향을 보이는 것 같습니다. 이외 alpha를 너무 작게 잡았을때보다, 수렴속도가 빠른 것을 확인 할 수 있습니다.

4. Basic gradient vs Mini-batch(Stochastic) gradient

이번에는 기본적인 gradient descent 와 데이터를 여러 배치로나눈 mini-batch gradient로 나누어 보겠습니다. (stochastic gradient는 batch size=1인 mini-batch gradient와 같습니다 ) 또한, stochastic의 강점을 조금 더 확실히 확인해보기 위해, 알기 쉽게 새롭게 데이터를 생성해 보았습니다. 데이터로 보기 편하도록 Y와 X가 선형관계에 있도록 만들어보았습니다.

\[Y=X+3\]

$X$는 0과 1사이의 난수를 생성해서 만들고, $Y$는 $X+3$에 오차를 조금 주어서 만들었습니다.

index=np.arange(0,101,1)
N = 100
theta = np.array([[1], [3]]) ##실제 theta값...
X = np.c_[np.random.rand(N,1), np.ones((N,1))]
y = np.dot(X,theta)+ 0.3*np.random.randn(N,1)
true_theta = np.linalg.lstsq(X, y, rcond=None)[0]
plt.scatter(X[:,0],y)

image

이 데이터를 바탕으로 일반적인 Gradient descent 와 SGD를 둘다 적용시켜 원하는 theta값을 찾아보겠습니다.

  • 코드 및 결과
y=np.squeeze(y)## making y from 2d to 1d
g_theta_hist,g_loss_hist,_ =gradDescent(X,y,alpha=0.01,num_iter=100,backtracking=False)

batchsize=[10,25,50]
plt.ylim(0,15)
index=np.arange(0,100,1)

plt.plot(index,g_loss_hist[0:100],label=r"$\alpha$={}".format(0.01))
b_loss_hist,b_theta_hist= minibatchgradDescent(X,y,batchsize=1,alpha=0.01,num_iter=100) ##setting batchsize=1-> SGD
plt.plot(index,np.mean(b_loss_hist,axis=1),label=r"sgd")
plt.legend()

image

데이터의 feature이 매우 적고, 예시로 만든 데이터에 불과하지만, SGD가 더 빠르게 수렴하는것을 확인할 수 있습니다. 더해서 수렴 속도의 차이를 확인 할 수 있도록, SGD와 GD의 Contour line 도 그려보았습니다.

##Draw a Contour line
Xmesh, Ymesh=np.meshgrid(np.linspace(0,5,100),np.linspace(0,5,100))
Z=np.zeros(Xmesh.shape)
for i in range(Xmesh.shape[0]):
    for j in range(Xmesh.shape[1]):
        theta1=Xmesh[i][j]
        theta2=Ymesh[i][j]
        Z[i,j]=SquareLossfunction(X,y,np.array([theta1,theta2]))

plt.figure(figsize=(15,10))
CS = plt.contour(Xmesh, Ymesh, Z, [0.03, 0.1, 0.3, 0.5, 1, 2, 3, 4, 5, 7])
plt.clabel(CS, inline=1, fontsize=10)
plt.plot(g_theta_hist[:,0],g_theta_hist[:,1],marker="o",label=r"$\alpha$={}".format(0.01))
plt.plot(b_theta_hist[:,99,0],b_theta_hist[:,99,1],marker="o",label="sgd")
plt.legend()

image

이처럼 SGD는 수렴만 할 수 있다면, 매우 빠르기 때문에 주로사용되는 기법 중 하나라고 합니다. 하지만, SGD는 하나의 데이터마다 update를 진행하기 때문에, alpha값에 매우 불안정하고, 데이터에 따라 수렴하지 않는 경우도 매우 많다고한다. 따라서 Step size를 조절하는 다양한 기법들이 존재하고, 기회가 된다면 더 다루어보기로 하겠습니다.

추가로, SGD를 구현할때는 Backtracking을 구현하지 않았는데, 이유는 다음과 같다. 다음은 위키피디아에서 설명하는 Stochastic gradient 의 일부입니다.

On the other hand, adaptive SGD does not guarantee the “descent property” which Backtracking line search enjoys .If the gradient of the cost function is globally Lipschitz continuous, with Lipschitz constant L, and learning rate is chosen of the order 1/L, then the standard version of SGD is a special case of backtracking line search.

Backtracking 부분을 해석해보면, Backtracking을 사용하는 이유는 올바른 stepsize를 매 iteration마다 찾기 위함인데, SGD에 backtracking을 사용하는것은 SGD가 가질 수 있는 속도라는 장점을 무효화시키는것이나 마찬가지이기 때문이라는겁니다. (cost function이 global하게 Lipschitz continuous 할 때 Stepsize를 1/L로 잡고 SGD를 하는 경우, Backtracking과 거의(Backtracking은 안전한 stepsize가 나올때까지 Stepsize를 줄이기 때문에)동일하고, L을 모르는 상태에서는 속도를위해 SGD를 적용할것인지, 조금 느리더라도 Backtracking을 적용할 것인지는 실험자가 선택해야하는 부분입니다. )

댓글남기기