読者です 読者をやめる 読者になる 読者になる

廿TT

譬如水怙牛過窓櫺 頭角四蹄都過了 因甚麼尾巴過不得

黄金分割法:一次元の最適化、黄金比による探索

R 最適化

一次元の最適化

この関数の最大値と、この関数を最大にするxがしりたいとする。

\displaystyle f(x)=\frac{\log(x + \log(x))}{\log(1+x)}

グラフはこう:


なので、最大値を与えるxは、区間(4,8)にいるだろう。

そこで、まず上記の関数を定義してやり、

f <- function(x) log(x + log(x))/log(1+x) 	# f を定義

optimizeを使う。
この関数の最大値を与えるxは、区間(4,8)にいるので、
lower = 4 , upper=8
最大値が知りたいので、
maximum = TRUE
にしてやる。
入力:

> optimize(f, lower = 4, upper = 8, maximum = TRUE)

出力:

$maximum
[1] 5.792299

$objective
[1] 1.055122

x=5.792299 で、最大値 f(x)=1.055122 が得られること

> f(5.792299)
[1] 1.055122

がわかった。超かんたん。

試しに点を打ってみる。

> points(5.792299 ,1.055122)


(あれ、そこ?…まあでもそんなもんか。とりあえずいいことにしよう)


optimize関数は、"golden section search(黄金分割法/黄金分割探索)"と"successive parabolic interpolation(断続放射状推定?)"というアルゴリズムの合わせ技でこれをやっているらしい。

Statistical Computing with R (Chapman & Hall/CRC The R Series)

Statistical Computing with R (Chapman & Hall/CRC The R Series)

↑ここまでがこの本の 11.5 章 Example11.11までに相当する。

黄金比による探索

↓ここからは

統計的データ解析のための数値計算法入門 (統計ライブラリー)

統計的データ解析のための数値計算法入門 (統計ライブラリー)

この本の4.1.1の手順4.1.1あたりに対応する。

successive parabolic interpolationはよく分からなかったので、golden section searchについて書くことにする。

黄金比による直接探索法の手順は以下の通り:

これは、さっきグラフでみたみたいな、ある区間[a,b]で一山形の関数を考えてるときに使う。f(x)をそのような関数とする。φは、
\phi = \frac {\sqrt{5}-1}{2}(ニアリーイコール0.618)。

(1)
 c = a + (b-a)\phi ^2, d=a+(b-a)\phi
とし、f(a) , f(b) , f(c) ,f(d) を計算する。
(2) f(c) ≧ f(d) ならば(3-1)へ、f(c) < f(d) ならば(3-2)へ。
(3-1) b <- d , d <- c として、新しい
 c = a + (b-a)\phi ^2
を計算して、関数値f(c)をもとめ、(2)へ。(aはそのまま)
(3-2) a <- c , c <- d として、新しい
 d = a + (b-a)\phi
を計算して、関数値f(d)をもとめ、(2)へ。(bはそのまま)

区間幅|b-a|はこの手順の一回の反復ごとにφ=0.618倍、m回反復したら元の幅の
 (0.6180)^m
となるので、必要な精度が得られる回数繰り返せばいい。

ふーん。ここで黄金比φに出会うとはおもわなかった。なんでφなの?←それを見ていこう。
それに、この手順がどういう意図でこれをやってるのか、これだけだとよくわからない。←そこも見てやろう。

上記の手順の考え方はこんなふうだ:
区間[a,b]のなかに a < c < d < b となる2点 c , d をとって、f(a) , f(b) , f(c) , f(d) を計算する。
(f は区間[a,b]で一山形の関数)

f(c) ≧ f(d) なら最大値は [d,b] にはない。[a,d]にある。
f(c) < f(d) なら最大値は [a,c] にはない。[c,b]にある。
このことは図を書いてながめてみるとなっとくがいく。

なので、

f(c) ≧ f(d) ならば、新[a,b] <- 旧[a,d]
f(c) < f(d) ならば、新[a,b] <- 旧[c,b]

として同じことを繰り返せば、区間の幅がだんだん小さくなっていって、最大値を与える点が分かってくる。


しかし
区間幅を更新する際,それまでに計算した関数値が次の探索でも使えるのが望ましい」
うん、たしかに。そこで、
 \phi < 1 , \phi ^2 = 1- \phi
となるような数、φを使う。

今、(簡単のため)[a,b]=[0,1]として、
 d = \phi , (\phi<1)
 c = \phi ^2,(\phi ^2 = 1 - \phi)
とする。
こうして、一回更新して[a,b]=[0,φ]になったとする。このとき
新dを、
 d = c (=\phi ^2)
として、
新c を
 c = \phi ^3
とすると、
 (a,b,c,d )=( 0, \phi ^2 ,\phi , 1)
だったのが
 (a,b,c,d )=( 0, \phi ^3 ,\phi ^2  ,\phi)
となった。
点の配置は元のやつをφ倍したものになっていて、比が保たれていることに注目。
(これは、更新後[a,b]=[1-φ,1]になった場合でも同じこと。)


でも、φってどんな数なの?
はい。ではこの式
 \phi ^2 = 1- \phi
をφについてとく。
 \phi ^2 = 1 -\phi \\ \phi ^2 + \phi -1 = 0 \\ \phi = \frac {-1 \pm \sqrt{1^2 -4 \cdot 1 \cdot (-1)}}{2 \cdot 1}
 \phi = \frac{\pm \sqrt{5} -1}{2}
だけど、区間[0,1]にあるように、
 \phi = \frac{\sqrt{5} -1}{2}
をとる。1:φ が黄金比(golden ratio)と呼ばれるやつである。


以上が、「黄金比による直接探索法の手順」の考え方である。おしまい。


でも、なんで区間[a,b]に2点とるの?一点じゃだめなの?
やってみよう。下の図を参照。

どちらも f(a) < f(c) かつ f(c) > f(b) だけど、最大値を与える点は、左は[c,b]にあり、右は[a,b]にある。
一点cしかとらないと最大値を与える点が、[a,c]にあるか、[c,b]にあるか、判断できない。最低でも、2点とらないといけないみたいだ。


今回はこのくらいにしとこう。
感謝しながら、参考になりそうなページを貼っておく。
ありがとうございます。
Rによる最適化、パラメータ推定入門 - Seeking for my unique color.

追記

一山形の関数書くのにベジエ曲線はあんま適当じゃなかったなー。まあいいや。

関連エントリ

abrahamcow.hatenablog.com