<jupyter_start><jupyter_text>4장 모델 훈련  핵심 주제 선형 회귀 모델의 훈련 작동과정의 상세 사항을 알아본다. 1. 선형 회귀    1. 정규방정식    1. 계산 복잡도1. 경사 하강법(GD)    1. 배치 경사 하강법    1. 확률적 경사 하강법(SGD)    1. 미니배치 경사 하강법1. 다항 회귀1. 학습 곡선    1. 과소적합    1. 과대적합    1. 편향/분산 트레이드오프1. 규제 선형 모델    1. 릿지 회귀    1. 라쏘 회귀    1. 엘라스틱넷    1. 조기 종료1. 로직스틱 회귀와 소프트맥스 회귀    1. 이진 분류: 로지스틱 회귀        1. 시그모이드 함수        1. 로그 손실 비용 함수        1. 결정 경계    1. 다중 클래스 분류: 소프트맥스 회귀        1. 소프트맥스 함수        1. 크로스 엔트로피 비용 함수  파이썬 기본 설정 * 파이선 버전 3.5 이상* 사이킷런 버전 0.20 이상* 모듈 불러오기    * numpy    * matplotlib    * matplotlib.pyplot* 그림 저장 위치 지정<jupyter_code># 파이썬 ≥3.5 필수
import sys
assert sys.version_info >= (3, 5)

# 사이킷런 ≥0.20 필수
import sklearn
assert sklearn.__version__ >= "0.20"

# 기본 모듈 임포트
import numpy as np
import os

# 랜덤 시드 설정. 동일한 결과 재생.
np.random.seed(42)

# 그래프 관련
%matplotlib inline
import matplotlib as mpl
import matplotlib.pyplot as plt
mpl.rc('axes', labelsize=14)
mpl.rc('xtick', labelsize=12)
mpl.rc('ytick', labelsize=12)

# 그림 저장 위치 지정
PROJECT_ROOT_DIR = "."
CHAPTER_ID = "training_linear_models"
IMAGES_PATH = os.path.join(PROJECT_ROOT_DIR, "images", CHAPTER_ID)
os.makedirs(IMAGES_PATH, exist_ok=True)

def save_fig(fig_id, tight_layout=True, fig_extension="png", resolution=300):
    path = os.path.join(IMAGES_PATH, fig_id + "." + fig_extension)
    print("그림 저장:", fig_id)
    if tight_layout:
        plt.tight_layout()
    plt.savefig(path, format=fig_extension, dpi=resolution)

# 불필요한 경고를 무시합니다 (사이파이 이슈 #5998 참조)
import warnings
warnings.filterwarnings(action="ignore", message="^internal gelsd")<jupyter_output><empty_output><jupyter_text>핵심 1: 선형 회귀 한 개의 특성 $x_1$을 사용하는 $i$번째 훈련 샘플에 대한 예측값을 구하는 선형 회귀 모델은 다음과 같다. $$\hat y^{(i)} = \theta_0 + \theta_1\, x_1^{(i)}$$$n\ge 1$ 개의 특성을 사용하는 $i$번째 훈련 샘플에 대한 예측값을 구하는 일반적인 선형 회귀 모델은 다음과 같다. $$\hat y^{(i)} = \theta_0 + \theta_1\, x_1^{(i)} + \cdots + \theta_n\, x_n^{(i)}$$사용된 기호들의 의미는 다음과 같다.* $\hat y^{(i)}$: $i$ 번째 훈련 샘플에 대한 예측값(타깃)* $n$: 특성 개수* $x_k^{(i)}$: $i$ 번째 훈련 샘플의 $k$ 번째 특성값* $\theta_i$    * $i > 0$: $i$ 번째 특성에 대한 **가중치** 파라미터    * $\theta_0$: **편향**을 의미함. **주의:** * 책에서 사용되는 인덱스가 약간 다르다.     하지만 여기서 사용하는 기호들이 보다 정확하며,     앞으로 계속해서 일반화될 수 있도록 준비된 기호와 인덱스 사용법이다. * 책에서 사용되는 수식은 모두 $\mathbf{x}$가 (n,1) 모양의 벡터로 가정한다.    하지만 실제 데이터셋에서는 각 샘플이 하나의 행으로 주어진다.     즉, (1, n) 모양의 벡터로 다루어진다.    따라서 행렬 연산의 통일성을 유지하기 위해     책에서 등장하는 많은 수식에서 사용되는 $\mathbf{x}$ 대신에 여기서는 $\mathbf{x}_b^T$를 사용한다.     여기서 $\mathbf x_b$는 편향을 다루기 위해 0번 인덱스에 1을 추가한 벡터이다.        실제로 161쪽에서 정규방정식을 이용하여 최적의 가중치 벡터를 계산할 때 $\mathbf{x}_b^T$가 사용된다.    사용되는 기호들의 모양에 대해 바로 이어서 구체적으로 설명한다.  행렬 연산과 선형 회귀 행렬 연산을 이용하면 전체 $m$개의 훈련 샘플을 대상으로 하는 선형 회귀 모델을 아래와 같이하나의 식으로 표현할 수 있다. $$\hat {\mathbf y}^T = h_\theta (\mathbf X) = \theta^{T}\, \mathbf X_b^T$$* $\hat {\mathbf y}$는 모든 훈련 샘플에 대한 예측값으로 구성된 $m\times 1$ 행렬이다.\begin{bmatrix}    y^{(0)} \\    y^{(1)} \\    \vdots \\    y^{(m-1)}\end{bmatrix}* $\theta$는 편형과 가중치로 이루어진 $(n+1)\times 1$ 행렬이다.\begin{bmatrix}    \theta_0 \\    \theta_1 \\    \vdots \\    \theta_n\end{bmatrix}    * $\mathbf X$는 모든 훈련 샘플에 대한 예측값으로 구성된 $m\times n$ 행렬이다.$$\begin{bmatrix}    \mathbf x^{(0)} \\    \mathbf x^{(1)} \\    \vdots \\    \mathbf x^{(m-1)}\end{bmatrix}=\begin{bmatrix}    x_1^{(0)} & \cdots & x_n^{(0)} \\    x_1^{(1)} & \cdots & x_n^{(1)} \\    \vdots \\    x_1^{(m-1)} & \cdots & x_n^{(m-1)}\end{bmatrix}$$* $\mathbf X_b$는 $\mathbf X$에 1로 이루어진 열을 맨 앞에 추가한 $m\times (n+1)$ 행렬이다.    이유는 각 샘플의 예측값을 계산할 때     편향($\theta_0$)을 더해주어야 하기 때문이다.    \begin{bmatrix}    1 & x_1^{(0)} & \cdots & x_n^{(0)} \\    1 & x_1^{(1)} & \cdots & x_n^{(1)} \\    \vdots \\    1 & x_1^{(m-1)} & \cdots & x_n^{(m-1)}\end{bmatrix}      평균 제곱 오차(MSE) 행렬을 이용한 평균 제곱 오차(MSE)의 정의는 다음과 같다.$$\mathrm{MSE}(\mathbf X, h_\theta) = \frac 1 m \sum_{i=0}^{m-1} (\theta^{T}\, \mathbf (\mathbf x_b^{(i)})^T - y^{(i)})^2$$단, $\mathbf{x}_b^{(i)} = \begin{bmatrix} 1 & x_1^{(i)} & \cdots & x_n^{(i)}\end{bmatrix}$ 이다.앞으로 $\mathrm{MSE}(\mathbf X, h_\theta)$ 대신에 좀 더 간략한 $\mathrm{MSE}(\theta)$를 사용한다.  넘파이 어레이 계산 행렬을 넘파이 어레이로 다루면 $\mathbf y, \hat{\mathbf y}, \theta, \mathbf X, \mathbf X_b$는 각각 아래 모양(shape)의 2차원 어레이로 구현된다.| 기호           |     수학 행렬의 모양 |  넘파이 어레이 모양||:-------------:|:---------------:|:------:|| $\mathbf y$, $\hat{\mathbf y}$  | $m \times 1$ | (m, 1) || $\theta$      | $(n+1)\times 1$ | (n+1, 1) || $\mathbf X$   | $m\times n$     | (m, n)   || $\mathbf X_b$ | $m\times (n+1)$ | (m, n+1) | 따라서 넘파이 함수인 `np.mean()`과 `np.square()` 함수를 이용하면 평균 제곱 오차를 아래와 같이 쉽게 구할 수 있다.$$\mathrm{MSE}(\mathbf X, h_\theta) = \texttt{np.mean}(\texttt{np.square}(\theta^{T}\, {\mathbf X}_b^T - {\mathbf y}))$$  정규방정식 선형 모델의 비용 함수를 최소화하는 $\theta$를 수학적으로 해결할 수 있다.선형회귀의 경우 아래 방정식이 만족되면 평균 제곱 오차(MSE)가 0이 되면서 당연히 최소값을 갖는다.$$\theta^{T}\, {\mathbf X}_b^T - {\mathbf y}^T = 0$$따라서 위 방정식을 만족시키는 가중치들의 리스트 $\hat\theta$를 찾으면 되는데,수학적으로 $\hat\theta$는 아래와 같다.$$\hat\theta = (\mathbf{X}_b^T \mathbf{X}_b)^{-1} \mathbf{X}_b^T \mathbf{y}$$위 방정식을 **정규방정식**(normal equation)이라 부른다.  증명 정규방정식을 유도하는 과정은 다음과 같다. 이부분은 이해가 어렵다면 건너 뛰어도 된다.먼저 $\hat\theta$에 의해 ${\hat\theta}^{T}\, {\mathbf X}_b^T = {\mathbf y}^T$가 성립한다고 가정하자.그러면 다음 등식이 성립한다.  $$\begin{align*}{\hat\theta}^{T}\, {\mathbf X}_b^T = {\mathbf y}^T \quad &\Longrightarrow\quad \mathbf X_b\, \hat\theta = \mathbf y\\&\Longrightarrow\quad \mathbf{X}_b^T\,\mathbf X_b\, \hat\theta = \mathbf{X}_b^T\,\mathbf y\\&\Longrightarrow\quad \hat\theta = (\mathbf{X}_b^T\,\mathbf X_b)^{-1}\, \mathbf{X}_b^T\,\mathbf y\\\end{align*}$$  예제 정규방정식을 이용하여 가중치를 계산하는 예제를 살펴보고자 한다.먼저 선형처럼 보이는 데이터를 인위적으로 생성한다.* 기본 데이터(`X`)    * (0,1) 사이의 값을 균등하게 무작위로 100개 선택하여 각각 두 배 해준 (100,1) 모양의 어레이 생성* 가우시안 잡음    * `np.random.randn(100,1)` 활용    * 정규분포(평균 0, 표준편차 1)를 따르는 100개의 숫자를 무작위로 선택한 후 (100,1) 모양의 어레이 생성    * 데이터에 잡읍을 추가하는 역할 수행<jupyter_code>import numpy as np

X = 2 * np.random.rand(100, 1)
y = 4 + 3 * X + np.random.randn(100, 1)<jupyter_output><empty_output><jupyter_text>`X`와 `y` 값들을 산포도로 확인하면 노이즈의 역할을 볼 수 있다.<jupyter_code>plt.plot(X, y, "b.")
plt.xlabel("$x_1$", fontsize=18)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.axis([0, 2, 0, 15])
save_fig("generated_data_plot")
plt.show()<jupyter_output>그림 저장: generated_data_plot<jupyter_text>`X`로 부터 `X_b`를 생성한다.앞서 설명한 $\mathbf X_b$를 생성하는 과정에 해당한다.<jupyter_code>X_b = np.c_[np.ones((100, 1)), X]  # 모든 샘플에 x0 = 1을 추가합니다.<jupyter_output><empty_output><jupyter_text>이제 정규방정식을 이용하여 최적의 가중치값들을 계산한다.* `np.linalg.inv()` 함수: 역행렬 계산* `T` 속성: 전치행렬을 가리킴.* `dot()` 메서드: 2차원 어레이에 대해서는 행렬의 곱에 해당하는 계산 수행.<jupyter_code>theta_best = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)

theta_best<jupyter_output><empty_output><jupyter_text>$\hat\theta_0 = 4.215$와 $\hat\theta_1 = 2.770$이 계산되었다.애초에 사용된 4, 3과는 약간 다르지만, 이것은 **가우시안 잡음** 때문에 어쩔 수 없다. 구해진 $\hat\theta$를 이용하여 예측을 할 수 있다.모든 샘플에 `x0=1`을 추가하여 `X_new_b`를 생성하는 이유는 앞서 설명한 그대로이다.<jupyter_code>X_new = np.array([[0], [2]])
X_new_b = np.c_[np.ones((2, 1)), X_new]  # 모든 샘플에 x0 = 1을 추가합니다.
y_predict = X_new_b.dot(theta_best)
y_predict<jupyter_output><empty_output><jupyter_text>예측된 선형회귀 그래프는 아래와 같다.<jupyter_code>plt.plot(X_new, y_predict, "r-")
plt.plot(X, y, "b.")
plt.axis([0, 2, 0, 15])
plt.show()<jupyter_output><empty_output><jupyter_text>책에 있는 그림은 범례와 축 레이블이 있는 그래프이다.<jupyter_code>plt.plot(X_new, y_predict, "r-", linewidth=2, label="Predictions")
plt.plot(X, y, "b.")
plt.xlabel("$x_1$", fontsize=18)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.legend(loc="upper left", fontsize=14)
plt.axis([0, 2, 0, 15])
save_fig("linear_model_predictions_plot")
plt.show()<jupyter_output>그림 저장: linear_model_predictions_plot<jupyter_text>사이키런의 선형 회귀 사이킷런의 `LinearRegression` 클래스가 선형 회귀 모델의 근사치를 계산해준다.사실 이 예제의 경우는 동일하다. 하지만 이는 특별한 경우이다.<jupyter_code>from sklearn.linear_model import LinearRegression

lin_reg = LinearRegression()
lin_reg.fit(X, y)
lin_reg.intercept_, lin_reg.coef_<jupyter_output><empty_output><jupyter_text>예측값도 동일하다. (물론 특별한 경우임.)<jupyter_code>lin_reg.predict(X_new)<jupyter_output><empty_output><jupyter_text>하지만 `LinearRegression.fit()` 메서드가 정규방정식을 이용하지는 않는다.이유는 행렬 계산 비용이 매우 비싸기 때문이다.반면에 `LinearRegression.fit()` 메서드는 **특잇값 분해**(singular value decomposition, SVD)라는표준 행렬 분해 기법을 사용하여 아래 성질을 만족하는 행렬 $\mathbf X^+$를 계산한다. $$\hat\theta \cong \mathbf X^+ \, \mathbf y$$$\mathbf X^+$는 $\mathbf X$의 **유사역행렬**(pseudoinverse)라고 불린다.**주의:** 위 예제의 경우 특별하게도 $\mathbf X^+ = \mathbf X^{-1}$이 성립한다. 하지만 일반적으로 $\mathbf X$의 역행렬은 존재하지 않는다. 반면에 $\mathbf X^+$는 항상 존재한다.  유사역행렬 계산 `LinearRegression` 클래스는 `scipy.linalg.lstsq()` 함수("least squares"의 약자)를 사용하는데,이 함수는 유사역행렬을 계산하여 이용한다.<jupyter_code>theta_best_svd, residuals, rank, s = np.linalg.lstsq(X_b, y, rcond=1e-6)
theta_best_svd<jupyter_output><empty_output><jupyter_text>유사역행렬 계산은 `np.linalg.pinv()` 함수를 이용한다.<jupyter_code>np.linalg.pinv(X_b).dot(y)<jupyter_output><empty_output><jupyter_text>계산 복잡도 정규방정식에 사용되는 $(n+1)\times (n+1)$ 모양의 $\mathbf{X}_b^T \mathbf{X}_b$의 역행렬을 계산해야 한다.그런데 역행렬 계산 복잡도는 $O(n^{2.4}) \sim O(n^3)$ 정도로 꽤 비용이 많이 든다.즉, 특성 수가 두 배로 늘면 계산 시간은 $2^{2.4}=5.3 \sim 2^3 = 8$ 배 정도 증가한다. 반면에 `LinearRegression` 에서 사용하는 SVD 방식의 시간 복잡도는 $O(n^2)$으로 상대적으로 저렴하다. 하지만 그래도 특성 수가 크면 역시 실용성이 떨어진다. 이렇게 특성 수가 매우 많거나 훈련 세트의 크기가 너무 큰 경우 하나의 행렬로 지정한 후 매우 큰 수의 계산 과정을 실행하는 것은 시간 뿐만 아니라 메모리 관련 문제도 발생한다.따라서 이럴 때는 정규방정식 방법 대신에 이어서 설명하는 경사 하강법을 사용해야 한다. **주의:** * 훈련 세트의 샘플 수에 대해서 계산복잡도가 $O(m)$이다.* 또한 이미 선형 모델 학습이 되었다면 예측을 위한 시간 복잡도는 $O(n)$, 즉, 특성 수에 선형적이다.   핵심 2: 경사 하강법 **경사 하강법**(gradient descent, GD)의 기본 아이디어는 훈련 세트를 이용한학습 과정을 통해 가중치 등의 파라미터를 조금씩 조정해 가는 것이며,사용되는 주요 개념은 다음과 같다.1. **최적 학습 모델**    * 경사 하강법의 목표는 **최적의 학습 모델**을 점진적으로 찾아가는 것이다.        1. **파라미터**    * 최적의 학습 모델을 찾으려면 학습 모델을 결정하는 최적의 파라미터를 찾아야 한다.    * 선형 회귀 모델의 경우, 파라미터는 편형과 가중치를 가리킨다.         이는 편향($\theta_0$)과 가중치($\theta_i$)가 결정되면 선형 회귀 모델이 결정되기 때문이다.     $$\hat {\mathbf y}^T = h_\theta (\mathbf X) = \theta^{T}\, \mathbf X_b^T$$        1. **비용 함수**    * **최적**의 학습 모델은 비용 함수가 최소가 되도록 하는 모델이다.    * 선형 회귀 모델의 경우 평균 제곱 오차(MSE)가 비용 함수이며, 파라미터 $\theta_i$에 의존한다.    $$\mathrm{MSE}(\mathbf X, h_\theta) = \frac 1 m \sum_{i=0}^{m-1} (\theta^{T}\, \mathbf (\mathbf x_b^{(i)})^T - y^{(i)})^2$$    1. **전역 최솟값**    * 비용 함수가 갖는 최소값을 전역 최솟값이라 하며, 전역 최솟값이 되도록 하는 파라미터를 찾아야 한다.    * 선형 회귀 모델의 경우 평균 제곱 오차(MSE) 함수가 최소가 되도록 하는         파라미터 $\theta_i$를 찾아야 한다.    1. **그레이디언트 벡터**    * 그레이디언트(gradient) 벡터는 다변수 함수의 미분값을 가리킨다.         다변수 함수의 미분값은 하나의 값이 아닌 여러 개의 값으로 이루어진 벡터를 이룬다.    * 선형 회귀의 경우 $n+1$개의 파라미터($\theta_i$)가 평균 제곱 오차(MSE) 함수의 매개변수로 사용되며,        이렇게 두 개 이상의 매개변수를 사용하는 함수를 **다변수 함수**라 부른다.         $$\nabla_\theta \textrm{MSE}(\theta) =\begin{bmatrix}    \frac{\partial}{\partial \theta_0} \textrm{MSE}(\theta) \\    \frac{\partial}{\partial \theta_1} \textrm{MSE}(\theta) \\    \vdots \\    \frac{\partial}{\partial \theta_n} \textrm{MSE}(\theta)\end{bmatrix} =\frac{2}{m}\, \mathbf{X}_b^T\, (\mathbf{X}_b\, \theta - \mathbf y)$$        * 다변수 함수의 미분값은 편미분을 이용하며, 미분보다 다루기가 좀 더 어렵다.        하지만 그레이디언트의 엄밀한 정의 보다는 의미를 이해하는 것이 훨씬 중요하다.        * 비용 함수의 전역 최솟값은 그레이디언트가 0 벡터가 되도록 하는 지점에서 비용 함수가 갖는 값이다.        * 선형 회귀의 경우 $\nabla_\theta \textrm{MSE}(\theta)=0$이 되게하는             파라미터 $\theta$를 찾아야 한다.        * 위 수식에서 보듯이 수학적으로 답은 정해져 있다.             하지만 행렬의 곱셈은 앞서 설명한데로 비용이 많이 든다.            따라서 임의의 파라미터($\theta$)로 시작해서 최적의 $\theta$를 찾아가는 방식을 사용한다.    * $\nabla_\theta \textrm{MSE}(\theta)$에 담긴 내용        * 지점 $\theta$에서의 편미분 값들로 이루어진 벡터        * 벡터는 **방향**과 **크기**에 대한 정보를 갖고 있음.        * 비용 함수의 전역 최솟값은 그레이디언트가 가리키는 방향의 **반대쪽**으로 움직여야 가장 빠르게 접근할 수 있음.        * 예제: 언덕을 가장 빠르게 내려가는 길은 가장 가파르게 내려가는 경로를 택하는 것이다.            이때, 가장 가파르게 내려가는 방향은 위치한 지점에서의 그레이디언트 벡터가 가리키는 값의 반대 방향이다.            1. **학습률**    * 비용 함수의 파라미터를 임의의 값으로 지정한 후 학습과정에서 **반복적으로** 그레이디언트를 계산하여        계산 결과에 따라 파라미터 값을 조정해야 하며,         이 조정 과정을 비용 함수가 특정 값에 수렴할 때까지 반복한다.     * 파라미터를 조정할 때마다 조정되는 값의 크기를 제한해야 한다.        * 너무 크게 움직이면 비용 함수가 수렴하지 않을 수 있으며,        * 너무 자게 움직이면 비용 함수가 전역 최소값에 너무 느리게 수렴하거나,            아니면 전역 최솟값이 아닌 지역 최솟값에 수렴할 수도 있다.             * 지역 최솟값: 특정 영역 내에서 최솟값을 가리킴.    * 파라미터 조정의 크기를 정하는 값이 **학습률**(learning rate)이며,         보통 그리스어 알파벳인 $\eta$(에타)로 표기한다.        학습률 크기에 따른 영향은 책 165쪽에서 167쪽 사이에 위치한 그림 4-3 ~ 그림 4.6을 참조한다.    * 앞서 설명한대로 파라미터 조정은 해당 지점의 그레이디언트 벡터가 가리키는 반대 방향으로 움직여야 한다.        예를 들어, 선형 회귀의 경우 새로운 파라미터 $\theta^{(\text{next step})}$을         학습률을 이용하여 다음과 같이 계산한다.        즉, 그레이디어트 벡터가 가리키는 반대 방향으로 학습률이 지정하는 만큼         파라미터를 조정한다.        $$\theta^{(\text{new})} = \theta^{(\text{old})}\, -\, \eta\cdot \nabla_\theta \textrm{MSE}(\theta^{(\text{old})})$$                    * 선형 회귀의 경우 학습률이 너무 크지만 않으면 언젠가는 전역 최솟값에 수렴한다는 것이 보장된다.        다만 학습률이 너무 크면 수렴하지 않고 무한으로 발산할 수도 있다.        * 이유는 비용함수가 항상 볼록 함수이기 때문이다.         * 볼록 함수는 지역 최솟값을 갖지 않는다.    * 예제: 언덕을 내려갈 때 보폭을 크게 해야 빠르게 내려온다. 반면에 보폭을 적게 하면 많은 내려오는 데에        보다 많은 시간이 걸린다. 학습률이 바로 보폭의 크기를 결정한다고 이해하면 된다.    1. **하이퍼파라미터**    * 학습률은 모델 학습 과정에서 변하는 값이 아니라, 학습 시작 전에 지정되어야 하는 값이다.    * 이런 값들을 **하이퍼파라미터**(hyperparameter)라고 한다.    * 학습률 이외에 대표적인 하이퍼파라미터는 전체 훈련 세트를 몇 번 학습할지를 결정하는 반복횟수이며,        **에포크**(epoch)라고 부른다.     * 학습률과 에포크는 모든 모델 학습에 기본으로 사용되는 하이퍼파라미터이다.    * **주의**: 아래에서 설명할 확률적 경사 하강법의 경우 학습 과정에서        지정된 비율로 학습률을 점점 작게 만드는 기법을 사용하기도 한다.         이를 **기본 학습 스케줄**이라 부른다.        하지만 언급한 대로 학습률을 작게하는 비율을 학습 이전에 미리 지정해야 한다는         의미에서 학습률을 하이퍼파라이터로 간주하는 것이 적절하다.    * 학습 모델에 따라 하이퍼파라미터가 추가될 수 있다.         예를 들어,몇 개의 훈련 샘플을 학습한 뒤에 그레이디언트 벡터를 계산할지를 결정하는        **배치 크기**(batch size)가 하이퍼파라미터로 종종 사용된다.        * 배치 경사 하강법            * 배치 크기: 훈련 세트의 크기            * 즉, 하나의 에포크가 지날 때마다 한 번 그레이디언트를 계산해서 파라미터를 조정한다.            * **주의**: 여기서 사용되는 '배치'의 의미가 '배치 크기'의 '배치'와 다른 의미이다.                책에서는 '배치 경사 하강법' 대신에 '전체 경사 하강법'이 적절하다고 언급하였다.        * 확률적 경사 하강법            * 배치 크기: 1            * 즉, 하나의 훈련 샘플을 학습할 때마다 그레이디언트를 계산해서 파라미터를 조정한다.        * 미니배치 경사 하강법            * 배치 크기: 2에서 수백 사이            * 최적 배치 크기: 경우에 따라 다름. 여러 논문이 32 이하 추천. **참조사항*** **스텝**(step)    * 지정된 배치 크기의 샘플을 학습한 후에 파라미터를 조정하는 단계를 가리킨다.     * 예를 들어, 훈련 세트의 크기가 1,000이고 배치 크기가 10이면,         하나의 에포크 기간동안 총 100번의 스텝이 지나간단.         즉, 하나의 에포크 기간동안 파라미터 조정이 100번 발생한다.     * 배치 경사하강법의 경우 스텝의 크기는 1이며        확률적 경사 하강법의 경우 스텝의 크기는 전체 훈련 세트의 크기이다.    * 스텝 크기가 하이퍼파라미터로 주어지는 경우도 있다.        그러면 배치 크기는 총 훈련 데이터 수를 스텝 크기로 나눈 값이 된다.* 특성 스케일링의 중요성    * 사용되는 모든 특성을 동일한 크기로 스케일링을 한 후에 경사 하강법을 적용해야 한다.     * 스케일링을 통해 모든 특성의 크기를 통일한 후의 학습이 보다 빠르게 이루어진다.* **차원**(dimension)과 파라미터 공간(parameter space)    * 사용되는 특성 수가 많아지면, 비용함수에 사용되는 파라미터 수도 그만큼 늘어난다.    * 즉, 조사해야 하는 파라미터들의 공간 또한 커지며, 따라서 검색도 그만큼 어려워진다.    * 또한 행렬 연산 등도 제곱 배 이상($O(n^2)$)으로 보다 많은 시간을 요구한다.    * 특성 수를 **차원**이라 부르며, 나중에 차원 축소 기법을 활용하여         보다 효율적인 모델 학습을 지원하는 방법을 알아볼 것이다.   배치 경사 하강법 앞서 설명한 대로 배치 크기에 따라 학습 과정 동안 파라미터 조정 횟수 및 시기가 달라진다.배치 경사 하강법은 전체 훈련 샘플을 모두 학습할 때마다 평균 제곱 오차(MSE)를 계산하고 파라미터를 조정한다. 즉, 하나의 에포크 마다 파라미터를 조정한다. 사이킷런은 배치 경사 하강법을 활용한 선형 회귀는 기본으로 지원하지는 않는다. 하지만 다음과 같이 쉽게 구현해서 사용해볼 수 있다. **주의사항*** 아래 예제는 아래 등식을 활용한다.$$\nabla_\theta \textrm{MSE}(\theta) =\frac{2}{m}\, \mathbf{X}_b^T\, (\mathbf{X}_b\, \theta - \mathbf y)$$        * 책 176쪽, 표 4-1에서 사이킷런의 SGDRegressor가 배치 경사 하강법을 지원한다고 명시되어 있지만    사실이 아니다.  아래 코드에서 사용된 하이퍼 파라미터는 다음과 같다.* `eta`: 학습률 0.1* `n_iterations`: 에포크 수 1000사용하는 데이터는 앞서 가우시안 잡음을 넣어 생성한 100개의 데이터이다.<jupyter_code>eta = 0.1  # 학습률
n_iterations = 1000 # 에포크 수

# 무작위로 생성한 100개의 샘플 수
m = 100

# theta0, theta1 파라미터 임의로 초기값 지정
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<jupyter_output><empty_output><jupyter_text>1000번의 에포크를 학습한 결과는 다음과 같으며, 앞서 정규방정식을 이용한 결과와 동일하다. 하지만 이는 우연에 불과하며 일반적으로는 오차가 발생한다.<jupyter_code>theta<jupyter_output><empty_output><jupyter_text>새로운 데이터에 대한 예측도 물론 동일하게 나온다.<jupyter_code>X_new_b.dot(theta)<jupyter_output><empty_output><jupyter_text>학습률을 달리할 때 발생하는 현상을 아래 코드로 확인할 수 있다.<jupyter_code># theta 기록 장치
theta_path_bgd = []

# theta가 조정될 때마다 학습모델 그래프를 그린다.
# theta_path가 옵션으로 사용될 경우 theta 기록 장치 활용
def plot_gradient_descent(theta, eta, theta_path=None):
    m = len(X_b)
    plt.plot(X, y, "b.")
    n_iterations = 1000
    for iteration in range(n_iterations):
        if iteration < 10:
            y_predict = X_new_b.dot(theta)
            style = "b-" if iteration > 0 else "r--"
            plt.plot(X_new, y_predict, style)
        gradients = 2/m * X_b.T.dot(X_b.dot(theta) - y)
        theta = theta - eta * gradients
        if theta_path is not None:
            theta_path.append(theta)
    plt.xlabel("$x_1$", fontsize=18)
    plt.axis([0, 2, 0, 15])
    plt.title(r"$\eta = {}$".format(eta), fontsize=16)<jupyter_output><empty_output><jupyter_text>사용된 학습률은 다음과 같다.* `eta=0.02`: 최적 선형 모델에 느리게 수렴* `eta=0.1`: 최적 선형 모델에 적절한 속도로 수렴* `eta=0.5`: 최적 선형 모델에 수렴하지 못하고 발산파라미터의 수렴 과정을 추적하기 위해 `eta=0.1`인 경우 `theta`가 조정될 때 마다 `theta_path_bgd` 리스트에 저장한다.잠시 뒤에 학습 방법에 따라 파라미터의 수렴 과정이 어떻게 다른지 설명할 때 사용될 것이다.<jupyter_code>np.random.seed(42)
theta = np.random.randn(2,1)  # random initialization

plt.figure(figsize=(10,4))
plt.subplot(131); plot_gradient_descent(theta, eta=0.02)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.subplot(132); plot_gradient_descent(theta, eta=0.1, theta_path=theta_path_bgd)
plt.subplot(133); plot_gradient_descent(theta, eta=0.5)

save_fig("gradient_descent_plot")
plt.show()<jupyter_output>그림 저장: gradient_descent_plot<jupyter_text>에포크와 허용오차 적절한 학습률은 사용되는 모델과 훈련 세트의 데이터에 따라 달라지기 때문에 정해진 기준이 없다.보통 그리드 탐색 등을 사용하여 적절한 학습률을 찾는다.하지만 예를 들어 그리드 탐색은 많은 학습을 시도하기 때문에 에포크(반복횟수)를 제한해야 한다.하지만 학습이 얼마나 잘 되는가를 미리 알 수 없기 때문에, 에포크(반복횟수)를 얼마로 제한해야 할지 알 수 없다.따라서 에포크는 보통 크게 설정하고 대신에 **허용오차**(tolerance, 보통 tol로 줄여서 사용)를활용한다.즉, 그레이디언트 벡터의 크기(노름, norm)가 허용오차 보다 작아지면 학습을 멈추도록 한다.이유는 파라미터가 최적점에 가까워지면 그레이디언트 벡터의 크기가 0에 가까워지기 때문이다.예를 들어, 포물선의 최소점에서의 접선의 기울기는 0이어야 하며, 최소점에 가까워질 수록 접선의 기울기가 0에 수렴한다.허용오차($\varepsilon$)와 에포크(반복횟수)는 서로 반비례의 관계이다.즉 에포크는 $(1/\varepsilon)$에 의존한다. 따라서 오차를 1/10로 줄이려 하면 에포크를 10배 늘려야 한다.  확률적 경사 하강법 배치 경사 하강법은 훈련 세트의 크면 그레이디언트를 계사한는 데에 많은 시간이 걸린다는 문제가 있다.또한 아주 많은 데이터를 저장해야 하는 메모리 문제도 발생할 수 있다. 확률적 경사 하강법은 훈련 샘플 하나를 학습할 때마다 평균 제곱 오차(MSE)를 계산하고 파라미터를 조정한다. 즉, 파라미터 조정 스텝이 샘플을 하나 학습할 때마다 이루어진다. 알고리즘 계산이 매우 빠르게 처리되어 학습 과정이 매우 빠르다.또한 하나의 샘플만 대상으로 하므로 매우 큰 훈련 세트로 처리할 수 있다. 반면에 파라미터가 조정되는 정도를 예상할 수 없고 요동을 많이 칠 수 있기 때문에 배치 경사 하강법보다 불안정하다.**확률적**이라고 부르는 이유가 여기에 있다. 확률적 경사 하강법의 장단점은 다음과 같다.* 매우 큰 훈련 세트를 다룰 수 있다. 하나의 샘플을 대상으로 학습이 이루어지기 때문이다.    * 예를 들어, 외부 메모리(out-of-core) 학습을 활용할 수 있다.    * 참조: 책 46쪽, 1장, 그림 1-14* 학습 과정이 매우 빠르다.* 학습 과정에서 파라미터의 동요가 심할 수 있다. * 파라미터 조정이 불안정 할 수 있기 때문에 오히려 지역 최솟값에 상대적으로 덜 민감하다.* 경우에 따라 전역 최솟값에 수렴하지 못하고 계속해서 요동칠 가능성도 높다.   학습 스케줄 요동치는 파라미터를 제어하기 위해 학습률을 학습 과정 동안천천히 줄어들게 만들 수 있다.이렇게 하면 파라미터가 요동치는 것을 어느 정도 제어할 수 있다.* 하지만 학습률이 너무 빨리 줄어들면, 지역 최솟값에 갇힐 수 있다.* 반면에 학습률이 너무 느리게 줄어들면 전역 최솟값에 제대로 수렴하지 못하고 주변만 맴돌 수도 있다.아래 코드는 **학습 스케줄**(learning schedule)을 활용한 확률적 경사 하강법을 직접 구현한 것이다. 사용되는 학습 스케줄 함수(`learning_schedule()`)는 에포크, 훈련 샘플 수, 학습되는 샘플의 인덱스에의존해서 학습이 이루어질 수록 작아지는 학습률을 계산한다.코드에 사용된 주요 요소들은 다음과 같다.* `theta_path_sgd`: 확률적 경사 하강법 학습 과정에서 변하는 파라미터 `theta` 값을 기록하기 위해 사용된다. * `np.random.randint()` 함수: 샘플 선택을 매번 무작위적으로 실행한다.<jupyter_code>theta_path_sgd = []
m = len(X_b)
np.random.seed(42)

# 에포크
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):
        
        # 학습 모델 그리기 (처음 20번만 그림)
        if epoch == 0 and i < 20:           
            y_predict = X_new_b.dot(theta) 
            style = "b-" if i > 0 else "r--"
            plt.plot(X_new, y_predict, style)

        # 샘플 무작위 선택
        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
        # 파라미터 기억해 두기
        theta_path_sgd.append(theta)               

plt.plot(X, y, "b.")                               
plt.xlabel("$x_1$", fontsize=18)                   
plt.ylabel("$y$", rotation=0, fontsize=18)         
plt.axis([0, 2, 0, 15])                            
save_fig("sgd_plot")                               
plt.show()<jupyter_output>그림 저장: sgd_plot<jupyter_text>학습된 파라미터의 최종 결과는 다음과 같다.<jupyter_code>theta<jupyter_output><empty_output><jupyter_text>50번의 에포크를 학습한 결과가 매우 우수하다. 앞서 살펴보았던 배치 경사하강법의 경우 에포크가 1000이었던 것에 비해 훨씬 좋은 성능을 확률적 경사 하강법이 발휘하였다. 하지만 이는 매우 특별한 경우이며보통은 더 많은 학습 반복을 요구한다.   사이킷런의 `SGDRegressor` 사이킷런의 `SGDRegressor` 예측기가 확률적 경사 하강법을 바로 지원한다. 사용되는 하이퍼파라미터는 다음과 같다.* `max_iter=1000`: 에포크 수 제한. 최대 1000번.* `tol=1e-3`: 허용 오차는 0.001.* `eta0=0,1`: SGDRegressor가 사용하는 학습 스케줄 함수에 사용되는 매개 변수임.    * 일종의 학습률 역할을 수행하지만 엄밀히 따지면 다르다.     * 사용되는 학습 스케줄 함수도 앞서 소개한 것과 다르다. 굳이 자세히 알 필요는 없음.* `penalty=None`: 규제를 사용하지 않음. 규제에 대한 설명은 잠시 뒤에 이루어짐.**주의:** `SGDRegressor는`의 `fit()` 추정기는 타깃 벡터 `y`가 1차원 어레이일 때만 다룬다. `y`는 원래 (100,1) 모양의 2차원 어레이이다. 따라서 `ravel()` 메서드를 이용하여 1차원 어레이로 변환시킨 후 `fit()` 추정기를 실행한다.<jupyter_code>from sklearn.linear_model import SGDRegressor

sgd_reg = SGDRegressor(max_iter=1000, tol=1e-3, penalty=None, eta0=0.1, random_state=42)
sgd_reg.fit(X, y.ravel())<jupyter_output><empty_output><jupyter_text>학습된 파라미터는 편향과 가중치로 구분되어 각각 아래 두 인스턴스 속성에 저장된다. * `intercept_`: 편향 저장* `coef_`: 가중치 저장결과는 기존과 거의 같다.<jupyter_code>sgd_reg.intercept_, sgd_reg.coef_<jupyter_output><empty_output><jupyter_text>미니배치 경사 하강법 각 스텝에서 학습하는 샘플의 개수를 지정하여 사용하는 경사 하강법을 **미니매치 경사 하강법**이라 부른다.* 배치 크기를 어느 정도 크게 하면 확률적 경사 하강법(SGD) 보다 파라미터의 움직임이 보다 덜 불규칙적이 된다.* 반면에 배치 경사 하강법보다 빠르게 학습한다. * 물론 지역 최솟값을 피하기가 덜 쉬워질 수 있다. 아래 코드는 배치 크기(`minibatch_size`)를 20으로 하여 경사 하강법을 적용한 코드이다. 사용된 하이퍼파라미터는 다음과 같다.* `n_iterations=50`: 에포크 수* `minibatch_size=20`: 배치 크기, 즉, 하나의 스텝에 사용되는 샘플 수또한 학습 스케줄을 사용하며,에포크가 달라질 때마다 전체 훈련 세트를 무작위로 섞는다. 각 학습법의 비교를 위해 여기서도 학습 과정에서 조정되는 파라미터(`theta`)의 변화를 추적한다.(`theta_path_mgd` 활용)<jupyter_code># 학습 과정에서 파라미터 변화 기록
theta_path_mgd = []

# 에포크
n_iterations = 50
# 배치 크기
minibatch_size = 20

np.random.seed(42)

# 파라미터 무작위 초기값
theta = np.random.randn(2,1)

# 학습 스케줄
t0, t1 = 200, 1000
def learning_schedule(t):
    return t0 / (t + t1)

t = 0

# 학습 시작
for epoch in range(n_iterations):
    # 매 에포크마다 훈련 세트 무작위로 섞기
    shuffled_indices = np.random.permutation(m)
    X_b_shuffled = X_b[shuffled_indices]
    y_shuffled = y[shuffled_indices]
    
    # 하나의 스텝에서 20개의 샘플을 학습한 후 파라미터를 조정한다.
    for i in range(0, m, minibatch_size):
        t += 1
        xi = X_b_shuffled[i:i+minibatch_size]
        yi = y_shuffled[i:i+minibatch_size]
        gradients = 2/minibatch_size * xi.T.dot(xi.dot(theta) - yi)
        eta = learning_schedule(t)
        theta = theta - eta * gradients
        # 조정된 파라미터 저장
        theta_path_mgd.append(theta)<jupyter_output><empty_output><jupyter_text>학습 결과는 이전과 매우 비슷하다.<jupyter_code>theta<jupyter_output><empty_output><jupyter_text>경사 하강법 비교 배치 경사 하강법, 확률적 경사 하강법, 미니배치 경사 하강법 실행 과정에서파라미터가 변하는 과정을 각각 아래 세 개의 리스트에 저장해 두었다.* `theta_path_bgd`: 배치 경사 하강법 파라미터 리스트* `theta_path_sgd`: 확률적 경사 하강법 파라미터 리스트* `theta_path_mgd`: 미니배치 경사 하강법 파라미터 리스트세 리스트를 어레이로 변환하여 사용한다.<jupyter_code>theta_path_bgd = np.array(theta_path_bgd)
theta_path_sgd = np.array(theta_path_sgd)
theta_path_mgd = np.array(theta_path_mgd)<jupyter_output><empty_output><jupyter_text>저장되는 파라리터는 편향과 가중치로 구성된다.* 편향($\theta_0$): 0번 열* 가중치($\theta_1$): 1번 열편향을 수평축, 가중치를 수직축으로 사용하여 그래프를 그리면 각각의 모델이학습되는 과정에서 파라미터가 어떻게 변하는가를 한 눈에 확인할 수 있다. * 배치 경사 하강법(파란색): 부드럽게 전역 최솟값으로 수렴한다.* 확률적 경사 하강법(빨간색): 전역 최솟값 근처에서 많이 요동친다.* 미니배치 경사 하강법(녹색): 확률적 경사 하강법 보다 덜 요동치지만, 여전히 전역 최솟값 근처에서 요동친다.<jupyter_code>plt.figure(figsize=(7,4))
plt.plot(theta_path_sgd[:, 0], theta_path_sgd[:, 1], "r-s", linewidth=1, label="Stochastic")
plt.plot(theta_path_mgd[:, 0], theta_path_mgd[:, 1], "g-+", linewidth=2, label="Mini-batch")
plt.plot(theta_path_bgd[:, 0], theta_path_bgd[:, 1], "b-o", linewidth=3, label="Batch")
plt.legend(loc="upper left", fontsize=16)
plt.xlabel(r"$\theta_0$", fontsize=20)
plt.ylabel(r"$\theta_1$   ", fontsize=20, rotation=0)
plt.axis([2.5, 4.5, 2.3, 3.9])
save_fig("gradient_descent_paths_plot")
plt.show()<jupyter_output>그림 저장: gradient_descent_paths_plot<jupyter_text>선형 회귀 알고리즘 비교 지금까지 소개한 선형 회귀 알고리즘을 간략하게 비교하면 다음과 같다.| 알고리즘   | 많은 샘플 수 | 외부 메모리 학습 | 많은 특성 수 | 하이퍼 파라미터 수 | 스케일 조정 | 사이킷런 지원 ||:--------:|:---------:|:---------:|:---------:|:---------:|:---------:|:---------:|| 정규방정식  | 빠름       | 지원 안됨      |  느림        | 0          | 불필요    | 지원 없음      || SVD      | 빠름       | 지원 안됨      |  느림        | 0          | 불필요     | LinearRegression     || 배치 GD   | 느림       | 지원 안됨      |  빠름        | 2          | 필요      | LogisticRegression      || SGD      | 빠름       | 지원          |  빠름        | >= 2       | 필요      | SGDRegressor || 미니배치 GD | 빠름       | 지원         |  빠름        | >=2        | 필요      | 지원 없음      |**주의**: 책 176쪽, 표 4-1에서는 SGGRegressor가 배치 경사 하강법과 미니배치 경사 하강법을 지원한다고 언급되었다.하지만 정확하지 않다. 이유는 다음과 같다.* 배치 경사 하강법 관련    * [SGDRegressor 문서](https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.SGDRegressor.html)에         배치에 대한 어떠한 언급도 없다.    * 또한 SGD와 배치 경사 하강법은 기본적으로 다른 아이디어를 사용한다.    * SGDRegressor에서 배치 경사 하강법을 지원할 이유가 없다.     * 반면에 잠시 뒤에 배울 로지스틱 회귀를 지원하는 사이킷런의 `LogisticRegression`이         배치 경사 하강법을 사용한다.* 미니배치 경사 하강법    * 책 176쪽, 20번 주석에 옮긴이가 언급한 내용을 보면,         SGDRegressor 클래스의 `partial_fit()` 메서드와 **미니배치와는 아무 상관이 없다**.    * `partial_fit()` 메서드를 실행하면 전체 샘플에 대해 한 번만 확률적 경사 하강법으로 학습한다.        즉, `max_iter=1`인 경우와 동일하게 작동한다.     * 참조: [SGDRegressor 문서](https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.SGDRegressor.html)  핵심 3: 다항 회귀 **다항 회귀**(polynomial regression)는 선형 회귀를 이용하여 비선형 데이터를 학습하는 기법이며,기본 아이디어는 다음과 같다.* 2장에서 배운 특성 조합을 이용한다.* 추가되는 특성은 기존의 특성값들의 거듭제곱, 특성값들 사이의 곱 등으로 이루어진다. * 즉, 특성 변수들의 다항식을 조합 특성으로 추가한다.  예제: 2차 방정식 활용 다항 회귀 2차 다항식 곡선에 대응하는 비선형 모델을 선형 회귀를 이용하여 학습 하는 방법을 소개한다.<jupyter_code>import numpy as np
import numpy.random as rnd

np.random.seed(42)<jupyter_output><empty_output><jupyter_text>먼저 구간 (0, 1)에서 균등하게 무작위로 100개의 샘플을 생성한 후,2차 함수 $f(x) = \frac 1 2 x^2 + x + 2$ 맞추어레이블 생성한다.제대로된 학습을 위해 가우시안 잡음도 추가한다.<jupyter_code>m = 100
X = 6 * np.random.rand(m, 1) - 3
y = 0.5 * X**2 + X + 2 + np.random.randn(m, 1)<jupyter_output><empty_output><jupyter_text>생성된 데이터의 분포를 확인하면 다음과 같이 비선형이다.<jupyter_code>plt.plot(X, y, "b.")
plt.xlabel("$x_1$", fontsize=18)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.axis([-3, 3, 0, 10])
save_fig("quadratic_data_plot")
plt.show()<jupyter_output>그림 저장: quadratic_data_plot<jupyter_text>비선형 문제를 선형 회귀를 단순히 앞서 배운대로 사용하면 아래에서처럼 학습이 제대로 이루어지지 않는다.<jupyter_code>lin_reg = LinearRegression()
lin_reg.fit(X, y)

X_new=np.linspace(-3, 3, 100).reshape(100, 1)
y_new = lin_reg.predict(X_new)
plt.plot(X, y, "b.")
plt.plot(X_new, y_new, "r-", linewidth=2, label="Predictions")
plt.xlabel("$x_1$", fontsize=18)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.legend(loc="upper left", fontsize=14)
plt.axis([-3, 3, 0, 10])
save_fig("linear_predictions_plot")
plt.show()<jupyter_output>그림 저장: linear_predictions_plot<jupyter_text>`PolynomialFeatures` 변환기 앞서 본 것처럼 주어진 특성들만 사용하여 선형 회귀를 적용하면 학습이 제대로 이루어지지 않는다.그런데 주어진 특성들의 거듭제곱을 특성으로 추가하여 선형 회귀를 적용하면 위 경우와 같은 문제는 해결할 수 있다.사이킷런의 `PolynomialFeatures` 변환기가주어진 특성들의 거듭제곱과 특성들 사이의 곱셈을 실행하여 특성을 추가하는 기능을 제공한다.* `degree=d`: $d\ge 2$에 대해 `PolynomialFeatures` 변환기가 하는 일은 다음과 같다.    * $x_1, \dots, x_n$가 n개의 특성 변수라 하자.        * 그러면 아래 다항식에 사용되는 모든 항들에 대한 특성 변수가 추가 된다.        $$\sum_{k=2}^d(x_1 + \cdots + x_n)^k$$                * $n, d \ge 2$인 경우 정확히 다음 숫자 만큼의 특성이 추가된다.            따라서 높은 차수를 선택하면 매우 많은 수의 특성이 추가됨에 주의해야 한다.            $$\frac{(n+d)!}{n! \, d!}-n$$             * 예를 들어 $n=2, d=3$인 경우         $(x_1+x_2)^2 + (x_1+x_2)^3$        에 사용된 아래 항들에 대한 특성 6개가 추가된다.         $$x_1^2,\,\, x_1 x_2,\,\, x_2^2,\,\, x_1^3,\,\, x_1^2 x_2,\,\, x_1 x_2^2,\,\, x_2^3$$        * 아래 예제에서는 $n=1, d=2$ 이기에 $x_1^2$에 대한 특성 변수만 추가된다.* `include_bias=False`: 기본값은 `True`임.    * `include_bias=True`인 경우 편향(bias) 특성이 0번 열에 추가된 후 변환 과정이 실행됨.    * 하지만 `LinearRegression`과 함께 사용될 경우 자동으로 0번 열에 1이 추가됨.    * 따라서 여기서는 `False`를 사용한 것임.    * 다른 모델과 함께 사용될 경우는 `True` 옵션을 그대로 사용하는 것이 추천됨. 이제 앞서 2차 함수를 활용하여 무작위로 생성한 데이터에 대한 선형 회귀 학습을 시작한다.먼저 2차 다항 회귀를 준비한 후 기존의 데이터를 변환한다.<jupyter_code>from sklearn.preprocessing import PolynomialFeatures
poly_features = PolynomialFeatures(degree=2, include_bias=False)
X_poly = poly_features.fit_transform(X)<jupyter_output><empty_output><jupyter_text>`X_poly`는 변환된 데이터를 갖고 있으며 0번 열은 기존의 특성을, 1번 열은 0번 열 값의 제곱을 항목으로 갖는다. 처음 10개의 샘플을 확인하면 다음과 같다.<jupyter_code>X_poly[:10]<jupyter_output><empty_output><jupyter_text>이제 선형 회귀 학습을 진행한다.<jupyter_code>lin_reg = LinearRegression()
lin_reg.fit(X_poly, y)<jupyter_output><empty_output><jupyter_text>학습된 파라미터 $\theta_0, \theta_1, \theta_2$는 다음과 같다.<jupyter_code>lin_reg.intercept_, lin_reg.coef_<jupyter_output><empty_output><jupyter_text>애초에 사용되었던 계수들은 차례대로 2, 1, 1/2 이었으며, 상당히 비슷한 결과를 얻었다.물론 잡음 때문에 완벽한 예측은 불가능하다. 이제 학습 결과를 이용하여 새로운 샘플 100개에 대한 예측을 실행하고실행 결과를 기존의 데이터와 함께 그래프로 그린다. 그러면 예측 결과가 멋진 2차원 그래프가 되었음을 확인할 수 있다.<jupyter_code># 새로운 샘플 100개 무작위 생성
X_new=np.linspace(-3, 3, 100).reshape(100, 1)
# 새로운 샘플도 변환은 필요
X_new_poly = poly_features.transform(X_new)
# 새로운 샘플에 대한 예측
y_new = lin_reg.predict(X_new_poly)

# 기존 훈련 세트와 새로운 샘플을 함께 그래프로 그리기
plt.plot(X, y, "b.")
plt.plot(X_new, y_new, "r-", linewidth=2, label="Predictions")
plt.xlabel("$x_1$", fontsize=18)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.legend(loc="upper left", fontsize=14)
plt.axis([-3, 3, 0, 10])
save_fig("quadratic_predictions_plot")
plt.show()<jupyter_output>그림 저장: quadratic_predictions_plot<jupyter_text>`include_bias=True` 인 경우 변환 과정으로 생성된 `X_poly`의 0번 열이 1로 채워진다.<jupyter_code>from sklearn.preprocessing import PolynomialFeatures
poly_features = PolynomialFeatures(degree=2, include_bias=True)
X_poly = poly_features.fit_transform(X)
X_poly[:10]<jupyter_output><empty_output><jupyter_text>이후 `LinearRegression` 이용한 학습 결과는 하지만 동일하다. 이유는 앞서 설명한 대로 `LinearRegression`가 자동으로 0번 열에 1을 추가하는 기능을 갖고 있기 때문이다.<jupyter_code>lin_reg = LinearRegression()
lin_reg.fit(X_poly, y)
lin_reg.intercept_, lin_reg.coef_
X_new=np.linspace(-3, 3, 100).reshape(100, 1)
X_new_poly = poly_features.transform(X_new)
y_new = lin_reg.predict(X_new_poly)
plt.plot(X, y, "b.")
plt.plot(X_new, y_new, "r-", linewidth=2, label="Predictions")
plt.xlabel("$x_1$", fontsize=18)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.legend(loc="upper left", fontsize=14)
plt.axis([-3, 3, 0, 10])
save_fig("quadratic_predictions_plot")
plt.show()<jupyter_output>그림 저장: quadratic_predictions_plot<jupyter_text>핵심 4: 학습 곡선 고차 다항 회귀를 적용하면 일반적으로 주어진 데이터에 더욱 잘 적응한다. 아래 코드는 기본 선형 회귀 모델과 2차, 300차 다항 회귀 모델의 학습 결과를 보여준다.빨간색, 파란색, 녹색 선은 새로운 데이터에 대한 예측 결과를 각각의 모델에 대한 그래프로 보여준다.<jupyter_code>from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import Pipeline

# 300차 다항 회귀(녹색), 2차 다항 회귀(파란색), 기본 다항 회귀(빨간색) 순서대로
# 새로운 데이터에 대한 예측 결과를 그래프로 그린다. 
# 기본 다항 회귀의 degree는 1을 사용한다. 
for style, width, degree in (("g-", 1, 300), ("b--", 2, 2), ("r-+", 2, 1)):
    polybig_features = PolynomialFeatures(degree=degree, include_bias=False)
    std_scaler = StandardScaler()
    lin_reg = LinearRegression()
    polynomial_regression = Pipeline([
            ("poly_features", polybig_features),
            ("std_scaler", std_scaler),
            ("lin_reg", lin_reg),
        ])
    polynomial_regression.fit(X, y)
    y_newbig = polynomial_regression.predict(X_new)
    plt.plot(X_new, y_newbig, style, label=str(degree), linewidth=width)

# 기존 훈련 세트 산점도 그리기
plt.plot(X, y, "b.", linewidth=3)

# 기타 그래프 관련 정보 설정
plt.legend(loc="upper left")
plt.xlabel("$x_1$", fontsize=18)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.axis([-3, 3, 0, 10])
save_fig("high_degree_polynomials_plot")
plt.show()<jupyter_output>그림 저장: high_degree_polynomials_plot<jupyter_text>과소적합 기본 선형 회귀의 경우 고차 항을 전혀 사용하지 않으며, 따라서 빨간색 직선으로 예측한다.이는 학습 과정에서 훈련 세트에 제대로 적응하지 않았음을 보이며, 이런 것을 **과소적합**(underfitting)이라 부른다.  과대 적합 반면에 300차 다항 회귀의 경우 주어진 훈련 세트에 과하게 적응하였으며, 이를 **과대적합**(overfitting)이라 부른다.과대적합 모델의 알고 있는 데이터에 너무 적응이 되어서 새로운 데이터에 대한 예측을 제대로 하지 못한다.   학습 곡선 학습된 모델의 과소적합 또는 과대적합 여부를 확인하려면 일반적으로 2장에서 배운 **교차 검증**을 사용한다.* 훈련 세트에 대한 검증은 우수하지만 교차 검증 점수가 낮으면 모델이 과대적합되었다.* 훈련 세트와 교차 검증 점수 모두 낮으면 모델이 과소적합되었다.또 다른 방법은 학습 곡선을 분석하는 것이다. 훈련 세트와 검증 세트에 대한 모델 성능을 모델이 학습하는 과정을 반복하여 추적한 결과를 그래프로 그린 결과가 **학습 곡선**이다.   훈련 세트, 검증 세트, 테스트 세트의 구분 머신러닝 모델 학습에 사용되는 데이터셋은 역할에 따라 다음 세 종류로 구분된다.* 훈련 세트    * 모델 훈련 과정에 사용되는 데이터 샘플들의 집합    * 즉, 비용함수를 계산하는 데에 직접 사용되는 데이터 샘플들의 집합    * 훈련 세트와 테스트 세트로 나뉘어질 때 생성된 전체 훈련 세트를 가리키거나,        아니면, 교차 검증 등에서처럼 반복 학습 과정에서 사용되는 (보다 작은) 훈련 세트를 가리키기도 함.* 검증 세트    * 교차 검증 기법 등에서 처럼, 전체 훈련 세트를 여러 개의 세트로 나뉘어진 부분 집합 중의 하나.    * 다른 훈련 세트를 통해 학습된 모델의 성능을 측정하기 위해 사용되는 데이터 샘플 집합    * 하지만, 다음 번 학습 과정에서는 훈련 세트에 포함될 수도 있는 전체 훈련 세트의 일부.* 테스트 세트    * 처음부터 전체 훈련 세트로부터 분리된 데이터셋.    * 모델 학습 과정에 전현 관여하지 않음.    * 검증 과정에도 관여하지 않음.    * 모든 훈련이 끝난 후에 구현된 모델의 실전 성능을 테스트 하기위해 사용되는 데이터셋을 가리킴.  학습 곡선 함수 아래 코드에서 정의된 `plot_learning_curves()`함수는훈련 세트의 부분 집합의 크기를 1부터 키워 가면서 훈련시킨 모델의 성능을훈련 과정에서 사용된 데이터에 대한 예측 결과와 검증 세트의 데이터, 즉 새로운 데이터에 대한 예측 결과에 대한 예측 결과를 비교하는 그래프를 그리는 함수이다. 사용되는 매개변수와 값들의 역할은 다음과 같다.* `model`: 학습 모델* `X`: 전체 데이터셋* `y`: 전체 레이블* 먼저, 8대 2로 훈련 세트(`X_train`, `y_train`)와 검증(`X_val`, `y_val`) 세트로 나눔.* 훈련 세트를 대상으로 지정된 인덱스를 기준으로 짤라서 아래 과정을 반복한다.     인덱스는 1부터 `len(X_train)`까지 움직인다.     * 지정된 인덱스까지의 훈련 샘플을 대상으로 학습    * 학습에 사용된 훈련 세트를 대상으로 예측 후, 평균 제곱근 오차 계산하여 `train_errors` 리스트에 추가.    * 검증 세트 `X_val`을 대상으로 예측 후, 평균 제곱근 오차 계산하여 `val_errors` 리스트에 추가* 최종적으로 `train_errors`와 `val_errors`를 이용하여 그래프를 그린다.<jupyter_code>from sklearn.metrics import mean_squared_error
from sklearn.model_selection import train_test_split

def plot_learning_curves(model, X, y):
    # 훈련 세트와 검증 세트로 분리 (8대 2)
    X_train, X_val, y_train, y_val = train_test_split(X, y, test_size=0.2, random_state=10)
    # 평균 제곱근 오차 기억 리스트
    train_errors, val_errors = [], []
    
    # 총 len(X_train) 번 학습 반복
    # 반복할 때마다 훈련 세트의 크기가 1씩 커진다.
    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="train")
    plt.plot(np.sqrt(val_errors), "b-", linewidth=3, label="val")
    plt.legend(loc="upper right", fontsize=14)  
    plt.xlabel("Training set size", fontsize=14)
    plt.ylabel("RMSE", fontsize=14)<jupyter_output><empty_output><jupyter_text>과소적합 모델의 학습 곡선 특징 기본 선형 회귀 모델을 대상으로 하면 과소적합이 발생하며, 과소적합에 대한 학습 곡선은 기본적으로 다음 특성을 보인다.* 훈련 데이터에 대한 성능    * 훈련 세트가 커지면서 RMSE(평균 제곱근 오차)가 커진다.        * 2차 다항 회귀가 아닌 기본 선형 회귀이기 때문에 어쩔 수 없다.    * 훈련 세트가 어느 정도 커지면 더 이상 RMSE가 변하지 않는다.* 검증 데이터에 대한 성능    * 검증 세트에 대한 성능이 훈련 세트에 대한 성능과 거의 비슷해진다.<jupyter_code>lin_reg = LinearRegression()
plot_learning_curves(lin_reg, X, y)
plt.axis([0, 80, 0, 3])                         
save_fig("underfitting_learning_curves_plot")   
plt.show()<jupyter_output>그림 저장: underfitting_learning_curves_plot<jupyter_text>과대적합 모델의 학습 곡선 특징 고차 다항 회귀 모델을 대상으로 하면 과대적합이 발생하며, 과대적합에 대한 학습 곡선은 기본적으로 다음 특성을 보인다.* 훈련 데이터에 대한 성능    * 훈련 데이터에 대한 평균 제곱근 오차가 다른 모델보다 훨씬 적다.* 검증 데이터에 대한 성능    * 훈련 데이터에 대한 성능과 차이가 크게 벌어진다. 즉, 훈련 데이터에 대한 성능이 훨씬 좋다.    아래 코드는 10차 다항 회귀를 이용하여 학습한 결과를 비교한다.<jupyter_code>from sklearn.pipeline import Pipeline

polynomial_regression = Pipeline([
        ("poly_features", PolynomialFeatures(degree=10, include_bias=False)),
        ("lin_reg", LinearRegression()),
    ])

plot_learning_curves(polynomial_regression, X, y)
plt.axis([0, 80, 0, 3])           # 책에는 없음
save_fig("learning_curves_plot")  # 책에는 없음
plt.show()                        # 책에는 없음<jupyter_output>그림 저장: learning_curves_plot<jupyter_text>과대적합 모델 개선법 검증 데이터에 대한 성능이 훈련 데이터에 대한 성능에 근접할 때가지 훈련 데이터의 크기를 늘린다.  편향/분산 트레이드오프 모델이 훈련 결과를 바탕으로 해서 새로운 데이터에 대한 예측을 할 때 발생하는 오차를**모델의 일반화 오차**(model's generalization error)라 부르며,다음 세 종류의 오차의 합으로 표현된다.* 편향    * 잘못된 가정으로 발생하는 오차    * 예를 들어, 2차 회귀 모델을 선형 모델로 구현하면 과소적합이 발생하는 데 이때 큰 편향 발생.    * 즉, 편향이 큰 모델은 과소적합 모델이 될 가능성 높음.    * **주의**: 선형 회귀 모델을 설정할 때 사용되는 편향($\theta_0)$과는 다른 의미로 사용됨.* 분산    * 훈련 샘플의 작은 변화에 민감하게 반응하는 경우 분산이 높음.    * 예를 들어, 고차 다항 회귀 모델은 분산이 크며, 그 결과 과대적합되는 경향 큼.* 줄일 수 없는 오차    * 데이터 자체에 포함되어 있는 잡음(가우시안 잡음 등)때문에 어쩔 수 없이 발생하는 오차    * 이전 예제들에서 보았듯이 가우시안 잡음으로 인해 원래 사용되었던 파라미터를 정확하게 찾을 수 없었음.    편향과 분산은 서로 반대로 움직이며, 이를 **편향/분산 트레이드오프**라 부른다. 예를 들어, 모델의 복잡도가 올라가면 분산은 늘고 편향은 줄어든다. 반대로 모델의 복잡도가 내려가면 분산은 줄고 편향은 늘어난다.  핵심 5: 규제 선형 모델 규제를 이용하여 과대적합을 감소하는 방법을 알아본다. **규제**(regularization)는 **자유도**(degree of freedom)를 제한한다는 의미로 사용된다.자유도는 학습 모델 결정에 영향을 주는 요소(특성)들의 개수를 가리킨다. 선형 회귀의 경우 특성 수로 , 다항 선형 회귀의 경우 차수로 자유도가 정해진다. 다항 선형 회귀 모델의 규제는 차수를 줄이는 방식으로 이루어지며,선형 회귀 모델에 대한 규제는 가중치의 역할을 제한하는 방식으로 이루어진다.여기서는 선형 회귀 모델의 가중치를 제한하는 세 가지 방식과 함께 검즈 세트에 대한 비용(오차)이 최소가 될 때바로 훈련을 중지시키는 **조기 종료**(early stopping) 기법을 소개한다. 가중치에 제한을 가하는 세 가지 방식은 다음과 같다.* 릿지 회귀(티코노프 규제)* 라쏘 회귀* 엘라스틱넷언급된 세 가지 방식은 모두 선형 회귀 학습 과정에서 사용되는 비용 함수인 MSE에 가중치 파라미터를제한하는 값을 추가하는 방식으로 이루어진다.   릿지 회귀 **릿지 회귀**(티코노프 규제)의 비용 함수 $J(\theta)$는 다음과 같이 정의된다.$$J(\theta) = \textrm{MSE}(\theta) + \alpha \, \frac{1}{2} \sum_{i=1}^{n}\theta_i^2$$평균 제곱 오차(MSE)에 추가된 규제항$\alpha \, \frac{1}{2} \sum_{i=1}^{n}\theta_i^2$의 역할은 다음과 같다. (**주의**: $i>0$임. 즉 $\theta_0$은 규제하지 않음)* $\alpha$(알파)    * 하이퍼파라미터로 지정됨.    * 규제 강도 지정    * $\alpha=0$이면 규제가 전혀 없는 기본 선형 회귀    * $\alpha$가 커질 수록 가중치의 역할이 줄어듦.         * 이유: 비용을 줄여야 하기 때문에 가중치를 줄이는 방향으로 학습할 것임.    * $\alpha$가 매우 크면 가중치는 거의 0에 가까울 것이며,        따라서 모델은 데이터의 평균 근처를 지나는 수평선이 된다.        * 이유: 가중치가 거의 0이기에 편향($\theta_0$)에 의해 모델이 결정됨.            따라서 $\theta_0$가 훈련 데이터의 평균값으로 수렴할 것임.            그래야 MSE가 최소화됨.            * $\theta_i$    * $\alpha>1$이면, $\theta_i^2$에 의해 MSE에 크게 영향을 줌.    * 따라서 $\theta_i$ 값을 제한하는 형태로 학습이 이루어지도록 유도됨.  참조 규제항은 훈련 과정에만 사용된다. 테스트 과정에는 다른 기준으로 성능을 평가한다.* 훈련 과정: 비용 최소화 목표* 테스트 과정: 최종 목표에 따른 성능 평가* 예제: 분류기 (잠시 뒤에 다룸)    * 훈련 과정: 로그 손실(log loss)을 비용함수로 활용    * 테스트 과정: 재현율/정밀도 기준으로 성능 평가  특성 스케일링 릿지 회귀를 포함하여 규제를 사용하는 모델은 기본적으로 특성들의 스케일에 민감하다.따라서 표준화를 지원하는 `StandardScaler` 등을 이용하여 스케일링 전처리를 실행해야 한다.  예제: 릿지 모델 훈련 릿지 모델 학습을 위해 사용되는 데이터셋은 다음과 같이 구한다.* `m`: 전체 훈련 샘플 수* `X`: 전체 훈련 세트    * (0,3) 구간에서 균등하게 무작위로 `m`개 선택* `y`: 전체 훈련 샘플들의 레이블. 가우시안 잡음 추가됨.* `X_new`: 테스트 세트    * (0,3) 구간을 100개의 작은 구간으로 균등하게 쪼갠 값들의 테스트 세트 (`np.linspace()` 함수 활용)<jupyter_code>np.random.seed(42)
m = 20
X = 3 * np.random.rand(m, 1)
y = 1 + 0.5 * X + np.random.randn(m, 1) / 1.5
X_new = np.linspace(0, 3, 100).reshape(100, 1)<jupyter_output><empty_output><jupyter_text>사이킷런에서 릿지 모델을 훈련시키려면 `Ridge` 예측기를 활용한다.`Ridge` 예측기에 사용되는 하이퍼파라미터 중에서 아래 코드에서 사용되는 값들은 다음과 같다.* `alpha`: 규제항에 사용되는 알파값* `solver`: 비용함수 최적화에 사용되는 알고리즘    * 기본값은 `auto`: 사용되는 데이터에 따라 자동으로 알고리즘 선택    * `cholesky`(숄레스키): 릿지 회귀의 비용함수에 대한 아래 모양의 정규방정식 활용. (자세한 내용은 중요하지 않음)        $$\hat\theta = (\mathbf{X}_b^T \mathbf{X}_b + \alpha\, \mathbf A)^{-1} \mathbf{X}_b^T \mathbf{y}$$            * `sag`: 확률적 평균 경사 하강법(Stochastic Average Gradient descent, SAG) 사용.        * SAG: 매 샘플에 대해 계산되는 그레이디언트 벡터의 구성원을            확률적으로 선택해서 조정하는 기법.        * 참조 논문: [Minimizing Finite Sums with the Stochastic Average Gradient](https://arxiv.org/abs/1309.2388)        * `random_state`: `solver`에 사용되는 알고리즘 중에서 `sag`와 `saga`의 경우 훈련 세트를 무작위로 섞음.    * `cholesky` 알고리즘을 사용할 경우 사용할 필요 없음.  `cholesky` 알고리즘 활용<jupyter_code>from sklearn.linear_model import Ridge
ridge_reg = Ridge(alpha=1, solver="cholesky")
ridge_reg.fit(X, y)

# 예측 예제
ridge_reg.predict([[1.5]])<jupyter_output><empty_output><jupyter_text>`random_state`를 지정해도 결과는 동일하다.<jupyter_code>ridge_reg = Ridge(alpha=1, solver="cholesky", random_state=42)
ridge_reg.fit(X, y)

# 예측 예제
ridge_reg.predict([[1.5]])<jupyter_output><empty_output><jupyter_text>릿지 회귀 비교 아래 코드는 선형 모델과 다중 회귀 모델을 **릿지 규제**를 추가하여 실행한 결과를 비교한다.* 왼쪽 이미지    * 선형 회귀 + 릿지 규제    * `alpha`: 0, 10, 100        * 0 인 경우: 일반 선형 회귀* 오른쪽 이미지    * 다중 회귀 + 릿지 규제        * `degree=10`: 10차 다중 회귀 사용    * `alpha`: 0, 1e-05, 1        * 0인 경우: 일반 다중 회귀 적용    * 전처리: `PolynomialFeatures` 변환기 + 표준화 스케일링`plot_model()` 함수에 사용되는 매개변수의 역할은 다음과 같다.* `model_classe`: 모델 지정    * `Ridge`        * `solver`는 `auto` 기본값 사용    * `Rasso`        * 바로 이어서 설명됨.* `polynomial`: 다중 회귀 사용 여부 결정* `alphas`: 알파 값들의 리스트* `**model_kargs`: 모델 학습에 사용될 추가 옵션.     * `random_state` 가 사용됨.<jupyter_code>from sklearn.linear_model import Ridge

def plot_model(model_class, polynomial, alphas, **model_kargs):
    # 알파값을 달리 하면서 그래프 그리기
    for alpha, style in zip(alphas, ("b-", "g--", "r:")):
        # 모델 선택
        model = model_class(alpha, **model_kargs) if alpha > 0 else LinearRegression()
        # 다항 회귀의 경우 따른 전처리 활용
        if polynomial:
            model = Pipeline([
                    ("poly_features", PolynomialFeatures(degree=10, include_bias=False)),
                    ("std_scaler", StandardScaler()),
                    ("regul_reg", model),
                ])
            
        # 모델 훈련
        model.fit(X, y)
        # 테스트 세트 예측
        y_new_regul = model.predict(X_new)
        # 선, 점선 두께 설정
        lw = 2 if alpha > 0 else 1
        # 그래프 그리기
        plt.plot(X_new, y_new_regul, style, linewidth=lw, label=r"$\alpha = {}$".format(alpha))
        
    # 기존 훈련 데이터 그래프 그리기
    plt.plot(X, y, "b.", linewidth=3)
    
    # 기타 설정
    plt.legend(loc="upper left", fontsize=15)
    plt.xlabel("$x_1$", fontsize=18)
    plt.axis([0, 3, 0, 4])

plt.figure(figsize=(8,4))
plt.subplot(121)
plot_model(Ridge, polynomial=False, alphas=(0, 10, 100), random_state=42)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.subplot(122)
plot_model(Ridge, polynomial=True, alphas=(0, 10**-5, 1), random_state=42)

save_fig("ridge_regression_plot")
plt.show()<jupyter_output>그림 저장: ridge_regression_plot<jupyter_text>위 그래프를 분석하면 다음과 같다.* 알파($\alpha$)가 커질 수록 모델이 직선에 가까워진다.     * 분산은 감소    * 편향은 증가* 왼쪽 그래프에서 $\alpha=100$ 인 경우    * 앞서 설명했듯이 가중치 $\theta_9$가 0에 수렴함.    * 따라서 평균값에 해당하는 수평선으로 모델 학습이 이루어짐.  규제 확률적 경사 하강법 확률적 경사 하강법에 규제를 가하는 방식으로 이루어진다. `SGDRegressor`의 경우 `penalty=l2` 로 지정해서 $\ell_2$ 노름을 지정한다. **주의**: 향후 버전이 바뀌더라도 동일한 결과를 만들기 위해 사이킷런 0.21 버전의 기본값인 `max_iter=1000`과 `tol=1e-3`으로 지정합니다.<jupyter_code>sgd_reg = SGDRegressor(penalty="l2", max_iter=1000, tol=1e-3, random_state=42)
sgd_reg.fit(X, y.ravel())

# 예측 예제
sgd_reg.predict([[1.5]])<jupyter_output><empty_output><jupyter_text>`sag` 알고리즘 활용 `Ridge` 예측기에 사용되는 알고리즘을 `sag`로 지정하면`SGDRegressor`를 좀 더 개선된 결과를 얻는다.<jupyter_code>ridge_reg = Ridge(alpha=1, solver="sag", random_state=42)
ridge_reg.fit(X, y)

# 예측 예제
ridge_reg.predict([[1.5]])<jupyter_output><empty_output><jupyter_text>라쏘 회귀 **라쏘 회귀**의 비용 함수 $J(\theta)$는 다음과 같이 정의된다.$$J(\theta) = \textrm{MSE}(\theta) + \alpha \, \sum_{i=1}^{n}\mid\theta_i\mid$$평균 제곱 오차(MSE)에 추가된 규제항$\alpha \, \sum_{i=1}^{n} \mid\theta_i\mid$의 역할은 다음과 같다. (**주의**: $i>0$임. 즉 $\theta_0$은 규제하지 않음)* $\alpha$(알파)    * 하이퍼파라미터로 지정됨.    * 규제 강도 지정    * $\alpha=0$이면 규제가 전혀 없는 기본 선형 회귀            * $\theta_i$: 덜 중요한 특성을 무시하기 위해 $\mid\theta_i\mid$가 0에 수렴하도록 학습 유도.**주의**: 절댓값을 사용한 라쏘 회귀의 비용함수는 $\theta_i = 0$ ($i = 1, 2, \dots, n$)에서 미분가능하지 않다.하지만 이 경우 다행히도 서브그레이디언트 벡터라는 것을 활용하여 경사 하강법을 적용할 수 있다.여기서는 이 정도만 알아두면 된다.다음은 `Lasso` 예측기의 간단한 활용 예제이다.<jupyter_code>from sklearn.linear_model import Lasso

lasso_reg = Lasso(alpha=0.1)
lasso_reg.fit(X, y)

# 예측 예제
lasso_reg.predict([[1.5]])<jupyter_output><empty_output><jupyter_text>라쏘 회귀 비교 아래 코드는 선형 모델과 다중 회귀 모델을 **라쏘 규제**를 추가하여 실행한 결과를 비교하기 위해앞서 정의된 `plot_model()` 함수 그대로 사용한다. * 왼쪽 이미지    * 선형 회귀 + 라쏘 규제    * `alpha`: 0, 10, 100* 오른쪽 이미지    * 다중 회귀 + 라쏘 규제        * `degree=10`: 10차 다중 회귀 사용    * `alpha`: 0, 1e-07, 1        * 좀 더 작은 알파($\alpha$) 사용<jupyter_code>from sklearn.linear_model import Lasso

plt.figure(figsize=(8,4))
plt.subplot(121)
plot_model(Lasso, polynomial=False, alphas=(0, 0.1, 1), random_state=42)
plt.ylabel("$y$", rotation=0, fontsize=18)
plt.subplot(122)
plot_model(Lasso, polynomial=True, alphas=(0, 10**-7, 1), random_state=42)

save_fig("lasso_regression_plot")
plt.show()<jupyter_output>/opt/conda/envs/tf2/lib/python3.7/site-packages/sklearn/linear_model/coordinate_descent.py:475: ConvergenceWarning: Objective did not converge. You might want to increase the number of iterations. Duality gap: 2.8028677038274545, tolerance: 0.0009294783355207351
  positive)<jupyter_text>오른쪽 그래프를 분석해보자.* $\alpha = 0$    * 10차 다중 회귀가 사용되어서인지 8차 다항식 그래프처럼 보인다.    * 과대적합된 것으로 보인다. * $\alpha =$ `1e-07`    * 10차 다중 회귀를 사용했음에도 불구하고 3차다항식 그래프처럼 보인다.     * 4차 이상의 특성값들에 대한 가중치가 모두 0이 된 것을 의미한다.     * 이렇게 가중치가 대부분 0이 되는 모델을 **희소 모델**(sparse model)이라 부름.* $\alpha = 1$    * 10차 다중 회귀를 사용했음에도 불구하고 데이터 평균값 근처를 수평선이 되었다.    * 모든 가중치가 0이 되었음을 의미한다.    * 전혀 의미 없는 모델이 생성되었다.    **주의:** 경고(Warning)가 뜨는 이유는 명확하지 않음. 허용오차(`tol`), 반복횟수(`max_iter`)를변경해도 수치가 달라질 뿐 동일한 경고가 계속 뜸. Lasso 예측기에 대한 보다 자세한 정보가 필요함.  엘라스틱넷 **엘라스틱넷**(elastic net)은 릿지 회귀와 라쏘 회귀를 절충한 모델이다.$$J(\theta) = \textrm{MSE}(\theta) + r\, \alpha \, \sum_{i=1}^{n}\mid\theta_i\mid + \,\frac{1-r}{2}\, \alpha\, \sum_{i=1}^{n}\theta_i^2$$$r$의 값에 따라 릿지 규제와 라쏘 규제를 적절하게 조절한다. 규제 사용 방법의 기본은 다음과 같다.* 대부분의 경우 약간이라도 규제 사용 추천* 릿지 규제가 기본* 유용한 속성이 많지 않다고 판단되는 경우     * 라쏘 규제나 엘라스틱넷 활용 추천    * 불필요한 속성의 가중치를 0으로 만들기 때문* 특성 수가 훈련 샘플 수보다 크거나 특성 몇 개가 강하게 연관되어 있는 경우    * 라쏘 규제는 적절치 않음.    * 엘라스틱넷 추천    아래 코드는 사이킷런의 `ElasticNet` 예측기의 간단한 활용 예제를 보여준다.* `l1_ratio`: 라쏘 규제의 홉합 비율 지정. 위 비용함수에서 사용된 $r$에 해당하는 값.<jupyter_code>from sklearn.linear_model import ElasticNet
elastic_net = ElasticNet(alpha=0.1, l1_ratio=0.5, random_state=42)
elastic_net.fit(X, y)

# 예측 예제
elastic_net.predict([[1.5]])<jupyter_output><empty_output><jupyter_text>라쏘 대 릿지 규제 비교 설명 아래 코드에서 그리는 그래프는 라쏘 회귀와 릿지 회귀 각각의 훈련 과정에서가중치 $\theta_1, \theta_2$와 비용 함수 $J(\theta) = J(\theta_1, \theta_2)$가 변화는 과정을 보여준다. * 왼편: 가중치들의 변화 그래프    * 상단: 라쏘 규제, 즉, $\ell_1$ 규제(페널티) 사용        * 등고선: $\ell_1$ 손실 ($\mid \theta_1 \mid + \mid \theta_2 \mid$) 나타냄        * $\theta$ 초기값: $\theta_1 = 2$, $\theta_2 = 0.5$        * 경사 하강법이 진행되면서  $\theta_2$가 먼저 0에 수렴. 이후 $\theta_1$도 0으로 수렴        * (0,0) 근처에서 진동. 이유는 미분이 불가능하기 때문임.    * 하단: 릿지 규제, 즉, $\ell_2$ 규제(페널티) 사용         * 등고선: $\ell_2$ 손실 ($\theta_1^2 + \theta_2^2$) 나타냄        * $\theta$ 초기값: $\theta_1 = 2$, $\theta_2 = 0.5$        * 경사 하강법이 진행되면서 (0,0)을 향해 직진하며 수렴.* 오른편: 비용함수 변화 그래프    * 상단: 라쏘 회귀 비용함수와 가중치 변화        * 등고선: 라쏘 회귀 비용함수. 색이 짙을 수록 낮은 값을 가리킴.        * $\theta$ 초기값: $\theta_1 = 0.25$, $\theta_2 = -1$        * 경사 하강법이 진행되면서  $\theta_2$가 먼저 0에 수렴.             이후 x축을 $\theta_1$이 전역  최솟점으로 수렴.             전역 최솟점 근처에서 진동함.        * $\alpha$가 감소하면 전역 최적점이 오른쪽으로 이동.             $\alpha$가 증가하면 전역 최적점이 왼쪽으로 이동.     * 하단: 릿지 회귀 비용함수와 가중치 변화        * 등고선: 릿지 회귀 비용함수. 색이 짙을 수록 낮은 값을 가리킴.        * $\theta$ 초기값: $\theta_1 = 0.25$, $\theta_2 = -1$        * 경사 하강법이 진행되면서  전역 최솟점으로 바로 수렴.            전역 최솟점에 가까워질 수록 경사 하강법이 느려지면서 수렴에 도움을 줌.        * $\alpha$가 감소하면 전역 최적점이 원점에서 멀어짐.            $\alpha$가 증가하면 전역 최적점이 원점으로 가까워짐.                 **참조**: 라쏘 규제를 사용할 경우 최적점 근처에서 진동을 줄이기 위해학습 스케줄을 사용하여 학습률을 점진적으로 줄이는 것을 추천한다.<jupyter_code>t1a, t1b = -1, 3
t2a, t2b = -1.5, 1.5

t1s = np.linspace(t1a, t1b, 500)
t2s = np.linspace(t2a, t2b, 500)
t1, t2 = np.meshgrid(t1s, t2s)

T = np.c_[t1.ravel(), t2.ravel()]
Xr = np.array([[1, 1], [1, -1], [1, 0.5]])
yr = 2 * Xr[:, :1] + 0.5 * Xr[:, 1:]

J = (1/len(Xr) * np.sum((T.dot(Xr.T) - yr.T)**2, axis=1)).reshape(t1.shape)

N1 = np.linalg.norm(T, ord=1, axis=1).reshape(t1.shape)
N2 = np.linalg.norm(T, ord=2, axis=1).reshape(t1.shape)

t_min_idx = np.unravel_index(np.argmin(J), J.shape)
t1_min, t2_min = t1[t_min_idx], t2[t_min_idx]

t_init = np.array([[0.25], [-1]])

def bgd_path(theta, X, y, l1, l2, core = 1, eta = 0.05, n_iterations = 200):
    path = [theta]
    for iteration in range(n_iterations):
        gradients = core * 2/len(X) * X.T.dot(X.dot(theta) - y) + \
                    l1 * np.sign(theta) + l2 * theta
        theta = theta - eta * gradients
        path.append(theta)
    return np.array(path)

fig, axes = plt.subplots(2, 2, sharex=True, sharey=True, figsize=(10.1, 8))

for i, N, l1, l2, title in ((0, N1, 2., 0, "Lasso"), (1, N2, 0,  2., "Ridge")):

    JR = J + l1 * N1 + l2 * 0.5 * N2**2
    
    tr_min_idx = np.unravel_index(np.argmin(JR), JR.shape)
    t1r_min, t2r_min = t1[tr_min_idx], t2[tr_min_idx]

    levelsJ=(np.exp(np.linspace(0, 1, 20)) - 1) * (np.max(J) - np.min(J)) + np.min(J)
    levelsJR=(np.exp(np.linspace(0, 1, 20)) - 1) * (np.max(JR) - np.min(JR)) + np.min(JR)
    levelsN=np.linspace(0, np.max(N), 10)
    
    path_J = bgd_path(t_init, Xr, yr, l1=0, l2=0)
    path_JR = bgd_path(t_init, Xr, yr, l1, l2)
    path_N = bgd_path(np.array([[2.0], [0.5]]), Xr, yr, np.sign(l1)/3, np.sign(l2), core=0)

    ax = axes[i, 0]
    ax.grid(True)
    ax.axhline(y=0, color='k')
    ax.axvline(x=0, color='k')
    ax.contourf(t1, t2, N / 2., levels=levelsN)
    ax.plot(path_N[:, 0], path_N[:, 1], "y--")
    ax.plot(0, 0, "ys")
    ax.plot(t1_min, t2_min, "ys")
    ax.set_title(r"$\ell_{}$ penalty".format(i + 1), fontsize=16)
    ax.axis([t1a, t1b, t2a, t2b])
    if i == 1:
        ax.set_xlabel(r"$\theta_1$", fontsize=16)
    ax.set_ylabel(r"$\theta_2$", fontsize=16, rotation=0)

    ax = axes[i, 1]
    ax.grid(True)
    ax.axhline(y=0, color='k')
    ax.axvline(x=0, color='k')
    ax.contourf(t1, t2, JR, levels=levelsJR, alpha=0.9)
    ax.plot(path_JR[:, 0], path_JR[:, 1], "w-o")
    ax.plot(path_N[:, 0], path_N[:, 1], "y--")
    ax.plot(0, 0, "ys")
    ax.plot(t1_min, t2_min, "ys")
    ax.plot(t1r_min, t2r_min, "rs")
    ax.set_title(title, fontsize=16)
    ax.axis([t1a, t1b, t2a, t2b])
    if i == 1:
        ax.set_xlabel(r"$\theta_1$", fontsize=16)

save_fig("lasso_vs_ridge_plot")
plt.show()<jupyter_output>그림 저장: lasso_vs_ridge_plot<jupyter_text>조기 종료 반복 훈련을 진행할 수록 훈련 세트에 대한 손실(RMSE 등)은 줄어든다. 모델이 훈련 데이터에 점점 더 익숙해지기 때문이다. 하지만 경우에 따라 과대적합이 발생하곤 한다.따라서 반복 훈련을 적절한 시기에 종료해야 하는데 그 기준을 검증 데이터에 대한 손실이 줄어 들다가 다시 커지는 순간으로 정할 수 있다. 이렇게 검증 오차가 최소에 다다렀을 때 반복 훈련을 멈추게 하는 기법을 **조기 종료**(early stoppy)이라 부른다. 다만, 확률적 경사 하강법, 미니배치 경사 하강법의 경우 손실 곡선이 매끄럽지 않고 진동이 발생할 수 있다.이런 경우에는 검증 오차가 한동안 최솟값보다 높게 유지될 때 반복 훈련을 멈추고검증 오차가 최소였을 때의 모델 파라미터를 확인하여 사용하면 된다.  예제: 조기 종료 구현 다음은 조기 종료를 구현하는 기본적인 방법을 보여주는 예제 코드이다.* 2차 다항 회귀 모델 적용을 위한 훈련 데이터 준비    * `m`: 전체 데이터셋 크기, 200    * 레이블 함수: $y = 2 + X + 1/2 X^2$ + 잡음    * 훈련 세트와 검증 세트를 반씩 나눔.<jupyter_code>np.random.seed(42)
m = 100
X = 6 * np.random.rand(m, 1) - 3
y = 2 + X + 0.5 * X**2 + np.random.randn(m, 1)

X_train, X_val, y_train, y_val = train_test_split(X[:50], y[:50].ravel(), test_size=0.5, random_state=10)<jupyter_output><empty_output><jupyter_text>먼저 훈련 데이터를 90차 다항 모델로 변환한다.<jupyter_code># 90차 다항 모델로 변환
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)<jupyter_output><empty_output><jupyter_text>규제 없는 확률적 경사 하강법 모델을 사용한다. 아래 코드는 최대 1000번의 에포크를 실행하면서 검증 오차가 최소가 되는 지점의 에포크와최상의 성능을 내는 학습 파라미터를 기억한다.즉, 최종적으로 최상의 검증 성능을 보이는 모델을 찾는다.* `warm_start=True`: 이전 에포크의 학습 결과를 이어 받아 학습함.     기본값은 `False`, 즉, 항상 새롭게 시작함.<jupyter_code>from sklearn.base import clone

# 확률적 경사 하강법 모델 설정. 규제 없음.
sgd_reg = SGDRegressor(max_iter=1, tol=-np.infty, warm_start=True,
                       penalty=None, learning_rate="constant", eta0=0.0005, random_state=42)

# 최대 1000번까지 전체 훈련 데이터 대상 학습을 반복하면서 평균 제곱근 오차(RMSE)가 작아질 때마다 
# 최소 RMSE 와 반복 횟수(에포크), 그리고 최고 성능의 모델 파라미터를 기억한다. 
# 최종적으로 가장 성능이 좋은 훈련 세트를 기억함. 
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)<jupyter_output><empty_output><jupyter_text>최상의 검증 모델을 찾기 위해 필요한 반복횟수(에포크)와 최상의 검증 모델은 다음과 같다.즉, 전체 훈련 세트를 대상으로 239번 반복 학습을 할 때 최상의 검증 성능이 나왔다.<jupyter_code>best_epoch, best_model<jupyter_output><empty_output><jupyter_text>이제, 반복 훈련을 최대 500번까지 할 때 훈련 세트 대상 오차(RMSE)와 검증 세트 대상 오차(RMSE)의 변화를그래프로 확인한다. 앞서 확인하였듯이 검증 성능이 가장 좋을 때는 239번 반복하였을 때이며,이 정보를 `best_epoch` 변수가 담고 있다.<jupyter_code># 확률적 경사 하강법 모델 설정. 규제 없음.
sgd_reg = SGDRegressor(max_iter=1, tol=-np.infty, warm_start=True,
                       penalty=None, learning_rate="constant", eta0=0.0005, random_state=42)

# 앞서 확인하였듯이 239번 반복하면 최상의 검증 성능이 나옴. 
# 따라서 500번 까지만 반복 훈련함.
n_epochs = 500

# 훈련 오차와, 검증 오차 기록 장치 
train_errors, val_errors = [], []

# 최대 500번 반복 훈련하면 훈련 오차, 검증 오차를 계산하여 저장
for epoch in range(n_epochs):
    sgd_reg.fit(X_train_poly_scaled, y_train)
    y_train_predict = sgd_reg.predict(X_train_poly_scaled)
    y_val_predict = sgd_reg.predict(X_val_poly_scaled)
    train_errors.append(mean_squared_error(y_train, y_train_predict))
    val_errors.append(mean_squared_error(y_val, y_val_predict))

# 검증 오차가 최소가 되는 지점 확인
best_epoch = np.argmin(val_errors)
# 검증 오차 최솟값 확인
best_val_rmse = np.sqrt(val_errors[best_epoch])

# 베스트 모델 지점 화살표 그리기
plt.annotate('Best model',
             xy=(best_epoch, best_val_rmse),
             xytext=(best_epoch, best_val_rmse + 1),
             ha="center",
             arrowprops=dict(facecolor='black', shrink=0.05),
             fontsize=16,
            )

# 최소 RMSE 값 약간 수정.
# 검은색 점선을 검증 오차 그래프에 접선으로 표시해서
# 보다 나은 그래프를 그리기 위한 것임.
best_val_rmse -= 0.03 

# 훈련 오차와 검증 오차의 변화 그래프 그리기
plt.plot([0, n_epochs], [best_val_rmse, best_val_rmse], "k:", linewidth=2)
plt.plot(np.sqrt(val_errors), "b-", linewidth=3, label="Validation set")
plt.plot(np.sqrt(train_errors), "r--", linewidth=2, label="Training set")
plt.legend(loc="upper right", fontsize=14)
plt.xlabel("Epoch", fontsize=14)
plt.ylabel("RMSE", fontsize=14)
save_fig("early_stopping_plot")
plt.show()<jupyter_output>그림 저장: early_stopping_plot<jupyter_text>핵심 6: 로지스트힉 회귀와 소프트맥스 회귀 로지스틱 회귀와 소프트맥스 회귀를 이용하여 분류를 위한 모델을 학습하는 과정을 살펴본다. * 이진 분류: 로지스틱 회귀 활용* 다중 클래스 분류: 소프트맥스 회귀 활용  로지스틱 회귀와 시그모이드 함수 로지스틱 회귀는 특성과 가중치의 곱한 값들을 더한 결과에 시그모이드 함수를 적용한 결과를 이용한다. 즉, 로지스틱 회귀 모델에서 샘플 $\mathbf x$에 대한 예측값은 다음과 같다.$$\hat p = h_\theta (\mathbf x) = \sigma(\theta^T \, \mathbf{x}_b^T) $$여기서 $\mathbf{x}_b = \begin{bmatrix} 1 & x_1 & \cdots & x_n\end{bmatrix}$이며, **시그모이드** 함수 $\sigma$는 정의는 다음과 같다.$$\sigma(t) = \frac{1}{1 + e^{-t}}$$시그모이드 함수의 그래프는 다음과 같다.<jupyter_code>t = np.linspace(-10, 10, 100)
sig = 1 / (1 + np.exp(-t))
plt.figure(figsize=(9, 3))
plt.plot([-10, 10], [0, 0], "k-")
plt.plot([-10, 10], [0.5, 0.5], "k:")
plt.plot([-10, 10], [1, 1], "k:")
plt.plot([0, 0], [-1.1, 1.1], "k-")
plt.plot(t, sig, "b-", linewidth=2, label=r"$\sigma(t) = \frac{1}{1 + e^{-t}}$")
plt.xlabel("t")
plt.legend(loc="upper left", fontsize=20)
plt.axis([-10, 10, -0.1, 1.1])
save_fig("logistic_function_plot")
plt.show()<jupyter_output>그림 저장: logistic_function_plot<jupyter_text>시그모이드 함수의 특징을 정리하면 다음과 같다. * $0< \sigma(t) < 1$* $\lim_{t \to \infty} \sigma(t) = 1$* $\lim_{t \to -\infty} \sigma(t) = -1$이제 결과 $\hat p$ 를 이용하여 로지스틱 회귀 모델의 예측값은 다음과 같이 정한다. $$\hat y = \begin{cases}0 & \text{if}\,\, \hat p < 0.5 \\1 & \text{if}\,\, \hat p \ge 0.5\end{cases}$$즉, 로지스틱 회귀는 $\mathbf{x}_b^T$가 0보다 같거나 클 때 1(양성 클래스), 작을 때 0(음성 클래스)으로 예측한다.   로그 손실 함수 로지스틱 회귀 모델의 학습 또한 경사하강법을 사용한다. 다만 사용되는 비용 함수가 선형 회귀의 그것과 많이 다르며, **로그 손실**(log loss) 함수라 부른다.$$J(\theta) = - \frac{1}{m}\, \sum_{i=0}^{m-1}\, [y^{(i)}\, \log(\,\hat p^{(i)}\,) + (1-y^{(i)})\, \log(\,1 - \hat p^{(i)}\,)]$$* $y^{(i)}$: $i$ 번째 샘플의 클래스 정보 레이블    * 양성 클래스에 속하면 1    * 음성 클래스에 속하면 0* $\hat{p}^{(i)}$: $i$ 번째 샘플이 양성 클래스 속할 확률 예측값로그 손실 함수를 이해하려면 하나의 샘플에 대한 아래의 값의 의미를 이해하면 된다.$$- [y^{(i)}\, \log(\,\hat p^{(i)}\,) + (1-y^{(i)})\, \log(\,1 - \hat p^{(i)}\,)]$$위 값은 틀린 예측을 하면 커지도록 설정되었으며, 이에 대한 설명은 다음과 같다. 먼저 로그 함수의 다음 성질에 주목한다. $$\lim_{t\to +0}\log(t) = -\infty$$이제 두 가지 틀린 예측을 하는 경우에 비용이 커지는 이유를 이해할 수 있다.* $y^{(i)} = 0$ 이지만 $\hat p^{(i)} \ge 0.5$ 인 경우:     $-(1-y^{(i)})\, \log(\,1 - \hat p^{(i)}\,)$ 가 커진다.* $y^{(i)} = 1$ 이지만 $\hat p^{(i)} < 0.5$ 인 경우:     $-y^{(i)}\, \log(\,\hat p^{(i)}\,)$ 가 커진다.**참고**: 로그 함수의 그래프는 다음과 같이 그려진다.<jupyter_code>t = np.linspace(1e-2, 4, 100)
log_t = np.log(t)
plt.figure(figsize=(9, 3))
plt.plot([0, 0], [-5, 2], "k-")
plt.plot([-0.1, 4], [0, 0], "k-")
plt.plot([1, 1], [-4.8, 2], "r:")

plt.plot(t, log_t, "b-", linewidth=2, label=r"$\log(t)$")
plt.legend(loc="lower right", fontsize=20)
plt.xlabel("t")
save_fig("logarithmic function")
plt.show()<jupyter_output>그림 저장: logarithmic function<jupyter_text>비용 함수의 최솟값 구하기 로지스틱 회귀 비용 함수의 최솟값을 위한 정규방정식은 알려지지 않았다.하지만 그레이디언트 벡터를 계산할 수 있으며 이를 이용하여 경사 하강법을 적용한다.$$\frac{\partial}{\partial \theta_j}\, J(\theta) =\frac{1}{m}\, \sum_{i=0}^{m-1} \, (\sigma(\theta^T\, (\mathbf{x}_b^{(i)})^T) - y^{(i)})\, x_j^{(i)}$$단, $x_0^{(i)} = 1$.  결정 경계 로지스틱 회귀의 작동법을 설명하기 위해 붓꽃 데이터셋을 사용한다. 붓꽃 데이터셋은 세 개의 품종과 관련된 총 150개의 샘플로 이루어진다.* Iris-Setosa* Iris-Versicolor* Iris-Virginica  붗꽃 데이터셋 가져오기 사이킷런에서 붓꽃 데이터셋을 제공한다.<jupyter_code>from sklearn import datasets
iris = datasets.load_iris()<jupyter_output><empty_output><jupyter_text>붓꽃 데이터셋은 사전 자료형으로 가져온다. 사용된 키를 확인하면 다음과 같다.<jupyter_code>list(iris.keys())<jupyter_output><empty_output><jupyter_text>데이터셋에 대한 설명은 `DESCR`의 키값을 확인하면 되며,내용을 정리하면 다음과 같다.* 데이터셋 크기: 150* 샘플 특성: 총 4개의 수치형 값으로 이루어짐. (`data`의 키값)    - 꽃받침(sepal) 길이(단위 cm)    - 꽃받침(sepal) 너비(단위 cm)    - 꽃잎(petal) 길이(단위 cm)    - 꽃잎(petal) 너비(단위 cm)* 샘플 타깃 (`target`의 키값)    - 0: Iris-Setosa 클래스    - 1: Iris-Versicolor 클래스    - 2: Iris-Virginica 클래스<jupyter_code>print(iris.DESCR)<jupyter_output>.. _iris_dataset:

Iris plants dataset
--------------------

**Data Set Characteristics:**

    :Number of Instances: 150 (50 in each of three classes)
    :Number of Attributes: 4 numeric, predictive attributes and the class
    :Attribute Information:
        - sepal length in cm
        - sepal width in cm
        - petal length in cm
        - petal width in cm
        - class:
                - Iris-Setosa
                - Iris-Versicolour
                - Iris-Virginica
                
    :Summary Statistics:

    ============== ==== ==== ======= ===== ====================
                    Min  Max   Mean    SD   Class Correlation
    ============== ==== ==== ======= ===== ====================
    sepal length:   4.3  7.9   5.84   0.83    0.7826
    sepal width:    2.0  4.4   3.05   0.43   -0.4194
    petal length:   1.0  6.9   3.76   1.76    0.9490  (high!)
    petal width:    0.1  2.5   1.20   0.76    0.9565  (high!)
    ============== ==== ==== ======= ===== ===========[...]<jupyter_text>예제 붓꽃 데이터와 로지스틱 회귀를 이용하여 이진 분류기를 학습시켜 보자.* 사용 특성: 꽃잎의 너비(`data` 값의 3번 열에 해당)* Iris-Versicolor 종인지 여부 판단(`target` 값의 1에 해당)훈련 세트를 아래와 같이 준비한다.<jupyter_code># 꽃잎 너비만 특성으로 사용
X = iris["data"][:, 3:]

# Iris virginica이면 1 아니면 0
y = (iris["target"] == 2).astype(np.int)<jupyter_output><empty_output><jupyter_text>이제 사이킷런의 `LogisticRegression` 예측기를 이용하여 모델을 학습한다. **참조**: 향후 버전이 바뀌더라도 동일한 결과를 만들기 위해 사이킷런 0.22 버전의 기본값인 `solver="lbfgs"`로 지정한다.<jupyter_code>from sklearn.linear_model import LogisticRegression
log_reg = LogisticRegression(solver="lbfgs", random_state=42)
log_reg.fit(X, y)<jupyter_output><empty_output><jupyter_text>`predict_proba()` 메서드는 음성과 양성에 대한 예측확률을 리턴한다. * 음성: Iris-Virginica 아님* 양성: Iris-Virginica 맞음아래 코드는 꽃잎의 너비가 0에서 3cm 사이인 붓꽃 1000개에 대한 모데릐 추정 확률을 보여주는 그래프를 그린다.<jupyter_code>X_new = np.linspace(0, 3, 1000).reshape(-1, 1)
y_proba = log_reg.predict_proba(X_new)

plt.plot(X_new, y_proba[:, 1], "g-", linewidth=2, label="Iris virginica")
plt.plot(X_new, y_proba[:, 0], "b--", linewidth=2, label="Not Iris virginica")<jupyter_output><empty_output><jupyter_text>위 그림을 책에서는 보다 멋지게 꾸몄다. 코드는 아래와 같으며, Iris-Virginica일 확률이 50%가 되도록 하는 **결정 경계**(`decision_bounday`)도 함께 계산하였다.<jupyter_code>X_new = np.linspace(0, 3, 1000).reshape(-1, 1)
y_proba = log_reg.predict_proba(X_new)
decision_boundary = X_new[y_proba[:, 1] >= 0.5][0]

plt.figure(figsize=(8, 3))
plt.plot(X[y==0], y[y==0], "bs")
plt.plot(X[y==1], y[y==1], "g^")
plt.plot([decision_boundary, decision_boundary], [-1, 2], "k:", linewidth=2)
plt.plot(X_new, y_proba[:, 1], "g-", linewidth=2, label="Iris virginica")
plt.plot(X_new, y_proba[:, 0], "b--", linewidth=2, 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.xlabel("Petal width (cm)", fontsize=14)
plt.ylabel("Probability", fontsize=14)
plt.legend(loc="center left", fontsize=14)
plt.axis([0, 3, -0.02, 1.02])
save_fig("logistic_regression_plot")
plt.show()<jupyter_output>그림 저장: logistic_regression_plot<jupyter_text>결정 경계는 약 1.66cm이다.<jupyter_code>decision_boundary<jupyter_output><empty_output><jupyter_text>실제로 Iris-Virginica의 꽃잎 너비는 1.4에서 2.5cm 정도이다. 반면에 다른 두 붓꽃 품종은 꽃잎 너비가 0.1에서 1.8cm 정도이며,따라서 약간 중첩되는 구간이 존재한다.하지만 `predict()` 메서드는 단순히 결정 경계를 이용하여 클래스를 분류한다.<jupyter_code>log_reg.predict([[1.67], [1.66]])<jupyter_output><empty_output><jupyter_text>하지만 꽃잎 너비가 1.66 근처인 경우 품종을 잘못 분류할 가능성이 높다. 아래 코드는 양성일 확률을 기준으로 결정 경계 근처에서 어떻게 품종이 분류되는지를 그래프로 보여준다. 그래프는 꽃잎 너비와 꽃잎 길이를 기준으로 그려졌다.확률 경계는 선형으로 그려지며, 사용된 양성일 확률은 15%, 30%, 45%. 60%, 75% 90%를 대상으로 한다.두꺼운 점선은 50% 경계, 즉 결정 경계를 나타낸다.* 15% 경계: 왼쪽 아래, 보라색 선    * 경계선 위쪽으로 15% 이상의 확률로 Iris-Virginica임을 확신한다는 의미.         즉, 믿음이 매우 약함.* 90% 경계: 오른쪽 위, 초록색 선    * 경계선 위쪽으로는 90% 이상의 확률로 Iris-Virginica임을 확신한다는 의미        즉, 믿음이 매우 강함.<jupyter_code>from sklearn.linear_model import LogisticRegression

X = iris["data"][:, (2, 3)]  # petal length, petal width
y = (iris["target"] == 2).astype(np.int)

log_reg = LogisticRegression(solver="lbfgs", C=10**10, random_state=42)
log_reg.fit(X, y)

x0, x1 = np.meshgrid(
        np.linspace(2.9, 7, 500).reshape(-1, 1),
        np.linspace(0.8, 2.7, 200).reshape(-1, 1),
    )
X_new = np.c_[x0.ravel(), x1.ravel()]

y_proba = log_reg.predict_proba(X_new)

plt.figure(figsize=(10, 4))
plt.plot(X[y==0, 0], X[y==0, 1], "bs")
plt.plot(X[y==1, 0], X[y==1, 1], "g^")

zz = y_proba[:, 1].reshape(x0.shape)
contour = plt.contour(x0, x1, zz, cmap=plt.cm.brg)


left_right = np.array([2.9, 7])
boundary = -(log_reg.coef_[0][0] * left_right + log_reg.intercept_[0]) / log_reg.coef_[0][1]

plt.clabel(contour, inline=1, fontsize=12)
plt.plot(left_right, boundary, "k--", linewidth=3)
plt.text(3.5, 1.5, "Not I<NAME>", fontsize=14, color="b", ha="center")
plt.text(6.5, 2.3, "<NAME>", fontsize=14, color="g", ha="center")
plt.xlabel("Petal length", fontsize=14)
plt.ylabel("Petal width", fontsize=14)
plt.axis([2.9, 7, 0.8, 2.7])
save_fig("logistic_regression_contour_plot")
plt.show()<jupyter_output>그림 저장: logistic_regression_contour_plot<jupyter_text>로지스틱 회귀 규제하기 로지스틱 회귀 모델에 `penalty`와 `C` 두 옵션 변수의 값을 지정하는 방식으로 규제를 가할 수 있다.* `penalty`: `l1`, `l2`, `elasticnet` 세 개중에 하나를 사용하며,    기본은 `l2`, 즉, $\ell_2$ 규제를 사용하는 릿지 규제이다.    * `elasticnet`을 선택한 경우 `l1_ration` 옵션 값을 지정해서 함께 사용한다.* `C`: 양의 실수값을 가지며, $\alpha$와 달리 0에 가까울 수록 강한 규제를 의미한다.  소프트맥스 회귀 로지스틱 회귀 모델을 일반화하여 다중 클래스 분류를 지원하도록 한 것이 **소프트맥스 회귀**(softmax regression)이다.소프트맥스 회귀는 **다항 로지스틱 회귀**(multinomial loogistic regression)라고도 부른다.  소프트맥스 회귀 개념 소프트맥스 회귀를 이용한 예측은 다음 세 단계를 거쳐서 진행된다.각 단계 설명에서 사용되는 변수들의 의미는 다음과 같다.* $\mathbf X$: 전체 훈련 세트* $\mathbf x$: 훈련 샘플* $m$: 전체 훈련 세트 크기* $n$: 특성 수* $k$: 클래스 $k$* $K$: 클래스 수 *단계 1: 소프트맥스 점수 계산* 샘플 $\mathbf x$에 대해 $s_k(\mathbf x)$를 계산한다.$s_k(\mathbf x)$의 정의는 다음과 같다.$$s_k(\mathbf x) = (\theta^{(k)})^T\, \mathbf{x}_b^T$$    * $\theta^{(k)}$: 클래스 $k$에 대한 가중치 행렬. 특성 수가 n 이라면 (n+1, 1) 모양의 행렬이다.     맨 처음에 설명했던 가중치 벡서 $\theta$와 동일한 개념이며,     여기서는 각 클래스별로 따로 가중치를 사용한다는 점만 다르다.    $$\begin{bmatrix}\theta_0^{(k)} \\\theta_1^{(k)} \\\vdots \\\theta_n^{(k)}\end{bmatrix}$$      * 전체 클래스를 대상으로 하는 가중치 벡터들의 모음을 $\Theta$로 표기한다.     $\Theta$는 (n+1, K) 모양의 행렬로 다루어진다.    $$\begin{bmatrix}\theta^{(0)} &\theta^{(1)} &\cdots &\theta^{(K-1)}\end{bmatrix}$$      * 이제 각 샘플에 대해 클래스별 점수를 담은 $\mathbf s(\mathbf x)$를 다음처럼 정의할 수 있다.    $$\mathbf s(\mathbf x) = \Theta^T\, \mathbf x_b^T$$          $\mathbf s(\mathbf x)$은 (K, 1) 모양의 행렬이다.     * 샘플 전체를 대상으로 하는 클래스별 점수를 담은 $\mathbf s(\mathbf X)$은 따라서 다음과 같다.    $$\mathbf s(\mathbf X) = \Theta^T\, \mathbf X_b^T$$          $\mathbf s(\mathbf X)$은 (K, m) 모양의 행렬이다.  *단계 2: 추정 확률 계산* 샘플 별로 소프트맥스 점수를 계산한 후에 **소프트맥스 함수**(softmax function) $\sigma()$를 이용하여각각의 클래스에 속할 확률 $\hat p_k$를 계산한다.$$\hat p_k = \sigma(\mathbf s(\mathbf x))_k =\frac{\exp(s_k(\mathbf x))}{\sum_{j=0}^{K-1}\exp(s_j(\mathbf x))}$$ *단계 3: 클래스 선택* 추정 확률이 가장 높은 클래스를 선택한다.즉, 샘플 $\mathbf x$에 대한 예측값 $\hat y$는 아래와 같이 정해진다.$$\hat y = \mathrm{argmax}_k \sigma(\mathbf s(\mathbf x))_k =\mathrm{argmax}_k s_k(\mathbf x) =\mathrm{argmax}_k ((\theta^{(k)})^T\, \mathbf{x}_b^T)$$`argmax`는 최대의 값을 갖는 인덱스를 확인한다.여기서는 추정 확률이 가장 높은 클래스를 선택한다.**참조**: 소프트맥스 회귀는 다중 출력 분류를 지원하지 못한다. 즉, 하나의 샘플이 주어졌을 때, 서로 배타적인 클래스로 이루어진 클래스들 중에서 하나의 클래스만 선택한다. 따라서 예를 들어, 하나의 사진에서 여러 사람의 얼굴을 인식하는 데에는 사용하지 못한다.  크로스 엔트로피 비용 함수 계산* 최적의 소프트맥스 회귀 모델을 경사 하강법을 이용하여 찾을 수 있다.학습 과정에 사용되는 비용 함수로는 다음과 같이 정의된 **크로스 엔트로피 비용 함수**를 사용한다. $$J(\Theta) = - \frac{1}{m}\, \sum_{i=0}^{m-1}\sum_{k=0}^{K-1} y^{(i)}_k\, \log(\hat{p}_k^{(i)})$$* $y^{(i)}_k$: $i$ 번째 샘플의 클래스 $k$ 관련 정보 레이블    * 클래스 $k$에 속하면 1    * 클래스 $k$에 속하지 않으면 0* $\hat{p}_k^{(i)}$: $i$ 번째 샘플이 클래스 $k$에 속할 확률 예측값**참조** * $K=2$이면 로지스틱 회귀의 로그 손실 함수와 정확하게 일치한다.* 크로스 엔트로피 개념은 정보 이론에서 유래하였다. (자세한 설명은 생략)  비용 함수의 그레이디어트 벡터 각각의 클래스에 대한 그레이디언트 벡터를 계산할 수 있다.예를 들어, 클래스 $k$에 대한 그레이디언트 벡터는 다음과 같다.$$\nabla_{\theta^{(k)}}\, J(\Theta) =\frac{1}{m}\, \sum_{i=0}^{m-1} \, (\hat p_k^{(i)} - y_k^{(i)})\, (\mathbf{x}_b^{(i)})^T$$따라서 최적화된 파라미터 행렬 $\Theta$를구하기 위해 경사 하강법 등 최적화 알고리즘을 사용할 수 있다.  예제: 멀티 클래스 분류 소프트맥스 회귀를 이용하여 붓꽃을 세 개의 클래스로 분류해보자.사이킷런의 `LogisticRegression` 예측기의 `multi_class` 옵션 변수의 값으로 `multinomial`를사용하면 소프트맥스 회귀를 사용할 수 있다.* `multi_class`: `auto`, `ovr`, `multinomial` 사용 가능. 기본값은 `auto`    * `ovr`: 로지스틱 회귀의 이진 분류 기능을 일대다 전략을 활용하여 다중 클래스 분류 지원    * `multinomial`: 다중 클래스 분류 지원        * `solver` 설정에서 `lbfgs`, `sag`, `saga` 등을 선택해야 함. 기본값은 `lbfgs`    * `auto`: 분류 방법에 따라 분류 모델 자동 선택* `penalty`와 `C`에 대한 설명은 앞서 언급한 규제 로지스틱 회귀에 대한 설명과 동일아래 코드는 꽃잎 길이와 너비를 이용하여 세 개의 클래스로 분류하는 소프트맥스 회귀 모델이다.즉, 특성 수 $n=2$이다.<jupyter_code># 꽃잎 길이, 꽃잎 너비 정보는 `data` 키의 값에서 2번, 3번 열에 들어 있다.
X = iris["data"][:, (2, 3)]
# 타깃
y = iris["target"]

# 소프트맥스 회귀 모델 학습
softmax_reg = LogisticRegression(multi_class="multinomial",
                                 solver="lbfgs", C=10, random_state=42)
softmax_reg.fit(X, y)<jupyter_output><empty_output><jupyter_text>학습을 통해 계산된 최적화된 특성들의 가중치는 `coef_` 속성에 들어 있다.**참조**: $\Theta^T$의 0번 열을 제외한 클래스별 가중치 벡터로 이루어진 (K,n) 모양의 어레이가 해당한다.여기서는 (3,2) 모양의 2차원 어레이로 표시된다.<jupyter_code>softmax_reg.coef_<jupyter_output><empty_output><jupyter_text>최적화된 편향 벡터는 `intercept_` 속성에 저장된다.**참조**: $\Theta$의 0번 행, $\Theta^T$의 0번 열에 해당한다.여기서는 (1,3) 모양의 행렬, 즉 길이가 3인 벡터로 표시된다.<jupyter_code>softmax_reg.intercept_<jupyter_output><empty_output><jupyter_text>붓꽃 데이터 전체를 대상으로 분류한 결과는 다음과 같다.(꽃잎 길이와 너비 기준으로 분류하였다.)<jupyter_code>x0, x1 = np.meshgrid(
        np.linspace(0, 8, 500).reshape(-1, 1),
        np.linspace(0, 3.5, 200).reshape(-1, 1),
    )
X_new = np.c_[x0.ravel(), x1.ravel()]


y_proba = softmax_reg.predict_proba(X_new)
y_predict = softmax_reg.predict(X_new)

zz1 = y_proba[:, 1].reshape(x0.shape)
zz = y_predict.reshape(x0.shape)

plt.figure(figsize=(10, 4))
plt.plot(X[y==2, 0], X[y==2, 1], "g^", label="Iris virginica")
plt.plot(X[y==1, 0], X[y==1, 1], "bs", label="Iris versicolor")
plt.plot(X[y==0, 0], X[y==0, 1], "yo", label="Iris setosa")

from matplotlib.colors import ListedColormap
custom_cmap = ListedColormap(['#fafab0','#9898ff','#a0faa0'])

plt.contourf(x0, x1, zz, cmap=custom_cmap)
contour = plt.contour(x0, x1, zz1, cmap=plt.cm.brg)
plt.clabel(contour, inline=1, fontsize=12)
plt.xlabel("Petal length", fontsize=14)
plt.ylabel("Petal width", fontsize=14)
plt.legend(loc="center left", fontsize=14)
plt.axis([0, 7, 0, 3.5])
save_fig("softmax_regression_contour_plot")
plt.show()<jupyter_output>그림 저장: softmax_regression_contour_plot<jupyter_text>위 그래프에서 사용된 결정 경계는 세 개의 클래스에 대해 선형으로 그려졌다.또한 Iris-Versicolor 클래스에 대한 확률을 곡선으로 나타냈다. 경우에 따라 모든 클래스에 대한 추정 확률이 50% 미만인 경우도 발생한다. 위그림에서 모든 결정 경계가 만나는 지점의 경우 모든 클래스에 대해 추정 확률이 33% 정도일 것이다. 그럼에도 불구하고 클래스는 선택된다. 예를 들어, 꽃잎의 길이가 2.6cm이고, 너비는 3cm인 붓꽃에 대한 예측은 Iris-Versicolor 이다.<jupyter_code>softmax_reg.predict([[2.6, 3]])<jupyter_output><empty_output><jupyter_text>이유는 `predict_proba()` 메서드가 예측하는 클래스별 추정 확률이 다음과 같기 때문인데,각 클래스에 대한 추정 확률이 모두 0.5보다 작다.<jupyter_code>softmax_reg.predict_proba([[2.6, 3]])<jupyter_output><empty_output><jupyter_text>연습문제 해답  1. to 11. 부록 A를 참고하세요.  12. 조기 종료를 사용한 배치 경사 하강법으로 소프트맥스 회귀 구현하기(사이킷런을 사용하지 않고) 먼저 데이터를 로드합니다. 앞서 사용했던 Iris 데이터셋을 재사용하겠습니다.<jupyter_code>X = iris["data"][:, (2, 3)]  # 꽃잎 길이, 꽃잎 넓이
y = iris["target"]<jupyter_output><empty_output><jupyter_text>모든 샘플에 편향을 추가합니다 ($x_0 = 1$):<jupyter_code>X_with_bias = np.c_[np.ones([len(X), 1]), X]<jupyter_output><empty_output><jupyter_text>결과를 일정하게 유지하기 위해 랜덤 시드를 지정합니다:<jupyter_code>np.random.seed(2042)<jupyter_output><empty_output><jupyter_text>데이터셋을 훈련 세트, 검증 세트, 테스트 세트로 나누는 가장 쉬운 방법은 사이킷런의 `train_test_split()` 함수를 사용하는 것입니다. 하지만 이 연습문제의 목적은 직접 만들어 보면서 알고리즘을 이해하는 것이므로 다음과 같이 수동으로 나누어 보겠습니다:<jupyter_code>test_ratio = 0.2
validation_ratio = 0.2
total_size = len(X_with_bias)

test_size = int(total_size * test_ratio)
validation_size = int(total_size * validation_ratio)
train_size = total_size - test_size - validation_size

rnd_indices = np.random.permutation(total_size)

X_train = X_with_bias[rnd_indices[:train_size]]
y_train = y[rnd_indices[:train_size]]
X_valid = X_with_bias[rnd_indices[train_size:-test_size]]
y_valid = y[rnd_indices[train_size:-test_size]]
X_test = X_with_bias[rnd_indices[-test_size:]]
y_test = y[rnd_indices[-test_size:]]<jupyter_output><empty_output><jupyter_text>타깃은 클래스 인덱스(0, 1 그리고 2)이지만 소프트맥스 회귀 모델을 훈련시키기 위해 필요한 것은 타깃 클래스의 확률입니다. 각 샘플에서 확률이 1인 타깃 클래스를 제외한 다른 클래스의 확률은 0입니다(다른 말로하면 주어진 샘플에 대한 클래스 확률이 원-핫 벡터입니다). 클래스 인덱스를 원-핫 벡터로 바꾸는 간단한 함수를 작성하겠습니다:<jupyter_code>def to_one_hot(y):
    n_classes = y.max() + 1
    m = len(y)
    Y_one_hot = np.zeros((m, n_classes))
    Y_one_hot[np.arange(m), y] = 1
    return Y_one_hot<jupyter_output><empty_output><jupyter_text>10개 샘플만 넣어 이 함수를 테스트해 보죠:<jupyter_code>y_train[:10]
to_one_hot(y_train[:10])<jupyter_output><empty_output><jupyter_text>잘 되네요, 이제 훈련 세트와 테스트 세트의 타깃 클래스 확률을 담은 행렬을 만들겠습니다:<jupyter_code>Y_train_one_hot = to_one_hot(y_train)
Y_valid_one_hot = to_one_hot(y_valid)
Y_test_one_hot = to_one_hot(y_test)<jupyter_output><empty_output><jupyter_text>이제 소프트맥스 함수를 만듭니다. 다음 공식을 참고하세요:$\sigma\left(\mathbf{s}(\mathbf{x})\right)_k = \dfrac{\exp\left(s_k(\mathbf{x})\right)}{\sum\limits_{j=1}^{K}{\exp\left(s_j(\mathbf{x})\right)}}$<jupyter_code>def softmax(logits):
    exps = np.exp(logits)
    exp_sums = np.sum(exps, axis=1, keepdims=True)
    return exps / exp_sums<jupyter_output><empty_output><jupyter_text>훈련을 위한 준비를 거의 마쳤습니다. 입력과 출력의 개수를 정의합니다:<jupyter_code>n_inputs = X_train.shape[1] # == 3 (특성 2개와 편향)
n_outputs = len(np.unique(y_train))   # == 3 (3개의 붓꽃 클래스)<jupyter_output><empty_output><jupyter_text>이제 좀 복잡한 훈련 파트입니다! 이론적으로는 간단합니다. 그냥 수학 공식을 파이썬 코드로 바꾸기만 하면 됩니다. 하지만 실제로는 꽤 까다로운 면이 있습니다. 특히, 항이나 인덱스의 순서가 뒤섞이기 쉽습니다. 제대로 작동할 것처럼 코드를 작성했더라도 실제 제대로 계산하지 못합니다. 확실하지 않을 때는 각 항의 크기를 기록하고 이에 상응하는 코드가 같은 크기를 만드는지 확인합니다. 각 항을 독립적으로 평가해서 출력해 보는 것도 좋습니다. 사실 사이킷런에 이미 잘 구현되어 있기 때문에 이렇게 할 필요는 없습니다. 하지만 직접 만들어 보면 어떻게 작동하는지 이해하는데 도움이 됩니다.구현할 공식은 비용함수입니다:$J(\mathbf{\Theta}) =- \dfrac{1}{m}\sum\limits_{i=1}^{m}\sum\limits_{k=1}^{K}{y_k^{(i)}\log\left(\hat{p}_k^{(i)}\right)}$그리고 그레이디언트 공식입니다:$\nabla_{\mathbf{\theta}^{(k)}} \, J(\mathbf{\Theta}) = \dfrac{1}{m} \sum\limits_{i=1}^{m}{ \left ( \hat{p}^{(i)}_k - y_k^{(i)} \right ) \mathbf{x}^{(i)}}$$\hat{p}_k^{(i)} = 0$이면 $\log\left(\hat{p}_k^{(i)}\right)$를 계산할 수 없습니다. `nan` 값을 피하기 위해 $\log\left(\hat{p}_k^{(i)}\right)$에 아주 작은 값 $\epsilon$을 추가하겠습니다.<jupyter_code>eta = 0.01
n_iterations = 5001
m = len(X_train)
epsilon = 1e-7

Theta = np.random.randn(n_inputs, n_outputs)

for iteration in range(n_iterations):
    logits = X_train.dot(Theta)
    Y_proba = softmax(logits)
    loss = -np.mean(np.sum(Y_train_one_hot * np.log(Y_proba + epsilon), axis=1))
    error = Y_proba - Y_train_one_hot
    if iteration % 500 == 0:
        print(iteration, loss)
    gradients = 1/m * X_train.T.dot(error)
    Theta = Theta - eta * gradients<jupyter_output>0 5.446205811872683
500 0.8350062641405651
1000 0.6878801447192402
1500 0.6012379137693313
2000 0.5444496861981872
2500 0.5038530181431525
3000 0.4729228972192248
3500 0.44824244188957774
4000 0.4278651093928793
4500 0.41060071429187134
5000 0.3956780375390374<jupyter_text>바로 이겁니다! 소프트맥스 모델을 훈련시켰습니다. 모델 파라미터를 확인해 보겠습니다:<jupyter_code>Theta<jupyter_output><empty_output><jupyter_text>검증 세트에 대한 예측과 정확도를 확인해 보겠습니다:<jupyter_code>logits = X_valid.dot(Theta)
Y_proba = softmax(logits)
y_predict = np.argmax(Y_proba, axis=1)

accuracy_score = np.mean(y_predict == y_valid)
accuracy_score<jupyter_output><empty_output><jupyter_text>와우, 이 모델이 매우 잘 작동하는 것 같습니다. 연습을 위해서 $\ell_2$ 규제를 조금 추가해 보겠습니다. 다음 코드는 위와 거의 동일하지만 손실에 $\ell_2$ 페널티가 추가되었고 그래디언트에도 항이 추가되었습니다(`Theta`의 첫 번째 원소는 편향이므로 규제하지 않습니다). 학습률 `eta`도 증가시켜 보겠습니다.<jupyter_code>eta = 0.1
n_iterations = 5001
m = len(X_train)
epsilon = 1e-7
alpha = 0.1  # 규제 하이퍼파라미터

Theta = np.random.randn(n_inputs, n_outputs)

for iteration in range(n_iterations):
    logits = X_train.dot(Theta)
    Y_proba = softmax(logits)
    xentropy_loss = -np.mean(np.sum(Y_train_one_hot * np.log(Y_proba + epsilon), axis=1))
    l2_loss = 1/2 * np.sum(np.square(Theta[1:]))
    loss = xentropy_loss + alpha * l2_loss
    error = Y_proba - Y_train_one_hot
    if iteration % 500 == 0:
        print(iteration, loss)
    gradients = 1/m * X_train.T.dot(error) + np.r_[np.zeros([1, n_outputs]), alpha * Theta[1:]]
    Theta = Theta - eta * gradients<jupyter_output>0 6.629842469083912
500 0.5339667976629507
1000 0.503640075014894
1500 0.4946891059460321
2000 0.4912968418075477
2500 0.48989924700933296
3000 0.48929905984511984
3500 0.4890351244397859
4000 0.4889173621830818
4500 0.48886433374493027
5000 0.48884031207388184<jupyter_text>추가된 $\ell_2$ 페널티 때문에 이전보다 손실이 조금 커보이지만 더 잘 작동하는 모델이 되었을까요? 확인해 보죠:<jupyter_code>logits = X_valid.dot(Theta)
Y_proba = softmax(logits)
y_predict = np.argmax(Y_proba, axis=1)

accuracy_score = np.mean(y_predict == y_valid)
accuracy_score<jupyter_output><empty_output><jupyter_text>와우, 완벽한 정확도네요! 운이 좋은 검증 세트일지 모르지만 잘 된 것은 맞습니다. 이제 조기 종료를 추가해 보죠. 이렇게 하려면 매 반복에서 검증 세트에 대한 손실을 계산해서 오차가 증가하기 시작할 때 멈춰야 합니다.<jupyter_code>eta = 0.1 
n_iterations = 5001
m = len(X_train)
epsilon = 1e-7
alpha = 0.1  # 규제 하이퍼파라미터
best_loss = np.infty

Theta = np.random.randn(n_inputs, n_outputs)

for iteration in range(n_iterations):
    logits = X_train.dot(Theta)
    Y_proba = softmax(logits)
    xentropy_loss = -np.mean(np.sum(Y_train_one_hot * np.log(Y_proba + epsilon), axis=1))
    l2_loss = 1/2 * np.sum(np.square(Theta[1:]))
    loss = xentropy_loss + alpha * l2_loss
    error = Y_proba - Y_train_one_hot
    gradients = 1/m * X_train.T.dot(error) + np.r_[np.zeros([1, n_outputs]), alpha * Theta[1:]]
    Theta = Theta - eta * gradients

    logits = X_valid.dot(Theta)
    Y_proba = softmax(logits)
    xentropy_loss = -np.mean(np.sum(Y_valid_one_hot * np.log(Y_proba + epsilon), axis=1))
    l2_loss = 1/2 * np.sum(np.square(Theta[1:]))
    loss = xentropy_loss + alpha * l2_loss
    if iteration % 500 == 0:
        print(iteration, loss)
    if loss < best_loss:
        best_loss = loss
    else:
        print(iteration - 1, best_loss)
        print(iteration, loss, "조기 종료!")
        break
logits = X_valid.dot(Theta)
Y_proba = softmax(logits)
y_predict = np.argmax(Y_proba, axis=1)

accuracy_score = np.mean(y_predict == y_valid)
accuracy_score<jupyter_output><empty_output><jupyter_text>여전히 완벽하지만 더 빠릅니다. 이제 전체 데이터셋에 대한 모델의 예측을 그래프로 나타내 보겠습니다:<jupyter_code>x0, x1 = np.meshgrid(
        np.linspace(0, 8, 500).reshape(-1, 1),
        np.linspace(0, 3.5, 200).reshape(-1, 1),
    )
X_new = np.c_[x0.ravel(), x1.ravel()]
X_new_with_bias = np.c_[np.ones([len(X_new), 1]), X_new]

logits = X_new_with_bias.dot(Theta)
Y_proba = softmax(logits)
y_predict = np.argmax(Y_proba, axis=1)

zz1 = Y_proba[:, 1].reshape(x0.shape)
zz = y_predict.reshape(x0.shape)

plt.figure(figsize=(10, 4))
plt.plot(X[y==2, 0], X[y==2, 1], "g^", label="Iris virginica")
plt.plot(X[y==1, 0], X[y==1, 1], "bs", label="Iris versicolor")
plt.plot(X[y==0, 0], X[y==0, 1], "yo", label="<NAME>")

from matplotlib.colors import ListedColormap
custom_cmap = ListedColormap(['#fafab0','#9898ff','#a0faa0'])

plt.contourf(x0, x1, zz, cmap=custom_cmap)
contour = plt.contour(x0, x1, zz1, cmap=plt.cm.brg)
plt.clabel(contour, inline=1, fontsize=12)
plt.xlabel("Petal length", fontsize=14)
plt.ylabel("Petal width", fontsize=14)
plt.legend(loc="upper left", fontsize=14)
plt.axis([0, 7, 0, 3.5])
plt.show()<jupyter_output><empty_output><jupyter_text>이제 테스트 세트에 대한 모델의 최종 정확도를 측정해 보겠습니다:<jupyter_code>logits = X_test.dot(Theta)
Y_proba = softmax(logits)
y_predict = np.argmax(Y_proba, axis=1)

accuracy_score = np.mean(y_predict == y_test)
accuracy_score<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.15080146049667245, 'fraction_numerical': 0.03344382007923971, 'mean_word_length': 3.392345313921747, 'pattern_counts': {'":': 0, '<': 348, '<?xml version=': 0, '>': 353, 'https://': 3, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '90', 'dir': 'jupyter-structured-clean-dedup', 'id': '563136', 'license': 'Apache-2.0', 'n_tokens_mistral': 52647, 'n_tokens_neox': 63161, 'n_words': 9722, 'path': 'notebooks/04_training_linear_models.ipynb', 'repo_name': 'Park-taenam/handson-ml2', 'size': '77254'}