Gauss–Newtons metod

Gauss-Newtons metod används för att lösa icke-linjära minsta kvadrat-problem. Dessa uppstår till exempel vid icke-linjär regression, där parametrar i en modell söks så att modellen stämmer väl överens med tillgängliga observationer.

Anpassning av en brusig kurva med en asymmetrisk modell för topparna med Gauss-Newton-algoritmen med variabel dämpningsfaktor α.

Överst: Rådata och modell.

Nederst: Utveckling av den normaliserade kvadratsumman av felen.

Det är en variant av Newtons metod för att hitta ett minimum av en funktion. Till skillnad från Newtons metod kan Gauss-Newton-algoritmen endast användas för att minimera summan av kvadrerade funktionsvärden, men den har fördelen att andraderivator, som kan vara svåra att beräkna, inte krävs.

Metoden är uppkallad efter matematikerna Carl Friedrich Gauss och Isaac Newton och presenterades först i Gauss verk från 1809 Theoria motus corporum coelestium in sectionibus conicis solem ambientum.

Beskrivning

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

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

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

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

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

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

Beräkningar

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

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

När $\mathbf {r}$ är komplex $\mathbf {r} :\mathbb {C} ^{n}\rightarrow \mathbb {C}$ den konjugerade formen ska användas: $\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

${\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 $\Delta$ . De kan lösas i ett steg, med hjälp av Choleskyuppdelning, eller, bättre, QR-faktorisering av $\mathbf {J_{r}}$. 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å $\mathbf {J_{r}} ^{T}\mathbf {J_{r}}$ blir singular.

Beräkningar för dataanpassning

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

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

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

${\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 $\left(\mathbf {J_{f}} ^{\mathsf {T}}\mathbf {J_{f}} \right)^{-1}\mathbf {J_{f}} ^{\mathsf {T}}$ är den vänstra pseudoinversen av $\mathbf {J_{f}}$ .

Exempel

Beräknad kurva erhållen med ${\hat {\beta }}_{1}=0,362$ och ${\hat {\beta }}_{2}=0,556$ (i blått) kontra observerade data (i rött).

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

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

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

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

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

minimeras.

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

${\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 $\beta _{1}=0,9$ och $\beta _{2}=0,2$ och efter fem iterationer av Gauss-Newton-metoden erhålls de optimala värdena ${\hat {\beta }}_{1}=0,362$ och ${\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.

Härledning från Newtons metod

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.

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

${\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 $S=\sum _{i=1}^{m}r_{i}^{2}$, ges gradienten av

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

Hessianens element beräknas genom att derivera gradientelementen, $g_{j}$, med avseende på $\beta _{k}$ :

$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

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

där $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

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

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

$\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:

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

Referenser

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

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