Gauss–Newtons metod

Givna ${\displaystyle m}$ funktioner ${\displaystyle {\textbf {r}}=(r_{1},\ldots ,r_{m})}$ (ofta kallade rester) av ${\displaystyle n}$ variabler ${\displaystyle {\boldsymbol {\beta }}=(\beta _{1},\ldots \beta _{n}),}$ med ${\displaystyle m\geq n,}$ ^[a] hittar Gauss-Newton-algoritmen iterativt värdet av variablerna som minimerar kvadratsumman ^[3]

${\displaystyle S({\boldsymbol {\beta }})=\sum _{i=1}^{m}r_{i}({\boldsymbol {\beta }})^{2}.}$

Man börjar med en första gissning ${\displaystyle {\boldsymbol {\beta }}^{(0)}}$ och fortsätter iterativt

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

där elementen i jakobianen är

${\displaystyle \left(\mathbf {J_{r}} \right)_{ij}={\frac {\partial r_{i}\left({\boldsymbol {\beta }}^{(s)}\right)}{\partial \beta _{j}}},}$

r och β är kolumnvektorer och symbolen ${\displaystyle ^{\mathsf {T}}}$ betecknar matristransponering.

Beräkningar

Vid varje iteration, kan uppdateringen ${\displaystyle \Delta ={\boldsymbol {\beta }}^{(s+1)}-{\boldsymbol {\beta }}^{(s)}}$ hittas genom att ordna om föregående ekvation i följande två steg:

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

${\textstyle \mathbf {J_{r}} ^{\mathsf {T}}\mathbf {J_{r}} \Delta =-\mathbf {J_{r}} ^{\mathsf {T}}\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right)}$

Med beteckningarna ${\textstyle A=\mathbf {J_{r}} ^{\mathsf {T}}\mathbf {J_{r}} }$ , ${\displaystyle \mathbf {b} =-\mathbf {J_{r}} ^{\mathsf {T}}\mathbf {r} \left({\boldsymbol {\beta }}^{(s)}\right)}$ , och ${\displaystyle {\mathbf {x}}=\Delta }$ , förvandlas detta till den vanliga matrisekvationen ${\displaystyle A{\mathbf {x}}={\mathbf {b}}}$ , som sedan kan lösas på en mängd olika metoder (se anmärkningar ).

När ${\displaystyle \mathbf {r} }$ är komplex ${\displaystyle \mathbf {r} :\mathbb {C} ^{n}\rightarrow \mathbb {C} }$ den konjugerade formen ska användas: ${\displaystyle \left({\overline {\mathbf {J_{r}} }}^{\mathsf {T}}\mathbf {J_{r}} \right)^{-1}{\overline {\mathbf {J_{r}} }}^{\mathsf {T}}}$ .Om m = n, kan iterationen förenklas till

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

vilket är en direkt generalisering av Newtons metod i en dimension.

Normalekvationerna är n samtidiga linjära ekvationer i okända steg ${\displaystyle \Delta }$ . De kan lösas i ett steg, med hjälp av Choleskyuppdelning, eller, bättre, QR-faktorisering av ${\displaystyle \mathbf {J_{r}} }$ .^[4] För stora system kan en iterativ metod, såsom konjugatgradientmetoden, vara mer effektiv. Om det finns ett linjärt beroende mellan kolumner i J _r kommer iterationerna att misslyckas, då ${\displaystyle \mathbf {J_{r}} ^{T}\mathbf {J_{r}} }$ blir singular.

Beräkningar för dataanpassning

Inom dataanpassning, där målet är att hitta parametrarna ${\displaystyle {\boldsymbol {\beta }}}$ så att en given modell fungerar ${\displaystyle \mathbf {f} (\mathbf {x} ,{\boldsymbol {\beta }})}$ passar bäst på vissa datapunkter ${\displaystyle (x_{i},y_{i})}$ , är funktionerna ${\displaystyle r_{i}}$ är residualerna :

${\displaystyle r_{i}({\boldsymbol {\beta }})=y_{i}-f\left(x_{i},{\boldsymbol {\beta }}\right).}$

Sedan kan Gauss-Newton-metoden uttryckas i termer av jakobianen ${\displaystyle \mathbf {J_{f}} }$ av funktionen ${\displaystyle \mathbf {f} }$ som

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

Observera att ${\displaystyle \left(\mathbf {J_{f}} ^{\mathsf {T}}\mathbf {J_{f}} \right)^{-1}\mathbf {J_{f}} ^{\mathsf {T}}}$ är den vänstra pseudoinversen av ${\displaystyle \mathbf {J_{f}} }$ .

I det här exemplet kommer Gauss-Newton-metoden att användas för att anpassa en modell till vissa data genom att minimera summan av kvadrater av fel mellan data och modellens förutsägelser.

I ett biologiskt experiment som studerade sambandet mellan substratkoncentration [S] och reaktionshastighet V i en enzymmedierad reaktion, erhölls data i följande tabell.

Det är önskvärt att hitta en kurva (modellfunktion) av formen

${\displaystyle V={\frac {V_{\text{max}}\cdot [S]}{K_{M}+[S]}}}$

som bäst passar data i minsta kvadrat-mening. Då bestäms parametrarna ${\displaystyle V_{\text{max}}}$ och ${\displaystyle K_{M}}$ .

Beteckna med ${\displaystyle x_{i}}$ och ${\displaystyle y_{i}}$ värdena för [S] (koncentration) och V (hastighet) för ${\displaystyle i=1,\dots ,7}$ . Låt ${\displaystyle \beta _{1}=V_{\text{max}}}$ och ${\displaystyle \beta _{2}=K_{M}}$ och hitta ${\displaystyle \beta _{1}}$ och ${\displaystyle \beta _{2}}$ så att summan av kvadraterna av residualerna

${\displaystyle r_{i}=y_{i}-{\frac {\beta _{1}x_{i}}{\beta _{2}+x_{i}}},\quad (i=1,\dots ,7)}$

minimeras.

Jakobianen ${\displaystyle \mathbf {J_{r}} }$ av vektorn av residualerna ${\displaystyle r_{i}}$ med hänsyn till de okända ${\displaystyle \beta _{j}}$ är en ${\displaystyle 7\times 2}$ -matrismed där den ${\displaystyle i}$ :te raden har elementen

${\displaystyle {\frac {\partial r_{i}}{\partial \beta _{1}}}=-{\frac {x_{i}}{\beta _{2}+x_{i}}};{\frac {\partial r_{i}}{\partial \beta _{2}}}={\frac {\beta _{1}\cdot x_{i}}{\left(\beta _{2}+x_{i}\right)^{2}}}.}$

Man börjar med de första uppskattningarna ${\displaystyle \beta _{1}=0,9}$ och ${\displaystyle \beta _{2}=0,2}$ och efter fem iterationer av Gauss-Newton-metoden erhålls de optimala värdena ${\displaystyle {\hat {\beta }}_{1}=0,362}$ och ${\displaystyle {\hat {\beta }}_{2}=0,556}$ erhålls. Summan av kvadraterna på residualerna minskade från initialvärdet 1,445 till 0,00784 efter den femte iterationen. Figuren till höger visar kurvan som bestäms av modellen för de optimala parametrarna med de observerade data.

I det följande kommer Gauss–Newton-metoden att härledas från Newtons metod för funktionsoptimering via en approximation. Som en konsekvens kan konvergenshastigheten för Gauss-Newton-metoden vara kvadratisk under vissa regularitetsförhållanden. I allmänhet (under svagare förhållanden) är konvergenshastigheten linjär.^[5]

Iterationsekvationen för Newtons metod för att minimera en funktion S av parametrarna ${\displaystyle {\boldsymbol {\beta }}}$ är

${\displaystyle {\boldsymbol {\beta }}^{(s+1)}={\boldsymbol {\beta }}^{(s)}-\mathbf {H} ^{-1}\mathbf {g} ,}$

där g betecknar gradientvektorn för S och H betecknar den hessianen för S .

Eftersom ${\displaystyle S=\sum _{i=1}^{m}r_{i}^{2}}$ , ges gradienten av

${\displaystyle g_{j}=2\sum _{i=1}^{m}r_{i}{\frac {\partial r_{i}}{\partial \beta _{j}}}.}$

Hessianens element beräknas genom att derivera gradientelementen, ${\displaystyle g_{j}}$ , med avseende på ${\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).}$

Gauss-Newton-metoden erhålls genom att försumma andra ordningens derivator (den andra termen i summanderna). Det vill säga, hessianen approximeras av

${\displaystyle H_{jk}\approx 2\sum _{i=1}^{m}J_{ij}J_{ik},}$

där ${\displaystyle J_{ij}={\frac {\partial r_{i}}{\partial \beta _{j}}}}$ är element i jakobianen J _r. Gradienten och den ungefärliga hessianen kan skrivas i matrisnotation som

${\displaystyle \mathbf {g} =2{\mathbf {J} _{\mathbf {r} }}^{\mathsf {T}}\mathbf {r} ,\quad \mathbf {H} \approx 2{\mathbf {J} _{\mathbf {r} }}^{\mathsf {T}}\mathbf {J_{r}} .}$

Dessa uttryck ersätts i iterationsekvationen ovan för att erhålla ekvationerna

${\displaystyle {\boldsymbol {\beta }}^{(s+1)}={\boldsymbol {\beta }}^{(s)}+\Delta ;\quad \Delta =-\left(\mathbf {J_{r}} ^{\mathsf {T}}\mathbf {J_{r}} \right)^{-1}\mathbf {J_{r}} ^{\mathsf {T}}\mathbf {r} .}$

Konvergens av Gauss-Newton-metoden garanteras inte i alla fall. Uppskattningen

${\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|}$

behöver gälla för att kunna försumma andra ordningens derivator. Det kan ske i två fall och då förväntas konvergens: ^[6]

1. Funktionsvärdena ${\displaystyle r_{i}}$ är små i storleksordningen, åtminstone runt minimum.
2. Funktionerna är bara "milt" olinjära, så att ${\displaystyle {\frac {\partial ^{2}r_{i}}{\partial \beta _{j}\partial \beta _{k}}}}$ är relativt liten i omfattning.

Den här artikeln är helt eller delvis baserad på material från engelskspråkiga Wikipedia.

Noter

Källor

• Björck, Åke (1996). Numerical methods for least squares problems. Philadelphia, Pa: SIAM, Society for Industrial and Applied Mathematics. ISBN 0898713609 
• Fletcher, Roger (1987). Practical methods of optimization (2nd). New York: John Wiley & Sons. ISBN 978-0-471-91547-8. https://archive.org/details/practicalmethods0000flet 
• Nocedal, Jorge; Wright, Stephen (1999). Numerical optimization. New York: Springer. ISBN 0-387-98793-2 
• Ramsin, Håkan (1976). Ickelinjär optimering. Lund: Liber läromedel. ISBN 91-40-04288-X