Linear Least-Squares Fitting with linear constraints

Solves linearly constrained linear least-squares problems.

lsqlincon(C, d,  A = NULL, b = NULL,
          Aeq = NULL, beq = NULL, lb = NULL,  ub = NULL)

Arguments

C

mxn-matrix defining the least-squares problem.

d

vector or a one colum matrix with m rows

A

pxn-matrix for the linear inequality constraints.

b

vector or px1-matrix, right hand side for the constraints.

Aeq

qxn-matrix for the linear equality constraints.

beq

vector or qx1-matrix, right hand side for the constraints.

lb

lower bounds, a scalar will be extended to length n.

ub

upper bounds, a scalar will be extended to length n.

Details

lsqlincon(C, d, A, b, Aeq, beq, lb, ub) minimizes ||C*x - d|| (i.e., in the least-squares sense) subject to the following constraints: A*x <= b, Aeq*x = beq, and lb <= x <= ub.

It applies the quadratic solver in quadprog with an active-set method for solving quadratic programming problems.

If some constraints are NULL (the default), they will not be taken into account. In case no constraints are given at all, it simply uses qr.solve.

Value

Returns the least-squares solution as a vector.

Note

Function lsqlin in pracma solves this for equality constraints only, by computing a base for the nullspace of Aeq. But for linear inequality constraints there is no simple linear algebra `trick', thus a real optimization solver is needed.

Author

HwB email: <hwborchers@googlemail.com>

References

Trefethen, L. N., and D. Bau III. (1997). Numerical Linear Algebra. SIAM, Society for Industrial and Applied Mathematics, Philadelphia.

See also

lsqlin, quadprog::solve.QP

Examples

##  MATLABs lsqlin example
C <- matrix(c(
    0.9501,   0.7620,   0.6153,   0.4057,
    0.2311,   0.4564,   0.7919,   0.9354,
    0.6068,   0.0185,   0.9218,   0.9169,
    0.4859,   0.8214,   0.7382,   0.4102,
    0.8912,   0.4447,   0.1762,   0.8936), 5, 4, byrow=TRUE)
d <- c(0.0578, 0.3528, 0.8131, 0.0098, 0.1388)
A <- matrix(c(
    0.2027,   0.2721,   0.7467,   0.4659,
    0.1987,   0.1988,   0.4450,   0.4186,
    0.6037,   0.0152,   0.9318,   0.8462), 3, 4, byrow=TRUE)
b <- c(0.5251, 0.2026, 0.6721)
Aeq <- matrix(c(3, 5, 7, 9), 1)
beq <- 4
lb <- rep(-0.1, 4)   # lower and upper bounds
ub <- rep( 2.0, 4)

x <- lsqlincon(C, d, A, b, Aeq, beq, lb, ub)
# -0.1000000 -0.1000000  0.1599088  0.4089598
# check A %*% x - b >= 0
# check Aeq %*% x - beq == 0
# check sum((C %*% x - d)^2)    # 0.1695104