Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Pairwise operations (distance) on two lists in numpy

Tags:

python

numpy

I have two lists of coordinates:

l1 = [[x,y,z],[x,y,z],[x,y,z],[x,y,z],[x,y,z]]
l2 = [[x,y,z],[x,y,z],[x,y,z]]

I want to find the shortest pairwise distance between l1 and l2. Distance between two coordinates is simply:

numpy.linalg.norm(l1_element - l2_element)

So how do I use numpy to efficiently apply this operation to each pair of elements?

like image 559
TDN169 Avatar asked Apr 13 '15 15:04

TDN169


1 Answers

Here is a quick performance analysis of the four methods presented so far:

import numpy
import scipy
from itertools import product
from scipy.spatial.distance import cdist
from scipy.spatial import cKDTree as KDTree

n = 100
l1 = numpy.random.randint(0, 100, size=(n,3))
l2 = numpy.random.randint(0, 100, size=(n,3))

# by @Phillip
def a(l1,l2):
    return min(numpy.linalg.norm(l1_element - l2_element) for l1_element,l2_element in product(l1,l2))

# by @Kasra
def b(l1,l2):
    return numpy.min(numpy.apply_along_axis(
        numpy.linalg.norm,
        2,
        l1[:, None, :] - l2[None, :, :]
    ))

# mine
def c(l1,l2):
    return numpy.min(scipy.spatial.distance.cdist(l1,l2))

# just checking that numpy.min is indeed faster.
def c2(l1,l2):
    return min(scipy.spatial.distance.cdist(l1,l2).reshape(-1))

# by @BrianLarsen
def d(l1,l2):
    # make KDTrees for both sets of points
    t1 = KDTree(l1)
    t2 = KDTree(l2)
    # we need a distance to not look beyond, if you have real knowledge use it, otherwise guess
    maxD = numpy.linalg.norm(l1[0] - l2[0]) # this could be closest but anyhting further is certainly not
    # get a sparce matrix of all the distances

    ans = t1.sparse_distance_matrix(t2, maxD)

    # get the minimum distance and points involved
    minD = min(ans.values())
    return minD

for x in (a,b,c,c2,d):
    print("Timing variant", x.__name__, ':', flush=True)
    print(x(l1,l2), flush=True)
    %timeit x(l1,l2)
    print(flush=True)

For n=100

Timing variant a :
2.2360679775
10 loops, best of 3: 90.3 ms per loop

Timing variant b :
2.2360679775
10 loops, best of 3: 151 ms per loop

Timing variant c :
2.2360679775
10000 loops, best of 3: 136 µs per loop

Timing variant c2 :
2.2360679775
1000 loops, best of 3: 844 µs per loop

Timing variant d :
2.2360679775
100 loops, best of 3: 3.62 ms per loop

For n=1000

Timing variant a :
0.0
1 loops, best of 3: 9.16 s per loop

Timing variant b :
0.0
1 loops, best of 3: 14.9 s per loop

Timing variant c :
0.0
100 loops, best of 3: 11 ms per loop

Timing variant c2 :
0.0
10 loops, best of 3: 80.3 ms per loop

Timing variant d :
0.0
1 loops, best of 3: 933 ms per loop
like image 142
PeterE Avatar answered Sep 19 '22 00:09

PeterE