Suppose there exists no $x$ such that $Ax=b$. How would you go about choosing $x$ to maximize $|\{i : A_i x = b_i \}|$?

  • 1
    $\begingroup$ Just to indicate some terms to search for: you want to find the sparsest vector in the affine subspace $\operatorname{ran} A - b$. $\endgroup$ – user147263 Dec 10 '15 at 20:39
  • $\begingroup$ Is this an abstract math problem? Or a real world, numerical optimization problem where, because of floating point error and imprecise data, constraint i will be considered satisfied if something like | A_i * x - b_i | < tolerance? $\endgroup$ – Matthew Gunn Dec 10 '15 at 23:06
  • 2
    $\begingroup$ To add on to Normal Human's comment: the sparsest vector is the one which minimizes the so-called "$\ell_0$ norm" $\|y\|_0=|\{i:y_i\ne0\}|$. In other words, you're looking for the $x$ which minimizes $\|Ax-b\|_0$. A common heuristic is to replace it with the $\ell_1$ norm $\|y\|_1=\sum |y_i|$, which leads to a linear programming problem that can be efficiently solved. $\endgroup$ – Rahul Dec 11 '15 at 1:03
  • $\begingroup$ @Rahul and minimizing $||Ax-b||_2$ would simply be least squares. $\mathcal{l}_0$ norm is a nice way to think about it! $\endgroup$ – Matthew Gunn Dec 11 '15 at 1:10
  • $\begingroup$ @MatthewGunn Fair question. I suppose it's an abstract math problem. If I can't satisfy a particular constraint exactly, I don't care how much I miss by. $\endgroup$ – Shane Dec 11 '15 at 1:54

Rahul in the comments made an excellent point that your problem is to minimize the $L_0$ norm (which should be mentioned doesn't satisfy all properties of a norm). The L0 norm is the number of non-zero elements. $$ \text{minimize (over $x$) } \|Ax - b\|_0 $$

This is not a convex problem, and may be difficult to solve. An alternative is the $L_1$ norm, which is a well defined norm and can be minimized efficiently on a computer. The $L_1$ norm may be a reasonable approximation of the $L_0$ norm.

$$ \text{minimize (over $x$) } \|Ax - b\|_1 $$

Moving perhaps further afield from your original objective, another norm to minimize is the 2-norm which would give the solution in the least squares sense.

$$ \text{minimize (over $x$) } \|Ax - b\|_2 $$

Minimizing the 2 norm has the well known analytic solution:

$$ x = (A^{\top}A)^{-1}A^{\top} b $$

where $(A^{\top}A)^{-1}A^{\top}$ is the pseudo inverse for a skinny matrix.


My first hunch would be to create a Mixed Integer Programming model: $$\begin{array}{ll} \min & \sum_i y_i \\ & -My_i \le \sum_j a_{i,j} x_j-b_i \le My_i \\ & y_i \in \{0,1\} \\ \end{array} $$ where $M$ is a large enough constant. There are some other tools developed to help diagnose infeasible Linear Programming models, such as IIS (Irreducable Infeasible Set).

See also the book John W. Chinneck: Feasibility and Infeasibility in Optimization.

In a different answer it is suggested to add a positive slack to each equation and then minimize the sum of the slacks. This is related but I think not the same. The first thing we should notice is that we need to allow positive and negative slacks. We then can either minimize the absolute values or the squared values. I.e.: $$\begin{array}{ll} \min & \sum_i |s_i| \\ & Ax = b + s \\ \end{array} $$

or $$\begin{array}{ll} \min & \sum_i s^2_i \\ & Ax = b + s \\ \end{array} $$ The latter is of course just ordinary least squares. Besides this technical issue there is a more fundamental difference. When minimizing slacks we minimize the sum of of the infeasibilities. This is very different from minimizing the number of infeasibilities. Minimizing the sum of the infeasibilities may lead to a large number of small infeasibilities. In my MIP formulation we explicitly minimize the number of the infeasibilities.

In a reaction we see two objections made to the MIP formulations:

  • How to choose $M$. This is valid point. We need in advance an upper bound on the maximum infeasibility we expect. This is a major drawback.
  • This method would not flag things exactly. The concept of exactness is somewhat difficult when solving a system of linear equations. At least in LP/MIP solvers there are a large number of tolerances involved. But even a straight equation solver will not find exact solutions unless very special techniques are used. I don't think the MIP approach adds much "inexactness".

Another option would be to define a positive penalty $p_i$ for every constraint $i$:

$$ a_ix_i+p_i =b_i $$

and then to minimize $$ \sum_i \omega_ip_i $$ subject to $$ a_ix_i+p_i =b_i\\ p_i\ge 0, $$

where $\omega_i$ is a positive weight associated to $p_i$. If the objective is 0, the solution is feasible. One advantage of this method is that you can control which constraint you want to satisfy in priority (by giving a large $\omega_i$). If you give the same $\omega_i$ to each constraint, you might violate every constraint, but by very little! up to you.

Note. As recommended in some of the comments below, it is a good idea to let $p_i$ be a variable in $\mathbb{R}$ and to minimize the squares, because if you impose $p_i\ge 0$, you will force $a_ix_i<b_i$, whereas it may be more interesting to have $a_ix_i>b_i$.

  • $\begingroup$ Does the minimizer choose the values of $p_i$ or does the person? I don't understand what you mean by, "One advantage of this method is that you can control which constraint you want to satisfy in priority (by giving a large penalty)." $\endgroup$ – Matthew Gunn Dec 11 '15 at 0:18
  • $\begingroup$ Yes you are right, this made no sense. I made the edit accordingly. $\endgroup$ – Kuifje Dec 11 '15 at 0:53
  • $\begingroup$ Minimization is over ($p$, $x$) where $p_i$ is a slack variable and objective is to minimize a weighted sum of the slack. Correct? Do you want to minimize $\sum_i \omega_i \mid p_i \mid$? Whether $a_ix_i$ is $<$ or $>$ than $b_i$ doesn't really matter? $\endgroup$ – Matthew Gunn Dec 11 '15 at 0:57
  • 1
    $\begingroup$ In this case, set $\omega_i=1$ for all $i$, and add the constraint $p_i \in \mathbb{N}$. It will work as, $\sum_i p_i$ then represents the number of violated constraints. $\endgroup$ – Kuifje Dec 11 '15 at 3:16
  • 1
    $\begingroup$ @Kuifje I do not think that that actually works. Your objective function evaluates to the same thing if I miss on two constraints (missing by one on each) as it would if I hit one constraint exactly and miss by two on the other, when I want it to favor the latter. $\endgroup$ – Shane Dec 11 '15 at 14:57

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.