Определение точек с наименьшим евклидовым расстоянием

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). Вы можете быть уверены, что две ближайшие точки будут ребром симплекса в триангуляции, который является гораздо меньшим подмножеством пар, чем любая комбинация.

Вот код (обновлен для общего N-D):

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

Принятый ответ подходит для небольших наборов данных, но время его выполнения масштабируется как 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.Однако обратите внимание, что cKDTree, в отличие от BallTree, плохо масштабируется для больших размеров

from scipy.spatial import cKDTree as tree

Answer 7

Насколько это быстро по сравнению с простым вложенным циклом и отслеживанием самой короткой пары?Я думаю, что создание огромного перекрестного массива - это то, что может повредить вам.Даже O (n ^ 2) все еще довольно быстр, если вы делаете только двухмерные точки.