Comment calculer la distance euclidienne avec NumPy?
-
05-07-2019 - |
Question
J'ai deux points en 3D:
(xa, ya, za)
(xb, yb, zb)
Et je veux calculer la distance:
dist = sqrt((xa-xb)^2 + (ya-yb)^2 + (za-zb)^2)
Quelle est la meilleure façon de faire cela avec NumPy ou avec Python en général? J'ai:
a = numpy.array((xa ,ya, za))
b = numpy.array((xb, yb, zb))
La solution
Utilisez numpy.linalg.norm
:
dist = numpy.linalg.norm(a-b)
Autres conseils
Il y a une fonction pour cela dans SciPy. C'est ce qu'on appelle Euclidean .
Exemple:
from scipy.spatial import distance
a = (1, 2, 3)
b = (4, 5, 6)
dst = distance.euclidean(a, b)
Pour ceux qui souhaitent calculer plusieurs distances à la fois, j'ai effectué une petite comparaison à l'aide de perfplot ( un petit projet à moi). Il s'avère que
a_min_b = a - b
numpy.sqrt(numpy.einsum('ij,ij->i', a_min_b, a_min_b))
calcule les distances des lignes dans a
et b
le plus rapide. C’est vrai aussi pour une seule ligne!
Code pour reproduire l'intrigue:
import matplotlib
import numpy
import perfplot
from scipy.spatial import distance
def linalg_norm(data):
a, b = data
return numpy.linalg.norm(a-b, axis=1)
def sqrt_sum(data):
a, b = data
return numpy.sqrt(numpy.sum((a-b)**2, axis=1))
def scipy_distance(data):
a, b = data
return list(map(distance.euclidean, a, b))
def mpl_dist(data):
a, b = data
return list(map(matplotlib.mlab.dist, a, b))
def sqrt_einsum(data):
a, b = data
a_min_b = a - b
return numpy.sqrt(numpy.einsum('ij,ij->i', a_min_b, a_min_b))
perfplot.show(
setup=lambda n: numpy.random.rand(2, n, 3),
n_range=[2**k for k in range(20)],
kernels=[linalg_norm, scipy_distance, mpl_dist, sqrt_sum, sqrt_einsum],
logx=True,
logy=True,
xlabel='len(x), len(y)'
)
Un autre exemple de cette méthode de résolution de problèmes :
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)
Je souhaite exposer la réponse simple avec diverses notes de performance. np.linalg.norm fera peut-être plus que ce dont vous avez besoin:
dist = numpy.linalg.norm(a-b)
Tout d’abord - cette fonction est conçue pour traiter une liste et renvoyer toutes les valeurs, par exemple. comparer la distance entre pA
et l'ensemble des points sP
:
sP = set(points)
pA = point
distances = np.linalg.norm(sP - pA, ord=2, axis=1.) # 'distances' is a list
Rappelez-vous plusieurs choses:
- Les appels de fonctions Python sont coûteux.
- [Normal] Python ne met pas en cache les recherches de noms.
Alors
def distance(pointA, pointB):
dist = np.linalg.norm(pointA - pointB)
return dist
n'est pas aussi innocent qu'il n'y paraît.
>>> 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
Tout d'abord, chaque fois que nous l'appelons, nous devons effectuer une recherche globale pour "np", une recherche ciblée pour "linalg". et une recherche ciblée pour "norme", et la surcharge de simplement appeler la fonction peuvent équivaloir à des dizaines d'instructions en python.
Enfin, nous avons perdu deux opérations pour stocker le résultat et le recharger pour le retour ...
Premier passage à l'amélioration: accélérez la recherche, ignorez le magasin
def distance(pointA, pointB, _norm=np.linalg.norm):
return _norm(pointA - pointB)
Nous obtenons le beaucoup plus simple:
>>> 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
Toutefois, le temps système d’appel de la fonction représente toujours du travail. Et vous voudrez faire des points de repère pour déterminer si vous feriez mieux de faire le calcul vous-même:
def distance(pointA, pointB):
return (
((pointA.x - pointB.x) ** 2) +
((pointA.y - pointB.y) ** 2) +
((pointA.z - pointB.z) ** 2)
) ** 0.5 # fast sqrt
Sur certaines plateformes, ** 0.5
est plus rapide que math.sqrt
. Votre kilométrage peut varier.
**** Notes de performances avancées.
Pourquoi calculez-vous la distance? Si le seul but est de l'afficher,
print("The target is %.2fm away" % (distance(a, b)))
avancez. Mais si vous comparez des distances, effectuez des contrôles de distance, etc., j'aimerais ajouter quelques observations utiles sur les performances.
Prenons deux cas: le tri par distance ou la sélection d'une liste d'éléments correspondant à une contrainte de plage.
# 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)
La première chose à retenir est que nous utilisons Pythagoras pour calculer le distance ( dist = sqrt (x ^ 2 + y ^ 2 + z ^ 2)
), nous faisons donc beaucoup d'appels sqrt
. Math 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
En bref: jusqu’à ce que nous demandions réellement la distance dans une unité de X plutôt que de X ^ 2, nous pouvons éliminer la partie la plus difficile des calculs.
# 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)
Génial, les deux fonctions ne font plus de racines carrées coûteuses. Ce sera beaucoup plus rapide. Nous pouvons également améliorer in_range en le convertissant en générateur:
def in_range(origin, range, things):
range_sq = range**2
yield from (thing for thing in things
if distance_sq(origin, thing) <= range_sq)
Cela présente des avantages si vous faites quelque chose comme:
if any(in_range(origin, max_dist, things)):
...
Mais si la prochaine chose que vous allez faire nécessite une distance,
for nearby in in_range(origin, walking_distance, hotdog_stands):
print("%s %.2fm" % (nearby.name, distance(origin, nearby)))
envisagez de céder des n-uplets:
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)
Cela peut être particulièrement utile si vous pouvez enchaîner les vérifications de plage ("rechercher des objets proches de X et à moins de Nm de Y", car vous n'avez pas à calculer de nouveau la distance).
Mais qu'en est-il si nous recherchons une très grande liste de éléments
et si nous anticipons qu'un grand nombre d'entre eux ne méritent pas d'être pris en compte?
Il existe en réalité une optimisation très simple:
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
Si cela est utile dépendra de la taille des "choses".
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 ...
Et encore une fois, envisagez de céder le dist_sq. Notre exemple de hot-dog devient alors:
# 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))
Je trouve une fonction 'dist' dans matplotlib.mlab, mais je ne pense pas que ce soit assez pratique.
Je le poste ici juste pour référence.
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)
Cela peut se faire comme suit. Je ne sais pas à quelle vitesse, mais il n’utilise pas 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)))
Vous pouvez simplement soustraire les vecteurs et ensuite le produit interne.
Suivant votre exemple,
a = numpy.array((xa, ya, za))
b = numpy.array((xb, yb, zb))
tmp = a - b
sum_squared = numpy.dot(tmp.T, tmp)
result sqrt(sum_squared)
Il s’agit d’un code simple et facile à comprendre.
J'aime np.dot
(produit scalaire):
a = numpy.array((xa,ya,za))
b = numpy.array((xb,yb,zb))
distance = (np.dot(a-b,a-b))**.5
Démarrage de Python 3.8
, le math < Le module / code>
fournit directement le dist
, qui renvoie la distance euclidienne entre deux points (exprimée en tuple de coordonnées):
from math import dist
dist((1, 2, 6), (-2, 3, 2)) # 5.0990195135927845
Si vous travaillez avec des listes au lieu de tuples:
dist(tuple([1, 2, 6]), tuple([-2, 3, 2]))
Avec a
et b
tels que vous les avez définis, vous pouvez également utiliser:
distance = np.sqrt(np.sum((a-b)**2))
Un beau one-liner:
dist = numpy.linalg.norm(a-b)
Cependant, si la vitesse est un problème, je vous recommande de faire des essais sur votre machine. J'ai constaté que l'utilisation de sqrt
de la bibliothèque math
avec l'opérateur **
pour le carré est beaucoup plus rapide sur ma machine que la solution à une seule ligne NumPy .
J'ai exécuté mes tests à l'aide de ce programme simple:
#!/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
Sur ma machine, math_calc_dist
est beaucoup plus rapide que numpy_calc_dist
: 1,5 seconde contre 23,5 secondes . . p>
Pour obtenir une différence mesurable entre greatest_calc_dist
et math_calc_dist
, je devais augmenter de TOTAL_LOCATIONS
à 6000. Ensuite, most_calc_dist
prend ~ 50 secondes tandis que math_calc_dist
prend ~ 60 secondes .
Vous pouvez également expérimenter avec numpy.sqrt
et numpy.square
bien qu'ils soient tous deux plus lents que les alternatives math
sur ma machine.
Mes tests ont été exécutés avec Python 2.6.6.
Voici un code concis pour la distance euclidienne en Python à partir de deux points représentés sous forme de listes en Python.
def distance(v1,v2):
return sum([(x-y)**2 for (x,y) in zip(v1,v2)])**(0.5)
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)
import math
dist = math.hypot(math.hypot(xa-xb, ya-yb), za-zb)
Vous pouvez facilement utiliser la formule
distance = np.sqrt(np.sum(np.square(a-b)))
qui ne fait en réalité rien d'autre que d'utiliser le théorème de Pythagore pour calculer la distance, en ajoutant les carrés de & x 916; x, & # 916; y et & # 916; z et en enracinant le résultat.
Calculez la distance euclidienne pour un espace multidimensionnel:
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
Trouvez la différence de deux matrices en premier. Ensuite, appliquez la multiplication élément par élément avec la commande multiply de numpy. Ensuite, trouvez la somme de l’élément multiplié par la nouvelle matrice. Enfin, recherchez la racine carrée de la somme.
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