Python SciPy : 多変数スカラー関数の制約付き局所的最適化

前回 に続いて SciPy の minimize を取りあげます。
今回は多変数スカラー関数の制約付き局所的最適化についてまとめます。

関数について

制約付き局所的最適化の method には表 1 のものがあります。

Table 1: 多変数スカラー関数の制約付き局所的最適化
関数 説明 勾配ベクトル、ヘッセ行列の要否
minimize(method="l-bfgs-b") 境界制約付き記憶制限 BFGS 法 勾配ベクトル
minimize(method="tnc") truncated Newton 法(Newton-CG 法) 勾配ベクトル、ヘッセ行列
minimize(method="cobyla") 線形近似法 不要
minimize(method="slsqp") 逐次最小 2 乗法 勾配ベクトル

L-BFGS-B は大規模問題において準 Newton 法を適用でるように計算容量を減らす工夫がされた方法です。
また、解の探索区間を指定できます。
偏導関数を与えなくても使用できますが、与えた方が高速です。

truncated Newton 法(切断 Newton 法? 短縮 Newton 法?)は直線探索の Newton-CG 法の別名です。
tnc は更に解の探索区間を指定できます。
特徴は newton-cg と同じですが、こちらは C 言語で実装されているため高速です。
L-BFGS-B より高速ですが、メモリの使用量が大きいため大規模な問題では使えないことがあります。

等式、不等式の制約条件がある問題の最適化には、COBYLA か SLSQP を使います。
COBYLA は偏導関数を使わないため遅いですが、比較的安定して解が得られます。
一方 SLSQP は偏導関数を与えられるため、収束が速いのが特徴です。

なお、多くの制約付き最適化問題は Lagrange の未定乗数法を用いて制約なしの最適化問題に置きかえれば、制約なし最適化手法で計算することができます。

SciPy の Tutorial の例を示します。

以下のように目的関数、その導関数、制約条件を定義します。

# coding: utf-8
import numpy as np
from scipy.optimize import minimize


def func(x, sign=1.0):
    return sign*(2*x[0]*x[1] + 2*x[0] - x[0]**2 - 2*x[1]**2)


def func_deriv(x, sign=1.0):
    dfdx0 = sign*(-2*x[0] + 2*x[1] + 2)
    dfdx1 = sign*(2*x[0] - 4*x[1])
    return np.array([dfdx0, dfdx1])

cons = ({'type': 'eq',
	 'fun' : lambda x: np.array([x[0]**3 - x[1]]),
	 'jac' : lambda x: np.array([3.0*(x[0]**2.0), -1.0])},
	{'type': 'ineq',
	 'fun' : lambda x: np.array([x[1] - 1]),
	 'jac' : lambda x: np.array([0.0, 1.0])})

制約条件を反映させない場合と

res1 = minimize(func, [-1.0, 1.0], args=(-1.0,), jac=func_deriv,
		  method='SLSQP')

res1.x
array([ 2.,  1.])

反映させた場合の結果です。

res2 = minimize(func, [-1.0, 1.0], args=(-1.0,), jac=func_deriv,
		  constraints=cons, method='SLSQP')

res2.x
array([ 1.00000009,  1.        ])

ついでに結果をプロットすると以下のようになります。

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

fig = plt.figure()
ax = fig.gca(projection='3d')

x = np.arange(-2, 5, 0.5)
y = np.arange(-2, 5, 0.5)
xx, yy = np.meshgrid(x, y, sparse=True)
z = 2*xx*yy + 2*xx - xx**2 - 2*yy**2

ax.plot_wireframe(xx, yy, z)
ax.scatter(res1.x[0], res1.x[1], func(res2.x), s=200, c="r")
ax.scatter(res2.x[0], res2.x[1], func(res2.x), s=200, c="y")
# plt.show()

15062301.png

コメント

Comments powered by Disqus
書籍更新情報
2017-01-16
Pythonによる科学技術計算 基礎編
1.3版への更新が可能になりました。
サポートページはこちら
電子書籍
Pythonによる科学技術計算 基礎編
Kindle ストアで販売中です
Pythonによる科学技術計算 基礎編
同人誌
技術書典(2016/6/25)
Emacs/org-modeのPDF作成術
電子版をBOOTHで販売中です
Emacs/org-modeのPDF作成術
Share