# What is the fastest way to solve Ax=b (subject to constraints and an absolute term)

I am trying to solve/optimize $$Ax=b$$ in the least squares sense subject to

1. box constraints;
2. a few (less than 5) equality/inequality constraints; and
3. an absolute function penalty (or some other piecewise linear convex function).

That is, min $$f(\mathbf{x})=|A\mathbf{x}-b|^2+\mathbf{c}\cdot|\mathbf{x}|$$ where $$\mathbf{c}\cdot|\mathbf{x}|=c_1|x_1|+...+c_n|x_n|$$

In particular, $$A\mathbf{x}=b$$ is massively underdetermined (i.e. $$A$$ has several thousand columns but hundreds of rows) and any valid solution will do if many exist. This code is performance critical. I've exhausted everything I can think of. Is there any obvious way to solve this efficiently (either a clever formulation or appropriate choice of algorithm)?

Language agnostic solutions would be appreciated.

tl;dr - Skip the rest. Here's what I have tried so far

### Approach 1 - Keep it simple

Use Bounded BFGS and evaluate $$f(\mathbf{x})$$ directly. I am trying to improve on this as it is too slow.

### Approach 2 - QP:

$$f(\mathbf{x})=\mathbf{x^T}A^TA\mathbf{x}-2\mathbf{x^T}A^Tb+\mathbf{c}\cdot|\mathbf{x}|=\frac{1}{2}\mathbf{x^T}Q\mathbf{x}+\mathbf{d}\cdot\mathbf{x}+\mathbf{c}\cdot|\mathbf{x}|$$

So this is basically standard quadratic form but $$Q$$ is now enormous and positive semi-definite making it difficult to work with. Slightly slower.

### Approach 2 - Use QR decomposition:

The idea is to reduce the size of the quadratic matrix using QR decomposition and a co-ordinate transformation.

$$A^T=QR$$ and $$\mathbf{x^T}Q=\mathbf{y^T}$$

$$f(\mathbf{x})=\mathbf{x^T}A^TA\mathbf{x}-2\mathbf{x^T}A^Tb+\mathbf{c}\cdot|\mathbf{x}|=\mathbf{y^T}RR^T\mathbf{y}+\mathbf{d}\cdot Q\mathbf{y}+\mathbf{c}\cdot|Q\mathbf{y}|$$

This makes $$A^TA$$ very small and fast to execute. However, the transformation forces the thousands of (quick) box constraints to become linear inequality constraints and this hammers the performance. A lot slower.

Any other ideas? Thanks.

• Is $A$ sparse or dense? Do you expect the solution to be extremely sparse or not sparse and/or have many variables at their upper/lower bounds? Have you considered scaling the problem so that the $c \cdot | x |$ term simply becomes $\| x \|_{1}$ and then using methods specialized for 1-norm regularized least squares? – Brian Borchers Mar 28 '18 at 1:44
• Are you going to be solving lots of these problems in sequence with small modifications each time? – Brian Borchers Mar 28 '18 at 1:47
• @BrianBorchers $A$ is dense. I expect the solution to have >90% of the variables at the lower bound. I have not considered that scaling as that is not something I am familiar with. If you could point me to any literature on the matter, I'd be more than happy to work through it. I do not expect to be solving lots of these problems with small modifications. – Big AL Mar 28 '18 at 6:58
• I'm assuming that every variable has a negative lower bound and a positive upper bound, right? – Brian Borchers Mar 28 '18 at 12:22
• @BrianBorchers That's correct, yes, although some lower bounds are fixed at zero. If I have understood your post below correctly, I would not be obliged to 'double up' those variables. – Big AL Mar 28 '18 at 20:31

This problem can be formulated as a standard positive semidefinite QP with bounded variables.

First, deal with the absolute value terms in the objective by letting $x=u-v$, where $u\geq 0$ and $v \geq 0$. Then write the problem as

$\min \;\; (A(u-v)-b)^{T}(A(u-v)-b)+c^{T}u+c^{T}v$

subject to

$B(u-v)=d$ (equality constraints)

Upper bounds on u and v.

$u,v \geq 0$.

Then use an active set method for QP to solve this QP. In this approach, most of the variables are fixed at their lower or upper bounds, which dramatically reduces the size of the QP. You solve a simple linear programming feasibility problem with the constraints $B(u-v)=d$ and the bounds constraints to find an initial feasible solution and an initial active set of variables.

In each iteration of the active set method you solve the reduced size QP over the current set of active variables, and then check optimality conditions to see if any of the fixed variables should be released from their bounds and whether any of the free variables should be pinned to their upper or lower bounds.

Active set methods for QP are widely available- don't reinvent the wheel.

• Wow. Thank you. I think I understood all of that but I might take some time to mull it over. Can I check, are $w$ and $z$ directly related to $u$ and $v$. Is it just a change in notation or have I fundamentally misunderstood? Sorry, this does not always come naturally to me. – Big AL Mar 28 '18 at 20:28
• Sorry, I switched from $w$ and $z$ to $u$ and $v$ half way through the answer. I've edited it to make it consistent. – Brian Borchers Mar 28 '18 at 22:16