The Gauss-Newton method (after Carl Friedrich Gauß and Isaac Newton ) is a numerical method for solving nonlinear minimization problems using the least squares method . The method is related to the Newton method for solving non-linear optimization problems, but has the advantage that the calculation of the 2nd derivative required for the Newton method is not necessary. The calculation of the 2nd derivative is often a limiting factor, especially for large problems with tens of thousands of parameters.
The optimization problem
The Gauss-Newton method solves problems in which the minimum of a sum of squares of continuously differentiable functions is sought, i.e.
![{\ displaystyle f_ {i} \ colon \ mathbb {R} ^ {n} \ to \ mathbb {R}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/15b8d5ffd08da06e3e721ba3637800ef97d19db3)
![{\ displaystyle \ min _ {x \ in \ mathbb {R} ^ {n}} \ left \ {{\ frac {1} {2}} \ sum _ {i = 1} ^ {m} \ left (f_ {i} (x) \ right) ^ {2} \ right \}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/6483ece2909a0eeff297744ef619033794576622)
with . With the Euclidean norm, this can also be written as
![\ | \ cdot \ |](https://wikimedia.org/api/rest_v1/media/math/render/svg/113f0d8fe6108fc1c5e9802f7c3f634f5480b3d1)
![{\ displaystyle \ min _ {x \ in \ mathbb {R} ^ {n}} \ left \ {{\ frac {1} {2}} \ | f (x) \ | ^ {2} \ right \} }](https://wikimedia.org/api/rest_v1/media/math/render/svg/e73a22b2fa282193ba3918d0a4a22b9695c2e2a6)
with . Problems of this form occur frequently in practice, in particular the nonlinear problem is equivalent to minimizing , provided that it has a zero . Is a linear map , there is the standard case, the method of least squares linear model function .
![{\ displaystyle f = (f_ {1}, \ dotsc, f_ {m}) \ colon \ mathbb {R} ^ {n} \ to \ mathbb {R} ^ {m}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/27a352b900bb828fb20a4b41f329e54c678b3314)
![f (x) = 0](https://wikimedia.org/api/rest_v1/media/math/render/svg/cf85883d74b75fe35ca8d3f2b44802df078e4fa1)
![{\ displaystyle {\ tfrac {1} {2}} \ | f (x) \ | ^ {2}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/117ac1d52b0ae73603f573d1858a6bbdc1ea2551)
![f](https://wikimedia.org/api/rest_v1/media/math/render/svg/132e57acb643253e7810ee9702d9581f159a1c61)
![f](https://wikimedia.org/api/rest_v1/media/math/render/svg/132e57acb643253e7810ee9702d9581f159a1c61)
The method
The basic idea of the Gauss-Newton method is to linearize the objective function and to optimize the linearization in the sense of the least squares. The linearization, i.e. the 1st order Taylor expansion , of in the point reads
![f](https://wikimedia.org/api/rest_v1/media/math/render/svg/132e57acb643253e7810ee9702d9581f159a1c61)
![f](https://wikimedia.org/api/rest_v1/media/math/render/svg/132e57acb643253e7810ee9702d9581f159a1c61)
![{\ displaystyle x ^ {0} \ in \ mathbb {R} ^ {n}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a017c1b48b3cbc52be9d856dc280b00cf89ce701)
-
.
The matrix is the Jacobi matrix and is often referred to as. The linear least squares problem is obtained
![{\ displaystyle \ nabla f (x ^ {0}) ^ {T}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/7de46f87b295e61204bb1eaa0cd884e889c51532)
![J](https://wikimedia.org/api/rest_v1/media/math/render/svg/359e4f407b49910e02c27c2f52e87a36cd74c053)
-
,
with gradient .
![{\ displaystyle \ nabla {\ frac {1} {2}} \ left \ Vert {\ tilde {f}} (x) \ right \ Vert ^ {2} = J ^ {T} \ left (J (xx ^ {0}) + f (x ^ {0}) \ right)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/11ea79126660c8fa4e18944b71a60fe9264d05bb)
Setting the gradient to zero provides the so-called normal equations
![{\ displaystyle J ^ {T} J (xx ^ {0}) = - J ^ {T} f (x ^ {0})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/dae3c78647d989328595dd5b2e493f360e064cca)
with the explicit solution
-
.
The Gauss-Newton iteration step results directly from this
,
|
making it clear that the Jacobi matrix is evaluated at this point and is a step size .
![{\ displaystyle J | _ {x ^ {k}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d0fe1d482e50acfdddf002b9e2109a1b70de6aef)
![x ^ k](https://wikimedia.org/api/rest_v1/media/math/render/svg/525e1133440e1565055dec6243aaf0f27d4d4e9b)
![{\ displaystyle \ alpha ^ {k} \ geq 0}](https://wikimedia.org/api/rest_v1/media/math/render/svg/6e739350b732eb5fec93f82fb585bfa6d537ccf6)
To solve the linear system of equations in the Gauss-Newton iteration step, there are different possibilities depending on the problem size and the structure:
- Small problems ( ) are best solved with QR decomposition
![{\ displaystyle n <1000, \ m <10000}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ad9d117fcbc1070f71a48b35abdff10c462671ea)
- The Cholesky decomposition is ideal for large problems , since the matrix is symmetrical by construction . There are specially adapted Cholesky variants for sparsely populated
![{\ displaystyle J ^ {T} J}](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3ee4504ba79b22818c8871d44729bce9a8339b4)
- The CG method can be used as a general option , although preconditioning is usually necessary here
convergence
The update vector in the Gauss-Newton iteration step has the form where . If full rank has, then and therefore also positive definite . On the other hand, the gradient of the quadratic problem , so a descent direction , i.e. H. it applies . From this follows (with a suitable choice of the step size ) the convergence of the Gauss-Newton method to a stationary point . From this representation it can also be seen that the Gauss-Newton method is essentially a scaled gradient method with the positively definite scaling matrix .
![{\ displaystyle d = -DJ ^ {T} f (x)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/97bc10b856ebc9ec5d295f2f3a3cb5fb13f02f57)
![{\ displaystyle D = \ left (J ^ {T} J \ right) ^ {- 1}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f116fae9ebee0de8fb57d4a00afea0a7e7a91dbc)
![{\ displaystyle J ^ {T} J}](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3ee4504ba79b22818c8871d44729bce9a8339b4)
![{\ displaystyle J ^ {T} f (x)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/de9d88c4d0f47d333b9de53a6308eaa904ce9cfc)
![{\ displaystyle \ min _ {x \ in \ mathbb {R} ^ {n}} {\ frac {1} {2}} \ | f (x) \ | ^ {2}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/0938b24395fcf39862512157a839b51fb85ef0d7)
![d](https://wikimedia.org/api/rest_v1/media/math/render/svg/e85ff03cbe0c7341af6b982e47e9f90d235c66ab)
![{\ displaystyle \ left (J ^ {T} f (x) \ right) ^ {T} d <0}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a6573e5e07c726357c239f00b166a9da8275b84d)
![D.](https://wikimedia.org/api/rest_v1/media/math/render/svg/f34a0c600395e5d4345287e21fb26efd386990e6)
In general, no statement can be made about the speed of convergence . If the starting point is very far from the optimum or the matrix is poorly conditioned , the Gauss-Newton method converges as slowly as desired. On the other hand, if the starting point is sufficiently close to the optimum, one can show that the Gauss-Newton method converges quadratically.
![x ^ {0}](https://wikimedia.org/api/rest_v1/media/math/render/svg/1871ffeb57c11624b375dbb7157d5887c706eb87)
![{\ displaystyle J ^ {T} J}](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3ee4504ba79b22818c8871d44729bce9a8339b4)
![x ^ {0}](https://wikimedia.org/api/rest_v1/media/math/render/svg/1871ffeb57c11624b375dbb7157d5887c706eb87)
extension
In order to improve the behavior in the case of badly conditioned or singular , the Gauss-Newton iteration step can be modified as follows
![{\ displaystyle J ^ {T} J}](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3ee4504ba79b22818c8871d44729bce9a8339b4)
-
,
where the diagonal matrix is chosen so that positive is definite. With the choice , i.e. a scalar multiple of the identity matrix , one obtains the Levenberg-Marquardt algorithm .
![{\ displaystyle \ Delta ^ {k}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ee42513cd347f0c36b666f632d1a8d515edd93d4)
![{\ displaystyle J | _ {x ^ {k}} ^ {T} J | _ {x ^ {k}} + \ Delta ^ {k}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/17bb6159245a00a65fcb541b74454b2327005212)
![{\ displaystyle \ Delta ^ {k} = \ beta ^ {k} I, \ \ beta \ geq 0}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d4a13916b42dad1246b9c88dfbf9ff4740fe3b5c)
example
The Rosenbrock function with .
![{\ displaystyle a = 1, \ b = 100}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d496b3e5f2df4a1ddbc780cea9185ddff46d65a7)
The Rosenbrock function
![{\ displaystyle g: \ mathbb {R} ^ {2} \ to \ mathbb {R}: x \ mapsto \ left (a-x_ {1} \ right) ^ {2} + b \ left (x_ {2} -x_ {1} ^ {2} \ right) ^ {2}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a3258ab7c236f1b2f3ac8233ce945d3b6fbdaff8)
is often used as a test for optimization methods because it is challenging because of the narrow and shallow valley in which iterative methods can only take small steps. The constants are usually chosen with , the global optimum in this case is with the function value .
![{\ displaystyle a = 1, \ b = 100}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d496b3e5f2df4a1ddbc780cea9185ddff46d65a7)
![{\ displaystyle x ^ {*} = (1,1)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/3a119753af09978d0ea1fe3a4f7690d1542d9837)
![{\ displaystyle g (x ^ {*}) = 0}](https://wikimedia.org/api/rest_v1/media/math/render/svg/07b1834847c045bf868b8caba66e9d1d5c359023)
In order to use the Gauss-Newton method, the Rosenbrock function must first be brought into the form "sum of squares of functions". Since the Rosenbrock function already consists of a sum of two terms, the approach is chosen
![{\ displaystyle {\ frac {1} {2}} \ left (f_ {1} (x) \ right) ^ {2} = \ left (a-x_ {1} \ right) ^ {2} \ \ \ Longleftrightarrow \ \ f_ {1} (x) = {\ sqrt {2}} (a-x_ {1})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/990f6da44f09d98617c02c87ab6ece41b95d5ba5)
and
-
.
The Gauss-Newton problem for the Rosenbrock function is thus
-
whereby .![{\ displaystyle f: \ mathbb {R} ^ {2} \ to \ mathbb {R} ^ {2}: x \ mapsto {\ begin {pmatrix} f_ {1} (x) \\ f_ {2} (x ) \ end {pmatrix}} = {\ begin {pmatrix} {\ sqrt {2}} (a-x_ {1}) \\ {\ sqrt {2b}} (x_ {2} -x_ {1} ^ { 2}) \ end {pmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/1496c62c998c780bf0af23076c0f6112fbb9d054)
The Jacobi matrix results as and thus is . Since has full rank , positive is definite and the inverse exists. The following simple line search is used to determine the step size :
![{\ displaystyle J = {\ begin {bmatrix} {\ tfrac {\ partial f_ {1}} {\ partial x_ {1}}} & {\ tfrac {\ partial f_ {1}} {\ partial x_ {2} }} \\ {\ tfrac {\ partial f_ {2}} {\ partial x_ {1}}} & {\ tfrac {\ partial f_ {2}} {\ partial x_ {2}}} \ end {bmatrix} } = {\ begin {bmatrix} - {\ sqrt {2}} & 0 \\ - 2x_ {1} {\ sqrt {2b}} & {\ sqrt {2b}} \ end {bmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/2f174f186a234e89d5b0cf0e6a8ccae0a9e286a7)
![{\ displaystyle D = J ^ {T} J = {\ begin {bmatrix} 8bx_ {1} ^ {2} + 2 & -4bx_ {1} \\ - 4bx_ {1} & 2b \ end {bmatrix}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/e858c48a744425741f456606daf11cdaa5e220da)
![J](https://wikimedia.org/api/rest_v1/media/math/render/svg/359e4f407b49910e02c27c2f52e87a36cd74c053)
![D.](https://wikimedia.org/api/rest_v1/media/math/render/svg/f34a0c600395e5d4345287e21fb26efd386990e6)
![D ^ {- 1}](https://wikimedia.org/api/rest_v1/media/math/render/svg/05cb013783b23fe6dc52b27a1f286f47c028bd52)
![\ alpha ^ k](https://wikimedia.org/api/rest_v1/media/math/render/svg/cbe2221400f734971e3aa0cd94f6703241e9fee1)
- Start with .
![{\ displaystyle \ alpha ^ {k} = 1}](https://wikimedia.org/api/rest_v1/media/math/render/svg/e5871f6a657499b6b1229e6ad1d5e8484a0ecbf9)
- Calculate the new point as well .
![{\ displaystyle {\ tilde {x}} = x ^ {k} + \ alpha ^ {k} d}](https://wikimedia.org/api/rest_v1/media/math/render/svg/aa95670f75b806dcc07ac99c47c344761e197a2e)
![{\ displaystyle d = - \ left ((J | _ {x ^ {k}}) ^ {T} J | _ {x ^ {k}} \ right) ^ {- 1} (J | _ {x ^ {k}}) ^ {T} f (x ^ {k})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/b0e9ccab5903254cb0692deceddd8d5588f78cb2)
- If sit and go to the next iteration.
![{\ displaystyle f ({\ tilde {x}}) <f (x ^ {k})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/5144e32bca9cbdcc039002a76b37cd011f6270fe)
![{\ displaystyle x ^ {k + 1} = {\ tilde {x}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/4d6b52df50a43c74d0c6aa9b951e367798903285)
- Otherwise halve and go to 2.
![\ alpha ^ k](https://wikimedia.org/api/rest_v1/media/math/render/svg/cbe2221400f734971e3aa0cd94f6703241e9fee1)
The line search forces the new function value to be less than the previous one; it is guaranteed to terminate (with possibly a very small one ), since it is a direction of descent.
![\ alpha ^ k](https://wikimedia.org/api/rest_v1/media/math/render/svg/cbe2221400f734971e3aa0cd94f6703241e9fee1)
![d](https://wikimedia.org/api/rest_v1/media/math/render/svg/e85ff03cbe0c7341af6b982e47e9f90d235c66ab)
The starting point is chosen. The Gauss-Newton method converges to the global optimum in a few iterations:
![{\ displaystyle x ^ {0} = (0, -0.1)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/38e8b91b43f46e69cf67b1cad2e8317710a5d423)
Optimization of the Rosenbrock function with the Gauss-Newton method
Optimization with the Gauss-Newton method
![k](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3c9a2c7b599b37105512c5d570edc034056dd40) |
![x ^ k](https://wikimedia.org/api/rest_v1/media/math/render/svg/525e1133440e1565055dec6243aaf0f27d4d4e9b) |
|
0 |
(0, -0.1) |
2
|
1 |
(0.1250, -0.0875) |
1.8291
|
2 |
(0.2344, -0.0473) |
1.6306
|
3 |
(0.4258, 0.0680) |
1.6131
|
4th |
(0.5693, 0.2186) |
1.3000
|
5 |
(0.7847, 0.5166) |
1.0300
|
6th |
(1.0, 0.9536) |
0.2150
|
7th
|
(1.0, 1.0)
|
1.1212e-27
|
The gradient method (with the same line search) provides the following result in comparison, it does not find the optimum even after 500 iterations:
Optimization of the Rosenbrock function with the gradient method
Optimization with gradient methods
|
|
|
0
|
(0, -0.1)
|
2
|
1
|
(0.0156, 0.0562)
|
1.2827
|
2
|
(0.0337, -0.0313)
|
1.0386
|
3
|
(0.0454, 0.0194)
|
0.9411
|
4th
|
(0.0628, -0.0077)
|
0.8918
|
5
|
(0.0875, 0.0286)
|
0.8765
|
|
|
|
500
|
(0.8513, 0.7233)
|
0.0223
|
literature
- Dimitri P. Bertsekas: '' Nonlinear Programming. '' Second Edition, Athena Scientific, 1995, ISBN 9781886529144 .
- Yurii Nesterov: " Introductory Lectures on Convex Optimization: A Basic Course. " Springer Science & Business Media, 2003, ISBN 978-1-4419-8853-9 .
- Jorge Nocedal, Stephen Wright: " Numerical Optimization. " Springer Science & Business Media, 2000, ISBN 9780387987934 .
- Amir Beck: " Introduction to Nonlinear Optimization. " SIAM, 2014, ISBN 978-1611973648 .
Individual evidence
-
↑ Ceres Solver documentation. Retrieved May 10, 2019 .