欧式距离的高效精确计算
Efficient and precise calculation of the euclidean distance
经过一些在线研究 (1, 2, numpy, scipy, scikit, math)，我发现了几种计算 Python 中的 欧氏距离的方法：
# 1
numpy.linalg.norm(a-b)
# 2
distance.euclidean(vector1, vector2)
# 3
sklearn.metrics.pairwise.euclidean_distances
# 4
sqrt((xa-xb)^2 + (ya-yb)^2 + (za-zb)^2)
# 5
dist = [(a - b)**2 for a, b in zip(vector1, vector2)]
dist = math.sqrt(sum(dist))
# 6
math.hypot(x, y)
我想知道是否有人可以提供关于以上哪一项（或我没有找到的任何其他项）被认为是最好的 效率和精度。如果有人知道讨论该主题的任何资源，那也很好。
我感兴趣的 context 是计算数字元组对之间的欧氏距离，例如(52, 106, 35, 12) 和 (33, 153, 75, 10) 之间的距离。
作为一般经验法则，尽可能坚持 scipy 和 numpy 实现，因为它们是矢量化的并且比原生 Python 代码快得多。 （主要原因是：C 中的实现，矢量化消除了循环所做的类型检查开销。）
（旁白：我的回答不包括这里的精度，但我认为同样的原则适用于精度和效率。）
作为额外奖励，我将提供一些有关如何分析代码以衡量效率的信息。如果您使用的是 IPython 解释器，秘诀就是使用 %prun 行魔法。
In [1]: import numpy
In [2]: from scipy.spatial import distance
In [3]: c1 = numpy.array((52, 106, 35, 12))
In [4]: c2 = numpy.array((33, 153, 75, 10))
In [5]: %prun distance.euclidean(c1, c2)
35 function calls in 0.000 seconds
Ordered by: internal time
ncalls tottime percall cumtime percall filename:lineno(function)
1 0.000 0.000 0.000 0.000 {built-in method builtins.exec}
1 0.000 0.000 0.000 0.000 linalg.py:1976(norm)
1 0.000 0.000 0.000 0.000 {built-in method numpy.core.multiarray.dot}
6 0.000 0.000 0.000 0.000 {built-in method numpy.core.multiarray.array}
4 0.000 0.000 0.000 0.000 numeric.py:406(asarray)
1 0.000 0.000 0.000 0.000 distance.py:232(euclidean)
2 0.000 0.000 0.000 0.000 distance.py:152(_validate_vector)
2 0.000 0.000 0.000 0.000 shape_base.py:9(atleast_1d)
1 0.000 0.000 0.000 0.000 misc.py:11(norm)
1 0.000 0.000 0.000 0.000 function_base.py:605(asarray_chkfinite)
2 0.000 0.000 0.000 0.000 numeric.py:476(asanyarray)
1 0.000 0.000 0.000 0.000 {method 'ravel' of 'numpy.ndarray' objects}
1 0.000 0.000 0.000 0.000 linalg.py:111(isComplexType)
1 0.000 0.000 0.000 0.000 <string>:1(<module>)
2 0.000 0.000 0.000 0.000 {method 'append' of 'list' objects}
1 0.000 0.000 0.000 0.000 {built-in method builtins.issubclass}
4 0.000 0.000 0.000 0.000 {built-in method builtins.len}
1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
2 0.000 0.000 0.000 0.000 {method 'squeeze' of 'numpy.ndarray' objects}
In [6]: %prun numpy.linalg.norm(c1 - c2)
10 function calls in 0.000 seconds
Ordered by: internal time
ncalls tottime percall cumtime percall filename:lineno(function)
1 0.000 0.000 0.000 0.000 {built-in method builtins.exec}
1 0.000 0.000 0.000 0.000 linalg.py:1976(norm)
1 0.000 0.000 0.000 0.000 {built-in method numpy.core.multiarray.dot}
1 0.000 0.000 0.000 0.000 <string>:1(<module>)
1 0.000 0.000 0.000 0.000 numeric.py:406(asarray)
1 0.000 0.000 0.000 0.000 {method 'ravel' of 'numpy.ndarray' objects}
1 0.000 0.000 0.000 0.000 linalg.py:111(isComplexType)
1 0.000 0.000 0.000 0.000 {built-in method builtins.issubclass}
1 0.000 0.000 0.000 0.000 {built-in method numpy.core.multiarray.array}
1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
%prun 所做的是告诉您函数调用到 运行 需要多长时间，包括一些跟踪以找出瓶颈可能在哪里。在这种情况下，scipy.spatial.distance.euclidean 和 numpy.linalg.norm 实现都非常快。假设您定义了一个函数 dist(vect1, vect2)，您可以使用相同的 IPython 魔法调用进行分析。作为另一个额外的好处，%prun 也可以在 Jupyter notebook 中运行，您可以 %%prun 分析整个代码单元，而不仅仅是一个函数，只需将 %%prun 作为第一个该单元格的行。
结论先行：
从使用timeit进行效率测试的测试结果，我们可以得出结论关于效率：
Method5 (zip, math.sqrt) > Method1 (numpy.linalg.norm) > Method2 (scipy.spatial.distance) > Method3 (sklearn.metrics.pairwise.euclidean_distances )
虽然我没有真正测试你的Method4，因为它不适合一般情况，它通常等同于Method5。
对于其余的，令人惊讶的是，Method5 是最快的。而对于使用 numpy 的 Method1，正如我们预期的那样，它在 C 中进行了大量优化，是第二快的。
对于scipy.spatial.distance，如果你直接进入函数定义，你会看到它实际上使用了numpy.linalg.norm，除了它会在实际[之前对两个输入向量执行验证=24=]。这就是为什么它比 numpy.linalg.norm.
稍慢的原因
最后 sklearn，根据文档：
This formulation has two advantages over other ways of computing distances. First, it is computationally efficient when dealing with sparse data. Second, if one argument varies but the other remains unchanged, then dot(x, x) and/or dot(y, y) can be pre-computed.
However, this is not the most precise way of doing this computation, and the distance matrix returned by this function may not be exactly symmetric as required
由于你的问题中你希望使用固定的一组数据，所以这种实现方式的优势并没有体现出来。并且由于性能和精度之间的权衡，它也给出了所有方法中最差的精度。
关于精度,Method5=Metho1=Method2>Method3
效率测试脚本：
import numpy as np
from scipy.spatial import distance
from sklearn.metrics.pairwise import euclidean_distances
import math
# 1
def eudis1(v1, v2):
return np.linalg.norm(v1-v2)
# 2
def eudis2(v1, v2):
return distance.euclidean(v1, v2)
# 3
def eudis3(v1, v2):
return euclidean_distances(v1, v2)
# 5
def eudis5(v1, v2):
dist = [(a - b)**2 for a, b in zip(v1, v2)]
dist = math.sqrt(sum(dist))
return dist
dis1 = (52, 106, 35, 12)
dis2 = (33, 153, 75, 10)
v1, v2 = np.array(dis1), np.array(dis2)
import timeit
def wrapper(func, *args, **kwargs):
def wrapped():
return func(*args, **kwargs)
return wrapped
wrappered1 = wrapper(eudis1, v1, v2)
wrappered2 = wrapper(eudis2, v1, v2)
wrappered3 = wrapper(eudis3, v1, v2)
wrappered5 = wrapper(eudis5, v1, v2)
t1 = timeit.repeat(wrappered1, repeat=3, number=100000)
t2 = timeit.repeat(wrappered2, repeat=3, number=100000)
t3 = timeit.repeat(wrappered3, repeat=3, number=100000)
t5 = timeit.repeat(wrappered5, repeat=3, number=100000)
print('\n')
print('t1: ', sum(t1)/len(t1))
print('t2: ', sum(t2)/len(t2))
print('t3: ', sum(t3)/len(t3))
print('t5: ', sum(t5)/len(t5))
效率测试输出：
t1: 0.654838958307
t2: 1.53977598714
t3: 6.7898791732
t5: 0.422228400305
精度测试脚本和结果：
In [8]: eudis1(v1,v2)
Out[8]: 64.60650122085238
In [9]: eudis2(v1,v2)
Out[9]: 64.60650122085238
In [10]: eudis3(v1,v2)
Out[10]: array([[ 64.60650122]])
In [11]: eudis5(v1,v2)
Out[11]: 64.60650122085238
我不知道精度和速度与您提到的其他库相比如何，但您可以使用内置的 math.hypot() 函数对 2D 向量执行此操作：
from math import hypot
def pairwise(iterable):
"s -> (s0, s1), (s1, s2), (s2, s3), ..."
a, b = iter(iterable), iter(iterable)
next(b, None)
return zip(a, b)
a = (52, 106, 35, 12)
b = (33, 153, 75, 10)
dist = [hypot(p2[0]-p1[0], p2[1]-p1[1]) for p1, p2 in pairwise(tuple(zip(a, b)))]
print(dist) # -> [131.59027319676787, 105.47511554864494, 68.94925670375281]
这并不能完全回答问题，但可能值得一提的是，如果您对实际的欧几里得距离不感兴趣，而只是想相互比较欧几里得距离，则平方根是单调函数，即x**(1/2) < y**(1/2) 当且仅当 x < y.
因此，如果您不想要显式距离，但例如只想知道 vector1 的欧氏距离是否更接近一个向量列表，称为 vectorlist，您可以避免昂贵的（就两者而言精度和时间）平方根，但可以用
min(vectorlist, key = lambda compare: sum([(a - b)**2 for a, b in zip(vector1, compare)])
这里有一个关于如何只使用 numpy 的例子。
import numpy as np
a = np.array([3, 0])
b = np.array([0, 4])
c = np.sqrt(np.sum(((a - b) ** 2)))
# c == 5.0
改进已接受的 的基准，我发现，假设您已经获得 numpy 数组格式的输入，方法 5 可以更好地写成：
import numpy as np
from numba import jit
@jit(nopython=True)
def euclidian_distance(y1, y2):
return np.sqrt(np.sum((y1-y2)**2)) # based on pythagorean
速度测试：
euclidian_distance(y1, y2)
# 2.03 µs ± 138 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
np.linalg.norm(y1-y2)
# 17.6 µs ± 5.08 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
有趣的是，您可以将 jit 添加到 numpy 函数：
@jit(nopython=True)
def jit_linalg(y1, y2):
return np.linalg.norm(y1-y2)
jit_linalg(y[i],y[j])
# 2.91 µs ± 261 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
经过一些在线研究 (1, 2, numpy, scipy, scikit, math)，我发现了几种计算 Python 中的 欧氏距离的方法：
# 1
numpy.linalg.norm(a-b)
# 2
distance.euclidean(vector1, vector2)
# 3
sklearn.metrics.pairwise.euclidean_distances
# 4
sqrt((xa-xb)^2 + (ya-yb)^2 + (za-zb)^2)
# 5
dist = [(a - b)**2 for a, b in zip(vector1, vector2)]
dist = math.sqrt(sum(dist))
# 6
math.hypot(x, y)
我想知道是否有人可以提供关于以上哪一项（或我没有找到的任何其他项）被认为是最好的 效率和精度。如果有人知道讨论该主题的任何资源，那也很好。
我感兴趣的 context 是计算数字元组对之间的欧氏距离，例如(52, 106, 35, 12) 和 (33, 153, 75, 10) 之间的距离。
作为一般经验法则，尽可能坚持 scipy 和 numpy 实现，因为它们是矢量化的并且比原生 Python 代码快得多。 （主要原因是：C 中的实现，矢量化消除了循环所做的类型检查开销。）
（旁白：我的回答不包括这里的精度，但我认为同样的原则适用于精度和效率。）
作为额外奖励，我将提供一些有关如何分析代码以衡量效率的信息。如果您使用的是 IPython 解释器，秘诀就是使用 %prun 行魔法。
In [1]: import numpy
In [2]: from scipy.spatial import distance
In [3]: c1 = numpy.array((52, 106, 35, 12))
In [4]: c2 = numpy.array((33, 153, 75, 10))
In [5]: %prun distance.euclidean(c1, c2)
35 function calls in 0.000 seconds
Ordered by: internal time
ncalls tottime percall cumtime percall filename:lineno(function)
1 0.000 0.000 0.000 0.000 {built-in method builtins.exec}
1 0.000 0.000 0.000 0.000 linalg.py:1976(norm)
1 0.000 0.000 0.000 0.000 {built-in method numpy.core.multiarray.dot}
6 0.000 0.000 0.000 0.000 {built-in method numpy.core.multiarray.array}
4 0.000 0.000 0.000 0.000 numeric.py:406(asarray)
1 0.000 0.000 0.000 0.000 distance.py:232(euclidean)
2 0.000 0.000 0.000 0.000 distance.py:152(_validate_vector)
2 0.000 0.000 0.000 0.000 shape_base.py:9(atleast_1d)
1 0.000 0.000 0.000 0.000 misc.py:11(norm)
1 0.000 0.000 0.000 0.000 function_base.py:605(asarray_chkfinite)
2 0.000 0.000 0.000 0.000 numeric.py:476(asanyarray)
1 0.000 0.000 0.000 0.000 {method 'ravel' of 'numpy.ndarray' objects}
1 0.000 0.000 0.000 0.000 linalg.py:111(isComplexType)
1 0.000 0.000 0.000 0.000 <string>:1(<module>)
2 0.000 0.000 0.000 0.000 {method 'append' of 'list' objects}
1 0.000 0.000 0.000 0.000 {built-in method builtins.issubclass}
4 0.000 0.000 0.000 0.000 {built-in method builtins.len}
1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
2 0.000 0.000 0.000 0.000 {method 'squeeze' of 'numpy.ndarray' objects}
In [6]: %prun numpy.linalg.norm(c1 - c2)
10 function calls in 0.000 seconds
Ordered by: internal time
ncalls tottime percall cumtime percall filename:lineno(function)
1 0.000 0.000 0.000 0.000 {built-in method builtins.exec}
1 0.000 0.000 0.000 0.000 linalg.py:1976(norm)
1 0.000 0.000 0.000 0.000 {built-in method numpy.core.multiarray.dot}
1 0.000 0.000 0.000 0.000 <string>:1(<module>)
1 0.000 0.000 0.000 0.000 numeric.py:406(asarray)
1 0.000 0.000 0.000 0.000 {method 'ravel' of 'numpy.ndarray' objects}
1 0.000 0.000 0.000 0.000 linalg.py:111(isComplexType)
1 0.000 0.000 0.000 0.000 {built-in method builtins.issubclass}
1 0.000 0.000 0.000 0.000 {built-in method numpy.core.multiarray.array}
1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}
%prun 所做的是告诉您函数调用到 运行 需要多长时间，包括一些跟踪以找出瓶颈可能在哪里。在这种情况下，scipy.spatial.distance.euclidean 和 numpy.linalg.norm 实现都非常快。假设您定义了一个函数 dist(vect1, vect2)，您可以使用相同的 IPython 魔法调用进行分析。作为另一个额外的好处，%prun 也可以在 Jupyter notebook 中运行，您可以 %%prun 分析整个代码单元，而不仅仅是一个函数，只需将 %%prun 作为第一个该单元格的行。
结论先行：
从使用timeit进行效率测试的测试结果，我们可以得出结论关于效率：
Method5 (zip, math.sqrt) > Method1 (numpy.linalg.norm) > Method2 (scipy.spatial.distance) > Method3 (sklearn.metrics.pairwise.euclidean_distances )
虽然我没有真正测试你的Method4，因为它不适合一般情况，它通常等同于Method5。
对于其余的，令人惊讶的是，Method5 是最快的。而对于使用 numpy 的 Method1，正如我们预期的那样，它在 C 中进行了大量优化，是第二快的。
对于scipy.spatial.distance，如果你直接进入函数定义，你会看到它实际上使用了numpy.linalg.norm，除了它会在实际[之前对两个输入向量执行验证=24=]。这就是为什么它比 numpy.linalg.norm.
最后 sklearn，根据文档：
This formulation has two advantages over other ways of computing distances. First, it is computationally efficient when dealing with sparse data. Second, if one argument varies but the other remains unchanged, then dot(x, x) and/or dot(y, y) can be pre-computed. However, this is not the most precise way of doing this computation, and the distance matrix returned by this function may not be exactly symmetric as required
由于你的问题中你希望使用固定的一组数据，所以这种实现方式的优势并没有体现出来。并且由于性能和精度之间的权衡，它也给出了所有方法中最差的精度。
关于精度,Method5=Metho1=Method2>Method3
效率测试脚本：
import numpy as np
from scipy.spatial import distance
from sklearn.metrics.pairwise import euclidean_distances
import math
# 1
def eudis1(v1, v2):
return np.linalg.norm(v1-v2)
# 2
def eudis2(v1, v2):
return distance.euclidean(v1, v2)
# 3
def eudis3(v1, v2):
return euclidean_distances(v1, v2)
# 5
def eudis5(v1, v2):
dist = [(a - b)**2 for a, b in zip(v1, v2)]
dist = math.sqrt(sum(dist))
return dist
dis1 = (52, 106, 35, 12)
dis2 = (33, 153, 75, 10)
v1, v2 = np.array(dis1), np.array(dis2)
import timeit
def wrapper(func, *args, **kwargs):
def wrapped():
return func(*args, **kwargs)
return wrapped
wrappered1 = wrapper(eudis1, v1, v2)
wrappered2 = wrapper(eudis2, v1, v2)
wrappered3 = wrapper(eudis3, v1, v2)
wrappered5 = wrapper(eudis5, v1, v2)
t1 = timeit.repeat(wrappered1, repeat=3, number=100000)
t2 = timeit.repeat(wrappered2, repeat=3, number=100000)
t3 = timeit.repeat(wrappered3, repeat=3, number=100000)
t5 = timeit.repeat(wrappered5, repeat=3, number=100000)
print('\n')
print('t1: ', sum(t1)/len(t1))
print('t2: ', sum(t2)/len(t2))
print('t3: ', sum(t3)/len(t3))
print('t5: ', sum(t5)/len(t5))
效率测试输出：
t1: 0.654838958307
t2: 1.53977598714
t3: 6.7898791732
t5: 0.422228400305
精度测试脚本和结果：
In [8]: eudis1(v1,v2)
Out[8]: 64.60650122085238
In [9]: eudis2(v1,v2)
Out[9]: 64.60650122085238
In [10]: eudis3(v1,v2)
Out[10]: array([[ 64.60650122]])
In [11]: eudis5(v1,v2)
Out[11]: 64.60650122085238
我不知道精度和速度与您提到的其他库相比如何，但您可以使用内置的 math.hypot() 函数对 2D 向量执行此操作：
from math import hypot
def pairwise(iterable):
"s -> (s0, s1), (s1, s2), (s2, s3), ..."
a, b = iter(iterable), iter(iterable)
next(b, None)
return zip(a, b)
a = (52, 106, 35, 12)
b = (33, 153, 75, 10)
dist = [hypot(p2[0]-p1[0], p2[1]-p1[1]) for p1, p2 in pairwise(tuple(zip(a, b)))]
print(dist) # -> [131.59027319676787, 105.47511554864494, 68.94925670375281]
这并不能完全回答问题，但可能值得一提的是，如果您对实际的欧几里得距离不感兴趣，而只是想相互比较欧几里得距离，则平方根是单调函数，即x**(1/2) < y**(1/2) 当且仅当 x < y.
因此，如果您不想要显式距离，但例如只想知道 vector1 的欧氏距离是否更接近一个向量列表，称为 vectorlist，您可以避免昂贵的（就两者而言精度和时间）平方根，但可以用
min(vectorlist, key = lambda compare: sum([(a - b)**2 for a, b in zip(vector1, compare)])
这里有一个关于如何只使用 numpy 的例子。
import numpy as np
a = np.array([3, 0])
b = np.array([0, 4])
c = np.sqrt(np.sum(((a - b) ** 2)))
# c == 5.0
改进已接受的 numpy 数组格式的输入，方法 5 可以更好地写成：
import numpy as np
from numba import jit
@jit(nopython=True)
def euclidian_distance(y1, y2):
return np.sqrt(np.sum((y1-y2)**2)) # based on pythagorean
速度测试：
euclidian_distance(y1, y2)
# 2.03 µs ± 138 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
np.linalg.norm(y1-y2)
# 17.6 µs ± 5.08 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
有趣的是，您可以将 jit 添加到 numpy 函数：
@jit(nopython=True)
def jit_linalg(y1, y2):
return np.linalg.norm(y1-y2)
jit_linalg(y[i],y[j])
# 2.91 µs ± 261 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)