Approximate subset sum with two-dimensional vectors

Consider the following optimization problem:

Given $n\leq 10^3$ vectors $v_i\in\mathbb{R}^2$, all of which are small, i.e., $\|v_i\| \leq 1$, find a subset $S$ of them that minimizes $\| w + \sum_{i\in S} v_i \|^2$, where $w$ is a fixed known vector.

I know that this can be reduced to an unconstrained 0-1 quadratic program (with unknowns $x_i\in\{0,1\}$, $x_i=1$ whenever $i\in S$), but that is a hard problem in general.

However, since the constraints are so specific, is there a way to get a reasonably simple efficient approximation algorithm for it? One that doesn't involve integer programming? It's not even clear to me how to efficiently implement a greedy algorithm for it.

Are there special cases of this problem that are easier? For example, would it help if all the vectors were unit length, $\|v_i\|=1$?

• – D.W. May 24 '17 at 23:33
• Your problem is the Closest Vector Problem (CVP) in an integer lattice. In general, the CVP is fairly hard. However, you are working in only 2 dimensions, so the problem might become tractable. It's known that the 2D SVP can be solved in polynomial time (algorithm due to Gauss) if you start with only 2 vectors. I don't know if it is solvable if you start with $n>2$ vectors, or if it extends to the CVP. – D.W. May 24 '17 at 23:50
• @D.W. I don't follow you: in a lattice, the coefficients are just arbitrary integers. How do you reduce the 0-1 problem to a lattice problem? – Kirill May 25 '17 at 0:22
• You are absolutely right; my mistake. I overlooked that we need the coefficients to be not just integers, but either zero or one. To try to save my idea: I think we still might be able to formulate as a lattice problem by replacing the vector $(x_i,y_i)$ with a vector of the form $(x_i,y_i,0,0,\dots,0,K,0,\dots,0)$ where $K$ is a sufficiently large constant, but then we're not in a lattice of dimension two. That might lead to an approximation algorithm using LLL lattice reduction, but I haven't thought about it carefully. – D.W. May 25 '17 at 16:42
• Why isn't the objective function to be minimized the following? $$\left\| \left( \sum_{i\in S} v_i \right) - w \right\|_2^2$$ – Rodrigo de Azevedo Jun 4 '17 at 2:45

Let $\rm V$ be the $2 \times n$ matrix whose $i$-th column is vector $\rm v_i$. We have the Boolean optimization problem in $\mathrm z \in \{0,1\}^n$

$$\min_{\mathrm z \in \{0,1\}^n} \| \mathrm V \mathrm z + \mathrm w \|_2^2$$

Let

$$\rm z = \frac 12 (x + 1_n)$$

where $\mathrm x \in \{\pm 1\}^n$. Hence, we have an instance of the Boolean least-squares (BLS) problem

$$\min_{\mathrm x \in \{\pm 1\}^n} \| \mathrm A \mathrm x - \mathrm b \|_2^2$$

where

$$\rm A = \frac 12 V \qquad \qquad \qquad b = - \left( \frac 12 V 1_n + w \right)$$

Since $\{\pm 1\}$ is the solution set of the quadratic equation $x_i^2 = 1$, we have the following (non-convex) quadratically constrained quadratic program (QCQP) in $\mathrm x \in \mathbb R^n$

$$\begin{array}{ll} \text{minimize} & \| \mathrm A \mathrm x - \mathrm b \|_2^2\\ \text{subject to} & x_i^2 = 1 \quad \forall i \in \{1,2,\dots,n\}\end{array}$$

Exhaustive evaluation of the objective function at all $2^n$ points is only feasible for small $n$.

SDP relaxation

Note that

$$\begin{array}{rl} \| \mathrm A \mathrm x - \mathrm b \|_2^2 &= (\mathrm A \mathrm x - \mathrm b)^{\top} (\mathrm A \mathrm x - \mathrm b)\\ &= \mbox{tr} \left( (\mathrm A \mathrm x - \mathrm b) (\mathrm A \mathrm x - \mathrm b)^{\top} \right)\\ &= \mbox{tr} \left( \begin{bmatrix} \mathrm A & -\mathrm b\end{bmatrix} \begin{bmatrix} \mathrm x\\ 1\end{bmatrix} \begin{bmatrix} \mathrm x\\ 1\end{bmatrix}^{\top} \begin{bmatrix} \mathrm A & -\mathrm b\end{bmatrix}^{\top} \right)\\ &= \mbox{tr} \left( \begin{bmatrix} \mathrm A^{\top}\\ -\mathrm b^{\top}\end{bmatrix} \begin{bmatrix} \mathrm A & -\mathrm b\end{bmatrix} \begin{bmatrix} \mathrm x\\ 1\end{bmatrix} \begin{bmatrix} \mathrm x^{\top} & 1\end{bmatrix} \right)\\ &= \mbox{tr} \left( \begin{bmatrix} \,\,\,\, \mathrm A^{\top} \mathrm A & -\mathrm A^{\top}\mathrm b\\ -\mathrm b^{\top}\mathrm A & \,\,\,\, \mathrm b^{\top}\mathrm b\end{bmatrix} \begin{bmatrix} \mathrm x \mathrm x^{\top} & \mathrm x\\ \mathrm x^{\top} & 1\end{bmatrix} \right)\end{array}$$

Since $\mathrm x \in \{\pm 1\}^n$, all $n$ entries on the main diagonal of $\mathrm x \mathrm x^{\top}$ are equal to $1$. Thus,

$$\begin{bmatrix} \mathrm x \mathrm x^{\top} & \mathrm x\\ \mathrm x^{\top} & 1\end{bmatrix}$$

is symmetric, positive semidefinite and has only ones on its main diagonal, i.e., it is a correlation matrix. Its rank is $1$. Let

$$\mathrm C := \begin{bmatrix} \,\,\,\, \mathrm A^{\top} \mathrm A & -\mathrm A^{\top}\mathrm b\\ -\mathrm b^{\top}\mathrm A & \,\,\,\, \mathrm b^{\top}\mathrm b\end{bmatrix}$$

Hence, the BLS problem can be written as the following rank-constrained optimization problem in $(n+1) \times (n+1)$ symmetric matrix $\rm Y$

$$\begin{array}{ll} \text{minimize} & \langle \mathrm C , \mathrm Y \rangle\\ \text{subject to} & y_{ii} = 1, \quad \forall i \in \{1,2,\dots,n+1\}\\ & \mathrm Y \succeq \mathrm O_{n+1}\\ & \mbox{rank} (\mathrm Y) = 1\end{array}$$

which is a hard problem due to the rank constraint. Relaxing the optimization problem above by discarding the rank constraint, we then obtain the following semidefinite program (SDP) in $\rm Y$

$$\boxed{\begin{array}{ll} \text{minimize} & \langle \mathrm C , \mathrm Y \rangle\\ \text{subject to} & y_{ii} = 1, \quad \forall i \in \{1,2,\dots,n+1\}\\ & \mathrm Y \succeq \mathrm O_{n+1}\end{array}}$$

which is convex and computationally tractable. This SDP relaxation provides a lower bound on the minimum of the original BLS problem. In the very, very fortunate case where the optimal solution of the SDP happens to be rank-$1$, we have solved the BLS problem.

References

Consider the following approach to approximate this problem.

First, construct an efficient data structure using your vectors that can perform Nearest Neighbor related searching, such as a KD Tree. Then, formulate your problem as an optimization that's a function of $N$, the number of nearest neighbors you will compute. To get there, we can do the following with the cost function:

\begin{align} J &= \left\lVert w + \sum_{i \in S} v_i\right\Vert^2 \\ &= \frac{N^2}{N^2}\left\lVert w + \sum_{i \in S_N} v_i\right\Vert^2 \\ &= N^2\left\lVert \frac{w}{N} + \frac{1}{N}\sum_{i \in S_N} v_i\right\Vert^2 \end{align}

Note that $S_N$ is then an approximate solution for $S$. We can then approximate this problem by choosing $N$ such that we generate a set of vectors, which can be represented by the index set $S_N$, that are the $N$ nearest neighbors to $\frac{-w}{N}$. The logic is that selecting $N$ nearest neighbors from $\lbrace v_i \rbrace$ about $\frac{-w}{N}$ will produce a value for $\frac{1}{N}\sum_{i \in S_N} v_i$ close to $\frac{-w}{N}$ and in turn help minimize the norm.

The resulting optimization then ends up being:

\begin{align} N^{*} = \arg \min_{N} N^2 \left\lVert \frac{w}{N} + \frac{1}{N}\sum_{i \in S_N} v_i\right\Vert^2 \end{align}

How you go about tackling the above optimization is then up to you.

• Thank you for the answer. Does this assume that the vectors are uniformly distributed? I don't want to assume that, so I'm not sure about the step where the average of $N$ nearest neighbours should be close to $\frac{-w}{N}$. – Kirill May 26 '17 at 20:46
• @Kirill This approximation scheme does not depend on that, though I expect it would perform better if a subset of the vectors were at least uniformly distributed near $\frac{-w}{N}$. The level of accuracy this might achieve certainly is dependent on the distribution of the vectors. – spektr May 26 '17 at 21:05