Let $\mathbf G$ be a given $m \times n$ matrix, $\mathbf y$ a given $m \times 1$ column vector and $\mathbf x$ an unknown $n \times 1$ column vector such that $\mathbf x \ge 0$.

1) How do you find $\displaystyle \min_{\mathbf x} (\mathbf y - \mathbf G\mathbf x)^T(\mathbf y - \mathbf G\mathbf x)$?

2) Whether there is a solution of that problem similar to the simple pseudo-inverse $x=(G^TG)^{−1}G^Ty$ used for solving the least-squares problem?

3) And what if we have $n \times n$ matrix $\mathbf C$ and the constraint $\mathbf C\mathbf x \ge 0$ instead of $\mathbf x \ge 0$?

  • $\begingroup$ Looks like least-squares with positivity constraint. I was about to write up on the three usual methods for solving least squares problems until I noticed the $\mathbf x \geq \mathbf 0$ part. I'll have to check on how the usual schemes are modified for finding positive solutions. $\endgroup$ – J. M. is a poor mathematician Dec 5 '10 at 17:31
  • $\begingroup$ This might have what you're looking for. I'm not posting this as an answer for now since I don't have access to the paper at the moment for me to substantially discuss this. $\endgroup$ – J. M. is a poor mathematician Dec 5 '10 at 17:46
  • 4
    $\begingroup$ This is just quadratic programming. There are lots of classical algorithms available. Your particular instance is nice because the matrix involved is positive definite, and the constraints are all bound constraints. Some recent algorithms focus on such special cases, such as the one J.M. pointed to. $\endgroup$ – Rahul Dec 5 '10 at 20:08
  • 2
    $\begingroup$ "the matrix involved is positive definite" - iff $\mathbf G$ has rank $n$; otherwise, it's positive semidefinite. $\endgroup$ – J. M. is a poor mathematician Dec 6 '10 at 0:12
  • 1
    $\begingroup$ @Max: You haven't said whether you want (a) an analytical solution (there isn't one), (b) a numerical algorithm (there are many), (c) an efficient numerical algorithm (). $\endgroup$ – Rahul Dec 8 '10 at 15:03

You can start by writing the expression as

\begin{equation} F(x) = (y^T - x^T G^T)(y - Gx) = y^T y - y^T Gx - x^TG^Ty + x^TG^TGx \end{equation}

You can then differentiate wrt $x$ and set the expression to zero.

\begin{equation} F'(x) = -2G^Ty + 2G^TG x = 0 \end{equation}

So you get $x = (G^TG)^{-1}G^T y$.

  • 3
    $\begingroup$ That's vanilla least squares. Did you notice the positivity constraint? $\endgroup$ – J. M. is a poor mathematician Dec 5 '10 at 17:44
  • $\begingroup$ Sorry I didn't notice the positivity constraint. Thanks for pointing that out. $\endgroup$ – svenkatr Dec 5 '10 at 18:04
  • $\begingroup$ In order for $G^TG$ to be invertible we need to know that $G$ has full rank. $\endgroup$ – Yiorgos S. Smyrlis Jan 7 '14 at 16:41

First, let's take a look at your objective function:

$$ \begin{align} f(x)&=(y-Gx)^{T}(y-Gx) \\ &=y^Ty-y^TGx-x^TGy+x^TG^TGx \\ &=y^Ty-y^TGx-y^TGx+x^TG^TGx \\ &=y^Ty-2y^TGx+x^TG^TGx \\ \end{align} $$

This is a quadratic function, and we can put it in a more standard form as follows: $$ \begin{align} a&=y^Ty \\ b^T&=-2y^TG \\ b&=-2G^Ty \\ C&=2G^TG \\ f(x)&=a + b^Tx + \frac{1}{2}x^TCx \\ \end{align} $$

Then you offer two possibilities for constraints: inequality constraints, or linear inequality constraints.

Calculus tells us that for a bounded optimization problem, the optimum value will be either a point in the interior of the permissible region at a spot where the gradient of the function is zero, or on the bounds of the permissible region.

This is pretty complicated to do explicitly, and not a viable strategy to code specifically for your problem unless $n\le2$.

All of this is really just to back up the following claim about your question (2): sadly, there is no nice simple formula like in the unconstrained case.

Fortunately, your problem fits into a common and well-researched class of problems: quadratic minimization with linear inequality constraints, which is a subset of quadratic programming. Thus, to answer (1) and (3), you simply need to look up and implement one of several well-known algorithms (a decent reference I googled on the topic: http://www.math.uh.edu/~rohop/fall_06/Chapter3.pdf), or use an existing implementation for your programming language/environment of choice. Examples include:

Googling "[your language name] quadratic programming" will likely give you some useful stuff for nearly any common language. Using an existing implementation is a much better idea than writing your own, it's pretty hard to get right.

Finally, if you can't find a specialized quadratic optimizer, you may find a general nonlinear optimizer that takes inequality constraints. It will be a performance hit relative to a specialized quadratic algorithm, but unless your problem is huge, or you have a very tight computation time budget, it will probably get the job done.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.