加速处理500万行坐标数据
speeding up processing 5 million rows of coordinate data
我有一个包含两列（纬度、经度）的 csv 文件，其中包含超过 500 万行的地理位置数据。
我需要确定不在列表中任何其他点 5 英里范围内的点，并将所有内容输出回另一个具有额外列 (CloseToAnotherPoint) 的 CSV，如果有另一个列，则为 True点在 5 英里以内，如果没有，则 False。
这是我目前使用 geopy 的解决方案（不进行任何网络调用，仅使用函数计算距离）：
from geopy.point import Point
from geopy.distance import vincenty
import csv
class CustomGeoPoint(object):
def __init__(self, latitude, longitude):
self.location = Point(latitude, longitude)
self.close_to_another_point = False
try:
output = open('output.csv','w')
writer = csv.writer(output, delimiter = ',', quoting=csv.QUOTE_ALL)
writer.writerow(['Latitude', 'Longitude', 'CloseToAnotherPoint'])
# 5 miles
close_limit = 5
geo_points = []
with open('geo_input.csv', newline='') as geo_csv:
reader = csv.reader(geo_csv)
next(reader, None) # skip the headers
for row in reader:
geo_points.append(CustomGeoPoint(row[0], row[1]))
# for every point, look at every point until one is found within 5 miles
for geo_point in geo_points:
for geo_point2 in geo_points:
dist = vincenty(geo_point.location, geo_point2.location).miles
if 0 < dist <= close_limit: # (0,close_limit]
geo_point.close_to_another_point = True
break
writer.writerow([geo_point.location.latitude, geo_point.location.longitude,
geo_point.close_to_another_point])
finally:
output.close()
正如您可能从中看到的那样，此解决方案非常慢。实际上太慢了，我让它 运行 3 天 但它仍然没有完成！
我考虑过尝试将数据分成块（多个 CSV 文件或其他），这样内部循环就不必查看其他所有点，但我必须弄清楚如何确保每个部分的边界与其相邻部分的边界相匹配，这看起来过于复杂，恐怕它会比它的价值更让人头疼。
关于如何加快速度的任何指示？
如果您按纬度 (n log(n)) 对列表进行排序，并且点分布大致均匀，则每个点将减少到 5 英里内约 1000 个点（餐巾纸数学，不精确） .通过仅查看纬度附近的点，运行时间从 n^2 变为 n*log(n)+.0004n^2。希望这足以加快速度。
我会分三步重做算法：
使用大圆距离，假设误差为 1%，因此使 limit 等于 1.01*limit。
将大圆距离编码为内联函数，这个测试应该很快
你会得到一些误报，你可以用 vincenty 进一步测试
我会 pandas 试一试。 Pandas 专为高效处理大量数据而设计。无论如何，这可能有助于提高 csv 部分的效率。但从它的声音来看，你自己解决了一个本质上效率低下的问题。您将第 1 点与其他 4,999,999 个点进行比较。然后取点 2 并将其与其他 4,999,998 个点进行比较，依此类推。算一算。那是您正在进行的 12.5 万亿 次比较。如果每秒可以进行 1,000,000 次比较，则计算需要 144 天。如果每秒可以进行 10,000,000 次比较，那就是 14 天。对于直接 python 中的加法，10,000,000 次操作大约需要 1.1 秒，但我怀疑您的比较是否与加法操作一样快。所以至少给它一两个星期。
或者，您可以想出另一种算法，但我没有想到任何特定的算法。
正如 的答案所指出的，这是 k-D 树的经典用例。
另一种方法是使用扫描线算法，如 How do I match similar coordinates using Python?
的答案所示
这是针对您的问题改编的扫描线算法。在我的笔记本电脑上，运行 通过 5M 个随机点需要 < 5 分钟。
import itertools as it
import operator as op
import sortedcontainers # handy library on Pypi
import time
from collections import namedtuple
from math import cos, degrees, pi, radians, sqrt
from random import sample, uniform
Point = namedtuple("Point", "lat long has_close_neighbor")
miles_per_degree = 69
number_of_points = 5000000
data = [Point(uniform( -88.0, 88.0), # lat
uniform(-180.0, 180.0), # long
True
)
for _ in range(number_of_points)
]
start = time.time()
# Note: lat is first in Point, so data is sorted by .lat then .long.
data.sort()
print(time.time() - start)
# Parameter that determines the size of a sliding lattitude window
# and therefore how close two points need to be to be to get flagged.
threshold = 5.0 # miles
lat_span = threshold / miles_per_degree
coarse_threshold = (.98 * threshold)**2
# Sliding lattitude window. Within the window, observations are
# ordered by longitude.
window = sortedcontainers.SortedListWithKey(key=op.attrgetter('long'))
# lag_pt is the 'southernmost' point within the sliding window.
point = iter(data)
lag_pt = next(point)
milepost = len(data)//10
# lead_pt is the 'northernmost' point in the sliding window.
for i, lead_pt in enumerate(data):
if i == milepost:
print('.', end=' ')
milepost += len(data)//10
# Dec of lead_obs represents the leading edge of window.
window.add(lead_pt)
# Remove observations further than the trailing edge of window.
while lead_pt.lat - lag_pt.lat > lat_span:
window.discard(lag_pt)
lag_pt = next(point)
# Calculate 'east-west' width of window_size at dec of lead_obs
long_span = lat_span / cos(radians(lead_pt.lat))
east_long = lead_pt.long + long_span
west_long = lead_pt.long - long_span
# Check all observations in the sliding window within
# long_span of lead_pt.
for other_pt in window.irange_key(west_long, east_long):
if other_pt != lead_pt:
# lead_pt is at the top center of a box 2 * long_span wide by
# 1 * long_span tall. other_pt is is in that box. If desired,
# put additional fine-grained 'closeness' tests here.
# coarse check if any pts within 80% of threshold distance
# then don't need to check distance to any more neighbors
average_lat = (other_pt.lat + lead_pt.lat) / 2
delta_lat = other_pt.lat - lead_pt.lat
delta_long = (other_pt.long - lead_pt.long)/cos(radians(average_lat))
if delta_lat**2 + delta_long**2 <= coarse_threshold:
break
# put vincenty test here
#if 0 < vincenty(lead_pt, other_pt).miles <= close_limit:
# break
else:
data[i] = data[i]._replace(has_close_neighbor=False)
print()
print(time.time() - start)
这只是第一遍，但到目前为止我已经通过使用 great_circle() instead of vincinty(), and cleaning up a couple of other things. The difference is explained here 将速度提高了一半，准确度损失约为 0.17%:
from geopy.point import Point
from geopy.distance import great_circle
import csv
class CustomGeoPoint(Point):
def __init__(self, latitude, longitude):
super(CustomGeoPoint, self).__init__(latitude, longitude)
self.close_to_another_point = False
def isCloseToAnother(pointA, points):
for pointB in points:
dist = great_circle(pointA, pointB).miles
if 0 < dist <= CLOSE_LIMIT: # (0, close_limit]
return True
return False
with open('geo_input.csv', 'r') as geo_csv:
reader = csv.reader(geo_csv)
next(reader, None) # skip the headers
geo_points = sorted(map(lambda x: CustomGeoPoint(x[0], x[1]), reader))
with open('output.csv', 'w') as output:
writer = csv.writer(output, delimiter=',', quoting=csv.QUOTE_ALL)
writer.writerow(['Latitude', 'Longitude', 'CloseToAnotherPoint'])
# for every point, look at every point until one is found within a mile
for point in geo_points:
point.close_to_another_point = isCloseToAnother(point, geo_points)
writer.writerow([point.latitude, point.longitude,
point.close_to_another_point])
我将进一步改进。
之前：
$ time python geo.py
real 0m5.765s
user 0m5.675s
sys 0m0.048s
之后：
$ time python geo.py
real 0m2.816s
user 0m2.716s
sys 0m0.041s
Oscar Smith 生成的更好的解决方案。你有一个 csv 文件，只是在 excel 中对它进行了排序，这是非常有效的）。然后在你的程序中使用二进制搜索来查找 5 英里以内的城市（你可以对二进制搜索方法进行小的改动，如果找到一个满足你条件的城市，它就会中断）。
另一项改进是设置一张地图，当您发现一个城市位于另一个城市内时，可以记住这对城市。例如，当您发现城市 A 在城市 B 的 5 英里范围内时，使用 Map 存储对（B 是键，A 是值）。所以下次遇到B，先去Map里找，如果有对应的值，就不用再查了。但它可能会使用更多的内存所以关心它。希望对你有帮助。
让我们看看你在做什么。
您将所有点读入名为 geo_points 的列表中。
现在，你能告诉我列表是否排序了吗？因为如果它已排序，我们 肯定 想知道。排序是有价值的信息，尤其是当您处理 500 万件东西时。
你遍历所有 geo_points。按你说是500万
在外部循环中，您再次 循环了全部 500 万 geo_points。
您计算两个循环项之间的距离（以英里为单位）。
如果距离小于阈值，则在第一个点上记录该信息，并停止内循环。
当内循环停止时，将有关外循环项的信息写入 CSV 文件。
注意几件事。首先，您在外循环中循环了 500 万次。然后你在内部循环中循环了 500 万次。
这就是 O(n²) 的意思。
下次您看到有人谈论 "Oh, this is O(log n) but that other thing is O(n log n)," 时请记住这次经历 - 您是 运行 一个 n² 算法，其中 n 在本例中为 5,000,000。糟透了，笨蛋？
总之，你遇到了一些问题。
问题 1：您最终会将每个点与自身进行比较。它的距离应该为零，这意味着它们都将被标记为在任何距离阈值内。如果您的程序完成，所有单元格都将标记为 True。
问题 2：当您将点 #1 与点 #12345 进行比较，并且它们彼此之间的距离在阈值范围内时，您正在记录关于点 #1 的信息。但是你没有记录关于另一点的相同信息。您 知道 #12345 (geo_point2) 点自反地位于点 #1 的阈值内，但您没有将其写下来。所以你错过了跳过超过 500 万次比较的机会。
问题 3：如果您比较点 #1 和点 #2，并且它们不在阈值距离内，那么当您比较点 #2 和点 #1 时会发生什么？您的内部循环每次都从列表的开头开始，但是您知道您已经比较了列表的开头和结尾。您可以将问题 space 减少一半，只需让您的外循环变为 i in range(0, 5million) 并且您的内循环变为 j in range(i+1, 5million).
答案?
考虑平面上的纬度和经度。您想知道 5 英里内是否有一个点。让我们考虑一个 10 平方英里的区域，以您的第 1 点为中心。这是一个以 (X1, Y1) 为中心的正方形，左上角位于 (X1 - 5miles, Y1 + 5miles)，右下角位于 (X1 + 5miles, Y1 - 5miles)。现在，如果一个点在那个方块内，它 可能不会 在你的点 #1 的 5 英里范围内。但是你可以打赌，如果它在那个广场之外，它就在5英里之外。
正如@SeverinPappadeaux 指出的那样，像地球这样的椭球体上的距离与平面上的距离并不完全相同。但那又怎样？将您的正方形设置得大一点以允许差异，然后继续！
排序列表
这就是排序很重要的原因。如果所有点都按 X 排序，然后是 Y（或 Y，然后是 X - 随便什么）你知道， 你真的可以加快速度。因为当 X（或 Y）坐标变得太大时，您可以简单地停止扫描，而不必经过 500 万个点。
那将如何运作？和以前一样，除了你的内部循环会有一些这样的检查：
five_miles = ... # Whatever math, plus an error allowance!
list_len = len(geo_points) # Don't call this 5 million times
for i, pi in enumerate(geo_points):
if pi.close_to_another_point:
continue # Remember if close to an earlier point
pi0max = pi[0] + five_miles
pi1min = pi[1] - five_miles
pi1max = pi[1] + five_miles
for j in range(i+1, list_len):
pj = geo_points[j]
# Assumes geo_points is sorted on [0] then [1]
if pj[0] > pi0max:
# Can't possibly be close enough, nor any later points
break
if pj[1] < pi1min or pj[1] > pi1max:
# Can't be close enough, but a later point might be
continue
# Now do "real" comparison using accurate functions.
if ...:
pi.close_to_another_point = True
pj.close_to_another_point = True
break
我在那里做什么？首先，我将一些数字放入局部变量中。然后我使用 enumerate 给我一个 i 值 和 对外部点的引用。 （你所谓的geo_point）。然后，我会快速查看我们是否已经知道这一点与另一点接近。
如果没有，我们将不得不扫描。所以我只扫描列表中的 "later" 点，因为我知道外循环扫描早期的点，我绝对不想将一个点与自身进行比较。我正在使用一些临时变量来缓存涉及外循环的计算结果。在内部循环中，我对临时对象做了一些愚蠢的比较。他们无法告诉我这两个点是否 彼此接近，但我可以检查它们是否确实 不 接近并向前跳.
最后，如果简单的检查通过，则继续进行昂贵的检查。如果检查确实通过了，请务必记录 both 点的结果，这样我们以后可以跳过第二点。
未排序列表
但如果列表未排序怎么办？
@RootTwo 将您指向一棵 kD 树（其中 D 代表 "dimensional"，在这种情况下 k 是“2”）。这个想法非常简单，如果您已经了解二叉搜索树：您循环遍历维度，比较树中偶数层的 X 和比较奇数层的 Y（反之亦然）。这个想法是这样的：
def insert_node(node, treenode, depth=0):
dimension = depth % 2 # even/odd -> lat/long
dn = node.coord[dimension]
dt = treenode.coord[dimension]
if dn < dt:
# go left
if treenode.left is None:
treenode.left = node
else:
insert_node(node, treenode.left, depth+1)
else:
# go right
if treenode.right is None:
treenode.right = node
else:
insert_node(node, treenode.right, depth+1)
这有什么用？这将为您提供一个可搜索的树，其中可以在 O(log n) 时间内插入点。这意味着整个列表的 O(n log n)，这比 n 平方好得多！ （500万的log底2基本就是23，所以n log n是500万乘以23，对比500万乘以500万！）
这也意味着您可以进行有针对性的搜索。由于树是有序的，因此查找 "close" 点非常简单（@RootTwo 的维基百科 link 提供了一种算法）。
建议
我的建议是，如果需要，只编写代码对列表进行排序。写起来更容易，手工检查也更容易，而且它是一个单独的通行证，您只需要制作一次。
对列表进行排序后，请尝试我上面展示的方法。它与您正在做的很接近，您应该很容易理解和编写代码。
这个问题可以用VP tree来解决。这些允许查询数据
距离是服从三角不等式的度量。
VP 树相对于 k-D 树的最大优势在于它们可以盲目地
无需担心应用于世界任何地方的地理数据
关于将其投影到合适的 2D space。另外一个真正的测地线
可以使用距离（无需担心之间的差异
测地线距离和投影中的距离）。
这是我的测试：在
世界。将这些放入 VP 树中。
遍历所有的点，查询VP树找到任意一个邻居a
距离在（0km，10km] 之外。（0km 不包括在这个集合中以避免
找到查询点。）计算没有这样的点的数量
邻居（在我的例子中是 229573）。
建立 VP 树的成本 = 5000000 * 20 次距离计算。
查询成本 = 5000000 * 23 次距离计算。
设置和查询时间为 5 分钟 7 秒。
我使用 C++ 和 GeographicLib 来计算距离，但是
该算法当然可以用任何语言实现，这里是
python version of GeographicLib.
ADDENDUM：给出了实现此方法的 C++ 代码 here。
我有一个包含两列（纬度、经度）的 csv 文件，其中包含超过 500 万行的地理位置数据。
我需要确定不在列表中任何其他点 5 英里范围内的点，并将所有内容输出回另一个具有额外列 (CloseToAnotherPoint) 的 CSV，如果有另一个列，则为 True点在 5 英里以内，如果没有，则 False。
这是我目前使用 geopy 的解决方案（不进行任何网络调用，仅使用函数计算距离）：
from geopy.point import Point
from geopy.distance import vincenty
import csv
class CustomGeoPoint(object):
def __init__(self, latitude, longitude):
self.location = Point(latitude, longitude)
self.close_to_another_point = False
try:
output = open('output.csv','w')
writer = csv.writer(output, delimiter = ',', quoting=csv.QUOTE_ALL)
writer.writerow(['Latitude', 'Longitude', 'CloseToAnotherPoint'])
# 5 miles
close_limit = 5
geo_points = []
with open('geo_input.csv', newline='') as geo_csv:
reader = csv.reader(geo_csv)
next(reader, None) # skip the headers
for row in reader:
geo_points.append(CustomGeoPoint(row[0], row[1]))
# for every point, look at every point until one is found within 5 miles
for geo_point in geo_points:
for geo_point2 in geo_points:
dist = vincenty(geo_point.location, geo_point2.location).miles
if 0 < dist <= close_limit: # (0,close_limit]
geo_point.close_to_another_point = True
break
writer.writerow([geo_point.location.latitude, geo_point.location.longitude,
geo_point.close_to_another_point])
finally:
output.close()
正如您可能从中看到的那样，此解决方案非常慢。实际上太慢了，我让它 运行 3 天 但它仍然没有完成！
我考虑过尝试将数据分成块（多个 CSV 文件或其他），这样内部循环就不必查看其他所有点，但我必须弄清楚如何确保每个部分的边界与其相邻部分的边界相匹配，这看起来过于复杂，恐怕它会比它的价值更让人头疼。
关于如何加快速度的任何指示？
如果您按纬度 (n log(n)) 对列表进行排序，并且点分布大致均匀，则每个点将减少到 5 英里内约 1000 个点（餐巾纸数学，不精确） .通过仅查看纬度附近的点，运行时间从 n^2 变为 n*log(n)+.0004n^2。希望这足以加快速度。
我会分三步重做算法：
使用大圆距离，假设误差为 1%，因此使 limit 等于 1.01*limit。
将大圆距离编码为内联函数，这个测试应该很快
你会得到一些误报，你可以用 vincenty 进一步测试
我会 pandas 试一试。 Pandas 专为高效处理大量数据而设计。无论如何，这可能有助于提高 csv 部分的效率。但从它的声音来看，你自己解决了一个本质上效率低下的问题。您将第 1 点与其他 4,999,999 个点进行比较。然后取点 2 并将其与其他 4,999,998 个点进行比较，依此类推。算一算。那是您正在进行的 12.5 万亿 次比较。如果每秒可以进行 1,000,000 次比较，则计算需要 144 天。如果每秒可以进行 10,000,000 次比较，那就是 14 天。对于直接 python 中的加法，10,000,000 次操作大约需要 1.1 秒，但我怀疑您的比较是否与加法操作一样快。所以至少给它一两个星期。
或者，您可以想出另一种算法，但我没有想到任何特定的算法。
正如
另一种方法是使用扫描线算法，如 How do I match similar coordinates using Python?
的答案所示这是针对您的问题改编的扫描线算法。在我的笔记本电脑上，运行 通过 5M 个随机点需要 < 5 分钟。
import itertools as it
import operator as op
import sortedcontainers # handy library on Pypi
import time
from collections import namedtuple
from math import cos, degrees, pi, radians, sqrt
from random import sample, uniform
Point = namedtuple("Point", "lat long has_close_neighbor")
miles_per_degree = 69
number_of_points = 5000000
data = [Point(uniform( -88.0, 88.0), # lat
uniform(-180.0, 180.0), # long
True
)
for _ in range(number_of_points)
]
start = time.time()
# Note: lat is first in Point, so data is sorted by .lat then .long.
data.sort()
print(time.time() - start)
# Parameter that determines the size of a sliding lattitude window
# and therefore how close two points need to be to be to get flagged.
threshold = 5.0 # miles
lat_span = threshold / miles_per_degree
coarse_threshold = (.98 * threshold)**2
# Sliding lattitude window. Within the window, observations are
# ordered by longitude.
window = sortedcontainers.SortedListWithKey(key=op.attrgetter('long'))
# lag_pt is the 'southernmost' point within the sliding window.
point = iter(data)
lag_pt = next(point)
milepost = len(data)//10
# lead_pt is the 'northernmost' point in the sliding window.
for i, lead_pt in enumerate(data):
if i == milepost:
print('.', end=' ')
milepost += len(data)//10
# Dec of lead_obs represents the leading edge of window.
window.add(lead_pt)
# Remove observations further than the trailing edge of window.
while lead_pt.lat - lag_pt.lat > lat_span:
window.discard(lag_pt)
lag_pt = next(point)
# Calculate 'east-west' width of window_size at dec of lead_obs
long_span = lat_span / cos(radians(lead_pt.lat))
east_long = lead_pt.long + long_span
west_long = lead_pt.long - long_span
# Check all observations in the sliding window within
# long_span of lead_pt.
for other_pt in window.irange_key(west_long, east_long):
if other_pt != lead_pt:
# lead_pt is at the top center of a box 2 * long_span wide by
# 1 * long_span tall. other_pt is is in that box. If desired,
# put additional fine-grained 'closeness' tests here.
# coarse check if any pts within 80% of threshold distance
# then don't need to check distance to any more neighbors
average_lat = (other_pt.lat + lead_pt.lat) / 2
delta_lat = other_pt.lat - lead_pt.lat
delta_long = (other_pt.long - lead_pt.long)/cos(radians(average_lat))
if delta_lat**2 + delta_long**2 <= coarse_threshold:
break
# put vincenty test here
#if 0 < vincenty(lead_pt, other_pt).miles <= close_limit:
# break
else:
data[i] = data[i]._replace(has_close_neighbor=False)
print()
print(time.time() - start)
这只是第一遍，但到目前为止我已经通过使用 great_circle() instead of vincinty(), and cleaning up a couple of other things. The difference is explained here 将速度提高了一半，准确度损失约为 0.17%:
from geopy.point import Point
from geopy.distance import great_circle
import csv
class CustomGeoPoint(Point):
def __init__(self, latitude, longitude):
super(CustomGeoPoint, self).__init__(latitude, longitude)
self.close_to_another_point = False
def isCloseToAnother(pointA, points):
for pointB in points:
dist = great_circle(pointA, pointB).miles
if 0 < dist <= CLOSE_LIMIT: # (0, close_limit]
return True
return False
with open('geo_input.csv', 'r') as geo_csv:
reader = csv.reader(geo_csv)
next(reader, None) # skip the headers
geo_points = sorted(map(lambda x: CustomGeoPoint(x[0], x[1]), reader))
with open('output.csv', 'w') as output:
writer = csv.writer(output, delimiter=',', quoting=csv.QUOTE_ALL)
writer.writerow(['Latitude', 'Longitude', 'CloseToAnotherPoint'])
# for every point, look at every point until one is found within a mile
for point in geo_points:
point.close_to_another_point = isCloseToAnother(point, geo_points)
writer.writerow([point.latitude, point.longitude,
point.close_to_another_point])
我将进一步改进。
之前：
$ time python geo.py
real 0m5.765s
user 0m5.675s
sys 0m0.048s
之后：
$ time python geo.py
real 0m2.816s
user 0m2.716s
sys 0m0.041s
Oscar Smith 生成的更好的解决方案。你有一个 csv 文件，只是在 excel 中对它进行了排序，这是非常有效的）。然后在你的程序中使用二进制搜索来查找 5 英里以内的城市（你可以对二进制搜索方法进行小的改动，如果找到一个满足你条件的城市，它就会中断）。 另一项改进是设置一张地图，当您发现一个城市位于另一个城市内时，可以记住这对城市。例如，当您发现城市 A 在城市 B 的 5 英里范围内时，使用 Map 存储对（B 是键，A 是值）。所以下次遇到B，先去Map里找，如果有对应的值，就不用再查了。但它可能会使用更多的内存所以关心它。希望对你有帮助。
让我们看看你在做什么。
您将所有点读入名为
geo_points的列表中。现在，你能告诉我列表是否排序了吗？因为如果它已排序，我们 肯定 想知道。排序是有价值的信息，尤其是当您处理 500 万件东西时。
你遍历所有
geo_points。按你说是500万在外部循环中，您再次 循环了全部 500 万
geo_points。您计算两个循环项之间的距离（以英里为单位）。
如果距离小于阈值，则在第一个点上记录该信息，并停止内循环。
当内循环停止时，将有关外循环项的信息写入 CSV 文件。
注意几件事。首先，您在外循环中循环了 500 万次。然后你在内部循环中循环了 500 万次。
这就是 O(n²) 的意思。
下次您看到有人谈论 "Oh, this is O(log n) but that other thing is O(n log n)," 时请记住这次经历 - 您是 运行 一个 n² 算法，其中 n 在本例中为 5,000,000。糟透了，笨蛋？
总之，你遇到了一些问题。
问题 1：您最终会将每个点与自身进行比较。它的距离应该为零，这意味着它们都将被标记为在任何距离阈值内。如果您的程序完成，所有单元格都将标记为 True。
问题 2：当您将点 #1 与点 #12345 进行比较，并且它们彼此之间的距离在阈值范围内时，您正在记录关于点 #1 的信息。但是你没有记录关于另一点的相同信息。您 知道 #12345 (geo_point2) 点自反地位于点 #1 的阈值内，但您没有将其写下来。所以你错过了跳过超过 500 万次比较的机会。
问题 3：如果您比较点 #1 和点 #2，并且它们不在阈值距离内，那么当您比较点 #2 和点 #1 时会发生什么？您的内部循环每次都从列表的开头开始，但是您知道您已经比较了列表的开头和结尾。您可以将问题 space 减少一半，只需让您的外循环变为 i in range(0, 5million) 并且您的内循环变为 j in range(i+1, 5million).
答案?
考虑平面上的纬度和经度。您想知道 5 英里内是否有一个点。让我们考虑一个 10 平方英里的区域，以您的第 1 点为中心。这是一个以 (X1, Y1) 为中心的正方形，左上角位于 (X1 - 5miles, Y1 + 5miles)，右下角位于 (X1 + 5miles, Y1 - 5miles)。现在，如果一个点在那个方块内，它 可能不会 在你的点 #1 的 5 英里范围内。但是你可以打赌，如果它在那个广场之外，它就在5英里之外。
正如@SeverinPappadeaux 指出的那样，像地球这样的椭球体上的距离与平面上的距离并不完全相同。但那又怎样？将您的正方形设置得大一点以允许差异，然后继续！
排序列表
这就是排序很重要的原因。如果所有点都按 X 排序，然后是 Y（或 Y，然后是 X - 随便什么）你知道， 你真的可以加快速度。因为当 X（或 Y）坐标变得太大时，您可以简单地停止扫描，而不必经过 500 万个点。
那将如何运作？和以前一样，除了你的内部循环会有一些这样的检查：
five_miles = ... # Whatever math, plus an error allowance!
list_len = len(geo_points) # Don't call this 5 million times
for i, pi in enumerate(geo_points):
if pi.close_to_another_point:
continue # Remember if close to an earlier point
pi0max = pi[0] + five_miles
pi1min = pi[1] - five_miles
pi1max = pi[1] + five_miles
for j in range(i+1, list_len):
pj = geo_points[j]
# Assumes geo_points is sorted on [0] then [1]
if pj[0] > pi0max:
# Can't possibly be close enough, nor any later points
break
if pj[1] < pi1min or pj[1] > pi1max:
# Can't be close enough, but a later point might be
continue
# Now do "real" comparison using accurate functions.
if ...:
pi.close_to_another_point = True
pj.close_to_another_point = True
break
我在那里做什么？首先，我将一些数字放入局部变量中。然后我使用 enumerate 给我一个 i 值 和 对外部点的引用。 （你所谓的geo_point）。然后，我会快速查看我们是否已经知道这一点与另一点接近。
如果没有，我们将不得不扫描。所以我只扫描列表中的 "later" 点，因为我知道外循环扫描早期的点，我绝对不想将一个点与自身进行比较。我正在使用一些临时变量来缓存涉及外循环的计算结果。在内部循环中，我对临时对象做了一些愚蠢的比较。他们无法告诉我这两个点是否 彼此接近，但我可以检查它们是否确实 不 接近并向前跳.
最后，如果简单的检查通过，则继续进行昂贵的检查。如果检查确实通过了，请务必记录 both 点的结果，这样我们以后可以跳过第二点。
未排序列表
但如果列表未排序怎么办？
@RootTwo 将您指向一棵 kD 树（其中 D 代表 "dimensional"，在这种情况下 k 是“2”）。这个想法非常简单，如果您已经了解二叉搜索树：您循环遍历维度，比较树中偶数层的 X 和比较奇数层的 Y（反之亦然）。这个想法是这样的：
def insert_node(node, treenode, depth=0):
dimension = depth % 2 # even/odd -> lat/long
dn = node.coord[dimension]
dt = treenode.coord[dimension]
if dn < dt:
# go left
if treenode.left is None:
treenode.left = node
else:
insert_node(node, treenode.left, depth+1)
else:
# go right
if treenode.right is None:
treenode.right = node
else:
insert_node(node, treenode.right, depth+1)
这有什么用？这将为您提供一个可搜索的树，其中可以在 O(log n) 时间内插入点。这意味着整个列表的 O(n log n)，这比 n 平方好得多！ （500万的log底2基本就是23，所以n log n是500万乘以23，对比500万乘以500万！）
这也意味着您可以进行有针对性的搜索。由于树是有序的，因此查找 "close" 点非常简单（@RootTwo 的维基百科 link 提供了一种算法）。
建议
我的建议是，如果需要，只编写代码对列表进行排序。写起来更容易，手工检查也更容易，而且它是一个单独的通行证，您只需要制作一次。
对列表进行排序后，请尝试我上面展示的方法。它与您正在做的很接近，您应该很容易理解和编写代码。
这个问题可以用VP tree来解决。这些允许查询数据 距离是服从三角不等式的度量。
VP 树相对于 k-D 树的最大优势在于它们可以盲目地 无需担心应用于世界任何地方的地理数据 关于将其投影到合适的 2D space。另外一个真正的测地线 可以使用距离（无需担心之间的差异 测地线距离和投影中的距离）。
这是我的测试：在 世界。将这些放入 VP 树中。
遍历所有的点，查询VP树找到任意一个邻居a 距离在（0km，10km] 之外。（0km 不包括在这个集合中以避免 找到查询点。）计算没有这样的点的数量 邻居（在我的例子中是 229573）。
建立 VP 树的成本 = 5000000 * 20 次距离计算。
查询成本 = 5000000 * 23 次距离计算。
设置和查询时间为 5 分钟 7 秒。
我使用 C++ 和 GeographicLib 来计算距离，但是 该算法当然可以用任何语言实现，这里是 python version of GeographicLib.
ADDENDUM：给出了实现此方法的 C++ 代码 here。