Broyden–Fletcher–Goldfarb–Shanno algorithm

In numerical optimization, the Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm is an iterative method for solving unconstrained nonlinear optimization problems. Like the related Davidon–Fletcher–Powell method, BFGS determines the descent direction by preconditioning the gradient with curvature information. It does so by gradually improving an approximation to the Hessian matrix of the loss function, obtained only from gradient evaluations (or approximate gradient evaluations) via a generalized secant method.
Since the updates of the BFGS curvature matrix do not require matrix inversion, its computational complexity is only formula_1, compared to formula_2 in Newton's method. Also in common use is L-BFGS, which is a limited-memory version of BFGS that is particularly suited to problems with very large numbers of variables (e.g., &gt;1000). The BFGS-B variant handles simple box constraints.
The algorithm is named after Charles George Broyden, Roger Fletcher, Donald Goldfarb and David Shanno.
Rationale.
The optimization problem is to minimize formula_3, where formula_4 is a vector in formula_5, and formula_6 is a differentiable scalar function. There are no constraints on the values that formula_4 can take.
The algorithm begins at an initial estimate for the optimal value formula_8 and proceeds iteratively to get a better estimate at each stage.
The search direction p"k" at stage "k" is given by the solution of the analogue of the Newton equation:
where formula_10 is an approximation to the Hessian matrix, which is updated iteratively at each stage, and formula_11 is the gradient of the function evaluated at x"k". A line search in the direction p"k" is then used to find the next point x"k"+1 by minimizing formula_12 over the scalar formula_13
The quasi-Newton condition imposed on the update of formula_10 is
Let formula_16 and formula_17, then formula_18 satisfies 
which is the secant equation. 
The curvature condition formula_20 should be satisfied for formula_18 to be positive definite, which can be verified by pre-multiplying the secant equation with formula_22. If the function is not strongly convex, then the condition has to be enforced explicitly e.g. by finding a point x"k"+1 satisfying the Wolfe conditions, which entail the curvature condition, using line search.
Instead of requiring the full Hessian matrix at the point formula_23 to be computed as formula_18, the approximate Hessian at stage "k" is updated by the addition of two matrices:
Both formula_26 and formula_27 are symmetric rank-one matrices, but their sum is a rank-two update matrix. BFGS and DFP updating matrix both differ from its predecessor by a rank-two matrix. Another simpler rank-one method is known as symmetric rank-one method, which does not guarantee the positive definiteness. In order to maintain the symmetry and positive definiteness of formula_18, the update form can be chosen as formula_29. Imposing the secant condition, formula_30. Choosing formula_31 and formula_32, we can obtain:
Finally, we substitute formula_35 and formula_36 into formula_29 and get the update equation of formula_18:
Algorithm.
From an initial guess formula_8 and an approximate Hessian matrix formula_41 the following steps are repeated as formula_42 converges to the solution:
formula_3 denotes the objective function to be minimized. Convergence can be checked by observing the norm of the gradient, formula_53. If formula_41 is initialized with formula_55, the first step will be equivalent to a gradient descent, but further steps are more and more refined by formula_56, the approximation to the Hessian.
The first step of the algorithm is carried out using the inverse of the matrix formula_10, which can be obtained efficiently by applying the Sherman–Morrison formula to the step 5 of the algorithm, giving
This can be computed efficiently without temporary matrices, recognizing that formula_59 is symmetric,
and that formula_60 and formula_61 are scalars, using an expansion such as
Therefore, in order to avoid any matrix inversion, the inverse of the Hessian can be approximated instead of the Hessian itself: formula_63
From an initial guess formula_8 and an approximate inverted Hessian matrix formula_65 the following steps are repeated as formula_42 converges to the solution:
In statistical estimation problems (such as maximum likelihood or Bayesian inference), credible intervals or confidence intervals for the solution can be estimated from the inverse of the final Hessian matrix . However, these quantities are technically defined by the true Hessian matrix, and the BFGS approximation may not converge to the true Hessian matrix.
Notable implementations.
Notable open source implementations are:
Notable proprietary implementations include: