在三维空间中有两个点
a = (ax, ay, az)
b = (bx, by, bz)
我想计算它们之间的距离:
dist = sqrt((ax-bx)^2 + (ay-by)^2 + (az-bz)^2)
我如何用NumPy做到这一点?我有:
import numpy
a = numpy.array((ax, ay, az))
b = numpy.array((bx, by, bz))
在三维空间中有两个点
a = (ax, ay, az)
b = (bx, by, bz)
我想计算它们之间的距离:
dist = sqrt((ax-bx)^2 + (ay-by)^2 + (az-bz)^2)
我如何用NumPy做到这一点?我有:
import numpy
a = numpy.array((ax, ay, az))
b = numpy.array((bx, by, bz))
这种解决问题方法的另一个例子:
def dist(x,y):
return numpy.sqrt(numpy.sum((x-y)**2))
a = numpy.array((xa,ya,za))
b = numpy.array((xb,yb,zb))
dist_a_b = dist(a,b)
使用numpy.linalg.norm:
dist = numpy.linalg.norm(a-b)
这是因为欧氏距离是l2范数,而numpy.linalg.norm中ord参数的默认值是2。 要了解更多理论,请参阅数据挖掘介绍:
一个很好的句子:
dist = numpy.linalg.norm(a-b)
但是,如果速度是一个问题,我建议在您的机器上进行试验。我发现,在我的机器上,使用数学库的√rt和**运算符对平方进行运算要比使用一行NumPy解决方案快得多。
我用这个简单的程序进行了测试:
#!/usr/bin/python
import math
import numpy
from random import uniform
def fastest_calc_dist(p1,p2):
return math.sqrt((p2[0] - p1[0]) ** 2 +
(p2[1] - p1[1]) ** 2 +
(p2[2] - p1[2]) ** 2)
def math_calc_dist(p1,p2):
return math.sqrt(math.pow((p2[0] - p1[0]), 2) +
math.pow((p2[1] - p1[1]), 2) +
math.pow((p2[2] - p1[2]), 2))
def numpy_calc_dist(p1,p2):
return numpy.linalg.norm(numpy.array(p1)-numpy.array(p2))
TOTAL_LOCATIONS = 1000
p1 = dict()
p2 = dict()
for i in range(0, TOTAL_LOCATIONS):
p1[i] = (uniform(0,1000),uniform(0,1000),uniform(0,1000))
p2[i] = (uniform(0,1000),uniform(0,1000),uniform(0,1000))
total_dist = 0
for i in range(0, TOTAL_LOCATIONS):
for j in range(0, TOTAL_LOCATIONS):
dist = fastest_calc_dist(p1[i], p2[j]) #change this line for testing
total_dist += dist
print total_dist
在我的机器上,math_calc_dist运行得比numpy_calc_dist快得多:1.5秒对23.5秒。
为了在fastst_calc_dist和math_calc_dist之间获得一个可测量的差异,我必须将TOTAL_LOCATIONS增加到6000。然后,fastst_calc_dist耗时约50秒,math_calc_dist耗时约60秒。
您也可以尝试使用numpy。SQRT和numpy。不过这两个运算都比我机器上的数学运算要慢。
我的测试使用Python 2.6.6运行。
你可以减去向量,然后内积。
以你为榜样,
a = numpy.array((xa, ya, za))
b = numpy.array((xb, yb, zb))
tmp = a - b
sum_squared = numpy.dot(tmp.T, tmp)
result = numpy.sqrt(sum_squared)
可以像下面这样做。我不知道它有多快,但它没有使用NumPy。
from math import sqrt
a = (1, 2, 3) # Data point 1
b = (4, 5, 6) # Data point 2
print sqrt(sum( (a - b)**2 for a, b in zip(a, b)))
我在matplotlib中找到了一个“dist”函数。mlab,但我认为它不够方便。
我把它贴在这里只是为了参考。
import numpy as np
import matplotlib as plt
a = np.array([1, 2, 3])
b = np.array([2, 3, 4])
# Distance between a and b
dis = plt.mlab.dist(a, b)
使用scipy.spatial.distance.euclidean:
from scipy.spatial import distance
a = (1, 2, 3)
b = (4, 5, 6)
dst = distance.euclidean(a, b)
这里有一些简洁的Python欧几里得距离代码,给出了Python中以列表表示的两个点。
def distance(v1,v2):
return sum([(x-y)**2 for (x,y) in zip(v1,v2)])**(0.5)
我喜欢np。点(点积):
a = numpy.array((xa,ya,za))
b = numpy.array((xb,yb,zb))
distance = (np.dot(a-b,a-b))**.5
计算多维空间的欧氏距离:
import math
x = [1, 2, 6]
y = [-2, 3, 2]
dist = math.sqrt(sum([(xi-yi)**2 for xi,yi in zip(x, y)]))
5.0990195135927845
我想用各种表演笔记来阐述这个简单的答案。Np.linalg.norm可能会做的比你需要的更多:
dist = numpy.linalg.norm(a-b)
首先,这个函数被设计用于处理一个列表并返回所有的值,例如比较pA到点集sP的距离:
sP = set(points)
pA = point
distances = np.linalg.norm(sP - pA, ord=2, axis=1.) # 'distances' is a list
记住几件事:
Python函数调用的开销很大。 [常规]Python不缓存名称查找。
So
def distance(pointA, pointB):
dist = np.linalg.norm(pointA - pointB)
return dist
并不像看上去那么无辜。
>>> dis.dis(distance)
2 0 LOAD_GLOBAL 0 (np)
2 LOAD_ATTR 1 (linalg)
4 LOAD_ATTR 2 (norm)
6 LOAD_FAST 0 (pointA)
8 LOAD_FAST 1 (pointB)
10 BINARY_SUBTRACT
12 CALL_FUNCTION 1
14 STORE_FAST 2 (dist)
3 16 LOAD_FAST 2 (dist)
18 RETURN_VALUE
首先,每次我们调用它时,我们都必须对“np”进行全局查找,对“linalg”进行范围查找,对“norm”进行范围查找,而仅仅调用这个函数的开销就相当于几十条python指令。
最后,我们浪费了两个操作来存储结果并重新加载它以返回…
改进的第一步:使查找更快,跳过存储
def distance(pointA, pointB, _norm=np.linalg.norm):
return _norm(pointA - pointB)
我们得到了更精简的:
>>> dis.dis(distance)
2 0 LOAD_FAST 2 (_norm)
2 LOAD_FAST 0 (pointA)
4 LOAD_FAST 1 (pointB)
6 BINARY_SUBTRACT
8 CALL_FUNCTION 1
10 RETURN_VALUE
不过,函数调用开销仍然需要一些工作。你会想要做基准测试,以确定你自己做数学是否会更好:
def distance(pointA, pointB):
return (
((pointA.x - pointB.x) ** 2) +
((pointA.y - pointB.y) ** 2) +
((pointA.z - pointB.z) ** 2)
) ** 0.5 # fast sqrt
在某些平台上,**0.5比math.sqrt快。你的里程可能会有所不同。
****高级性能说明。
你为什么要计算距离?如果唯一的目的是展示它,
print("The target is %.2fm away" % (distance(a, b)))
沿着。但是如果你在比较距离,进行范围检查等等,我想添加一些有用的性能观察。
让我们以两种情况为例:按距离排序或将列表剔除到满足范围约束的项。
# Ultra naive implementations. Hold onto your hat.
def sort_things_by_distance(origin, things):
return things.sort(key=lambda thing: distance(origin, thing))
def in_range(origin, range, things):
things_in_range = []
for thing in things:
if distance(origin, thing) <= range:
things_in_range.append(thing)
我们需要记住的第一件事是,我们使用毕达哥拉斯来计算距离(dist =根号(x²+ y²+ z²)),所以我们做了很多根号调用。数学101:
dist = root ( x^2 + y^2 + z^2 )
:.
dist^2 = x^2 + y^2 + z^2
and
sq(N) < sq(M) iff M > N
and
sq(N) > sq(M) iff N > M
and
sq(N) = sq(M) iff N == M
简而言之:直到我们真正需要以X为单位的距离,而不是X^2,我们才能消除计算中最难的部分。
# Still naive, but much faster.
def distance_sq(left, right):
""" Returns the square of the distance between left and right. """
return (
((left.x - right.x) ** 2) +
((left.y - right.y) ** 2) +
((left.z - right.z) ** 2)
)
def sort_things_by_distance(origin, things):
return things.sort(key=lambda thing: distance_sq(origin, thing))
def in_range(origin, range, things):
things_in_range = []
# Remember that sqrt(N)**2 == N, so if we square
# range, we don't need to root the distances.
range_sq = range**2
for thing in things:
if distance_sq(origin, thing) <= range_sq:
things_in_range.append(thing)
很好,这两个函数都不再做昂贵的平方根了。这样会快得多,但在进一步讨论之前,请检查自己:为什么sort_things_by_distance两次都需要一个“天真”的免责声明?在最下面回答(*a1)。
我们可以通过将in_range转换为生成器来改进它:
def in_range(origin, range, things):
range_sq = range**2
yield from (thing for thing in things
if distance_sq(origin, thing) <= range_sq)
如果你在做以下事情,这尤其有好处:
if any(in_range(origin, max_dist, things)):
...
但如果你接下来要做的事需要一段距离,
for nearby in in_range(origin, walking_distance, hotdog_stands):
print("%s %.2fm" % (nearby.name, distance(origin, nearby)))
考虑生成元组:
def in_range_with_dist_sq(origin, range, things):
range_sq = range**2
for thing in things:
dist_sq = distance_sq(origin, thing)
if dist_sq <= range_sq: yield (thing, dist_sq)
如果你可能要进行连锁范围检查(“找到X附近和Y的Nm范围内的东西”,因为你不需要再次计算距离),这可能特别有用。
但如果我们在搜索一个很大的列表,我们预计其中有很多不值得考虑呢?
其实有一个很简单的优化:
def in_range_all_the_things(origin, range, things):
range_sq = range**2
for thing in things:
dist_sq = (origin.x - thing.x) ** 2
if dist_sq <= range_sq:
dist_sq += (origin.y - thing.y) ** 2
if dist_sq <= range_sq:
dist_sq += (origin.z - thing.z) ** 2
if dist_sq <= range_sq:
yield thing
这是否有用取决于“事物”的大小。
def in_range_all_the_things(origin, range, things):
range_sq = range**2
if len(things) >= 4096:
for thing in things:
dist_sq = (origin.x - thing.x) ** 2
if dist_sq <= range_sq:
dist_sq += (origin.y - thing.y) ** 2
if dist_sq <= range_sq:
dist_sq += (origin.z - thing.z) ** 2
if dist_sq <= range_sq:
yield thing
elif len(things) > 32:
for things in things:
dist_sq = (origin.x - thing.x) ** 2
if dist_sq <= range_sq:
dist_sq += (origin.y - thing.y) ** 2 + (origin.z - thing.z) ** 2
if dist_sq <= range_sq:
yield thing
else:
... just calculate distance and range-check it ...
同样,考虑生成dist_sq。热狗的例子就变成了:
# Chaining generators
info = in_range_with_dist_sq(origin, walking_distance, hotdog_stands)
info = (stand, dist_sq**0.5 for stand, dist_sq in info)
for stand, dist in info:
print("%s %.2fm" % (stand, dist))
(*a1: sort_things_by_distance的排序键为每一项调用distance_sq,而那个看起来无辜的键是一个lambda,这是必须调用的第二个函数…)
对于那些对一次计算多个距离感兴趣的人来说,我已经使用perfplot(我的一个小项目)做了一些比较。
第一个建议是组织数据,使数组具有维数(3,n)(显然是c连续的)。如果加法发生在连续的第一维中,事情会更快,如果你使用带有axis=0的sqrt-sum, linalg,它也不会太重要。轴=0的范数,或
a_min_b = a - b
numpy.sqrt(numpy.einsum('ij,ij->j', a_min_b, a_min_b))
这是,以微弱优势,最快的变种。(这实际上也只适用于一行。)
在第二个轴上求和的变量,轴=1,都要慢得多。
代码重现情节:
import numpy
import perfplot
from scipy.spatial import distance
def linalg_norm(data):
a, b = data[0]
return numpy.linalg.norm(a - b, axis=1)
def linalg_norm_T(data):
a, b = data[1]
return numpy.linalg.norm(a - b, axis=0)
def sqrt_sum(data):
a, b = data[0]
return numpy.sqrt(numpy.sum((a - b) ** 2, axis=1))
def sqrt_sum_T(data):
a, b = data[1]
return numpy.sqrt(numpy.sum((a - b) ** 2, axis=0))
def scipy_distance(data):
a, b = data[0]
return list(map(distance.euclidean, a, b))
def sqrt_einsum(data):
a, b = data[0]
a_min_b = a - b
return numpy.sqrt(numpy.einsum("ij,ij->i", a_min_b, a_min_b))
def sqrt_einsum_T(data):
a, b = data[1]
a_min_b = a - b
return numpy.sqrt(numpy.einsum("ij,ij->j", a_min_b, a_min_b))
def setup(n):
a = numpy.random.rand(n, 3)
b = numpy.random.rand(n, 3)
out0 = numpy.array([a, b])
out1 = numpy.array([a.T, b.T])
return out0, out1
b = perfplot.bench(
setup=setup,
n_range=[2 ** k for k in range(22)],
kernels=[
linalg_norm,
linalg_norm_T,
scipy_distance,
sqrt_sum,
sqrt_sum_T,
sqrt_einsum,
sqrt_einsum_T,
],
xlabel="len(x), len(y)",
)
b.save("norm.png")
import numpy as np
from scipy.spatial import distance
input_arr = np.array([[0,3,0],[2,0,0],[0,1,3],[0,1,2],[-1,0,1],[1,1,1]])
test_case = np.array([0,0,0])
dst=[]
for i in range(0,6):
temp = distance.euclidean(test_case,input_arr[i])
dst.append(temp)
print(dst)
这个公式很容易用
distance = np.sqrt(np.sum(np.square(a-b)))
它实际上只是使用毕达哥拉斯定理来计算距离,通过将Δx, Δy和Δz的平方相加,并对结果进行根运算。
首先求两个矩阵的差。然后,使用numpy的multiply命令应用元素乘法。然后,求元素与新矩阵相乘的和。最后,求求和的平方根。
def findEuclideanDistance(a, b):
euclidean_distance = a - b
euclidean_distance = np.sum(np.multiply(euclidean_distance, euclidean_distance))
euclidean_distance = np.sqrt(euclidean_distance)
return euclidean_distance
从Python 3.8开始,math模块直接提供dist函数,它返回两点之间的欧几里得距离(以元组或坐标列表的形式给出):
from math import dist
dist((1, 2, 6), (-2, 3, 2)) # 5.0990195135927845
如果你使用列表:
dist([1, 2, 6], [-2, 3, 2]) # 5.0990195135927845
从Python 3.8开始
从Python 3.8开始,数学模块包含了math.dist()函数。 请看这里https://docs.python.org/3.8/library/math.html#math.dist。
数学。dist (p1, p2) 返回两点p1和p2之间的欧氏距离, 每一个都以坐标序列(或可迭代对象)给出。
import math
print( math.dist( (0,0), (1,1) )) # sqrt(2) -> 1.4142
print( math.dist( (0,0,0), (1,1,1) )) # sqrt(3) -> 1.7321
使用Python 3.8,这非常简单。
https://docs.python.org/3/library/math.html#math.dist
math.dist(p, q)
返回给定两点p和q之间的欧氏距离 作为坐标序列(或可迭代对象)。这两点一定有 相同的维度。 大致相当于: √(sum((px - qx) ** 2.0 for px, qx in zip(p, q)))
import numpy as np
# any two python array as two points
a = [0, 0]
b = [3, 4]
首先将list更改为numpy数组,并像这样做:print(np.linalg.norm(np.array(a) - np.array(b)))。第二种方法直接从python列表as: print(np.linalg.norm(np.subtract(a,b)))
其他答案适用于浮点数,但不能正确计算整数dtype的距离,因为整数dtype容易溢出和下溢。注意,即使scipy.distance.euclidean也有这个问题:
>>> a1 = np.array([1], dtype='uint8')
>>> a2 = np.array([2], dtype='uint8')
>>> a1 - a2
array([255], dtype=uint8)
>>> np.linalg.norm(a1 - a2)
255.0
>>> from scipy.spatial import distance
>>> distance.euclidean(a1, a2)
255.0
这是很常见的,因为许多图像库将图像表示为dtype="uint8"的ndarray。这意味着,如果你有一个由深灰色像素组成的灰度图像(比如所有像素的颜色都是#000001),你将它与黑色图像(#000000)进行区分,你最终可以在所有单元格中得到x-y为255的图像,这表明两张图像彼此相距很远。对于无符号整数类型(例如uint8),你可以安全地在numpy中计算距离:
np.linalg.norm(np.maximum(x, y) - np.minimum(x, y))
对于有符号整型,可以先转换为浮点类型:
np.linalg.norm(x.astype("float") - y.astype("float"))
对于图像数据,你可以使用opencv的norm方法:
import cv2
cv2.norm(x, y, cv2.NORM_L2)
用NumPy或一般的Python做这件事的最好方法是什么?我有:
最好的方法是最安全的,也是最快的
我建议使用低流量来获得可靠的结果,因为与编写自己的平方根计算器相比,下溢和溢出的几率非常小
我们来看看数学。函数,np。Hypot vs vanilla np.sqrt(np.sum(np.数组([i, j, k])) ** 2,轴=1))
i, j, k = 1e+200, 1e+200, 1e+200
math.hypot(i, j, k)
# 1.7320508075688773e+200
np.sqrt(np.sum((np.array([i, j, k])) ** 2))
# RuntimeWarning: overflow encountered in square
速度智慧的数学。Hypot看起来更好
%%timeit
math.hypot(i, j, k)
# 100 ns ± 1.05 ns per loop (mean ± std. dev. of 7 runs, 10000000 loops each)
%%timeit
np.sqrt(np.sum((np.array([i, j, k])) ** 2))
# 6.41 µs ± 33.3 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
下溢
i, j = 1e-200, 1e-200
np.sqrt(i**2+j**2)
# 0.0
溢出
i, j = 1e+200, 1e+200
np.sqrt(i**2+j**2)
# inf
没有下溢
i, j = 1e-200, 1e-200
np.hypot(i, j)
# 1.414213562373095e-200
没有溢出
i, j = 1e+200, 1e+200
np.hypot(i, j)
# 1.414213562373095e+200
请参考
对于大量的距离,我能想到的最快的解决方案是使用numexpr。在我的机器上,它比使用numpy einsum更快:
import numexpr as ne
import numpy as np
np.sqrt(ne.evaluate("sum((a_min_b)**2,axis=1)"))