$$(\frac{\partial y}{\partial x_0}, \frac{\partial y}{\partial x_1})$$
勾配とは、すべての変数で偏微分したときの値をベクトルとしてまとめたものです。
$$(\frac{\partial y}{\partial x_0}, \dots, \frac{\partial y}{\partial x_n})$$
偏微分の実装を始めよう
ではもとの関数でおさらいです。
$$f(x_0, x_1)=x_0^2+x_1^2$$
def function_2(x):
return x[0]**2 + x[1]**2この実装ではxを配列としてfunction_2に渡していますね。
変数はx[0]とx[1]ではなくて、xとyの2変数でもいいんですが、配列にしておいた方が、あとが楽なんです。
配列の要素数が変数の数と一致するからです（3変数なら、x[0], x[1], x[2]の3要素ってことです）。
先に実装した数値微分の利用
まずは先に実装した数値微分のおさらいです。
def numerical_diff(function, x):
h = 1e-4
return (function(x + h) - function(x - h)) / (2*h)さて重要なのは、function(x + h)のところですよね。微小なhを足すってことです。これを\(f(x_0, x_1)=x_0^2+x_1^2\)でも同じことをすればいいんです。
\(f(x_0, x_1)=x_0^2+x_1^2\)に\(h\)を足す
さっきとの違いは、変数が複数あることです。
でも基本的には同じです。function(x + h)を求めればいいんです。
まずはx[0]について考えます。
def calc_partial(x):
# 微小変化量をdeltaとします。
delta = 1e-4
# もとのx[0]をoriginとします。
origin = x[0]
plus_d = origin + delta
minus_d = origin - delta
# 微小変化量を加えた配列を作ります。
# x[1]は固定なのでそのまま。
x_p = [plus_d, x[1]]
x_m = [minus_d, x[1]]
# 新しい配列に対する関数の値
f_p = function_2(x_p)
f_m = function_2(x_m)
# x_0についての偏微分 partial_x_0
partial_x_0 = (f_p - f_m) / (2*delta)
return partial_x_0
print(calc_partial([3.0, 4.0]))
# 6.00000000000378かなり丁寧に書いたのでくどいですが、理解はしやすいと思います。
x[1]についての偏微分を考えるなら、x_p = [plus_d, x[1]]のところをx_p = [x[0], plus_d]とすればいいですね（他にも変えなきゃいけないけど、x[0]の方を固定すればOKってことです）。
勾配を出力させる
ではこれを勾配が出力されるように書き換えましょう。
def calc_grad(x):
# 出力用の勾配（配列）を用意します。
grad = [0, 0] # 初期値は0にしました。
# ここは同じ
delta = 1e-4
for i in (0, 1):
# ここは同じ
origin = x[i]
plus_d = origin + delta
minus_d = origin - delta
# x[0]だけを変えたい
# もとの配列xの0番目にだけ代入
x[i] = plus_d
f_p = function_2(x)
# いまの配列xは、(plus_d, x[1])だよ
# 今度は配列xを(minus_d, x[1])にしたい
x[i] = minus_d
f_m = function_2(x)
# i番目の要素についての偏微分
partial_x_i = (f_p - f_m) / (2*delta)
# gradに代入
grad[i] = partial_x_i
# 次の周回のためにx[i]の値を戻す
x[i] = origin
return grad
print(calc_grad([3.0, 4.0]))
# [6.00000000000378, 7.999999999999119]どうですか?読めますか?
はっきり言って無駄が多いですが、論理の流れはなんとなくわかると思います。
先ほどとの違いは、for i in (0,1) で配列を回すところですね。
historoidゆっくりでいいから、理解できるまで読んでください!!
変数は2つだけとは限らない
まだまだ書き直しますよ!!
今度は変数が何個か分からない場合にも対応しましょう!! つまり先ほどのfor i in (0, 1)が使えないってことです。
# 3変数関数を用意
def func_3(x):
return x[0]**2 + x[1]**2 + x[2]**2
def calc_grad_mod(x):
# 変数の個数は分からない!!
# でも配列xの要素数と同じなのはわかっている。
grad = np.zeros_like(x) # xと同じ要素数でゼロ初期化
# ここは同じ
delta = 1e-4
# 変数と同じ回数だけ回す
for i in range(x.size): # これを使うためにはxはNumPy配列でないとだめ
# ここは同じ
origin = x[i]
plus_d = origin + delta
minus_d = origin - delta
# ここは同じ
x[i] = plus_d
f_p = func_3(x)
# ここは同じ
x[i] = minus_d
f_m = func_3(x)
# ここは同じ
partial_x_i = (f_p - f_m) / (2*delta)
# ここは同じ
grad[i] = partial_x_i
x[i] = origin
return grad
print(calc_grad_mod(np.array([3.0, 4.0, 5.0])))
# [ 6. 0. 20.]3変数でもちゃんと動きましたね。
くどい書き方をシンプルに
これで書き直しは最後にします。
関数自体も引数としてとれるように書き換えつつ、くどい部分を直していきましょう。
def numerical_grad(f, x):
grad = np.zeros_like(x)
delta = 1e-4
for i in range(x.size):
origin = x[i]
x[i] = origin + delta
f_p = f(x)
x[i] = origin - delta
f_m = f(x)
partial_x_i = (f_p - f_m) / (2*delta)
grad[i] = partial_x_i
x[i] = origin
return grad
# 4変数関数を用意
def n4(x):
return x[0]**3 + x[1]**3 + x[2]**3 + x[3]**3
numerical_grad(n4, np.array([1.0, 2.0, 3.0, 4.0]))
# array([ 3.00000001, 12.00000001, 27.00000001, 48.00000001])いやー。長かったですね。
これが教科書に書いてあるコードです。理解できたでしょうか。
これは教科書に任せて図示だけします。
以下が元の関数です。
この関数の勾配を可視化してみましょう。
可視化のコードはGitHub上で公開されています。
z軸方向（真上）から見ているわけですね。グラフの底へ落ちって行く感じが伝わって来ますね。
矢印の長さが勾配の大きさです。周辺の方が急勾配になっているのが、グラフと一致しています。
\(f(x, y) = \exp(-\frac{x^2+y^2}{2}) \frac{x^2+y^2}{2\pi}\)の勾配も可視化
せっかくなので他の関数でも試してみました。
$$f(x, y) = \exp(-\frac{x^2+y^2}{2}) \frac{x^2+y^2}{2\pi}$$
def _numerical_gradient_no_batch(f, x):
h = 1e-4
grad = np.zeros_like(x)
for idx in range(x.size):
tmp_val = x[idx]
x[idx] = float(tmp_val) + h
fxh1 = f(x)
x[idx] = tmp_val - h
fxh2 = f(x)
grad[idx] = (fxh1 - fxh2) / (2*h)
x[idx] = tmp_val
return grad
def numerical_gradient(f, X):
if X.ndim == 1:
return _numerical_gradient_no_batch(f, X)
else:
grad = np.zeros_like(X)
for idx, x in enumerate(X):
grad[idx] = _numerical_gradient_no_batch(f, x)
return grad
def make_plot(x):
a = x[0]**2 + x[1]**2
return np.exp(-a/2) * (a/(np.pi*2))
def tangent_line(f, x):
d = numerical_gradient(f, x)
print(d)
y = f(x) - d*x
return lambda t: d*t + y
if __name__ == '__main__':
x0 = np.arange(-4, 4, 0.3)
x1 = np.arange(-4, 4, 0.3)
X, Y = np.meshgrid(x0, x1)
X = X.flatten()
Y = Y.flatten()
grad = numerical_gradient(make_plot, np.array([X, Y]).T).T
plt.figure(figsize=(10, 10))
plt.quiver(X, Y, -grad[0], -grad[1], angles="xy",color="#4ECDC4")
plt.xlim([-4, 4])
plt.ylim([-4, 4])
plt.xlabel('x0')
plt.ylabel('x1')
plt.grid()
plt.draw()
plt.show()勾配が急激に変化している様子がわかりますね。
\(f(x, y) = \frac{\sin(\sqrt{x^2+y^2}}{\sqrt{x^2+y^2}}\)の勾配も可視化
メキシカンハットといわれている関数の勾配も見てみましょう。
繰り返し、勾配が0になるところが見られますね（多重円が見える）。
こういうのを見ると勾配消失問題がよく理解できますよね。
勾配法
勾配とはすなわち関数の値を最も減らす方向（変数ごとの偏微分値からなるベクトル）です。
したがって矢印の方向へ進むと関数の値は小さくなります。
ニューラルネットワークでは、損失関数の勾配によって重みを調整していきます（第5章から解説）。
しかし前述の例で示したとおり、勾配を下っていけば必ず最小値になるわけではありません。
機械学習では、最適化関数や損失関数の種類を変えることで、より適切な重みを探す方法が模索されていますが、「最小値になっているかはわからない」という認識を忘れないでください。
数式で見る勾配法
では実際に勾配降下法を数式で見てみましょう。
関数を\(f(x,y)\)としたとき、
$$\begin{align}
x_{n+1}=x_n-\eta\frac{\partial f}{\partial x_n}\\
y_{n+1}=y_n-\eta\frac{\partial f}{\partial y_n}
\end{align}$$
というように表すことができます。
\(\eta\)は、学習率と呼ばれるハイパーパラメータです。自分で設定することができます。
勾配の矢印の大きさを増やしている感じですね。
\(eta\)が大きいほど、極小値に向かいやすくなります。
機械学習のパラメータは、とくにニューラルネットワークでは「重み」と「バイアス」です。
パラメータは学習の結果獲得されるものであり、自分で設定することはできません。
これに対しハイパーパラメータは、学習モデル内部ではなく、学習モデルを構築する際に設定するパラメータです。
例えば、層の数、ニューロンの数、学習率、ドロップアウト率などさまざなな値があります。
とくに難しくないので、さっと実装してしまいましょう。
すでに偏微分の関数はできていますので、学習率をかける部分を書けば良いですね。
def gradient_descent(f, init_x, lr=0.01, step_num=100):
x = init_x
for i in range(step_num):
grad = numerical_gradient(f, x)
x = x - lr * grad
return x勾配降下法の可視化
以下の関数で勾配降下法を試してみます。
勾配はこんなふうに見えていましたね。
# ライブラリのインポート等は省略
init_x = np.array([-3.0, 4.0])
lr = 0.1
step_num = 20
x, x_history = gradient_descent(function_2, init_x, lr=lr, step_num=step_num)
plt.plot( [-5, 5], [0,0], '--b')
plt.plot( [0,0], [-5, 5], '--b')
plt.plot(x_history[:,0], x_history[:,1], 'o')
plt.xlim(-3.5, 3.5)
plt.ylim(-4.5, 4.5)
plt.xlabel("X0")
plt.ylabel("X1")
plt.show()
味気ないですけど、こんな風に点が移動していきます。これが勾配降下法です。
初期値をinit_x = np.array([-3.0, 4.0])と設定したので、そこからプロット始まっていますね。
step_numが20なので、点が20個打たれています。早期に原点\((0,0)\)に辿り着いているので20個も見えませんが、実装上は20個あります。
学習率を変えたときの挙動
ではここで学習率を\(1/10\)にしてみましょう。つまりlr=0.01です。
最小値まで辿り着いていませんね!!
でもステップ数を増やせば極小値まで辿り着けそうです。
では逆に学習率を増やして、\(1.025\)にしてみます。
わかりますか? 最小値に向かうどころか拡散しています（軸のスケールが変わっているのに注意）。
ちなみにlr=1.0にすると、原点対称を繰り返し、見かけ上2点しかないように見えます。
このように、学習率の設定次第では極小値を得ることはできないのです。
例えばこんなNNがあるとします。
それぞれの矢印に重みが設定されているので、矢印の数だけ重みがあります。
矢印の数は、ニューロンの組み合わせの数なので\(2\times3\)ですね。
そのうち「パラメータの数を計算して」という場面が来ると思いますが、このかけ算をすればOKです。
数式で書くと以下のようになります。
$$ W = \left( \begin{matrix}
w_{11} & w_{12} & w_{13} \\
w_{21} & w_{22} & w_{23}
\end{matrix} \right) $$
ここで損失関数を\(L(w)\)としたとき、勾配は以下のように示されます。
$$ \frac{\partial L}{\partial W} =
\left( \begin{matrix}
\frac{\partial L}{\partial w_{11}} & \frac{\partial L}{\partial w_{12}} & \frac{\partial L}{\partial w_{13}} \\
\frac{\partial L}{\partial w_{21}} & \frac{\partial L}{\partial w_{22}} & \frac{\partial L}{\partial w_{23}}
\end{matrix} \right) $$
では、NNで勾配降下法を実装しましょう。これまでやってきたことの総集編です。
def softmax(x):
x = x - np.max(x, axis=-1, keepdims=True)
return np.exp(x) / np.sum(np.exp(x), axis=-1, keepdims=True)# これから分類タスクを解くのでクロスエントロピーを使います。
def cross_entropy_error(y, t):
if y.ndim == 1:
t = t.reshape(1, t.size)
y = y.reshape(1, y.size)
# 教師データがone-hot-vectorの場合、正解ラベルのインデックスに変換
if t.size == y.size:
t = t.argmax(axis=1)
batch_size = y.shape[0]
return -np.sum(np.log(y[np.arange(batch_size), t] + 1e-7)) / batch_sizedef _numerical_gradient_1d(f, x):
h = 1e-4 # 0.0001
grad = np.zeros_like(x)
for idx in range(x.size):
tmp_val = x[idx]
x[idx] = float(tmp_val) + h
fxh1 = f(x) # f(x+h)
x[idx] = tmp_val - h
fxh2 = f(x) # f(x-h)
grad[idx] = (fxh1 - fxh2) / (2*h)
x[idx] = tmp_val # 値を元に戻す
return gradここまではこれまでにやった内容なので、忘れていたら教科書を戻ってください。
つぎにNNをクラスとして定義します。
class simpleNet:
def __init__(self):
# 重みの初期値をランダムに設定
self.W = np.random.randn(2,3)
def predict(self, x):
return np.dot(x, self.W)
def loss(self, x, t):
z = self.predict(x)
y = softmax(z)
loss = cross_entropy_error(y, t)
return lossこのクラスでは、\(2\times3\)のNNが定義されています。そして属性として重みを持っているわけですね。
そしてメソッドとして、predictとlossがあります。
ではまず、初期の重みを見てみましょう。
net = simpleNet()
net.W
# array([[ 0.17686805, 0.97862957, -1.31785384],
# [-1.00737797, -0.76290503, 0.20865116]])$$ \begin{align}
W &=
\left( \begin{matrix}
w_{11} & w_{12} & w_{13} \\
w_{21} & w_{22} & w_{23}
\end{matrix} \right) \\
&=
\left( \begin{matrix}
0.17686805 & 0.97862957 & -1.31785384 \\
-1.00737797 & -0.76290503 & 0.20865116
\end{matrix} \right)
\end{align}$$
これはランダムな初期値なので、とくに意味はありません。これが学習によって変化するのを確認しましょう。
# 入力層にある2つのニューロンにそれぞれ値を渡します。
x = np.array([0.6, 0.9])
# 初期値の重みで予測します。
net.predict(x)
# array([-0.80051934, -0.09943679, -0.60292626])出力層の3つのニューロンに入る値がこの3つです。
では推測してみましょう。
# 最大値のインデックスを表示させます。
np.argmax(net.predict(x))
# 1つまり、1番目のカテゴリであると判断されました。
これは初期値での予測なので、全く意味はないですよ!!
では損失値も計算しましょう。
# 正解ラベルを設定
t = np.array([0, 0, 1])
# 損失関数で評価
net.loss(x, t)
# 1.2456482976828025この\(1.2456482976828025\)は、惜しくもなんともない値です。「ぜんぜん違うじゃん」と言ってもいいでしょう。
では、これに対して勾配を求めてみましょう。
# numerical_gradientは、関数を引数にとるので関数をつくってあげます。
def f(W):
return net.loss(x, t)
# 偏微分します
dW = numerical_gradient(f, net.W)
# array([[ 0.14169652, 0.28565081, -0.42734733],
# [ 0.21254478, 0.42847622, -0.641021 ]])$$ \begin{align}
\frac{\partial L}{\partial W} &=
\left( \begin{matrix}
\frac{\partial L}{\partial w_{11}} & \frac{\partial L}{\partial w_{12}} & \frac{\partial L}{\partial w_{13}} \\
\frac{\partial L}{\partial w_{21}} & \frac{\partial L}{\partial w_{22}} & \frac{\partial L}{\partial w_{23}}
\end{matrix} \right) \\
&=
\left( \begin{matrix}
0.14169652 & 0.28565081 & -0.42734733 \\
0.21254478 &0.42847622 & -0.641021
\end{matrix} \right)
\end{align}$$
「これが何なの?」と思うかもしれません。
この行列の要素はそれぞれの重みの偏微分でした。\(\partial L / \partial w_{11} = 0.14169652\)は、「重み\(w_{11}\)に関していえば、\(h=0.001\)だけ増やすと、\(w_{11}\)はおよそ\(0.141h\)だけ増えます」ということです。
ということは勾配を下りたいのだから、\(w_{11}\)は小さくした方がいいということです。
これが勾配降下法の基本的な考え方です。
実際に重みを更新する箇所は、これ以降の部分になります。
損失を前のニューロンに伝える手法ですが、これはあと20ページ以上あとに始まります。
そして、それを実装するのはさらに30ページあとです。
まだまだ先は長いですが頑張っていきましょう。
historoidいやー。長かったし、大変だった。
これで私の担当箇所はおしまいです。めちゃくちゃ長かったですが、理解が深まるかと思います。
もしE資格受験に役立ったら教えて下さい。
