算法实现LM算法

A Brief Description of theLevenberg-Marquardt Algorithm Implemenedby levmarManolis I. A. LourakisInstitute of Computer ScienceFoundation for Research and Technology - Hellas (FORTH)Vassilika Vouton, P.O. Box 1385, GR711 10Heraklion, Crete, GREECEFebruary 11, 2005AbstractThe Levenberg-Maiquardt (LM) algorithm is an iterative technique that locates the minimum of a function that is expressed as the sum of squares of nonlinear functions. It has become a standard technique for nonlinear least-squares problems and can be thought of as a combination of steepest descent and the Gauss-Newton method. This dDcumcnt briefly describes the nialiicmalius behind levmar, a free LM C/C+ iinplciiiciilalioii that uun be found at http: /www. ics . forth. gr/ lourakis/levmar.IntroductionThe Levenberg-Marquardt (LM) algorithm is an iterative technique that locates the minimum of a multivariate function that is expressed as the sum of squares of non-linear real-valued functions |4, 6|. It has become a standard technique for non-linear least-squares problems 7, widely adopted in a broad spectrum of disciplines. LM can be thought of as a combination of steepest descent and the Gauss-Newton method When the current solution is far from the correct one, the algorithm behaves like a steepest descent method: slow, but guaranteed to converge When the current solution is close to the correct solution, it becomes a Gauss-Newton method. Next, a short description of the LM algorithm based on the material in 5 is supplied Note, however, that a detailed analysis of the LM algorithm is beyond the scope of this report and the interested reader is referred to 5, & 9, 2, 10 for more comprehensive treatments.The Levenberg-Marquardt AlgorithmIn the following, vectors and arrays appear in boldface and T is used to denote transposition. Also, |.| and denote (he 2 and infinity norms respectively. Let / be an assumed functional relation which maps a parameter vector p e 冗 to an estimated measurement vector x = /(p)5 x e 冗.An initial parameter estimate p0 and a measured vector x are provided and it is desired to find the vector p+ that best satisfies the functional relation /, i.e. minimizes the squared distance f1( with e = x x. The basis of the LM algorithm is a linear approximation to f in the neighborhood of p. For a small | 0. The strategy of altering the diagonal elements of J7 J is called damping and “ is referred to as the damping term. If the updated parameter vector p + dp with dp computed from Eq. (3) leads to a reduction in the error c the update is accepted and the process repeats with a decreased damping term Otherwise, the damping term is increased, the augmented normal equations are solved again and the process iterates until a value of 6P that decreases error is found The process of repeatedly solving Eq. (3) for different values of the damping term until an acceptable update to the parameter vector is found corresponds to one iteration of the LM algorithmIn LM. the damping term is adjusted at each iteration to assure a reduction in the error & If the damping is set to a large value, matrix N in Eq(3) is nearly diagonal and the LM update step 6P is near the steepest descent direction. Moreover, the magnitude of 6P is reduced in this case Damping also handles situations where the Jacobian is rank deficient and J7 J is therefore singular 3. In this way, LM can defensively navigate a region of the parameter space in which the model is highly nonlinear. If the damping is small, the LM step approximates the exact quadratic step appropriate for a fully linear problem LM is adaptive because it controls its own damping: it raises the damping if a step fails to reduce e； otherwise it reduces the damping In this way LM is capable to alternate between a slow descent approach when being far from the minimum and a fast convergence when being at the minimums neighborhood 3 The LM algorithm terminates when at least one of the following conditions is met: The magnitude of the gradient of i.e. J re in the right hand side of Eq. (2). drops below a threshold “ The relative change in the magnitude of 6p drops below a threshold f.2 The error 6re drops below a threshold e3 A maximum number of iterations kmax is completedIf a covariance matrix Sx for the measured vector x is available, it can be incorporated into the LM algorithm by minimizing the squared Snorm c7 instead of the Euclidean e7 e. Accordingly, the minimum is found by solving a weighted least squares problem defined by the weighted normal equationsThe rest of the algorithm remains unchanged. The complete LM algorithm is shown in pseudocode in Fig. 1. It is derived by slight modification of algorithm3.16 in page 27 of |5; more details regarding the LM algorithm can be found there. Indicative values for the user-defined parameters are r = 103,= e2 =e3 = 1015, kmax = 100. levmar is a free C/C+ implementation of this LM algorithm that can be found at http: /www. ics forth . gr/ * lourakis/levmar.References1 G.H. Golub and C.F. Van Loan. Matrix Computations. Johns Hopkins University Press, 3rd editioir 1996.2 C.T. Kelley. Iterative Methods for Optimization. SIAM Press, Philadelphia. 1999.3 M Lainpton.Damping-Undamping Strategies for the Levenberg-Marquardt Nonlinear Least-Squares Method. Computers in Physics Journal. 11(1):110-115, Jan./Feb. 1997.4 K. Levenberg. A Method for the Solution of Certain Non-linear Problems in Least Squares. Quarterly of Applied Mathematics. 2(2):164-16 & Jul. 1944.5 K. Madsen, H.B. Nielsen, and O. Tingleff.Methodsfor Non-Linear Least Squares Problems Technical University of Denmark, 2004Lecture notes, available athttp:/wwwimm.dtudk/courses/02611/nllsq.pdf.|6| D.W. Marquardt. An Algorithm for the Least-Squares Estimation of Nonlinear Parameters. SIAM Journal of Applied Mathematics, 11 (2):431141, Jim. 1963.7 H.D Mittelmann The Least Squares Problemweb pagehttp:/platoasuedu/topics/problems/nlolsq.html, Jul. 2004. Accessed on 4 Aug. 2004.8 H.B. Nielsen. Damping Parameter in Marquardts Method. Technical Report IMM-REP-1999-05, Technical University of Denmark.1999. Available at http:/wwwimm.dtu.dk/hbn.J. Noceclal and S.J. Wright. Numerical Optimization. Springer New York, 1999Input: A vector function f : fJVn t TV with n m, a measurement vector x e TV1 and an initial parameters estimate p0 EOutput: A vector p+ e R minimizing |x - /(p)|2.Algorithm:R := 0; “ := 2; p := p0；A ：= J7 J；p ：= X - /(p)； g ：= JOp；stop：=(|g|oo 6)；“ ：= 丁* ma&iwSii)； while (not stop) and (k 0P = PnewiA := J7 J；6p := x - /(p)； g ：= Jrep;StOP：=(l|g|oo i)or(|p|2 0) or (stop)endwhileFigure 1: Levenberg-Marquardt non-linear least squares algorithm; see text and 5, 8 for details.|10 W.H. Press, S.A. Teukolsky, A.W.T. Vetterling, and B.P. Flannery. Numerical Recipes in C: The Art of Scientific Computing Cambridge University Press. New York, 1992.