多変数関数の微分

前回1変数関数の微分をやった.

今回は多変数関数の微分をやってみた. 多変数関数を微分して得られる導関数は勾配(gradient)と呼ばれる.

Source

# -*- coding:utf-8 -*-
"""
多変数関数の微分を計算機で行う.
今回は円を扱う
"""

import numpy as np
import matplotlib.pyplot as plt

# n変数関数fを微分して導関数を返す.
def gradient(f, n, h=1./10**2):
    I = np.identity(n)
    def partial(x, i):
        return (f(x+I[i]*h) - f(x)) / h
    def grad(x):
        return np.array([partial(x, i) for i in np.arange(n)], dtype=float)
    return grad

if __name__ == '__main__':
    # target function
    c = np.array([2, 1])
    R = 2
    f = lambda x: (x-c).dot(x-c) - R**2

    # for plot data
    d = 2
    n = 100
    x = np.linspace(-5, 5, n)
    X, Y = np.meshgrid(x, x)
    PX = np.array([np.array([xx,yy]) for x, y in zip(X, Y) for xx, yy in zip(x, y)]) # like expand.grid on R
    Z = np.array([f(px) for px in PX]).reshape((n, n))

    fx = gradient(f, d, h=1./10**2) # fの導関数
    p = np.array([-4., 1.]) 
    gp = fx(p) # pにおける勾配
    plt.contour(X, Y, Z, levels=[0])
    plt.plot(p[0], p[1], "ro")
    plt.plot(gp[0], gp[1], "bo")
    # pにおける勾配gpを直線で表す.直線の方程式t*gp + (1-t)*p
    t = np.linspace(-1.2, 1.1, n)
    pl = np.array([tt*gp + (1-tt)*p for tt in t])
    plt.plot(pl[:,0], pl[:,1], "b-")

    plt.show()

Result

点pにおける勾配を求めて,その勾配方向を直線で描いてみた. 赤い点が点pで,青い点がpの勾配,青い直線が勾配方向を示している. 勾配法を使えばこの直線上を下って行くことになるのかな? wpid-multvalscalc.png

コメントを残す

以下に詳細を記入するか、アイコンをクリックしてログインしてください。

WordPress.com ロゴ

WordPress.com アカウントを使ってコメントしています。 ログアウト / 変更 )

Twitter 画像

Twitter アカウントを使ってコメントしています。 ログアウト / 変更 )

Facebook の写真

Facebook アカウントを使ってコメントしています。 ログアウト / 変更 )

Google+ フォト

Google+ アカウントを使ってコメントしています。 ログアウト / 変更 )

%s と連携中