ガウス・ニュートン法

ガウス・ニュートン法は、非線形最小二乗問題を解くために使用され、これは二乗関数値の合計を最小化することと同等です。これは、非線形関数の最小値を見つけるためのニュートン法の拡張です。二乗和は非負でなければならないため、このアルゴリズムは、ニュートン法を使用して和の成分のゼロを反復的に近似し、和を最小化するものと見なすことができます。この意味で、このアルゴリズムは、過剰決定方程式を解くための効果的な方法でもあります。計算が難しい場合がある2次導関数が必要ないという利点があります。^[1]

非線形最小二乗問題は、例えば、モデルが利用可能な観測値とよく一致するようにモデル内のパラメータが求められる 非線形回帰で発生します。

この方法は数学者のカール フリードリッヒ ガウスとアイザック ニュートンにちなんで名付けられ、ガウスの 1809 年の著作『Theoria motus corporum coelestium in Sectionibus conicis solem ambientum 』に初めて登場しました。^[2]

ガウス・ニュートン法は、変数の関数（残差と呼ばれることが多い）を反復的に求め、平方和を最小にする値を求める^[3]${\displaystyle m}$${\displaystyle {\textbf {r}}=(r_{1},\ldots,r_{m})}$${\displaystyle n}$${\displaystyle {\boldsymbol {\beta }}=(\beta _{1},\ldots \beta _{n}),}$${\displaystyle m\geq n,}$${\displaystyle \beta}$ $${\displaystyle S({\boldsymbol {\beta }})=\sum _{i=1}^{m}r_{i}({\boldsymbol {\beta }})^{2}.}$$

最小値の初期推定から始めて、この方法は反復によって進行する。 ${\displaystyle {\boldsymbol {\beta }}^{(0)}}$$${\displaystyle {\boldsymbol {\beta }}^{(s+1)}={\boldsymbol {\beta }}^{(s)}-\left(\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {J_{r}} \right)^{-1}\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right),}$$

ここで、rとβが列ベクトルである場合、ヤコビ行列の要素は $${\displaystyle \left(\mathbf {J_{r}} \right)_{ij}={\frac {\partial r_{i}\left({\boldsymbol {\beta }}^{(s)}\right)}{\partial \beta _{j}}},}$$

記号は行列の転置を表します。 ${\displaystyle ^{\operatorname {T} }}$

各反復において、次の 2 つの手順で前の式を並べ替えることで 更新を見つけることができます。${\displaystyle \Delta ={\boldsymbol {\beta }}^{(s+1)}-{\boldsymbol {\beta }}^{(s)}}$

• ${\displaystyle \Delta =-\left(\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {J_{r}} \right)^{-1}\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right)}$
• ${\displaystyle \mathbf {J_{r}} ^{\operatorname {T} }\mathbf {J_{r}} \Delta =-\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right)}$

、、を代入すると、これは という形式の従来の行列方程式になり、さまざまな方法で解くことができます (注を参照)。 ${\textstyle A=\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {J_{r}} }$${\displaystyle \mathbf {b} =-\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right)}$${\displaystyle \mathbf {x} =\Delta }$${\displaystyle A\mathbf {x} =\mathbf {b} }$

m = nの場合、反復は次のように単純化される。

$${\displaystyle {\boldsymbol {\beta }}^{(s+1)}={\boldsymbol {\beta }}^{(s)}-\left(\mathbf {J_{r}} \right)^{-1}\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right),}$$

これはニュートン法を 1 次元に直接一般化したものです。

データフィッティングでは、与えられたモデル関数がいくつかのデータポイントに最もよく適合するようなパラメータを見つけることが目的であり、関数は残差です。 ${\displaystyle {\boldsymbol {\beta }}}$${\displaystyle \mathbf {f} (\mathbf {x} ,{\boldsymbol {\beta }})}$${\displaystyle (x_{i},y_{i})}$${\displaystyle r_{i}}$$${\displaystyle r_{i}({\boldsymbol {\beta }})=y_{i}-f\left(x_{i},{\boldsymbol {\beta }}\right).}$$

そして、ガウス・ニュートン法は、関数のヤコビアンを使って次のよう に表すことができます。${\displaystyle \mathbf {J_{f}} =-\mathbf {J_{r}} }$${\displaystyle \mathbf {f} }$$${\displaystyle {\boldsymbol {\beta }}^{(s+1)}={\boldsymbol {\beta }}^{(s)}+\left(\mathbf {J_{f}} ^{\operatorname {T} }\mathbf {J_{f}} \right)^{-1}\mathbf {J_{f}} ^{\operatorname {T} }\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right).}$$

は の左擬似逆であることに注意してください。 ${\displaystyle \left(\mathbf {J_{f}} ^{\operatorname {T} }\mathbf {J_{f}} \right)^{-1}\mathbf {J_{f}} ^{\operatorname {T} }}$${\displaystyle \mathbf {J_{f}} }$

アルゴリズムステートメントでは、m ≥ nという仮定が必要です。そうでないと、行列は逆行列にならず、正規方程式を解くことができません (少なくとも一意には)。 ${\displaystyle \mathbf {J_{r}} ^{T}\mathbf {J_{r}} }$

ガウス・ニュートン法は、関数r _iのベクトルを線形近似することによって導出できます。テイラーの定理を使用すると、各反復で次のように記述できます。 $${\displaystyle \mathbf {r} ({\boldsymbol {\beta }})\approx \mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right)+\mathbf {J_{r}} \left({\boldsymbol {\beta }}^{(s)}\right)\Delta }$$

とします。右辺の二乗和を最小化する 問題、すなわち、${\displaystyle \Delta ={\boldsymbol {\beta }}-{\boldsymbol {\beta }}^{(s)}}$${\displaystyle \Delta }$$${\displaystyle \min \left\|\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right)+\mathbf {J_{r}} \left({\boldsymbol {\beta }}^{(s)}\right)\Delta \right\|_{2}^{2},}$$

線形最小二乗問題であり、明示的に解くことができ、アルゴリズム内の正規方程式を生成します。

正規方程式は、未知の増分 におけるn個の同時線形方程式です。これらは、コレスキー分解、またはより適切な のQR分解を使用して 1 ステップで解くことができます。大規模なシステムでは、共役勾配法などの反復法 の方が効率的な場合があります。J _rの列間に線形依存関係がある場合、 が特異になるため、反復は失敗します 。 ${\displaystyle \Delta }$${\displaystyle \mathbf {J_{r}} }$${\displaystyle \mathbf {J_{r}} ^{T}\mathbf {J_{r}} }$

が複素数の場合は共役形式を使用する必要があります。 ${\displaystyle \mathbf {r} }$${\displaystyle \mathbf {r} :\mathbb {C} ^{n}\to \mathbb {C} }$${\displaystyle \left({\overline {\mathbf {J_{r}} }}^{\operatorname {T} }\mathbf {J_{r}} \right)^{-1}{\overline {\mathbf {J_{r}} }}^{\operatorname {T} }}$

この例では、ガウス・ニュートンアルゴリズムを使用して、データとモデルの予測値間の誤差の二乗和を最小化することにより、モデルをデータに適合させます。

酵素媒介反応における基質濃度[ S ]と反応速度の関係を調べる生物学実験で、次の表のデータが得られました。

次のような曲線（モデル関数）を求めることが望ましい。 $${\displaystyle {\text{rate}}={\frac {V_{\text{max}}\cdot [S]}{K_{M}+[S]}}}$$

最小二乗法の意味でデータに最もよく適合するもので、パラメータとが決定されます。 ${\displaystyle V_{\text{max}}}$${\displaystyle K_{M}}$

とをそれぞれ[ S ]と率の値を、と表記する。ととする。残差の二乗和 が${\displaystyle x_{i}}$${\displaystyle y_{i}}$${\displaystyle i=1,\dots ,7}$${\displaystyle \beta _{1}=V_{\text{max}}}$${\displaystyle \beta _{2}=K_{M}}$${\displaystyle \beta _{1}}$${\displaystyle \beta _{2}}$$${\displaystyle r_{i}=y_{i}-{\frac {\beta _{1}x_{i}}{\beta _{2}+x_{i}}},\quad (i=1,\dots ,7)}$$

最小化されます。

未知数に関する残差ベクトルのヤコビアンは、次の要素を持つ行を持つ 行列である。${\displaystyle \mathbf {J_{r}} }$${\displaystyle r_{i}}$${\displaystyle \beta _{j}}$${\displaystyle 7\times 2}$${\displaystyle i}$$${\displaystyle {\frac {\partial r_{i}}{\partial \beta _{1}}}=-{\frac {x_{i}}{\beta _{2}+x_{i}}};\quad {\frac {\partial r_{i}}{\partial \beta _{2}}}={\frac {\beta _{1}\cdot x_{i}}{\left(\beta _{2}+x_{i}\right)^{2}}}.}$$

およびの初期推定値から始めて、ガウス・ニュートン法を 5 回繰り返した後、最適値およびが得られます。残差の二乗和は、5 回目の繰り返し後に初期値の 1.445 から 0.00784 に減少しました。右の図のプロットは、観測データを使用してモデルによって決定された最適パラメータの曲線を示しています。 ${\displaystyle \beta _{1}=0.9}$${\displaystyle \beta _{2}=0.2}$${\displaystyle {\hat {\beta }}_{1}=0.362}$${\displaystyle {\hat {\beta }}_{2}=0.556}$

ガウス・ニュートン反復法は、次の4つの条件下で局所的最小値に向かって収束することが保証されています。^[4]関数は開凸集合 で2回連続微分可能、ヤコビアンはフル列ランク、初期反復は 付近、局所的最小値は小さい。 の場合には収束は2次式です。 ${\displaystyle {\hat {\beta }}}$${\displaystyle r_{1},\ldots ,r_{m}}$${\displaystyle D\ni {\hat {\beta }}}$${\displaystyle \mathbf {J} _{\mathbf {r} }({\hat {\beta }})}$${\displaystyle \beta ^{(0)}}$${\displaystyle {\hat {\beta }}}$${\displaystyle |S({\hat {\beta }})|}$${\displaystyle |S({\hat {\beta }})|=0}$

^[5]によれば、増分 Δ はSの下降方向であり、アルゴリズムが収束する場合、極限はSの停留点となる。しかし、最小値が大きい場合、収束は保証されず、ニュートン法のような局所収束や、通常のウルフ条件下での収束さえも保証されない。^[6]${\displaystyle |S({\hat {\beta }})|}$

ガウス・ニュートン法の収束率は2乗に近づくことがある。^[7]初期推定値が最小値から遠い場合や行列の条件が悪ければ、収束が遅くなったり、まったく収束しないこともある。例えば、次式で表される方程式と変数を持つ 問題を考える。${\displaystyle \mathbf {J_{r}^{\operatorname {T} }J_{r}} }$${\displaystyle m=2}$${\displaystyle n=1}$$${\displaystyle {\begin{aligned}r_{1}(\beta )&=\beta +1,\\r_{2}(\beta )&=\lambda \beta ^{2}+\beta -1.\end{aligned}}}$$

最適値は です。(実際はの場合、ですが であるため、で最適値は です。) の場合、問題は実際に線形であり、この方法は 1 回の反復で最適値を見つけます。|λ| < 1 の場合、この方法は線形に収束し、反復ごとに誤差は係数 |λ| で漸近的に減少します。しかし、|λ| > 1 の場合、この方法は局所的にも収束しません。^[8]${\displaystyle \beta =0}$${\displaystyle \beta =-1}$${\displaystyle \lambda =2}$${\displaystyle S(0)=1^{2}+(-1)^{2}=2}$${\displaystyle S(-1)=0}$${\displaystyle \lambda =0}$

ガウス・ニュートン反復法は 、および の形の過剰決定方程式系を 解くための効果的な方法です 。ここで、は のヤコビ行列のムーア・ペンローズ逆行列（擬似逆行列とも呼ばれる）です。これはニュートン法の拡張とみなすことができ、孤立した正則解に向かって 同様に局所的な二次収束^{[4]を示します。}$${\displaystyle \mathbf {x} ^{(k+1)}=\mathbf {x} ^{(k)}-J(\mathbf {x} ^{(k)})^{\dagger }\mathbf {f} (\mathbf {x} ^{(k)})\,,\quad k=0,1,\ldots }$$${\displaystyle \mathbf {f} (\mathbf {x} )=\mathbf {0} }$$${\displaystyle \mathbf {f} (\mathbf {x} )={\begin{bmatrix}f_{1}(x_{1},\ldots ,x_{n})\\\vdots \\f_{m}(x_{1},\ldots ,x_{n})\end{bmatrix}}}$$${\displaystyle m>n}$${\displaystyle J(\mathbf {x} )^{\dagger }}$ ${\displaystyle J(\mathbf {x} )}$${\displaystyle \mathbf {f} (\mathbf {x} )}$

解が存在しないが、初期反復が平方和が小さな局所的最小値に達する点に近い場合、ガウス・ニュートン反復は に線形収束します。この点は、過剰決定システムの最小二乗解と呼ばれることがよくあります。 ${\displaystyle \mathbf {x} ^{(0)}}$${\displaystyle {\hat {\mathbf {x} }}=({\hat {x}}_{1},\ldots ,{\hat {x}}_{n})}$${\textstyle \sum _{i=1}^{m}|f_{i}(x_{1},\ldots ,x_{n})|^{2}\equiv \|\mathbf {f} (\mathbf {x} )\|_{2}^{2}}$${\displaystyle {\hat {\mathbf {x} }}}$${\displaystyle {\hat {\mathbf {x} }}}$

以下では、近似法によって関数最適化のためのニュートン法からガウス・ニュートン法を導出します。その結果、ガウス・ニュートン法の収束率は、一定の規則性条件下では2次式になることがあります。一般的には（より弱い条件下では）、収束率は線形です。^[9]

パラメータ関数Sを最小化するニュートン法の再帰関係は、 ${\displaystyle {\boldsymbol {\beta }}}$$${\displaystyle {\boldsymbol {\beta }}^{(s+1)}={\boldsymbol {\beta }}^{(s)}-\mathbf {H} ^{-1}\mathbf {g} ,}$$

ここでg はSの勾配ベクトルを表し、H はSのヘッセ行列を表します。

なので、勾配は次のように与えられる。 ${\textstyle S=\sum _{i=1}^{m}r_{i}^{2}}$$${\displaystyle g_{j}=2\sum _{i=1}^{m}r_{i}{\frac {\partial r_{i}}{\partial \beta _{j}}}.}$$

ヘッセ行列の要素は、勾配要素 をについて微分することによって計算されます。 ${\displaystyle g_{j}}$${\displaystyle \beta _{k}}$$${\displaystyle H_{jk}=2\sum _{i=1}^{m}\left({\frac {\partial r_{i}}{\partial \beta _{j}}}{\frac {\partial r_{i}}{\partial \beta _{k}}}+r_{i}{\frac {\partial ^{2}r_{i}}{\partial \beta _{j}\partial \beta _{k}}}\right).}$$

ガウス・ニュートン法は、2次微分項（この式の2番目の項）を無視することで得られる。つまり、ヘッセ行列は次のように近似される。 $${\displaystyle H_{jk}\approx 2\sum _{i=1}^{m}J_{ij}J_{ik},}$$

ここで、 はヤコビ行列J _rの要素です。正確なヘッセ行列が正確な近似値の近くで評価されると、 はほぼゼロになるため、2 番目の項もほぼゼロになり、近似が正当化されます。勾配と近似ヘッセ行列は、行列表記で次のように表すことができます。 ${\textstyle J_{ij}={\partial r_{i}}/{\partial \beta _{j}}}$${\displaystyle r_{i}}$$${\displaystyle \mathbf {g} =2{\mathbf {J} _{\mathbf {r} }}^{\operatorname {T} }\mathbf {r} ,\quad \mathbf {H} \approx 2{\mathbf {J} _{\mathbf {r} }}^{\operatorname {T} }\mathbf {J_{r}} .}$$

これらの式を上記の再帰関係に代入すると、演算方程式が得られる。 $${\displaystyle {\boldsymbol {\beta }}^{(s+1)}={\boldsymbol {\beta }}^{(s)}+\Delta ;\quad \Delta =-\left(\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {J_{r}} \right)^{-1}\mathbf {J_{r}} ^{\operatorname {T} }\mathbf {r} .}$$

ガウス・ニュートン法の収束は、すべての場合に保証されるわけではない。近似 $${\displaystyle \left|r_{i}{\frac {\partial ^{2}r_{i}}{\partial \beta _{j}\partial \beta _{k}}}\right|\ll \left|{\frac {\partial r_{i}}{\partial \beta _{j}}}{\frac {\partial r_{i}}{\partial \beta _{k}}}\right|}$$

2次微分項を無視するために必要な式は、収束が期待される2つのケースで有効である可能性がある：^[10]

1. 関数値は、少なくとも最小値付近では大きさが小さいです。${\displaystyle r_{i}}$
2. 関数は「軽度」非線形であるため、その大きさは比較的小さくなります。${\textstyle {\frac {\partial ^{2}r_{i}}{\partial \beta _{j}\partial \beta _{k}}}}$

ガウス・ニュートン法では、残差Sの二乗和は反復ごとに減少しない可能性があります。ただし、 Δ は下降方向であるため、 が定常点でない限り、十分に小さい に対してが成り立ちます。したがって、発散が発生した場合、1 つの解決策は、更新式で増分ベクトル Δ の一部を使用することです。 ${\displaystyle S\left({\boldsymbol {\beta }}^{s}\right)}$${\displaystyle S\left({\boldsymbol {\beta }}^{s}+\alpha \Delta \right)<S\left({\boldsymbol {\beta }}^{s}\right)}$${\displaystyle \alpha >0}$${\displaystyle \alpha }$$${\displaystyle {\boldsymbol {\beta }}^{s+1}={\boldsymbol {\beta }}^{s}+\alpha \Delta .}$$

言い換えれば、増分ベクトルは長すぎるが、それでも「下り坂」を指しているので、途中まで進むだけで目的関数S は減少する。 の最適値は直線探索アルゴリズムを使用して見つけることができる。つまり、 の大きさは、通常、区間内の直接探索法またはアルミヨ直線探索などのバックトラッキング直線探索を使用して、 S を最小化する値を見つけることで決定される。通常、 は、ウルフ条件またはゴールドスタイン条件を満たすように選択する必要がある。^[11]${\displaystyle \alpha }$${\displaystyle \alpha }$${\displaystyle 0<\alpha <1}$${\displaystyle \alpha }$

シフトベクトルの方向が最適分数αがゼロに近い場合、発散を処理するための代替方法は、信頼領域法であるレベンバーグ・マルカートアルゴリズムを使用することです。^{[3]正規方程式は、増分ベクトルが}最急降下の方向に回転するように変更されます。 $${\displaystyle \left(\mathbf {J^{\operatorname {T} }J+\lambda D} \right)\Delta =-\mathbf {J} ^{\operatorname {T} }\mathbf {r} ,}$$

ここで、D は正の対角行列です。D が単位行列 I で のとき、となり、 Δの方向は負の勾配 の方向に近づくことに注意してください。 ${\displaystyle \lambda \to +\infty }$${\displaystyle \lambda \Delta =\lambda \left(\mathbf {J^{\operatorname {T} }J} +\lambda \mathbf {I} \right)^{-1}\left(-\mathbf {J} ^{\operatorname {T} }\mathbf {r} \right)=\left(\mathbf {I} -\mathbf {J^{\operatorname {T} }J} /\lambda +\cdots \right)\left(-\mathbf {J} ^{\operatorname {T} }\mathbf {r} \right)\to -\mathbf {J} ^{\operatorname {T} }\mathbf {r} }$${\displaystyle -\mathbf {J} ^{\operatorname {T} }\mathbf {r} }$

いわゆる Marquardt パラメータも直線探索によって最適化できますが、が変更されるたびにシフト ベクトルを再計算する必要があるため、これは非効率的です。より効率的な戦略は次のとおりです。発散が発生した場合、Sが減少するまで Marquardt パラメータを増やします。次に、1 つの反復から次の反復まで値を保持しますが、可能であればカットオフ値に達するまで値を減らします。カットオフ値に達すると、Marquardt パラメータを 0 に設定できます。これにより、 Sの最小化は標準的なガウス - ニュートン最小化になります。 ${\displaystyle \lambda }$${\displaystyle \lambda }$

大規模な最適化では、ガウス・ニュートン法が特に重要です。なぜなら、行列は近似ヘッセ行列よりも疎であることが多いからです（ただし、常にそうであるとは限りません）。そのような場合、ステップ計算自体は通常、共役勾配法などの大規模で疎な問題に適した近似反復法を使用して実行する必要があります。 ${\displaystyle \mathbf {J} _{\mathbf {r} }}$${\displaystyle \mathbf {J} _{\mathbf {r} }^{\operatorname {T} }\mathbf {J_{r}} }$

このようなアプローチを機能させるには、少なくとも積を計算する効率的な方法が必要である。 $${\displaystyle {\mathbf {J} _{\mathbf {r} }}^{\operatorname {T} }\mathbf {J_{r}} \mathbf {p} }$$

あるベクトルpに対して、 が成り立ちます。疎行列のストレージでは、 の行を圧縮した形式 (ゼロ要素なしなど) で保存することが一般的に実用的ですが、転置のため上記の積を直接計算するのは困難です。ただし、c _{i を}行列 の行iとして定義すると、次の単純な関係が成り立ちます。 ${\displaystyle \mathbf {J} _{\mathbf {r} }}$${\displaystyle \mathbf {J} _{\mathbf {r} }}$$${\displaystyle {\mathbf {J} _{\mathbf {r} }}^{\operatorname {T} }\mathbf {J_{r}} \mathbf {p} =\sum _{i}\mathbf {c} _{i}\left(\mathbf {c} _{i}\cdot \mathbf {p} \right),}$$

すべての行が積に独立して加算的に寄与するようになります。実用的なスパース ストレージ構造を尊重することに加えて、この式は並列計算にも適しています。すべての行c _{i は}対応する残差r _iの勾配であることに注意してください。これを念頭に置くと、上記の式は残差が互いに独立して問題に寄与するという事実を強調しています。

Davidon、Fletcher、Powellによる準ニュートン法や Broyden–Fletcher–Goldfarb–Shanno ( BFGS 法)などの準ニュートン法では、完全なヘッセ行列の推定値は、 n 回の改良サイクルの後、パフォーマンスの点でニュートン法に近似するように、1 次導関数のみを使用して数値的に構築されます。準ニュートン法は一般的な実数値関数を最小化できますが、ガウス–ニュートン法、レーベンバーグ–マルカート法などは非線形最小二乗問題にのみ適合することに注意してください。 ${\textstyle {\frac {\partial ^{2}S}{\partial \beta _{j}\partial \beta _{k}}}}$${\textstyle {\frac {\partial r_{i}}{\partial \beta _{j}}}}$

1 次導関数のみを使用して最小化問題を解く別の方法は、勾配降下法です。ただし、この方法では 2 次導関数が近似的にも考慮されません。その結果、多くの関数、特にパラメータに強い相互作用がある場合には非常に非効率的です。

Julia での次の実装では、提供されたヤコビアンを使用する 1 つの方法と、自動微分を使用して計算する別の方法が提供されます。

"""
    ガウスニュートン(r,J,β₀,maxiter,tol)

ガウス・ニュートン最適化を実行して、ヤコビアン `J` を `β₀` から開始して残差関数 `r` を最小化します。ステップのノルムが `tol` 未満になるか、`maxiter` 回の反復後にアルゴリズムは終了します。
""" 
function gaussnewton ( r , J , β₀ , maxiter , tol ) β = copy ( β₀ ) for _ in 1 : maxiter Jβ = J ( β ); Δ = - ( Jβ '* Jβ ) \ ( Jβ '* r ( β )) β += Δ if sqrt ( sum ( abs2 , Δ )) < tol break end end return β end 
      
       
          
              
          
           
            
        
    
     


import AbstractDifferentiation as AD , Zygote backend = AD . ZygoteBackend () # 他のバックエンドも利用可能    
   

"""
    ガウスニュートン(r,β₀,maxiter,tol)

ガウス・ニュートン最適化を実行して、`β₀` から始めて残差関数 `r` を最小化します。関連するヤコビアンは自動微分を使用して計算されます。ステップのノルムが `tol` 未満になるか、`maxiter` 回の反復後にアルゴリズムは終了します。
""" 
function gaussnewton ( r , β₀ , maxiter , tol ) β = copy ( β₀ ) for _ in 1 : maxiter rβ , Jβ = AD . value_and_jacobian ( backend , r , β ) Δ = - ( Jβ [ 1 ] '* Jβ [ 1 ]) \ ( Jβ [ 1 ] '* rβ ) β += Δ if sqrt ( sum ( abs2 , Δ )) < tol break end end return β end 
      
       
           
              
          
           
            
        
    
     

• Björck, A. (1996).最小二乗問題の数値解析法. SIAM, フィラデルフィア. ISBN 0-89871-360-9。
• フレッチャー、ロジャー（1987年）。最適化の実践的方法（第2版）。ニューヨーク：ジョン・ワイリー・アンド・サンズ。ISBN 978-0-471-91547-8。。
• Nocedal, Jorge; Wright, Stephen (1999).数値最適化. ニューヨーク: Springer. ISBN 0-387-98793-2。

• 確率、統計、推定 アルゴリズムの詳細が説明され、この記事で例として説明されている生物学実験に適用されます (推定値の不確実性については 84 ページを参照)。

• Artelys Knitro は、ガウス・ニュートン法を実装した非線形ソルバーです。C で記述されており、C++/C#/Java/Python/MATLAB/R へのインターフェースを備えています。