An efficient approximation for point-set diameter
in higher dimensions

Given a finite set $\mathcal{S}$ of $n$ points, the diameter of $\mathcal{S}$, denoted by $D(\mathcal{S})$ is the maximum distance between two points of $\mathcal{S}$. Namely, we want to find a diametrical pair $p$ and $q$ such that $D(\mathcal{S})=\max_{\forall{p,q\in\mathcal{S}}}(||p-q||)$. Computing the diameter of a set of points has a large history, and it may be required in various fields such as database, data mining, and vision. A trivial brute-force algorithm for this problem is as follows: compute the distance between each pair of points and then choose the maximum distance. Since computing the distance takes constant time, this algorithm takes $O(n^{2})$ time, but this is too slow for large-scale datasets that occur in the fields. Hence, we need a faster algorithm which may be exact or is an approximation.

By reducting from the set disjointness problem, it can be shown that computing the diameter of $n$ points in $\mathbb{R}^{d}$ requires $\Omega(n\log n)$ operations in the algebraic computation-tree model [1]. It is shown by Yao that it is possible to compute the diameter in sub-quadratic time in each dimension [2]. There are well-known solutions in two and three dimensions. In the plane, this problem can be computed in optimal time $O(n\log n)$, but in three dimensions, it is more difficult. Clarkson and Shor [3] present an $O(n\log n)$-time randomized algorithm. Their algorithm needs to compute the intersection of $n$ balls (with the same radius) in $\mathbb{R}^{3}$. It may be slower than the brute-force algorithm for the most practical datasets. Moreover, it is not an efficient method for higher dimensions because the intersection of $n$ balls with the same radius has a large size. Some deterministic algorithms with running time $O(n\log^{3}n)$ [4, 5] and $O(n\log^{2}n)$ [6, 7] are found for solving this problem. Finally, Ramos [8, 9] introduced an optimal deterministic $O(n\log n)$-time algorithm in $\mathbb{R}^{3}$. Cheong et al. [10] present an $O(n\log^{2}n)$ randomized algorithm that computes the all-pairs farthest neighbors for $n$ points on the convex position in $\mathbb{R}^{3}$.
   In the absence of fast algorithms, many attempts have been done to approximate the diameter in low and high dimensions. A 2-approximation algorithm with $O(dn)$ time in $d$ dimensions can be found easily by selecting a point $x$ of $\mathcal{S}$ and then finding the farthest point $y$ from it by brute-force manner. The first non-trivial approximation algorithm for the diameter is presented by Egecioglu and Kalantari [11] that approximates the diameter with factor $\sqrt{3}$ and operations cost $O(dn)$ in $d$ dimensions. They also present an iterative algorithm with $t\leq n$ iterations and the cost $O(dn)$ for each iteration that has approximation factor $\sqrt{5-2\sqrt{3}}$. Agarwal et al. [12] present a $(1+\varepsilon)$-approximation algorithm in $\mathbb{R}^{d}$ with $O(n/\varepsilon^{(d-1)/2})$ running time by projection to directions. Barequet and Har Peled [13] present a $\sqrt{d}$-approximate diameter method with $O(dn)$ time. They also describe a $(1+\varepsilon)$-approximate approach for computing the diameter with $O(n+1/\varepsilon^{2d})$ time in $\mathbb{R}^{d}$. They show that the running time can be improved to $O(n+1/\varepsilon^{2(d-1)})$. Similarly, Har Peled [14] presents an approach which for the most inputs is able to compute very fast the exact diameter, or an approximation. Although, in the worst case, the algorithm running time is still quadratic, and it is sensitive to the hardness of the input. His algorithm is able to return a pair of points $p$ and $q$ such that $||p-q||\geq(1-\varepsilon)D$, for each value $\varepsilon>0$ in each dimension with $O((n+1/\varepsilon^{2d})\log 1/\varepsilon)$ running time. He shows that with a complicated analysis, this running time can be reduced to $O((n+1/\varepsilon^{3(d-1)/2})\log 1/\varepsilon)$. Simultaneously, Maladain and Boissonnat [15] present an exact algorithm for the diameter by computing the double normals in each dimension, but their algorithm is not worst-case optimal. They also show that with having double normals, a $\sqrt{3}$-approximation of the diameter in each dimension is provided. Moreover, Finocchiaro and Pellegrini [16] describe an algorithm that finds in $O(dn\log n+n^{2})$ time with high probability a $(1+\varepsilon)$-approximation for the diameter of a set of $n$ points in $d$-dimensional Euclidean space. Chan [17] observes that a combination of two approaches in [12] and [13] yields a $(1+\varepsilon)$-approximation algorithm with $O(n+1/\varepsilon^{3(d-1)/2})$ time and a $(1+O(\varepsilon))$-approximation algorithm with $O(n+1/\varepsilon^{d-\frac{1}{2}})$ time. He also introduces a core-set theorem, and shows that using this theorem, a $(1+O(\varepsilon))$-approximation for the diameter in $O(n+1/\varepsilon^{d-\frac{3}{2}})$ time can be found [18]. Recently, Chan [19] has proposed an approximation algorithm with $O((n/\sqrt{\varepsilon}+1/\varepsilon^{\frac{d}{2}+1})(\log\frac{1}{\varepsilon})^{O(1)})$ time by applying the Chebyshev polynomials for the diameter in low constant dimensions, and Arya et al. [20] show that by applying an efficient decomposition of a convex body using a hierarchy of Macbeath regions, it is possible to compute an approximation for the diameter of a point set in $O(n\log\frac{1}{\varepsilon}+1/\varepsilon^{\frac{(d-1)}{2}+\alpha})$ time, where $\alpha$ is an arbitrarily small positive constant. Table 1 provides a summary on some non-constant approximation algorithm for the point-set diameter.

‎

In this section, we describe our new approximation algorithm to compute the diameter of a set $\mathcal{S}$ of $n$ points in $\mathbb{R}^{d}$. In order to follow our algorithm, we first find extreme points in each coordinate and compute the axis-parallel bounding box of $\mathcal{S}$, which is denoted by $B(\mathcal{S})$. We use the largest length side $\ell$ of $B(\mathcal{S})$ to impose grids on the point set. In fact, we first decompose $B(\mathcal{S})$ to a grid of regular hypercubes with side length $\xi$, where $\xi=\varepsilon\ell/2\sqrt{d}$. We call each hypercube a cell. Then, each point of $\mathcal{S}$ is rounded to its corresponding central cell-point. Figure 1 shows an example of the rounding process for a point set in $\mathbb{R}^{2}$.

In the following, we impose again a $\xi_{1}$-grid to $B(\mathcal{S})$ for $\xi_{1}=\sqrt{\varepsilon}\ell/2\sqrt{d}$. Then, we round each point of the rounded point set $\hat{\mathcal{S}}$ to its nearest grid-point in this new grid that results in a point set $\hat{\mathcal{S}_{1}}$. Let, $\mathcal{B}_{\delta}(p)$ be a hypercube with side length $\delta$ and central-point $p$. We restrict our search for finding diametrical pairs of the first rounded point set $\hat{\mathcal{S}}$ into two hypercubes $\mathcal{B}_{2\xi_{1}}(\hat{p}_{1})$ and $\mathcal{B}_{2\xi_{1}}(\hat{q}_{1})$ corresponding to two diametrical pair $\hat{p}_{1}$ and $\hat{q}_{1}$ in the point set $\hat{\mathcal{S}_{1}}$. Let us use two point sets $\mathcal{B}_{1}$ and $\mathcal{B}_{2}$ for maintaining points of the rounded point set $\hat{\mathcal{S}}$, which are inside two hypercubes $\mathcal{B}_{2\xi_{1}}(\hat{p}_{1})$ and $\mathcal{B}_{2\xi_{1}}(\hat{q}_{1})$, respectively. See Figure 2. Then, it is sufficient to find a diameter between points of $\hat{\mathcal{S}}$, which are inside two point sets $\mathcal{B}_{1}$ and $\mathcal{B}_{2}$. We use notation $Diam(\mathcal{B}_{1},\mathcal{B}_{2})$ for the process of computing the diameter of point set $\mathcal{B}_{1}\cup\mathcal{B}_{2}$. Altogether, we can present the following algorithm.

‎

‎
Algorithm 1‎: APPROXIMATE DIAMETER $(\mathcal{S},\varepsilon)$

‎
‎ ‎ ‎Input‎: a set $\mathcal{S}$ of $n$ points in $\mathbb{R}^{d}$ and an error parameter $\varepsilon$.

‎Output‎: approximate diameter $\tilde{D}$.

1:   Compute the axis-parallel bounding box $B(\mathcal{S})$ for a point set $\mathcal{S}$.

2:   $\ell\leftarrow$ Find the length of the largest side in $B(\mathcal{S})$.

3:   Set $\xi\leftarrow\varepsilon\ell/2\sqrt{d}$ and $\xi_{1}\leftarrow\sqrt{\varepsilon}\ell/2\sqrt{d}$.

4:   $\hat{\mathcal{S}}\leftarrow$ Round each point of $\mathcal{S}$ to its central-cell point in a $\xi$-grid.

5:   $\hat{\mathcal{S}_{1}}\leftarrow$ Round each point of $\hat{\mathcal{S}}$ to its nearest grid-point in a $\xi_{1}$-grid.

6:   $\hat{D}_{1}\leftarrow$ Compute the diameter of the point set $\hat{\mathcal{S}_{1}}$ by brute-force manner, and
        simultaneously, a list of the diametrical pair $(\hat{p}_{1},\hat{q}_{1})$, such that $\hat{D}_{1}=||\hat{p}_{1}-\hat{q}_{1}||$.

7:   Find points of $\hat{\mathcal{S}}$ which are in two hypercubes $\mathcal{B}_{1}=\mathcal{B}_{2\xi_{1}}(\hat{p}_{1})$ and $\mathcal{B}_{2}=\mathcal{B}_{2\xi_{1}}(\hat{q}_{1})$
      for each diametrical pair $(\hat{p}_{1},\hat{q}_{1})$.

8:   $\hat{D}\leftarrow$ Compute $Diam(\mathcal{B}_{1},\mathcal{B}_{2})$, corresponding to each diametrical pair $(\hat{p}_{1},\hat{q}_{1})$
       by brute-force manner and return the maximum value between them.

9:   $\tilde{D}\leftarrow\hat{D}+\varepsilon\ell/2$.

10: Output $\tilde{D}$. ‎

‎ ‎ ‎

Note that we only need $O(dn)$ time to round points to their central-cell points. We need $O(d)$ time for computing the central cell-point for each point. Thus, we can round all point of a set of $n$ points to their central-cell points in $O(dn)$ time. Similarly, rounding $n$ points to their nearest grid-point can be done in $O(dn)$ time.