leonfrank
leonfrank

Reputation: 201

calculate distance of 2 list of points in numpy

I have 2 lists of points as numpy.ndarray, each row is the coordinate of a point, like:

a = np.array([[1,0,0],[0,1,0],[0,0,1]])
b = np.array([[1,1,0],[0,1,1],[1,0,1]])

Here I want to calculate the euclidean distance between all pairs of points in the 2 lists, for each point p_a in a, I want to calculate the distance between it and every point p_b in b. So the result is

d = np.array([[1,sqrt(3),1],[1,1,sqrt(3)],[sqrt(3),1,1]])

How to use matrix multiplication in numpy to compute the distance matrix?

Upvotes: 7

Views: 15333

Answers (4)

haofeng
haofeng

Reputation: 651

L2 distance = (a^2 + b^2 - 2ab)^0.5

a = np.random.randn(5, 3)
b = np.random.randn(2, 3)
a2 = np.sum(np.square(a), axis = 1)[..., None]
b2 = np.sum(np.square(b), axis = 1)[None, ...]
ab = -2*np.dot(a, b.T)

dist = np.sqrt(a2 + b2 + ab)

Upvotes: 0

bfree67
bfree67

Reputation: 725

If you have 2 each 1-dimensional arrays, x and y, you can convert the arrays into matrices with repeating columns, transpose, and apply the distance formula. This assumes that x and y are coordinated pairs. The result is a symmetrical distance matrix.

x = [1, 2, 3]
y = [4, 5, 6]
xx = np.repeat(x,3,axis = 0).reshape(3,3)
yy = np.repeat(y,3,axis = 0).reshape(3,3)
dist = np.sqrt((xx-xx.T)**2 + (yy-yy.T)**2)


dist
Out[135]: 
array([[0.        , 1.41421356, 2.82842712],
       [1.41421356, 0.        , 1.41421356],
       [2.82842712, 1.41421356, 0.        ]])

Upvotes: 0

Divakar
Divakar

Reputation: 221504

To compute the squared euclidean distance for each pair of elements off them - x and y, we need to find :

(Xik-Yjk)**2 = Xik**2 + Yjk**2 - 2*Xik*Yjk

and then sum along k to get the distance at coressponding point as dist(Xi,Yj).

Using associativity, it reduces to :

dist(Xi,Yj) = sum_k(Xik**2) + sum_k(Yjk**2) - 2*sum_k(Xik*Yjk)

Bringing in matrix-multiplication for the last part, we would have all the distances, like so -

dist = sum_rows(X^2), sum_rows(Y^2), -2*matrix_multiplication(X, Y.T)

Hence, putting into NumPy terms, we would end up with the euclidean distances for our case with a and b as the inputs, like so -

np.sqrt((a**2).sum(1)[:,None] + (b**2).sum(1) - 2*a.dot(b.T))

Leveraging np.einsum, we could replace the first two summation-reductions with -

np.einsum('ij,ij->i',a,a)[:,None] + np.einsum('ij,ij->i',b,b) 

More info could be found on eucl_dist package's wiki page (disclaimer: I am its author).

Upvotes: 5

jakevdp
jakevdp

Reputation: 86310

Using direct numpy broadcasting, you can do this:

dist = np.sqrt(((a[:, None] - b[:, :, None]) ** 2).sum(0))

Alternatively, scipy has a routine that will compute this slightly more efficiently (particularly for large matrices)

from scipy.spatial.distance import cdist
dist = cdist(a, b)

I would avoid solutions that depend on factoring-out matrix products (of the form A^2 + B^2 - 2AB), because they can be numerically unstable due to floating point roundoff errors.

Upvotes: 13

Related Questions