# A quick algorithm for calculating the $\ell_1$-distance between two finite sets on the real line?

For two non-empty finite sets $A,B$ in the real line define the $\ell_1$-distance $d_1(A,B)$ between $A$ and $B$ as the smallest Lebesgue measure of a closed subset $\Gamma\subset \mathbb R$ such that $A\cup B\subset\Gamma$ and each connected component of $\Gamma$ intersects both sets $A$ and $B$. It is clear that such minimal set $\Gamma$ is a disjoint union of closed intervals with ends in the set $A\cup B$. As was noted by @fedja in his comment, the number of such unions is $\le 2^{n-1}$ where $n$ is the cardinality of the set $A\cup B$.

So, a brute force algorithm for calculating the distance $d_1(A,B)$ has exponential complexity.

Problem. Is there a polynomial-time algorithm for calculating $d_1(A,B)$? Or this problem is NP-hard? NP-complete?

Remark. A similar problem for the plane seems to be NP-hard, see Remark 2 here.

• The number of unions of disjoint intervals with endpoints in $S$ is $2^{|S|-1}$ (each interval between two consecutive points is either in the union or not), so the factorial bound is certainly not optimal. – fedja Aug 2 '17 at 6:26
• This is surely an inefficient way of going about things, but you can at least solve this as a linear program. As you note, the optimal solution is a disjoint union of closed intervals with ends in $A\cup B$, and each such interval contains at least one element from each set. There are $M=O(n^2)$ such intervals. – John Gunnar Carlsson Aug 2 '17 at 6:38
• (continued) Now, let $c_j$ be the length of interval $j$ (according to some enumeration), and let $Q$ be a matrix such that $q_{ij}=1$ if point $i$ belongs to interval $j$ and $0$ otherwise. You then want to solve the linear program of choosing $x_j$'s to minimize $\sum_{j=1}^M c_j x_j$ subject to the constraint that $\sum_{j} q_{ij}x_j \geq 1$ for all $i$ and that $x_j\geq 0$ for all $j$. The matrix $Q$ is totally unimodular because it has the consecutive ones property, so there is guaranteed to be a solution with all $x_j$'s being either $0$ or $1$. – John Gunnar Carlsson Aug 2 '17 at 6:38
• ...and those entries $j$ for which the optimal solution $x_j^*=1$ are the optimal intervals for your problem. – John Gunnar Carlsson Aug 2 '17 at 6:52
• I assume that you want $A \cup B \subset \Gamma.$ What you wrote seems to require just the shortest interval with one endpoint in $A$ and the other in $B.$ – Aaron Meyerowitz Aug 2 '17 at 9:07

While it doesn’t seem to be possible to beat the $O(n\log n)$ bound for general input, the problem can be solved in time $O(n)$ if $A\cup B$ is given sorted.

This builds on the ideas in fedja’s and Aaron Meyerowitz’s answers. As implicit in the other posts, I am assuming unit-cost real RAM or a similar model, so that arithmetic operations on real numbers are exact and take constant time.

Let $A\cup B=\{x_1,\dots,x_n\}$ with $x_1<x_2<\dots<x_n$. In order to compute $d_1(A,B)$, we will compute the sequence of numbers \begin{align*} u_i&=d_1\bigl(A\cap[x_1,x_i],B\cap[x_1,x_i]\bigr),\\ v_i&=d_1\bigl((A\cup\{x_i\})\cap[x_1,x_i],(B\cup\{x_i\})\cap[x_1,x_i]\bigr) \end{align*} for $i=1,\dots,n$. This can be done in linear time using the following recurrence: \begin{align*} v_1&=0,\\ u_1&=\begin{cases}0&x_1\in A\cap B,\\\infty&\text{otherwise,}\end{cases}\\ v_{i+1}&=\min\bigl\{u_i,v_i+(x_{i+1}-x_i)\bigr\},\\ u_{i+1}&=\begin{cases} v_{i+1}&x_{i+1}\in A\cap B,\\ v_i+(x_{i+1}-x_i)&x_{i+1}\in A\smallsetminus B\text{ and }x_i\in B,\text{ or v.v.,}\\ u_i+(x_{i+1}-x_i)&x_{i+1},x_i\in A\smallsetminus B,\text{ or v.v.} \end{cases} \end{align*} (Some optimization is possible, treating blocks of elements belonging to the same set in bulk. However, this does not improve the worst-case asymptotic complexity.) Note that the two cases in the minimum in the computation of $v_{i+1}$ correspond to the choice whether $[x_i,x_{i+1}]$ is in $\Gamma$ or not. It is straightforward to extend the algorithm so that it also outputs a witnessing sequence of intervals.

Hausdorff distance can be also computed in time $O(n)$ on sorted input, using the following pseudocode:

/* for each i, find distance to the maximal y \le x_i that belongs to the other set: */
lastA := lastB := -\infty
for i := 1 to n do:
if x_i in A then lastA := x_i
if x_i in B then lastB := x_i
if x_i in A then d_i := x_i - lastB
else d_i := x_i - lastA
/* do the same thing backwards, and compute the result along the way */
lastA := lastB := \infty
d := 0
for i := n downto 1 do:
if x_i in A then lastA := x_i
if x_i in B then lastB := x_i
if x_i in A then d_i := min (d_i, lastB - x_i)
else d_i := min (d_i, lastA - x_i)
d := max (d, d_i)
return d

• By the way, what the complexity of calculation of Hausdorff distance between two finite sets in $\mathbb R^n$ (endowed with the $\ell_p$-metric for $p\in\{1,2,\infty\}$)? Is it still $O(n\ln n)$ (or the answer depends on $p$)? The brute force algorith has compexity $O(n^2)$. – Taras Banakh Aug 4 '17 at 12:36
• I don’t know. However, I suspect that you can’t significantly beat the brute force algorithm, at least if you want to compute the distance exactly, because of the following heuristic: on the one hand, if all pairwise distances are not too different from each other, you need to more-or-less know them all to compute the result. Now, in dimension $d=1$, all the $n^2$ pairwise distances are determined by just $n$ of them (say, distances to a fixed point $x_1$) plus $O(n)$ bits of discrete information (say, the relative order of the points wrt $x_1$). This is not true in higher dimensions. ... – Emil Jeřábek Aug 4 '17 at 13:38
• ... While the triangle inequality gives you bounds, it doesn’t let you to compute the distance between two points from distances to some other points. Thus, I assume that all the pairwise distances are not determines by any number $o(n^2)$ of them. – Emil Jeřábek Aug 4 '17 at 13:40
• At least for the plane endowed with the $\ell_p$-metric for $p\in\{1,\infty\}$ there exists an algorithm of complexity $O(n\ln n)$ calculating the Hausdorff distance $d_\infty(A,B)$ between two non-empty finite sets $A,B$. The algorithm uses the standard sweep line'' approach, see Exercise 2 here (wuecampus2.uni-wuerzburg.de/moodle/pluginfile.php/250296/…). Maybe there exists a better reference than this homework task from Wuerzburg University. I do not known what happen in the Euclidean metric (maybe seepline works too)? – Taras Banakh Aug 4 '17 at 14:11
• Oh, I see. This should give something like $O(dn\log n)$ for $\ell_\infty$, and $O(2^dn\log n)$ for $\ell_1$, in dimension $d$. However, I don’t see how this could be adapted to other $\ell_p$; balls of distinct diameter are not piecewise translations of each other. – Emil Jeřábek Aug 4 '17 at 15:26

Of course one does not want brute force, but the case count is not as bad as suggested. I think a reasonable upper bound on the number of cases might be (a shift of the) the Padovan Sequence $1,1,1,2,2,3,4,5,7,12,\cdots$ with $p_n=p_{n-2}+p_{n-3}$ which is $O(\theta^n)$ for $\theta\approx 1.3245\cdots$ the unique real root of $t^3-t-1$. This is based on the idea that the hardest case is if the points are alternately in $A,B,A,B,\cdots.$

I will sketch this and give similar bounds. Consider the $n-1$ intervals $I_j=[x_j,x_{j+1}]$ (with the points $x_1 \lt x_2\lt \cdots)$ Then $\Gamma$ is one of the $2^{n-1}$ ways to pick some of those sub-intervals. If $A,B$ are disjoint that falls to $2^{n-3}$ since we can't omit $I_1$ nor $I_{n-1}$ and in fact to $F_{n-2} \lt 2^{n-3}$ since we can't omit two in a row. Of course there is no advantage to keeping all the intervals. We haven't taken into account $A$ and $B.$ Here is why strick alternation is harder than the other cases. Given a run $x_{i+1},\cdots, x_j$ all in $A$ but not $B,$ We can't omit more than one of $I_i,\cdots I_j$ and we should omit one. If there is one of $I_2,\cdots,I_{j-1}$ which is of the maximal length we should omit it. Otherwise the choices are the longest of those (which is second or third longest) or $I_1$ or $I_j$. We may need to see the other choices first before deciding on the second or third longest.) A point in $A \cap B$ splits the problem into two almost disjoint ones.

Assuming the $n$ points $x_1,x_2,x_3,\cdots,x_n$ do alternate $A,B,A,B,\cdots,$ let $c_n$ be the number of possible cases for $\Gamma$ where we assume nothing about the actual values of the $x_i$ except that they increase. Then $c_2=c_3=c_4=1.$ I will now observe that $c_n=c_{n-2}+c_{n-3}$ for $n \ge 5.$

We will have to use the interval $I_1$ . Record this, then shift to considering $x_2,x_3,\cdots,x_{n}$ where the $x_2$ is now considered universal , i.e. in $A \cap B,$ If we do not use [x_2,x_3] then we have the $c_{n-2}$ cases for $x+3,x_4,\cdots,x_n$ If we do use $[x_2,x_3],$ it is for the benefit of $x_3$ in which case we do not use $[x_3,x_4].$ This gives the $c_{n-3}$ cases for $x_4,x_5,\cdots, x_n.$

Once the members of $A \cup B$ are listed in increasing order we can in $O(n)$ time list for each $x_i$ the closest $x_j \le x_i \le x_k$ with $x_j,x_k \in B$ ($A$ resp.) if $x_i \in A$ ($B$ resp.) I feel that $O(n)$ time might be enough to finish, but I could be wrong.

LATER A few more observations:

Every complete interval in $\Gamma$ is of type $AB,ABA,BA$ or $BAB$ where type $AB$ means one or more points, all from $A$ (in increasing order) followed by one or more all from $B.$ However type $ABA$ means that with one exception every point is from $A$ and the excepional $B$ is not an endpoint.

In particular, if the points happen to be alternately from $A$ and $B$ (the hardest case in my opinion) then $\Gamma$ consists of a subset of the intervals $I_1,\cdots,I_{n-1}$ which never uses three in a row nor skips two subintervals in a row. In other words, looking at how the points are grouped, the number of possibilities is the number of ways to write $n$ as an ordered sum of $2$s and $3$s

If it happens that $I_k$ is longer than $I_{k-1}$ and $I_{k+1}$ combined then we must not use $I_k$ and hence are forced to use both $I_{k-1}$ and $I_{k+1}$

All this could be adapted to the general case.

• The greedy method does not guarantee a minimum. I can imagine you being successful on one point set with a given order type and certain interval weights, and fail on a similar example with the same order type and the same relative weight sizes, but with the weights tweaked slightly. Gerhard "Works For Upper Bounds Though" Paseman, 2017.08.03. – Gerhard Paseman Aug 3 '17 at 23:48
• @Gerhard You are right. I rolled it back. – Aaron Meyerowitz Aug 4 '17 at 0:09
• It seems that there exists an algorithm of complexity $O(n\ln n)$ for calculating $d_1(A,B)$. Now I will try to write down the proof. – Taras Banakh Aug 4 '17 at 5:16

Somebody should check the argument below because I'm still terribly sleepy but it looks like you can get away with simple dynamic programming. Indeed, let's go from the left and see what's going on. First we are forced to connect the leftmost point with the next point in the other set. That creates certain length and a leftmost "universal point" for the remaining configuration so that any interval containing that universal point satisfies the requirement that it has both types of points in it. You can now join a few next points backward to that universal point but once you make a gap, you start all over and are forced to create the next universal point at some particular position.

This calls for going from the right and solving the optimization problem with a universal point to the right of that point for each position of the universal point in $A\cup B$ . We need just to record the total length and the last point joined to the universal point in the optimal configuration to keep track of everything. At each step consider all possibilities to make a gap and the next universal point ($n$ choices at most) and minimize the sum of the required length between the universal points (finding that length takes time $n$ at most, but you can utilize the information about the previous gap choice, so, probably, constant on average) and the minimum for the next universal point (already computed). So, the total time is at most $n^3$ for sure, probably $n^2$.

Again, sorry if I wrote some nonsense.

• In the simplest situation, when you just delete the lefmost point of $A\cup B$ and reduce the problem to smaller number of points, you obtain exponential (more precisely Fibonaccy) complexity: the complexity of calculation of the distance for $n$ points is $\le$ the complexity for $n-1$ points (with one lefmost point deleted) + complexity for $n-2$ points (with two lefmost points removed). I did not catch the trick which lowers the complexity to $O(n^3)$. – Taras Banakh Aug 2 '17 at 7:37
• @TarasBanakh What exactly is unclear? – fedja Aug 2 '17 at 8:17
• I have understood your idea. Indeed, it can work but some more justifications should be made. The essence that at each next step we should not return back very far to make a gap (at most to the previous return point) and this gives you a polynomial complexity. Anyway this should be elaborated with more precision. – Taras Banakh Aug 2 '17 at 8:25
• I think the description is unnecessarily complicated. The only information you need to record for each point is the $\ell_1$-distance of the parts of $A$ and $B$ to the right of that point (inclusive). Then, in order to compute this for the next point $x$, try all points $y\ge x$ such that $[x,y]$ intersects both $A$ and $B$; for each of them, the cost is $Y-X$ plus the already recorded distance to the right of the first point to the right of $y$. So, you have $\le n$ choices for each of the $n$ points, hence you only need $O(n^2)$ arithmetic operations. – Emil Jeřábek Aug 2 '17 at 14:27
• @TarasBanakh Please, read properly what I wrote. When considering $x$, I am computing $d_1(A\cap[x,\infty),B\cap[x,\infty))$. So the only way I could have an interval $[z,x]$ is in the degenerate case $z=x$ (which is allowed iff $x\in A\cap B$), but then it is of the form $[x,y]$ for $y=x$. – Emil Jeřábek Aug 2 '17 at 14:54

Solution (developing ideas of @fedja, @EmilJerabek, and @AaronMeyerowitz):
There exists an algorithm of complexity $O(n\ln n)$ for calculating $d_1(A,B)$
(and simultaneously presenting the witness set $\Gamma$ of Lebesgue measure equal to $d_1(A,B)$).

The algorithm consists of 8 steps.

1. Given two non-empty finite sets $A,B\subset\mathbb R$ write the union $A\cup B$ as $\{x_1,\dots,x_n\}$ where $x_1<\cdots<x_n$. Also keep the sets $A,B$ sorted in binary search trees in order to make quick search.

2. Find the smallest number $s\in\{1,\dots,n\}$ such that the segment $[x_1,x_s]$ intersects $A$ and $B$ (so, $x_s=\max\{\min A,\min B\}$).

3. For every $i\in\{s,\dots,n\}$ find the largest number $l(i)\in\{1,\dots,i\}$ such that the segment $[x_{l(i)},x_i]$ intersects both sets $A$ and $B$;

4. Find the largest number $\lambda\in\{1,\dots,n\}$ such that the segment $[x_\lambda,x_n]$ intersects $A$ and $B$;

5. For every $i\in\{1,\dots,\lambda\}$ find the smallest number $r(i)\in\{i,\dots,n\}$ such that the segment $[x_i,x_{r(i)}]$ intersects both sets $A$ and $B$;

6. Let $D(s):=|x_s-x_1|$ and $g(s)=1$ (here $D$ stands for "distance" and $g$ for "gap");

7. For every $k\in\{s+1,\dots,n\}$ define the real number $D(k)$ and an integer number $g(k)$ as follows.
$\bullet$ If $r(g(k-1))\ge l(k)$, then put $g(k):=g(k-1)$ and $D(k):=D(k-1)+x_k-x_{k-1}$.
$\bullet$ If $r(g(k-1))<l(k)$, then put $D(k):=\min_{i}D(i-1)+|x_k-x_{i}|$ where the minimum is taken over all integer numbers $i$ such that $r(g(k-1))<i\le l(k)$. Put $g(k)$ be the largest integer number $i$ such that $r(g(k-1))<i\le l(k)$ and $D(k)=D(i-1)+|x_k-x_i|$.

8. Print the number $D(n)$. It is equal to $d_1(A,B)$. More generally, $D(k)$ is equal to $d_1(A\cap{\downarrow}x_k],B\cap{\downarrow}x_k)$ where ${\downarrow}x_k:=(-\infty,x_k]$.

Simultaneously we can also find a closed set $\Gamma\subset\mathbb R$ of the smallest Lebesgue measure such that $A\cup B\subset\Gamma$ and each connected component intersects both sets $A$ and $B$. For this put $\Gamma(0)=\emptyset$ and for every $k\in\{s,\dots, n\}$ let $\Gamma(k):=\Gamma(g(k)-1)\cup [x_{g(k)},x_k]$. Then the set $\Gamma:=\Gamma(n)$ will have the required property.

Now let us evaluate the complexity of this algorithm.

The step 1 (sorting) requires $O(n\ln n)$ operations.
The step 2 (finding $s$) requires $O(\ln n)$ operations.
The step 3 (calculating the function $l(\cdot)$) requires $O(n\ln n)$ operations using two the binary searches through the sets $A,B$).
By analogy, the steps 4,5 require $O(n\ln n)$ operations
The steps 6--8 require $O(n)$ operations (since the search of $i$ goes through disjoint sets).

So, in total the algorithm has complexity $O(n\ln n)$.

Remark. The presented algorithm for calculating $d_1(A,B)$ has the same complexity $O(n\ln n)$ as the (evident) algorithm for calculating the Hausdorff distance $$d_\infty(A,B)=\max\{\max_{a\in A}\min_{b\in B}|a-b|,\max_{b\in B}\min_{a\in A}|a-b|\}$$ (using the binary search for calculating $\min_{b\in B}|b-a|$ in time $O(\ln n)$).

Acknowledgements: I would like to thank the MathOverflow users @fedja, @EmilJerabek, @AaronMeyerowitz and all others participating in the discussion for presenting their valuable ideas.