Jak obliczyć odległość euklidesową za pomocą NumPy?

Mam dwa punkty w 3D:

(xa, ya, za)
(xb, yb, zb)

I chcę obliczyć odległość:

dist = sqrt((xa-xb)^2 + (ya-yb)^2 + (za-zb)^2)

Jaki jest najlepszy sposób, aby to zrobić za pomocą NumPy lub ogólnie Pythona? Mam:

import numpy
a = numpy.array((xa ,ya, za))
b = numpy.array((xb, yb, zb))

Użyj numpy.linalg.norm:

dist = numpy.linalg.norm(a-b)

Teorię tę można znaleźć we wstępie do eksploracji danych

Działa to, ponieważ odległość euklidesowa jest normą l2, a domyślną wartością parametru ord w numpy.linalg.norm jest 2.

Jest taka funkcja w SciPy. To się nazywa euklidesowe .

Przykład:

from scipy.spatial import distance
a = (1, 2, 3)
b = (4, 5, 6)
dst = distance.euclidean(a, b)

Dla wszystkich zainteresowanych obliczeniem wielu odległości naraz zrobiłem małe porównanie przy użyciu perfplot ( mały projekt).

Pierwsza rada to takie uporządkowanie danych, aby tablice miały wymiar (3, n)(i oczywiście są ciągłe w C). Jeśli dodawanie odbywa się w ciągłym pierwszym wymiarze, rzeczy są szybsze i nie ma to większego znaczenia, jeśli używasz sqrt-sumz axis=0, linalg.normz axis=0lub

a_min_b = a - b
numpy.sqrt(numpy.einsum('ij,ij->j', a_min_b, a_min_b))

który jest, z niewielkim marginesem, najszybszym wariantem. (Tak naprawdę dotyczy to tylko jednego rzędu).

Warianty, w których sumuje się na drugiej osi axis=1, są znacznie wolniejsze.

Kod do odtworzenia fabuły:

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


perfplot.save(
    "norm.png",
    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,
    ],
    logx=True,
    logy=True,
    xlabel="len(x), len(y)",
)

Chcę wyjaśnić prostą odpowiedź z różnymi notatkami dotyczącymi wydajności. np.linalg.norm zrobi może więcej niż potrzebujesz:

dist = numpy.linalg.norm(a-b)

Po pierwsze - ta funkcja została zaprojektowana do pracy nad listą i zwracania wszystkich wartości, np. W celu porównania odległości od pAzestawu punktów sP:

sP = set(points)
pA = point
distances = np.linalg.norm(sP - pA, ord=2, axis=1.)  # 'distances' is a list

Pamiętaj o kilku rzeczach:

• Wywołania funkcji w języku Python są drogie.
• [Zwykły] Python nie buforuje wyszukiwania nazw.

Więc

def distance(pointA, pointB):
    dist = np.linalg.norm(pointA - pointB)
    return dist

nie jest tak niewinny jak się wydaje.

>>> 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

Po pierwsze - za każdym razem, gdy go nazywamy, musimy przeprowadzić globalne wyszukiwanie dla „np”, wyszukiwanie w zakresie dla „linalg” i wyszukiwanie w zakresie dla „normy”, a narzut związany z samym wywoływaniem funkcji może być równy kilkudziesięciu pytonom instrukcje.

Na koniec zmarnowaliśmy dwie operacje, aby zapisać wynik i załadować go ponownie w celu zwrotu ...

Pierwszy krok do poprawy: przyspiesz wyszukiwanie, pomiń sklep

def distance(pointA, pointB, _norm=np.linalg.norm):
    return _norm(pointA - pointB)

Otrzymujemy znacznie bardziej usprawnione:

>>> 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

Narzut związany z wywoływaniem funkcji nadal jednak wymaga pewnej pracy. I będziesz chciał przeprowadzić testy porównawcze, aby ustalić, czy lepiej byłoby, gdybyś sam wykonał matematykę:

def distance(pointA, pointB):
    return (
        ((pointA.x - pointB.x) ** 2) +
        ((pointA.y - pointB.y) ** 2) +
        ((pointA.z - pointB.z) ** 2)
    ) ** 0.5  # fast sqrt

Na niektórych platformach **0.5jest szybszy niż math.sqrt. Twój przebieg może się różnić.

**** Zaawansowane uwagi dotyczące wydajności.

Dlaczego obliczasz odległość? Jeśli jedynym celem jest wyświetlenie go,

 print("The target is %.2fm away" % (distance(a, b)))

poruszać się. Ale jeśli porównujesz odległości, sprawdzasz zasięg itp., Chciałbym dodać kilka przydatnych obserwacji wydajności.

Weźmy dwa przypadki: sortowanie według odległości lub ubijanie listy do elementów spełniających ograniczenie zakresu.

# 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)

Pierwszą rzeczą, o której musimy pamiętać, jest to, że używamy Pitagorasa do obliczania odległości ( dist = sqrt(x^2 + y^2 + z^2)), więc wykonujemy wiele sqrtpołączeń. Matematyka 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

Krótko mówiąc: dopóki faktycznie nie będziemy potrzebować odległości w jednostce X zamiast X ^ 2, możemy wyeliminować najtrudniejszą część obliczeń.

# 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)

Świetnie, obie funkcje nie powodują już żadnych drogich pierwiastków kwadratowych. To będzie znacznie szybsze. Możemy również ulepszyć in_range, przekształcając go w generator:

def in_range(origin, range, things):
    range_sq = range**2
    yield from (thing for thing in things
                if distance_sq(origin, thing) <= range_sq)

Ma to szczególne zalety, jeśli robisz coś takiego:

if any(in_range(origin, max_dist, things)):
    ...

Ale jeśli następna rzecz, którą zamierzasz zrobić, wymaga dystansu,

for nearby in in_range(origin, walking_distance, hotdog_stands):
    print("%s %.2fm" % (nearby.name, distance(origin, nearby)))

rozważ uzyskanie krotek:

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)

Może to być szczególnie przydatne, jeśli możesz połączyć sprawdzanie zasięgu („znajdź rzeczy, które są w pobliżu X i w Nm od Y”, ponieważ nie musisz ponownie obliczać odległości).

Ale co, jeśli szukamy naprawdę dużej listy thingsi spodziewamy się, że wiele z nich nie będzie wartych rozważenia?

W rzeczywistości istnieje bardzo prosta optymalizacja:

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

To, czy jest to przydatne, będzie zależeć od wielkości „rzeczy”.

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 ...

I jeszcze raz zastanów się nad uzyskaniem dist_sq. Nasz przykład Hot-Dog to:

# 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))

Kolejne wystąpienie tej metody rozwiązywania problemów :

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)

Zaczynając Python 3.8, mathmoduł bezpośrednio udostępnia distfunkcję, która zwraca odległość euklidesową między dwoma punktami (podanymi jako krotki lub listy współrzędnych):

from math import dist

dist((1, 2, 6), (-2, 3, 2)) # 5.0990195135927845

A jeśli pracujesz z listami:

dist([1, 2, 6], [-2, 3, 2]) # 5.0990195135927845

Można to zrobić w następujący sposób. Nie wiem, jak szybko to jest, ale nie używa 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)))

Znajduję funkcję „dist” w matplotlib.mlab, ale nie sądzę, aby była wystarczająco przydatna.

Zamieszczam go tutaj tylko w celach informacyjnych.

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)

Lubię np.dot(produkt kropka):

a = numpy.array((xa,ya,za))
b = numpy.array((xb,yb,zb))

distance = (np.dot(a-b,a-b))**.5

Miły jednowarstwowy:

dist = numpy.linalg.norm(a-b)

Jeśli jednak chodzi o szybkość, zalecam eksperymentowanie na komputerze. Przekonałem się, że używając mathbiblioteki sqrtz** operatorem kwadratu jest znacznie szybsze na moim komputerze niż jedno-liniowe rozwiązanie NumPy.

Przeprowadziłem testy przy użyciu tego prostego programu:

#!/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

Na mojej maszynie math_calc_distdziała znacznie szybciej niżnumpy_calc_dist : 1,5 sekundy w porównaniu do 23,5 sekundy.

Aby uzyskać mierzalną różnicę między fastest_calc_disti math_calc_distmusiałem TOTAL_LOCATIONSdo 6000. Potem fastest_calc_distzajmuje ~ 50 sekundmath_calc_dist zajmuje ~ 60 sekund.

Możesz także eksperymentować z numpy.sqrti numpy.squarechoć oba były wolniejsze niżmath alternatywy na mojej maszynie.

Moje testy zostały uruchomione w języku Python 2.6.6.

Możesz po prostu odjąć wektory, a następnie produkt wewnętrzny.

Idąc za twoim przykładem,

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)

Posiadanie ai bjak je określono, można użyć również:

distance = np.sqrt(np.sum((a-b)**2))

Z Python 3.8 jest to bardzo łatwe.

https://docs.python.org/3/library/math.html#math.dist

math.dist(p, q)

Zwraca odległość euklidesową między dwoma punktami p i q, każdy podany jako sekwencja (lub iterowalna) współrzędnych. Dwa punkty muszą mieć ten sam wymiar.

Z grubsza odpowiada:

sqrt(sum((px - qx) ** 2.0 for px, qx in zip(p, q)))

Oto zwięzły kod odległości euklidesowej w Pythonie, biorąc pod uwagę dwa punkty reprezentowane jako listy w Pythonie.

def distance(v1,v2): 
    return sum([(x-y)**2 for (x,y) in zip(v1,v2)])**(0.5)

Od wersji Python 3.8

Od wersji Python 3.8 mathmoduł zawiera funkcję math.dist().
Zobacz tutaj https://docs.python.org/3.8/library/math.html#math.dist .

math.dist (p1, p2)
Zwraca odległość euklidesową między dwoma punktami p1 i p2, każdy podany jako ciąg (lub iterowalny) współrzędnych.

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

Oblicz odległość euklidesową dla przestrzeni wielowymiarowej:

 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

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)

Możesz łatwo użyć formuły

distance = np.sqrt(np.sum(np.square(a-b)))

co właściwie nie robi nic więcej niż użycie twierdzenia Pitagorasa do obliczenia odległości, poprzez dodanie kwadratów ofx, yy i zz i zrootowanie wyniku.

Najpierw znajdź różnicę dwóch macierzy. Następnie zastosuj mnożenie elementu za pomocą polecenia mnożenia numpy. Następnie znajdź sumę elementu pod kątem pomnożonej nowej macierzy. Na koniec znajdź pierwiastek kwadratowy sumy.

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

import numpy as np
# any two python array as two points
a = [0, 0]
b = [3, 4]

Państwo listę pierwsza zmiana na numpy tablicy i zrobić tak: print(np.linalg.norm(np.array(a) - np.array(b))). Druga metoda bezpośrednio z listy python jako:print(np.linalg.norm(np.subtract(a,b)))