SciPyによる最適化(scipy.optimize)③｜SciPy入門 #7

f:id:lib-arts:20190225182115g:plain

SciPyについて色々と話題になり面白そうだったので公式チュートリアルを元にまとめています。

SciPy Tutorial — SciPy v1.2.1 Reference Guide
#5ではscipy.optimizeから制約条件のない際の最適化、#6では制約条件がある場合の最適化や最小二乗法などに関して取り扱いました。

#7では1変数の関数の最適化や求根などについて取り扱っていければと思います。

Optimization (scipy.optimize) — SciPy v1.2.1 Reference Guide
以下目次になります。

1. Univariate function minimizers (minimize_scalar)
2. Root finding
3. まとめ

1. Univariate function minimizers (minimize_scalar)
まずは1変数関数(Univariate function)の取り扱いについてまとめていきます。スカラー値のような1変数の入力においてはより速く動く最適化アルゴリズムが開発されており、それはminimize_scalar関数を用いることで実行することができます。

・制約なしの最適化(Unconstrained minimization)
brentとgoldenの二つの手法があるが、goldenはほぼアカデミックでしか使われないので、brentが中心に解説されています。使用サンプルとして下記のコードを実行してみてください。

from scipy.optimize import minimize_scalar
f = lambda x: (x + 1)**2
res = minimize_scalar(f, method='brent')
print(res.x)

上記を実行することで、 $f(x)=(x+1)^2$ の最小値問題を解くことができます。答えとしては、 $f'(x)=2x=0$ の答えとなる $-1$ が出力されます。

・制約ありの最適化(Bounded minimization)
制約ありの最適化では先ほどmethod='brent'と引数を与えていたのをmethod='bounded'に変更し、boundsで定義域を与えることで実行することができます。

from scipy.optimize import minimize_scalar
f = lambda x: (x + 1)**2
res = minimize_scalar(f, bounds=(4, 7), method='bounded')
print(res.x)

上記を実行することで、 $f(x)=(x+1)^2$ において制約条件として $4 \leq x \leq 7$ を設定した際の最適化問題を解くことができます。答えは今まで通りres.xで得ることができ、4に近い値となります。(近似的に解いているため、4の周辺の実数となる可能性もあります。)

2. Root finding
次に求根(方程式の解を求める)問題についてです。まずは簡単な関数での実装例から見ていきます。

・方程式が一つの場合

import numpy as np
from scipy.optimize import root
def func(x):
return (x + 1)**2 - 4
sol = root(func, 0.3)
print(sol.x)

上記はoptimize.rootメソッドを用いて $(x+1)^2-4 = 0$ の解を求めています。手元での実行結果は $1$ になりました。ここで注目したいのがrootの引数として、funcと0.3を与えていることです。0.3はドキュメントではInitial guessとされており、初期値と解釈して良さそうです。

for init_v in np.arange(-2,0,0.1):
sol = root(func2, init_v)
print(str(init_v)+":"+str(sol.x))

上記を実行すると、初期値に近い解に収束していることがわかります。 $(x+1)^2-4 = 0$ の解は-3と1なので、-1を境にどちらに収束するかが決まります。

f:id:lib-arts:20190305132324p:plain
コードの実行結果は上記のようになり、仮説通りの結果が得られています。

・連立方程式
 チュートリアルを参考に実装例として下記をまとめました。

def func2(x):
f = [x[0] **2 - x[1] +1,
2*x[0] - x[1] + 1]
df = np.array([[2*x[0], -1],
[2, -1]])
return f, df
sol = root(func2, [3, 1], jac=True, method='lm')
print(sol.x)

これは $x^2 - y = -1$ と $2x - y = -1$ の連立方程式を解いており、解としてx[0]=2、x[1]=5を得ています。また、こちらも初期値の検証を下記のコードで行ってみました。

for init_x in np.arange(0,2,0.5):
for init_y in np.arange(1,5,1):
sol = root(func2, [init_x, init_y], jac=True, method='lm')
print("("+str(init_x)+","+str(init_y)+") :"+str(sol.x))

実際の解はx[0]が0と2なので、1を境に結果が違うようになるのではと推測できます。

f:id:lib-arts:20190305134914p:plain
実行結果は上記のようになり、概ね仮説通りの結果が得られているということがわかります。