一次元の解探索、f(x)=0を解く、uniroot()、ブレント法（Brent's method）

（前回）
いま読んでいるMaria L. Rizzo "Statistical Computing with R"では，
「ブレント法（Brent's method）は，root bracketingと二分法（bisection），逆二次補間の組み合わせである」
と書いてあるんですが，ウィキペディアでは
「ブレント法は、二分法、割線法、逆2次補間を組み合わせた、複雑ではあるが広く用いられている、数値解析における求根アルゴリズムの1つである」
となっている：http://ja.wikipedia.org/wiki/%E3%83%96%E3%83%AC%E3%83%B3%E3%83%88%E6%B3%95

逆二次補間と虐殺器官って似てるなあ

虐殺器官 (ハヤカワ文庫JA)

割線法（secant mathod；わりせんほう、と読むらしい）とroot bracketingというのが同じものなんでしょうか

• 割線法について解説している：http://next1.cc.it-hiroshima.ac.jp/MULTIMEDIA/numeanal1/node8.html
• root bracketingについて解説している：http://www.scribd.com/doc/27726572/8/Root-bracketing-methods

root bracketing method
f(x)は区間[a,b]で連続f(a)f(b)<0と仮定する．中間値の定理により，f(p)=0となるような(a,b)の要素，ある数pの存在を保証される．

root bracketing methodは入れ子になっている区間[a(n),b(n)]の数列から構成される．区間[a(n),b(n)]はそれぞれf(x)=0の解を含む．それぞれのステップで，[a(n),b(n)]の要素p(n)の計算と進行は以下の通り：

もしf(p(n))=0ならば，反復をストップしてa(n+1)=a(n)，b(n+1)=p(n)
さもなくば，
もしf(a(n))f(p(n))<0ならば，a(n+1)=a(n)，b(n+1)=p(n)．
さもなくば，a(n+1)=p(n)，b(n+1)=b(n)．

http://www.scribd.com/doc/27726572/8/Root-bracketing-methods より

どうもちがうような気がする
ちがうよなあ
コード読むか…きついなあ…

さて、Rでは、uniroot関数がブレント法を用いているらしい
前回と同じ方程式

a=1/2，n=20をunirootを使って解いてみる

a <- 0.5
n <- 20
out <- uniroot(function(y) {
	a^2 + y^2 + 2*a*y/(n-1)-(n-2)},  #関数を入れて
	lower=0,upper=n*5)  #[0,5n]の範囲で解く

結果は、

$root
[1] 4.186870

$f.root
[1] 0.0002381408

$iter
[1] 14

$estim.prec
[1] 6.103516e-05

解は 4.186870
また区間をひっくり返して[-5n,0]の範囲で解くと、

uniroot(function(y) {a^2 + y^2 + 2*a*y/(n-1)-(n-2)},
			interval = c(-n*5, 0))$root
[1] -4.239501

解は -4.239501

正確な解は，y=4.186841， -4.239473だった
前回（二分法）と比べてみる
前回：真値-近似値 = 4.186841-4.186845= -0.000004
今回：真値-近似値 = 4.186841-4.186870= -0.000029
あれ遠くなった…

また、

$f.root
[1] 0.0002381408

は
0.0002381408 > 0.0001220703=.Machine$double.eps^0.25
なんだけど，いいのか？これ？

でも

$iter
[1] 14

ということは繰り返し回数が14回なんだよね？
前回の21回より早い
ここは勝ってる

…なんか全体的に釈然としないな

Rには多項式の解を求めるpolyroot関数というのもあって，複素数の範囲での解も求められるらしいのでチェックしておきたい：
http://cse.naro.affrc.go.jp/takezawa/r-tips/r/38.html

↑しかしこちらではuniroot()はニュートン法を使っているって書いてるんだよね…また新説が…
ところで、このR-tipsは素晴らしいページだ
いつも感謝してます

以上がMaria L. Rizzo "Statistical Computing with R"のExample11.7の内容なんだけど…
また戻って読み直さないと駄目だろうな

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