Pythonで基礎から機械学習 「勾配法」

はじめに

　この「Pythonで基礎から機械学習」シリーズの目的や、環境構築方法、シリーズの他の記事などは以下まとめページを最初にご覧下さい。

Pythonで基礎から機械学習まとめ

　今回は、Pythonで基礎から機械学習 「単回帰分析」、Pythonで基礎から機械学習 「重回帰分析」を読んだことが前提の内容となっております。

　また、本記事は、初学者が自分の勉強のために個人的なまとめを公開している記事になります。そのため、記事中に誤記・間違いがある可能性が大いにあります。あらかじめご了承下さい。

　より良いものにしていきたいので、もし間違いに気づいた方は、編集リクエストやコメントをいただけましたら幸いです。

勾配法とは

　前回の記事で取り上げたロジスティック回帰では、単回帰・重回帰のように最小二乗法をもとに行列計算で解析解を求めることはできず、勾配法等の探索アルゴリズムで、任意の値から逐次計算をしていくことで解空間を探索して近似解を求める必要があるという話をしました。今回は、その勾配法に関して紹介します。

　といっても、勾配法自体を説明したサイトは、数多くあるので、そちらに説明を譲ります。以下記事等参照下さい。

最適化のための勾配法

数学 - 勾配法について可視化して理解する。

　この記事では、今までそれぞれの解き方で計算する必要があった単回帰・重回帰の解が、勾配法という一つのアルゴリズムで解けることと、必ずしも正しい解を求められるわけではないことを確認したいと思います。

単回帰分析を勾配法で計算

　まずは、単回帰分析を勾配法で行います。

単回帰分析の解析解の確認

　単回帰分析の記事の復習です。

　使用するデータは、以下コマンドでダウンロードします。

!wget https://raw.githubusercontent.com/karaage0703/machine-learning-study/master/data/karaage_data.csv

　データを可視化した結果は以下です。

　最小二乗法を元に解いた、単回帰分析の結果は以下となります。

　可視化と計算の詳細は、単回帰分析の記事を参照下さい。

勾配法を用いて単回帰分析

　単回帰分析を勾配法を用いて行います。

　まずは、必要なライブラリをインポートします。

import numpy as np
import pandas as pd

　続いて、データを読み込んで、計算できる形に前処理しておきます（データのダウンロードに関しては、前の節の冒頭を参照下さい）。

df = pd.read_csv('karaage_data.csv')
df.head()

x = df[['x']].values
x = np.insert(x, 0, 1.0, axis=1)
y = df[['y']].values.T[0]

　勾配法の設定値（ハイパーパラメータと呼ばれます）を設定します。

# 繰り返し回数（10より大きい数を入れる）
EPOCHS = 50

# 学習率
LEARNING_RATE = 0.01

　さて、ここでハイパーパラメータは、どうやって決めるのでしょうか？実は基本的に人間が経験則で決めます。もちろん決め方の指針に関する研究や自動化に関する研究もあるのですが、問題設定によって大きく変わるので、とりあえず色々試して良かったものを選ぶというのが一般的です。

　そんなんでいいの？と思うかもしれませんが、実はディープラーニングでも状況は同じで、ハイパーパラメータは試行錯誤しながらいまだに人間が調整することが多かったりします。AutoML等、自動化も進んできましたが、膨大な計算資源が必要となります。

　続いて勾配降下法を実装します。

def pred(x, w):
    return x @ w

def gradient_method(x, y, epochs=EPOCHS):
    # 重みベクトル初期化
    w = np.ones(x.shape[1])
    # 結果履歴初期化
    history = np.zeros((0,2))

    # 逐次計算（学習）
    for i in range(epochs):
        yp = pred(x, w)
        yd = yp - y
        w = w - LEARNING_RATE * (x.T @ yd) / x.shape[0]
        if ( i % ((int)(EPOCHS / 10)) == 0):
            loss = np.mean(yd ** 2) / 2
            history = np.vstack((history, np.array([i, loss])))

    return w, history

　やっていることは、損失関数の空間で勾配を計算して、損失関数を最小化する方向に重みを更新していくだけですね。

　勾配法による計算を行います。

weights, history = gradient_method(x, y, epochs=EPOCHS)

　計算結果を表示します。

import matplotlib.pyplot as plt
import seaborn as sns

xall = x[:,1].ravel()
xl = np.array([[1, xall.min()],[1, xall.max()]])
yl = pred(xl, weights)

plt.figure(figsize=(6,6))
plt.plot(x[:,1], y, 'o')
plt.plot(xl[:,1], yl, linestyle="solid")
plt.show()

　むむ、なんか変ですね。

　以下コマンドで学習曲線を表示させます。学習曲線は、損失関数の値（loss）の試行回数ごとの変化具合で、lossが小さいほど正しい解が得られているという指標になります。

plt.plot(history[1:,0], history[1:,1])
plt.show()

　lossがめちゃめちゃ発散していますね。

　学習率が大きすぎたかもしれないので、以下のように学習率を変えてみます。

LEARNING_RATE = 0.001

　再度計算します。

weights, history = gradient_method(x, y, epochs=EPOCHS)

　結果を先ほどと同じ要領で表示します。

　今度は良い感じですね。値自体も、ほぼ解析解と一致しています。

　学習曲線を表示させます。

　良い感じにlossが下がっていますね。

重回帰分析を勾配法で計算

重回帰分析の解析解の確認

　重回帰分析の記事の復習です。利用するデータは以下コマンドでダウンロードします。

!wget http://pythondatascience.plavox.info/wp-content/uploads/2016/07/winequality-red.csv

　データを可視化します。

　重回帰分析の解析解を可視化すると以下のようになりました。

　詳細は、重回帰分析の記事を参照下さい。

勾配法を用いて重回帰分析

　単回帰と同様に計算していきます。

　以下でデータを読み込み前処理します。

import pandas as pd

df = pd.read_csv('winequality-red.csv', sep=';')

x = df[['density', 'volatile acidity']]
y = df[['alcohol']]

x = x.values
x = np.insert(x, 0, 1.0, axis=1)
y = y.values.T[0]

　ハイパーパラメータの設定です。

# 繰り返し回数（10より大きい数を入れる）
EPOCHS = 500

# 学習率
LEARNING_RATE = 0.01

　勾配法を計算します。ここで勾配法のアルゴリズム（関数）自体は、単回帰と全く同じものを使用します。

weights, history = gradient_method(x, y, epochs=EPOCHS)

　可視化します。

from mpl_toolkits.mplot3d import Axes3D  #3Dplot
import matplotlib.pyplot as plt
import seaborn as sns

fig=plt.figure()
ax=Axes3D(fig)

y = df[['alcohol']]
x1 = df[['density']]
x2 = df[['volatile acidity']]

ax.scatter3D(x1, x2, y)
ax.set_xlabel("x1")
ax.set_ylabel("x2")
ax.set_zlabel("y")

mesh_x1 = np.arange(x1.min()[0], x1.max()[0], (x1.max()[0]-x1.min()[0])/20)
mesh_x2 = np.arange(x2.min()[0], x2.max()[0], (x2.max()[0]-x2.min()[0])/20)
mesh_x1, mesh_x2 = np.meshgrid(mesh_x1, mesh_x2)
mesh_y = weights[1] * mesh_x1 + weights[2] * mesh_x2 + weights[0]
ax.plot_wireframe(mesh_x1, mesh_x2, mesh_y)
plt.show()

　3次元なので、分かりづらいですが、あんまりフィットしていない感じがしますね。

　学習曲線を可視化します。

plt.plot(history[:,0], history[:,1])
plt.show()

　良い感じにlossが下がっていますが、まだまだ下がっていきそうな雰囲気がありますね。

　繰り返し回数が少なかったのかもしれません。回数を増やしてみましょう。

EPOCHS = 5000

　再度計算します。

weights, history = gradient_method(x, y, epochs=EPOCHS)

　それっぽい感じになりました。

　学習曲線も良い感じにlossが下がっています。

　ただし、解析解とはかなり異なる値になっていることがわかります。変数が増えるほど、損失関数が最小に見える解（局所解）が増えるので、解析解が得られないことは増えていきます。

ロジスティック回帰分析を勾配法で計算

　ロジスティック回帰も、もちろん勾配法で同様に解けるはずなのですが、今回は省略します。もし実施したら追記します。

TensorFlowを使って勾配法を実装

　TensorFlowは、ディープラーニングのフレームワークとして有名ですが、その実は自動微分のライブラリなので、勾配法はお手の物です。この記事では省略しますが、興味ある方は以下サイトなど参照下さい。

まとめ

　勾配法に関して学習しました。単回帰をはじめとして、重回帰、ロジスティック回帰も同じ方法で近似解を計算できてしまいます。一方、ハイパーパラメータを適切に設定しないと、適切な解を得られない場合があることも分かりました。

　実は、ディープラーニングでもニューラルネットワークの膨大な変数の損失関数の解を探索する方法のベースとして勾配法が使われています。そして、ハイパーパラメータによって、適切な解（モデル）が得られないのは、ディープラーニングでも同じだったりします。

　大局的に見れば、結局パラメータが増えただけで、単回帰もディープラーニングもやっていることは一緒と言えないこともないですね。