pracma (version 1.0.5)

gaussNewton: Gauss-Newton Function Minimization

Description

Gauss-Newton method of minimizing a term $f_1(x)^2 + \ldots + f_m(x)^2$ or $F' F$ where $F = (f_1, \ldots, f_m)$ is a multivariate function of $n$ variables, not necessarily $n = m$.

Usage

gaussNewton(x0, Ffun, Jfun = NULL,
                        maxiter =100, tol = .Machine$double.eps^(1/2), ...)

Arguments

Ffun

m functions of n variables.

Jfun

function returning the Jacobian matrix of Ffun; if NULL, the default, the Jacobian will be computed numerically. The gradient of f will be computed internally from the Jacobian (i.

x0

Numeric vector of length n.

maxiter

Maximum number of iterations.

tol

Tolerance, relative accuracy.

...

Additional parameters to be passed to f.

Value

• List with components: xs the minimum or root found so far, fs the square root of sum of squares of the values of f, iter the number of iterations needed, and relerr the absoulte distance between the last two solutions.

Details

Solves the system of equations applying the Gauss-Newton's method. It is especially designed for minimizing a sum-of-squares of functions and can be used to find a common zero of several function. This algorithm is described in detail in the textbook by Antoniou and Lu, incl. different ways to modify and remedy the Hessian if not being positive definite. Here, the approach by Goldfeld, Quandt and Trotter is used, and the hessian modified by the Matthews and Davies algorithm if still not invertible.

To accelerate the iteration, an inexact linesearch is applied.

References

Antoniou, A., and W.-S. Lu (2007). Practical Optimization: Algorithms and Engineering Applications. Springer Business+Science, New York.

See Also

newtonsys, softline

Examples

f1 <- function(x) c(x[1]^2 + x[2]^2 - 1, x[1] + x[2] - 1)
gaussNewton(c(4, 4), f1)

f2 <- function(x) c( x[1] + 10*x[2], sqrt(5)*(x[] - x[4]),
                    (x[2] - 2*x[3])^2, 10*(x[1] - x[4])^2)
gaussNewton(c(-2, -1, 1, 2), f2)

##  Examples (from Matlab)
F1 <- function(x) c(2*x[1]-x[2]-exp(-x[1]), -x[1]+2*x[2]-exp(-x[2]))
gaussNewton(c(-5, -5), F1)

# Find a matrix X such that X %*% X %*% X = [1 2; 3 4]
F2 <- function(x) {
    X <- matrix(x, 2, 2)
    D <- X %*% X %*% X - matrix(c(1,3,2,4), 2, 2)
    return(c(D))
}
sol <- gaussNewton(ones(2,2), F2)
(X  <- matrix(sol$xs, 2, 2))
#            [,1]      [,2]
# [1,] -0.1291489 0.8602157
# [2,]  1.2903236 1.1611747
X %*% X %*% X