UNISIA-SE Tech Blog

1. 前提条件

▼ 下記ページを理解していること。
[Python] [10] 損失関数と数値微分の実装サンプル

▼ Python3.6、NumPyがインストールされていること。
このページでは、venvの仮想環境(Python3.6)上にNumPyをインストールした環境で、Python対話モード(Pythonインタプリタ)にて実装サンプルを記載している。

※ Python対話モード、NumPyについては下記を参考。
Python対話モード：[Python] 対話モード (インタプリタ) の使用方法
NumPy：[Python] [NumPy] インストールとnumpy.ndarrayの使用方法

2. 偏微分のおさらい

3. 偏微分のPython実装サンプル

$ python
 >>> from mpl_toolkits.mplot3d import Axes3D
 >>> import matplotlib.pyplot as plt
 >>> import numpy as np
 >>>
 >>> def func_ex(x0, x1):
 ...     return x0**2 + x1**2
 ...
 >>> x = np.arange(-3.0, 3.0, 0.1)
 >>> x0 = np.arange(-3.0, 3.0, 0.1)
 >>> x1 = np.arange(-3.0, 3.0, 0.1)
 >>> X0, X1 = np.meshgrid(x0, x1)
 >>> Z = func_ex(X0, X1)
 >>> fig = plt.figure()
 >>> ax = Axes3D(fig)
 >>> ax.set_xlabel("x0")
 Text(0.5, 0, 'x0')
 >>> ax.set_ylabel("x1")
 Text(0.5, 0, 'x1')
 >>> ax.set_zlabel("f(x0, x1)")
 Text(0.5, 0, 'f(x0, x1)')
 >>> ax.set_title("f(x0, x1) = x0^2+x1^2 # arange:-3.0, 3.0, 0.1, label:f(x0, x1), x0, x1")
 Text(0.5, 0.92, 'f(x0, x1) = x0^2+x1^2 # arange:-3.0, 3.0, 0.1, label:f(x0, x1), x0, x1')
 >>> ax.plot_wireframe(X0, X1, Z)
 [<mpl_toolkits.mplot3d.art3d.Line3DCollection object at 0x7f6ae0cec6d8>]
 >>> plt.savefig('/var/www/vops/ops/macuos/static/macuos/img/b_id48_1.png')
 >>>

$ python
 >>> def num_dif(f, x):    # 下記 (※ num_dif) 参照
 ...     h = 1e-4
 ...     return (f(x+h) - f(x-h)) / (2 * h)
 ...
 >>> def func_partial_dif(x0):
 ...     return x0*x0 + 2.0**10.0
 ...
 >>> num_dif(func_partial_dif, 5.0)
 9.999999999976694
 >>>

$ python
 >>> def num_dif(f, x):
 ...     h = 1e-4
 ...     return (f(x+h) - f(x-h)) / (2 * h)
 ...
 >>> def func_partial_dif(x1):
 ...     return 2.0**5.0 + x1*x1
 ...
 >>> num_dif(func_partial_dif, 10.0)
 19.99999999995339
 >>>

4. 勾配のPython実装サンプル

$ python
 >>> import numpy as np
 >>>
 >>> def func_ex(x):    # このページで使用しているサンプル関数の定義
 ...     return x[0]**2 + x[1]**2
 ...
 >>> def num_gradient(f,x):    # 勾配関数
 ...     h = 1e-4
 ...     grad = np.zeros_like(x)    # xと同じ形状の配列で値がすべて 0
 ...
 ...     for idx in range(x.size):    # x の次元分ループする。 (下記例は、5.0, 10.0 の 2 周ループ)
 ...         idx_val = x[idx]
 ...         x[idx] = idx_val + h
 ...         fxh1 = f(x)    # f(x + h)の算出
 ...
 ...         x[idx] = idx_val - h
 ...         fxh2 = f(x)    # f(x - h)の算出
 ...
 ...         grad[idx] = (fxh1 - fxh2) / (2 * h)
 ...         x[idx] = idx_val    # 値をループ先頭の状態に戻す。
 ...     return grad
 ...
 >>>

$ python
 >>> num_gradient(func_ex, np.array([5.0, 10.0]))
 array([10., 20.])
 >>> num_gradient(func_ex, np.array([0.0, 10.0]))
 array([ 0., 20.])
 >>> num_gradient(func_ex, np.array([5.0, 0.0]))
 array([10.,  0.])
 >>> num_gradient(func_ex, np.array([2.0, -4.0]))
 array([4., -8.])
 >>> num_gradient(func_ex, np.array([-3.0, 4.0]))
 array([-6., 8.])
 >>>