识别具有最小欧氏距离的点

Question

我有一个 n 维点的集合，我想找到哪 2 个最接近。我能想到的最好的二维是：

from numpy import *
myArr = array( [[1, 2],
                [3, 4],
                [5, 6],
                [7, 8]] )

n = myArr.shape[0]
cross = [[sum( ( myArr[i] - myArr[j] ) ** 2 ), i, j]
         for i in xrange( n )
         for j in xrange( n )
         if i != j
         ]

print min( cross )

这给出了

[8, 0, 1]

但这对于大型数组来说太慢了。我可以对其应用什么样的优化？

相关：

---

两个不同 Numpy 数组中的点之间的欧几里得距离，不在范围内

Answer 1

尝试 scipy.spatial.distance.pdist(myArr) 。这将为您提供一个压缩距离矩阵。您可以对其使用 argmin 并找到最小值的索引。这可以被转换成配对信息。

Answer 2

关于这个问题有一个完整的维基百科页面，请参阅：http://en.wikipedia.org/wiki/Closest_pair_of_points

摘要：您可以使用递归分治算法实现 O(n log n)（在上面的 Wiki 页面上概述）。

Answer 3

您可以利用最新版本的 SciPy (v0.9) Delaunay 三角测量工具。您可以确定最近的两个点将是三角剖分中单纯形的边，这是比执行每个组合要小得多的对的子集。

这是代码（针对一般 ND 进行了更新）：

import numpy
from scipy import spatial

def closest_pts(pts):
    # set up the triangluataion
    # let Delaunay do the heavy lifting
    mesh = spatial.Delaunay(pts)

    # TODO: eliminate reduncant edges (numpy.unique?)
    edges = numpy.vstack((mesh.vertices[:,:dim], mesh.vertices[:,-dim:]))

    # the rest is easy
    x = mesh.points[edges[:,0]]
    y = mesh.points[edges[:,1]]

    dists = numpy.sum((x-y)**2, 1)
    idx = numpy.argmin(dists)

    return edges[idx]
    #print 'distance: ', dists[idx]
    #print 'coords:\n', pts[closest_verts]

dim = 3
N = 1000*dim
pts = numpy.random.random(N).reshape(N/dim, dim)

看起来很接近 O(n):

Answer 4

有一个 scipy 函数 pdist 可以以相当有效的方式获取数组中点之间的成对距离：

http://docs.scipy.org/doc/scipy/reference/spatial.distance.html

输出 N*(N-1)/ 2 个唯一对（因为 r_ij == r_ji）。然后，您可以搜索最小值并避免代码中的整个循环混乱。

Answer 5

也许您可以沿着这些思路继续进行：

In []: from scipy.spatial.distance import pdist as pd, squareform as sf
In []: m= 1234
In []: n= 123
In []: p= randn(m, n)
In []: d= sf(pd(p))
In []: a= arange(m)
In []: d[a, a]= d.max()
In []: where(d< d.min()+ 1e-9)
Out[]: (array([701, 730]), array([730, 701]))

有了更多的点，您需要能够以某种方式利用集群的层次结构。

Answer 6

与仅执行嵌套循环并跟踪最短对相比，它有多快？我认为创建一个巨大的交叉数组可能会伤害你。如果你只处理二维点，即使 O(n^2) 仍然相当快。

Answer 7

对于小型数据集，接受的答案是可以的，但其执行时间为 n**2。然而，正如 @payne 所指出的，最佳解决方案可以实现 n*log(n) 计算时间缩放。

这个最佳解决方案可以使用 sklearn.neighbors.BallTree如下。

import matplotlib.pyplot as plt
import numpy as np
from sklearn.neighbors import BallTree as tree

n = 10
dim = 2
xy = np.random.uniform(size=[n, dim])

# This solution is optimal when xy is very large
res = tree(xy)
dist, ids = res.query(xy, 2)
mindist = dist[:, 1]  # second nearest neighbour
minid = np.argmin(mindist)

plt.plot(*xy.T, 'o')
plt.plot(*xy[ids[minid]].T, '-o')

对于非常大的 xy 值集，甚至对于大尺寸 dim（尽管示例说明了 dim=2 的情况），此过程可以很好地扩展。结果输出如下所示

使用 scipy.spatial.cKDTree，通过将 sklearn 导入替换为以下 Scipy 导入。但请注意，与 BallTree 不同，cKDTree 不能很好地适应高维度

from scipy.spatial import cKDTree as tree