핸즈온 머신러닝 2 복습하기(챕터 4: 모델 훈련)

 

CH 4. 모델 훈련

정규방정식(normal equation)

비용 함수를 최소화하는 $\theta$ 값을 찾기 위한 해석적인 방법이 있다. 이를 정규방정식이라고 한다. 

$\Theta=(X^TX)^{-1}X^Ty$

이 공식을 테스트하기 위해 선형처럼 보이는 데이터를 생성해본다. 선형 함수는 $y = 4 + 3x_1 + gaussian_noise$ 이다. 

import numpy as np

# np.random.randint(n, m): n과 m-1 사이에서 랜덤 숫자 뽑기
X = 2 * np.random.rand(100, 1) # rand: 0~1 사이의 균일분포 표준정규분포 난수 생성
y = 4 + 3 * X + np.random.randn(100, 1) # 기댓값이 0이고 표준편차가 1인 가우시안 표준정규분포 난수 생성
X_b = np.c_[np.ones((100, 1)), X] # 모든 샘플에 x0 = 1을 추가
theta_best = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)
# 선형대수 모듈(np.linalg)에 있는 inv() 함수를 사용해 역행렬 계산 
# dot() 메서드를 통해 행렬 곱셈 

본래 $\theta_0$=4, $\theta_1=3$이 가장 정확하지만 아래 출력값을 보면 잡음 때문에 완벽하게 재현하지 못했다. 

theta_best
>>> array([[4.21509616],
           [2.77011339]])

어쨌거나 $\hat{\theta}$을 가지고 예측을 해보자. 

X_new = np.array([[0], [2]])
X_new_b = np.c_[np.ones((2, 1)), X_new]
y_predict = X_new_b.dot(theta_best)
y_predict

사이킷런에서 선형 회귀를 수행하는 것은 간단하다. 특성의 가중치(coef_)와 편향(intercept_)을 분리하여 저장한다. 

from sklearn.linear_model import LinearRegression
lin_reg = LinearRegression()
lin_reg.fit(X, y)
lin_reg.intercept_, lin_reg.coef_ # coef=coefficient(계수), 즉 파라미터를 가리킴
lin_reg.predict(X_new)

LinearRegression 클래스는 scipy.linalg.lstsq() 함수를 기반으로 하며 이 함수를 직접 호출할 수 있다. 

theta_best_svd, residuals, rank, s = np.linalg.lstsq(X_b, y, rcond=1e-6)
theta_best_svd

위 관련 내용은 유사역행렬이란 개념에 대한 이해가 필요하다. 현재는 다루지 않고 이후에 공부해서 살펴볼 예정이다. 이후 나오는 계산 복잡도(computational complexity) 역시 책 설명만 보고 이해하지 못했기 때문에 다시 공부할 필요가 있다. 일단 이해한만큼 적어보겠다. 

정규방정식은 (n+1) x (n+1) 크기가 되는 $X^TX$의 역행렬을 계산한다. 이를 계산하는 계산 복잡도는 일반적으로 $O(n^.24)$에서 $O(n^3)$ 사이이다. 반면, 사이킷런의 LinearRegression 클래스가 사용하는 SVD 방법은 약 $O(n^2)$이다. 주의할 점은 두 방식 모두 훈련 세트의 샘플 수에 대해서는 선형적으로 증가, 즉 계산 복잡도가 $O(m)$이라는 것이다. 이유는 $X^TX$는 (n$\times$m)(m$\times$n) = (n$\times$n) 크기의 행렬로 샘플 수(m)는 복잡도를 증가시키지 않고 점곱의 양만 선형적으로 증가시키기 때문이다. 

경사 하강법(gradient descent)

경사 하강법은 여러 종류의 문제에서 최적의 해법을 찾을 수 있는 일반적인 최적화 알고리즘이다. 기본 아이디어는 비용 함수를 최소화하기 위해 반복해서 파라미터를 조정해가는 것이다. 비용 함수는 기본적으로 볼록(convex) 형태이고, 최소화하는 지점은 바로 미분값이 0이 되는 곳이다. 구체적으로 말하면, $\theta$를 임의의 값으로 시작하여(무작위 초기화, random initialization) 한 번에 조금씩 비용 함수(예: MSE)가 감소되는 방향으로 진행하여 알고리즘이 최솟값에 수렴할 때까지 점진적으로 향상시키는 것이다. 

또 하나 중요한 파라미터는 스텝(step)의 크기로, 학습률(learning rate) 하이퍼파라미터로 결정하면 된다. 학습률이 너무 작으면 알고리즘이 수렴하기 위해 반복을 많이 진행해야 하므로 시간이 오래 걸리는 단점이 있다. 반대로 학습률이 너무 크면 알고리즘을 더 큰 값으로 발산하게 만들어 적절한 해법을 찾지 못하게 한다. 

중요한 것은 비용 함수가 항상 볼록 형태가 아니라는 것이다. 여러 개의 볼록 형태가 나타날 경우, 지역 최솟값(local minimum)에 수렴하기 때문에 진정한 전역 최솟값(global minimum)을 찾지 못한다. 선형 회귀를 위한 MSE 함수는 언제나 볼록 함수이다. 이는 지역 최솟값이 없고 하나의 전역 최솟값만 있다는 뜻이며, 연속된 함수이고 기울기가 갑자기 변하지 않는다.

(참고) 이런 함수의 도함수를 립시츠 연속(Lipschitz continous)이라 하는데, 도함수가 일정한 범위 안에서 변하는 성질을 가진다. 예를 들어, sin(x)는 립시츠 연속 함수지만, $\sqrt{x}$의 경우 x = 0일 때, 기울기가 무한대가 되므로 립시츠 연속 함수가 아니다. MSE는 x가 무한대일 때 기울기가 무한대가 되므로 국부적인(locally) 립시츠 함수라고 한다. 

경사 하강법을 사용할 때는 반드시 모든 특성이 같은 스케일을 갖도록 만들어야 한다. 사이킷런의 StandardScaler를 사용하면 데이터의 각 특성에서 평균을 빼고 표준편차로 나눠 평균을 0, 분산을 1로 만든다. 모델 훈련이 비용 함수를 최소화하는 모델 파라미터의 조합을 찾는 일을 모델의 파라미터 공간(space)에서 찾는다고 말하기도 한다. 

배치 경사 하강법(Batch gradient descent)

경사 하강법을 구현하려면 각 모델 파라미터 $\theta_j$에 대해 비용 함수의 그레이디언트를 계산해야 한다. 즉 $\theta_j$가 조금 변경될 때 비용 함수가 얼마나 바뀌는지 계산해야 하는 것이며 이를 편도함수(partial derivative)라고 한다. 다음은 비용 함수의 편도함수다. 

$\partial\over\partial \theta_j$ $MSE(\theta)$ = $2\over{m}$ $\sum_{i=1}^m (\theta^T x^{(i)} - y^{(i)}) x^{(i)}_j$

이를 벡터식으로 표현하면 다음과 같다. 매 경사 하강법 스텝에서 전체 훈련 세트 X에 대해 계산하고 이것이 바로 배치 경사 하강법이라고 한다. 아래 식의 최종 계산 결과는 파라미터 수와 같은 n x 1 벡터가 된다. 

$2\over{m}$ $X^T(X\theta - y)$

다음은 이제 $\theta$에서 $\Delta$ $MSE(\theta)$ 값을 빼야 한다. 여기서 학습률 $\eta$가 사용된다. 여기서 학습률은 우리가 흔히 보던 $\alpha$인 듯하다. 아래는 경사 하강법 스텝을 업데이트하는 알고리즘을 구현한 식이다. 

eta = 0.1
n_iterations = 1000
m = 100

theta = np.random.randn(2, 1)

for iteration in range(n_iterations):
    gradients = 2/m * X_b.T.dot(X_b.dot(theta) - y)
    theta = theta - eta * gradients

여기서 반복 횟수를 적절하게 지정해야 한다. 간단한 해결책은 반복 횟수를 아주 크게 지정하고 그레이디언트 벡터가 아주 작아지면, 즉 벡터의 노름이 어떤 값 $\epsilon$(허용오차)보다 작아지면 경사 하강법이 거의 최솟값에 도달한 것으로 알고리즘을 중지하는 것이다. 

확률적 경사 하강법(Stocastic gradient descent)

매 스텝에서 한 개의 샘플을 무작위로 선택하고 그 하나의 샘플에 대한 그레이디언트를 계산한다. 비용 함수가 최솟값에 다다를 때까지 부드럽게 감소하지 않고 위아래로 요동치며 평균적으로 감소한다. 이는 결국 최솟값에 가까워지나 안착하지는 못하는 것을 의미한다. 이를 해결하는 방법은 학습률을 점진적으로 감소시키는 것이다. 매 반복에서 학습률을 결정하는 함수를 학습 스케줄(learning schedule)이라 부른다. 아래 학습 스케줄을 사용한 확률적 경사 하강법의 예시이다.

n_epochs = 50
t0, t1 = 5, 50 # 학습 스케줄 하이퍼파라미터

def learning_schedule(t):
    return t0 / (t + t1)
    
theta = np.random.randn(2, 1) # 무작위 초기화

for epoch in range(n_epochs):
    for i in range(m):
        random_index = np.random.randint(m)
        xi = X_b[random_index:random_index+1] # 샘플 하나하나
        yi = y[random_index:random_index+1]
        gradients = 2 * xi.T.dot(xi.dot(theta) - yi)
        eta = learning_schedule(epoch * m + i)
        theta = theta - eta * gradients

확률적 경사 하강법은 무작위성의 위험이 있다. 알고리즘이 모든 샘플을 사용하게 하려면 훈련 세트를 섞은 후 (입력 특성과 레이블을 동일하게 섞어야 함) 차례대로 하나씩 선택하고 다음 에포크에서 다시 섞는 식의 방법을 사용할 수 있다. 사이킷런에서는 SGD 방식으로 선형 회귀를 사용하려면 기본값으로 제곱 오차 비용 함수를 최적화하는 SGDRegressor 클래스를 사용한다. 다음 예시를 살펴보자.

from sklearn.linear_model import SGDRegressor
sgd_reg = SGDRegressor(max_iter=1000, tol=1e-3, penalty=None, eta0=0.1)
sgd_reg.fit(X, y.ravel())

sgd_reg.intercept_, sgd_reg.coef_

max_iter은 최대 반복 횟수를 의미하고, tol=1e-3은 한 에포크에서 0.001보다 적게 손실이 줄어들 때까지 실행하라는 뜻이다. eta0은 학습률을 의미한다. 

미니배치 경사 하강법(Mini-batch gradient descent)

미니배치라 부르는 임의의 작은 샘플 세트에 대해 그레이디언트를 계산한다. 미니배치를 어느 정도 크게 하면 알고리즘은 파라미터 공간에서 SGD보다 덜 불규칙하게 움직이며 이는 결국 SGD보다 최솟값에 더 가까이 도달하게 됨을 의미한다. 하지만 지역 최솟값에서 빠져나오기는 더 힘들지도 모른다. 

배치 경사 하강법, 확률적 경사 하강법, 미니배치 경사 하강법을 간단히 비교하자면 다음과 같다. 보통 배치 경사 하강법은 실제로 최솟값에 도달할 수 있는 반면, 확률적 경사 하강법과 미니배치 경사 하강법은 근처에서 맴돌게 된다. 하지만 배치 경사 하강법에는 매 스텝에서 많은 시간이 소요되고, 확률적 경사 하강법과 미니배치 경사 하강법은 적절한 학습 스케줄을 사용하면 최솟값에 도달할 수 있다. 

 

알고리즘	m이 클 때	외부 메모리 학습 지원	n이 클 때	하이퍼파라미터 수	스케일 조정 필요	사이킷런
정규방정식	빠름	No	느림	0	No	N/A
SVD	빠름	No	느림	0	No	LinearRegression
배치경사하강법	느림	No	빠름	2	Yes	SGDRegressor
확률적 경사 하강법	빠름	Yes	빠름	$\ge$2	Yes	SGDRegressor
미니배치 경사 하강법	빠름	Yes	빠름	$\ge$2	Yes	SGDRegressor

다항 회귀

다항 회귀는 간단히 말해서 비선형 함수다. 직선으로 표현되지 않고 거듭제곱, 세제곱 등 다항이 포함된 것이다. 비선형으로 판단되는 데이터가 있다면 훈련 데이터 역시 새로운 특성이(예: 제곱항) 추가되어야 한다. 이는 사이킷런의 PolynomialFeatures를 사용하여 훈련 데이터를 변환하면 된다.

from sklearn.preprocessing import PolynomialFeatures
poly_features = PolynomialFeatures(degree=2, include_bias=False) # 기본값: degree=2 / include_bias=True
# include_bias=True -> 편향을 위한 특성(x0)인 1이 추가
X_poly = poly_features.fit_transform(X)
X_poly[0]

특성이 여러 개일 때 다항 회귀는 이 특성 사이의 관계를 찾을 수 있다. PolynomialFeatures가 주어진 차수까지 특성 간의 모든 교차항을 추가하기 때문이다. 예를 들어, 두 개의 특성 $a$, $b$가 있을 때 degree=3이라고 한다면 $a^2$, $a^3$, $b^2$, $b^3$뿐 아니라 $ab$, $a^2b$, $ab^2$도 특성으로 추가한다.

(참고) PolynomialFeatures(degree=d)는 특성이 n개인 배열을 특성이 $(n+d)!\over d!n!$개인 배열로 변환한다. 

학습 곡선

모델의 일반화 성능을 추정하기 위해선 교차 검증이 필요하다. 하지만 이외에도 사용할 수 있는 방법은 학습 곡선을 살펴보는 것이다. 이는 훈련 세트와 검증 세트의 모델 성능을 훈련 세트 크기(또는 훈련 반복)의 함수로 나타낸다. 이 그래프를 생성하기 위해선 단순히 훈련 세트에서 크기가 다른 서브 세트를 만들어 모델을 여러 번 훈련시키면 된다. 아래는 모델의 학습 곡선을 그리는 함수이다.

from sklearn.metrics import mean_squared_error
from sklearn.model_selection import train_test_split

def plot_learning_curves(model, X, y):
    X_train, X_val, y_train, y_val = train_test_split(X, y, test_size=0.2)
    train_errors, val_errors = [], []
    for m in range(1, len(X_train)):
        model.fit(X_train[:m], y_train[:m])
        y_train_predict = model.predict(X_train[:m])
        y_val_predict = model.predict(X_val)
        train_errors.append(mean_squared_error(y_train[:m], y_train_predict))
        val_errors.append(mean_squared_error(y_val, y_val_predict))
    plt.plot(np.sqrt(train_errors), "r-+", linewidth=2, label="훈련 세트") #RMSE 측정
    plt.plot(np.sqrt(val_errors), "b-", linewidth=3, label="검증 세트")

학습 곡선이 나온걸 볼 때 훈련 세트와 검증 세트 각각 그래프가 일정 훈련 세트 크기를 지나면서 거의 붙어있는 형태가 나오면 과소적합한 것이다. 반대로 과대적합한 모델은 일정 훈련 세트 크기를 지나면 훈련 세트, 검증 세트 그래프 간의 공간이 비교적 넓다. 더 큰 훈련 세트를 사용하면 두 그래프가 점점 가까워질 것이다.

여기서 잠깐! 편향/분산 트레이드오프에 대해 살펴보고 가자.

모델 복잡도가 커지면(다항 회귀 등) 통상적으로 분산이 늘어나고 편향이 줄어든다. 반대로 복잡도가 줄어들면 편향이 커지고 분산이 작아지게 되어 트레이드오프라고 부른다. 

규제가 있는 선형 모델

과대적합을 감소시키는 좋은 방법은 모델을 규재하는 것이다. 일반적으로 자유도(DF, degree of freedom)를 줄이면 데이터에 과대적합하기 어려워지는데, 다항 회귀 모델의 경우 다항식의 차수를 감소시키면 된다. 여기서 주의할 점은 10차원을 3차원으로 아예 줄인다는 것이 아니라, 변수는 그대로 두고 계수 값을 0에 가깝게 하여 사실상 사용하지 않는 것처럼 만든다는 것이다. 선형 회귀 모델은 가중치를 제한함으로써 규제를 가한다. 각기 다른 방법으로 라쏘 회귀(L1 Norm), 릿지 회귀(L2 Norm), 엘라스틱넷(L1+L2 Norm)이 있다.

<출처: Zou, H., & Hastie, T. (2005)>

라쏘(Lasso, least absolute shrinkage and selection operator) 회귀

릿지 회귀처럼 비용 함수에 규제항을 더하지만, l2 노름의 제곱을 2로 나눈 것 대신 가중치 벡터의 l1 노름을 사용한다. 아래 라쏘 회귀 비용 함수를 보면 알 수 있다. 

$J(\theta)=MSE(\theta) + \alpha$$\sum_{i=1}^n|\theta_i|$

라쏘 휘귀의 중요한 특징은 덜 중요한 특성의 가중치를 제거하려고 한다는 점이다. 다시 말해 라쏘 회귀는 자동으로 특성 선택(variable selection)을 하고 희소 모델(sparse model)을 만든다(변수 대부분이 0이고, 일부 변수만 0이 아님). 라쏘와 다른 릿지의 특징은 첫째, 파라미터가 전역 최적점(global minimum)에 가까워질수록 그레이디언트가 작아진다. 따라서 경사 하강법이 자동으로 느려지고 수렴에 도움이 된다. 둘째, $\alpha$를 증가시킬수록 최적의 파라미터가 원점에 더 가까워진다. 

라쏘의 비용 함수는 $\theta_i$=0(i=1, 2, ... n)에서 미분 가능하지 않다. 하지만 서브그레이디언트 벡터(subgradient vector) g를 사용하면 문제가 없다. 다음은 서브그레이디언트 벡터 공식이다.

$g(\theta, J)=MSE(\theta)+\alpha(sign(\theta_i))$ , 여기서 $sign(\theta_i)=-1(\theta_i<0일 때), 0(\theta_i=0일 때), 1(\theta_i>0일 때)$

다음은 Lasso 클래스를 사용한 간단한 사이킷런 예제다. Lasso 대신 SGDRegressor(penalty="l1")을 사용할 수도 있다.

from sklearn.linear_model import Lasso
lasso_reg = Lasso(alpha=0.1)
lasso_reg.fit(X, y)
lasso_reg.predict([[1.5]])

릿지(ridge) 회귀 = 티호노프(Tikhonov) 규제

규제항이 $\alpha \sum_{i=1}^n \theta^{2}_i$ 비용 함수에 추가된 선형 회귀 버전으로 학습 알고리즘을 데이터에 맞추는 것뿐만 아니라 모델의 가중치가 가능한 한 적게 유지되도록 노력한다. 규제항은 훈련하는 동안에만 비용 함수에 추가되고, 훈련이 끝나면 규제가 없는 성능 지표로 모델의 성능을 평가한다.  

하이퍼파라미터 $\alpha$는 모델을 얼마나 많이 규제할지 조절하는 것으로, 0이면 릿지 회귀는 선형 회귀와 같아진다. $\alpha$가 아주 크면 모든 가중치가 거의 0에 가까워지고 결국 데이터의 평균을 지나는 수평선이 된다. 아래는 릿지 회귀의 비용 함수이다. 

$J(\theta)=MSE(\theta) + \alpha$$1\over{2}$$\sum_{i=1}^n\theta^{2}_i$

규제항에 $1\over{2}$을 곱한 것은 미분 결과를 간단하게 만들기 위함이다. 또 편향 $\theta_0은 규제되지 않아 합 기호가 i=0이 아닌 i=1에서 시작한다. 가중치 벡터를 w로 정의하면 결국 규제항은 w의 l2 노름이다. 한 가지 주의할 점은 릿지 회귀는 입력 특성의 스케일에 민감하기 때문에 데이터 스케일을 우선적으로 맞춰야 한다. 이는 규제가 있는 모델 대부분이 그렇다.

결과적으로 $\alpha$를 증가시킬수록 직선에 가까워지고, 이는 모델의 분산이 줄어들지만 편향은 커지는 것이다. 그리고 l2 노름의 특성상, 변수 대부분은 0보다는 큰 값에 분포하게 된다.

아래는 사이킷런에서 정규방정식을 사용한 릿지 회귀를 적용하는 예이다(안드레 루이 숄레스키가 발견한 행렬 분해를 사용하여 기존의 정규방정식을 변형한 버전이다).

from sklearn.linear_model import Ridge
ridge_reg = Ridge(alpha=1, solver="cholesky")
ridge_reg.fit(X, y)
ridge_reg.predict([[1.5]])

# 확률적 경사 하강법을 사용한 버전
sgd_reg = SGDRegressor(penalty="l2")
sgd_reg.fit(X, y.ravel())
sgd_reg.predict([[1.5]])

엘라스틱넷(Elastic Net)

엘라스틱넷은 릿지 회귀와 라쏘 회귀를 절충한 모델이다. 규제항은 릿지와 회귀의 규제항을 단순히 더해서 사용하며, 혼합 정도는 혼합 비율 r을 사용해 조절한다. r=0이면 릿지 회귀와 같아지고, r=1이면 라쏘 회귀와 같다.

$J(\theta)=MSE(\theta)+r\alpha\sum_{i=1}^n|\theta_i|$+$1-r\over{2}$$\alpha\sum_{i=1}^n\theta_i^2$

보통 규제가 없는 선형 회귀 모델보단 적어도 규제가 있는 것을 사용한다. 이때 릿지가 기본이 되지만 쓰이는 특성이 몇 개뿐이라고 의심되면 라쏘나 엘라스틱넷이 낫다. 라쏘나 엘라스틱넷은 불필요한 특성의 가중치를 0으로 만들어주기 때문이다. 특성 수가 훈련 샘플 수보다 많거나 특성 몇 개가 강하게 연관되어 있을 때는 엘라스틱넷이 낫다. ^[각주:1]

다음은 사이킷런의 ElasticNet을 사용한 예제이다.(l1_ratio가 혼합 비율 r이다)

from sklearn.linear_model import ElasticNet
elastic_net = ElasticNet(alpha=0.1, l1_ratio=0.5)
elastic_net.fit(X, y)
elastic_net.predict([[1.5]])

조기 종료(Early Stopping)

경사 하강법과 같은 반복적인 학습 알고리즘을 규제하는 색다른 방식은 검증 에러가 최솟값에 도달하면 바로 훈련을 중지시키는 것이다. 에포크가 진행되면서 알고리즘은 점차 학습되어 훈련 세트에 대한 예측 에러(RMSE)와 검증 세트에 대한 예측 에러가 줄어든다. 하지만 감소하던 검증 에러가 잠시 멈추었다 상승하는 지점이 있는데, 모델이 훈련 데이터에 과대적합되기 시작했다는 뜻이다. 조기 종료는 검증 에러가 최소에 도달하는 즉시 훈련을 멈추는 것이다. 

다음은 조기 종료를 위한 코드이다.

from sklearn.base import clone

poly_scaler = Pipeline([
        ("poly_features", PolynomialFeatures(degree=90, include_bias=False)),
        ("std_scaler", StandardScaler())
    ])
X_train_poly_scaled = poly_scaler.fit_transform(X_train)
X_val_poly_scaled = poly_scaler.transform(X_val)

# warm_start=True: fit() 메서드가 호출될 때 처음부터 다시 시작하지 않고 이전 모델 파라미터에서 훈련 이어감
sgd_reg = SGDRegressor(max_iter=1, tol=-np.infty, warm_start=True, 
                       penalty=None, learning_rate="constant", eta0=0.0005)
                       
minimum_val_error = float("inf")
best_epoch = None
best_model = None
for epoch in range(1000):
    sgd_reg.fit(X_train_poly_scaled, y_train)
    y_val_predict = sgd_reg.predict(X_val_poly_scaled)
    val_error = mean_squared_error(y_val, y_val_predict)
    if val_error < minimum_val_error:
        minimum_val_error = val_error
        best_epoch = epoch
        best_model = clone(sgd_reg)

로지스틱 회귀(Logistic Regression)

로지스틱 회귀는 샘플이 특정 클래스에 속할 확률을 추정하는 데 널리 사용된다. 추정 확률이 50%가 넘으면 모델은 그 샘플이 해당 클래스에 속한다고 예측한다.(레이블이 '1'인 양성 클래스) 50%가 넘지 않으면 클래스에 속하지 않는다고 예측한다.(레이블이 '0'인 음성 클래스)

로지스틱 회귀 모델의 확률 추정 식은 다음과 같다.

$\hat{p}=h_\theta(x)=\sigma(\theta^Tx)$

로지스틱($\sigma$로 표시)은 0과 1 사이의 값을 출력하는 시그모이드(sigmoid) 함수이며 S자 형태의 곡선이다. 시그모이드 함수는 $\sigma(t)$=$1\over{1+exp(-t)}$로 표현된다. 시그모이드 함수는 t=0일 때, 0.5가 되고, t가 0보다 크면 > 0.5 , 0보다 작으면 < 0.5가 된다.

로지스틱 회귀 모델이 샘플 x가 양성 클래스에 속할 확률 $\hat{p}=h_\theta(x)$를 추정하면 이에 대한 예측 $\hat{y}$을 쉽게 구할 수 있다. 로지스틱 회귀 모델 예측은 다음과 같은 식을 사용한다.

$\hat{y}=0$($\hat{p}<0.5$일 때), $1(\hat{p}\ge0$일 때)

따라서 로지스틱 회귀 모델은 $\theta^Tx$가 양수일 때 1이라고 예측하고, 음수일 때 0이라고 예측한다. 훈련의 목적은 양성 샘플(y=1)에 대해선 높은 확률을 추정하고 음성 샘플(y=0)에 대해선 낮은 확률을 추정하는 모델의 파라미터 벡터 $\theta$를 찾는 것이다. 로지스틱 회귀 비용 함수는 아래와 같다.

$J(\theta)=$$-1\over{m}$$\sum_{i=1}^m$$[y^{(i)}log(\hat{p}^{(i)})+(1-y^{(i)})log(1-\hat{p}^{(i)})]$

이 비용 함수는 최솟값을 계산하는 정규방정식 같은 해가 없다. 하지만 볼록(convex) 함수이므로 경사 하강법이 전역 최솟값을 찾는 것을 보장한다. 이때 필요한 편도함수는 아래와 같다. 비용 함수의 j번째 모델 파라미터 $\theta_j$에 대해 편미분을 하면 

$d\over{d\theta_j}$=$1\over{m}$$sum_{i=1}^m$$(\sigma(\theta^Tx^{(i)})-y^{(i)})x_j^{(i)}$

로 각 샘플에 대해 예측 오차를 계산하고 j번째 특성값을 곱해서 모든 훈련 샘플에 대해 평균을 낸다. 

결정 경계(Decision Boundary)

다른 클래스로 분류한 그 기준값을 찾아보기 위해 붓꽃 데이터를 사용해 코드를 구현해본다.

from sklearn import datasets
import numpy as np
import matplotlib.pyplot as plt

from sklearn import datasets
iris = datasets.load_iris()
list(iris.keys())
X = iris['data'][:, 3:] # 꽃잎 너비
y = (iris['target'] == 2).astype(np.int) # Iris-Virginica면 1, 그렇지 않으면 0

# 로지스틱 회귀 모델 훈련
from sklearn.linear_model import LogisticRegression

log_reg = LogisticRegression()
log_reg.fit(X, y)

X_new = np.linspace(0, 3, 1000).reshape(-1, 1) 
y_proba = log_reg.predict_proba(X_new)
# 샘플 1에 속할 확률을 나열한 리스트 중 가장 첫번째로(=[0]) 샘플 1에 속하게 된 확률 값이 결정 경계
decision_boundary = X_new[y_proba[:, 1] >= 0.5][0] 
 
plt.figure(figsize=(8, 3))

# 예측값이 아닌 실제값 X,y 표시
plt.plot(X[y==0], y[y==0], "bs")
plt.plot(X[y==1], y[y==1], "g^")

# 선형 그래프로 앞의 [,]: x축 데이터 / 뒤의 [,]: y축 데이터
plt.plot([decision_boundary, decision_boundary], [-1, 2], "k:", linewidth=2) 

plt.plot(X_new, y_proba[:, 1], "g--", label="Iris virginica")
plt.plot(X_new, y_proba[:, 0], "b--", label="Not Iris virginica")

plt.text(decision_boundary+0.02, 0.15, "Decision boundary", fontsize=14, color="k", ha="center")

# 화살표 표시
plt.arrow(decision_boundary, 0.08, -0.3, 0, head_width=0.05, head_length=0.1, fc='b', ec='b') 
plt.arrow(decision_boundary, 0.92, 0.3, 0, head_width=0.05, head_length=0.1, fc='g', ec='g')

plt.legend(loc="center left", fontsize=14); plt.xlabel("petal width", fontsize=14); plt.ylabel("proba", fontsize=14)
plt.axis([0, 3, -0.02, 1.02])
plt.show()

-------------------------------------------------참고-------------------------------------------------

# linspace(n, m, k): n과 m 사이에서 k개의 숫자를 리스트 형태로 출력, 각 숫자는 서로 동일한 차이값을 가짐 
# reshpae(-1, 1): 앞의 -1은 k개의 column을 row로 변환함을 의미, 뒤의 1은 column의 개수를 정함, column 수에 따라 앞의 -1에 따른 row 수는 상대적으로 변함
(np.linspace(0, 1, 6).reshape(-1, 1)).shape()
>>> (6, 1)
(np.linspace(0, 1, 6).reshape(-1, 2)).shape()
>>> (3, 2)

소프트맥스 회귀(Softmax Regression)

다중 클래스를 지원하는 모델로 샘플 x가 주어지면 각 클래스 k에 대한 점수 $s_k(x)$를 계산하고, 그 점수에 소프트맥스 함수를 적용하여 각 클래스의 확률을 추정한다. 클래스 k에 대한 소프트맥스 점수를 구하는 식은 아래와 같다.

$s_k(\mathbf{x})=(\theta^{(k)})^T\mathbf{x}$

각 클래스는 자신만의 파라미터 벡터 $\theta^{(k)}$가 있고, 이들은 파라미터 행렬 $\Theta$에 행으로 저장된다. 소프트맥스 함수는 아래와 같다.

$\hat{p}_k=\sigma(s(\mathbf{x}))_k$=$exp(s_k(\mathbf{x}))\over{\sum_{j=1}^Kexp(s_j(\mathbf{x}))}$

K는 클래스 수, s(x)는 샘플 x에 대한 각 클래스의 점수를 담은 벡터, $\sigma(s(\mathbf{x}))_k$는 샘플 x에 대한 각 클래스의 점수가 주어졌을 때 이 샘플이 클래스 k에 속할 추정 확률이다.

이렇게 추정 확률을 구한 뒤에는 추정 확률이 가장 높은 클래스(즉 가장 높은 점수를 가진 클래스)를 선택한다. 소프트맥스 회귀 분류기는 한 번에 하나의 클래스만 예측하는 것으로 다중 클래스(multiclass)지만 다중 출력(multioutput)은 아니다. 

소프트맥스 회귀의 훈련 방법은 크로스 엔트로피(cross entropy) 비용 함수를 사용하는 것이다. 크로스 엔트로피는 추정된 클래스의 확률이 타깃 클래스에 얼마나 잘 맞는지 측정하는 용도로 종종 사용된다. 

$J(\Theta)$=$-1\over{m}$$\sum_{i=1}^m\sum_{k=1}^K$$y_k^{(i)}log(\hat{p}_k^{(i)})$

$y_k^{(i)}$는 i번째 샘플이 클래스 k에 속할 타깃 확률이다. 일반적으로 샘플이 클래스에 속하는지 아닌지에 따라 1 또는 0이 된다. 

이 비용 함수의 $\theta^{(k)}$에 대한 그레이디언트 벡터는 다음과 같다. $\Delta_\theta^{(k)} J(\Theta)$=$1\over{m}$$\sum_{i=1}^m$$(\hat{p}_k^{(i)}-y_k^{(i)})x^{(i)}$

사이킷런의 LogisticRegression을 사용하고, multi_class 매개변수를 "multinomial"로 바꾸면 소프트맥스 회귀를 사용할 수 있다.

X = iris["data"][:, (2, 3)] # 꽃잎 길이, 꽃잎 너비
y = iris["target"]

# 소프트맥스 회귀를 위해선 solver 매개변수에 "sag", "saga", "newton-cg"와 같은 알고리즘 지정
# 하이퍼파라미터 C를 사용하여 조절할 수 있는 l2 규제 적용
softmax_reg = LogisticRegression(multi_class="multinomial", solver="lbfgs", C=10)
softmax_reg.fit(X, y)

softmax_reg.predict([[5, 2]])
softmax_reg.predict_proba([[5, 2]])

1. 라쏘는 특성 수가 샘플 수(n)보다 많으면 최대 n개의 특성을 선택한다. 또한 여러 특성이 강하게 연관되어 있으면 이들 중 임의의 특성 하나를 선택한다. [본문으로]