从 Python 中的集合中高效地找到最近的坐标对

Efficiently finding the closest coordinate pair from a set in Python

问题

想象一下我站在机场。给定一对地理坐标,如何有效地确定我站在哪个机场?

输入

期望输出

一组机场坐标对中的一个坐标对 (a,b),表示距离点 (x,y) 最近的机场。

低效解决方案

这是我解决这个问题的低效尝试。机场集的长度显然是线性的。

shortest_distance = None
shortest_distance_coordinates = None

point = (50.776435, -0.146834)

for airport in airports:
    distance = compute_distance(point, airport)
    if distance < shortest_distance or shortest_distance is None:
        shortest_distance = distance
        shortest_distance_coordinates = airport

问题

如何改进这个解决方案?这可能涉及到根据我们当前所处位置的坐标预先过滤机场列表的某种方式,或者预先按特定顺序对它们进行排序。

从这个SO question:

import numpy as np
def closest_node(node, nodes):
    nodes = np.asarray(nodes)
    deltas = nodes - node
    dist_2 = np.einsum('ij,ij->i', deltas, deltas)
    return np.argmin(dist_2)

其中 node 是具有两个值 (x, y) 的元组,nodes 是具有两个值 ([(x_1, y_1), (x_2, y_2),])

的元组数组

如果您的坐标未排序,假设它是 (latitude,longitude),您的搜索只能略微改进,方法是像地球一样先过滤纬度

1 degree of latitude on the sphere is 111.2 km or 69 miles

但这不会带来很大的加速。

如果您首先按纬度对机场进行排序,那么您可以使用二进制搜索来查找第一个 可以 匹配 (airport_lat >= point_lat-tolerance) 的机场,然后只进行比较最后一个 可以 匹配 (airport_lat <= point_lat+tolerance) - 但注意 0 度等于 360。虽然你不能直接使用那个库,bisect 的来源是实施二进制搜索的良好开端。

虽然从技术上讲,这种搜索方式仍然是 O(n),但实际距离计算(取决于公差)和纬度比较要少得多。所以你会有一个巨大的加速。

使用k维树:

>>> from scipy import spatial
>>> airports = [(10,10),(20,20),(30,30),(40,40)]
>>> tree = spatial.KDTree(airports)
>>> tree.query([(21,21)])
(array([ 1.41421356]), array([1]))

其中1.41421356是查询点到最近邻的距离,1是近邻的索引

参见:http://docs.scipy.org/doc/scipy/reference/generated/scipy.spatial.KDTree.query.html#scipy.spatial.KDTree.query

@Juddling的回答很好,但是KDTree不支持haversine distance,比较适合latitude/longitude坐标。 对于 haversine 距离,您可以使用 BallTree。请注意,您需要先将坐标转换为弧度。

from math import radians
from sklearn.neighbors import BallTree
import numpy as np

airports = [(10,10),(20,20),(30,30),(40,40)]
airports_rad = np.array([[radians(x[0]), radians(x[1])] for x in airports ])
tree = BallTree(airports_rad , metric = 'haversine')
result = tree.query([(radians(21),radians(21))])
print(result)

给予

(array([[0.02391369]]), array([[1]], dtype=int64))

要将距离转换为米,您需要乘以地球半径(以米为单位)。

earth_radius = 6371000 # meters in earth
print(result[0][0] * earth_radius)
[152354.11114795]