lsei {limSolve} | R Documentation |
Solves an lsei inverse problem (Least Squares with Equality and Inequality Constraints)
\min(||Ax-b||^2)
subject to
Ex=f
Gx>=h
Uses either subroutine lsei (FORTRAN) from the LINPACK package, or
solve.QP
from R-package quadprog
.
In case the equality constraints Ex=f cannot be satisfied, a generalized inverse solution residual vector length is obtained for f-Ex.
This is the minimal length possible for ||f-Ex||^2.
lsei (A = NULL, B = NULL, E = NULL, F = NULL, G = NULL, H = NULL, Wx = NULL, Wa = NULL, type = 1, tol = sqrt(.Machine$double.eps), tolrank = NULL, fulloutput = FALSE, verbose = TRUE)
A |
numeric matrix containing the coefficients of the quadratic
function to be minimised, ||Ax-B||^2; if the columns of |
B |
numeric vector containing the right-hand side of the quadratic function to be minimised. |
E |
numeric matrix containing the coefficients of the equality
constraints, Ex=F; if the columns of |
F |
numeric vector containing the right-hand side of the equality constraints. |
G |
numeric matrix containing the coefficients of the inequality
constraints, Gx>=H; if the columns of |
H |
numeric vector containing the right-hand side of the inequality constraints. |
Wx |
numeric vector with weighting coefficients of unknowns (length = number of unknowns). |
Wa |
numeric vector with weighting coefficients of the quadratic function (Ax-B) to be minimised (length = number of number of rows of A). |
type |
integer code determining algorithm to use 1= |
tol |
tolerance (for singular value decomposition, equality and inequality constraints). |
tolrank |
only used if |
fulloutput |
if |
verbose |
logical to print error messages. |
a list containing:
X |
vector containing the solution of the least squares problem. |
residualNorm |
scalar, the sum of absolute values of residuals of equalities and violated inequalities. |
solutionNorm |
scalar, the value of the minimised quadratic function at the solution, i.e. the value of ||Ax-b||^2. |
IsError |
logical, |
type |
the string "lsei", such that how the solution was obtained can be traced. |
covar |
covariance matrix of the solution; only returned if
|
RankEq |
rank of the equality constraint matrix.; only returned if
|
RankApp |
rank of the reduced least squares problem (approximate
equations); only returned if |
See comments in the original code for more details; these comments are included in the ‘docs’ subroutine of the package.
Karline Soetaert <karline.soetaert@nioz.nl>
K. H. Haskell and R. J. Hanson, An algorithm for linear least squares problems with equality and nonnegativity constraints, Report SAND77-0552, Sandia Laboratories, June 1978.
K. H. Haskell and R. J. Hanson, Selected algorithms for the linearly constrained least squares problem - a users guide, Report SAND78-1290, Sandia Laboratories,August 1979.
K. H. Haskell and R. J. Hanson, An algorithm for linear least squares problems with equality and nonnegativity constraints, Mathematical Programming 21 (1981), pp. 98-118.
R. J. Hanson and K. H. Haskell, Two algorithms for the linearly constrained least squares problem, ACM Transactions on Mathematical Software, September 1982.
Berwin A. Turlach R and Andreas Weingessel (2007). quadprog: Functions to solve Quadratic Programming Problems. R package version 1.4-11. S original by Berwin A. Turlach R port by Andreas Weingessel.
solve.QR
the original function from package quadprog
.
# ------------------------------------------------------------------------------ # example 1: polynomial fitting # ------------------------------------------------------------------------------ x <- 1:5 y <- c(9, 8, 6, 7, 5) plot(x, y, main = "Polynomial fitting, using lsei", cex = 1.5, pch = 16, ylim = c(4, 10)) # 1-st order A <- cbind(rep(1, 5), x) B <- y cf <- lsei(A, B)$X abline(coef = cf) # 2-nd order A <- cbind(A, x^2) cf <- lsei(A, B)$X curve(cf[1] + cf[2]*x + cf[3]*x^2, add = TRUE, lty = 2) # 3-rd order A <- cbind(A, x^3) cf <- lsei(A, B)$X curve(cf[1] + cf[2]*x + cf[3]*x^2 + cf[4]*x^3, add = TRUE, lty = 3) # 4-th order A <- cbind(A, x^4) cf <- lsei(A, B)$X curve(cf[1] + cf[2]*x + cf[3]*x^2 + cf[4]*x^3 + cf[5]*x^4, add = TRUE, lty = 4) legend("bottomleft", c("1st-order", "2nd-order","3rd-order","4th-order"), lty = 1:4) # ------------------------------------------------------------------------------ # example 2: equalities, approximate equalities and inequalities # ------------------------------------------------------------------------------ A <- matrix(nrow = 4, ncol = 3, data = c(3, 1, 2, 0, 2, 0, 0, 1, 1, 0, 2, 0)) B <- c(2, 1, 8, 3) E <- c(0, 1, 0) F <- 3 G <- matrix(nrow = 2, ncol = 3, byrow = TRUE, data = c(-1, 2, 0, 1, 0, -1)) H <- c(-3, 2) lsei(E = E, F = F, A = A, B = B, G = G, H = H)