Эффективное вычисление расстояний от точек до линии в векторной форме

Question

Вычисление расстояния от точки до линии в векторной форме является тривиальным.

Однако я реализовал его следующим образом, который очень медленный:

def compute_point_distance_to_line(point):
    dist = np.linalg.norm((a - point) - np.vdot((a - point), n) * n)
    return dist
np.apply_along_axis(compute_point_distance_to_line, 1, xyz)

Я использовал обозначения из википедии, форма xyz (2521909, 3), форма a, n и точка, следовательно, (3,)

Я пробовал это следующим образом:

def compute_point_distance_to_line2(points):
    _a = np.tile(a, (points.shape[0], 1))
    _n = np.tile(n, (points.shape[0], 1))
    _n_t = np.ascontiguousarray(np.swapaxes(_n, 0, 1))

    diffs = _a - points
    vdots_scaled = np.dot(diffs, _n_t) * n
    diffs = diffs - vdots_scaled

    return np.linalg.norm(diffs, axis=1)

К сожалению, для меня это приводит к ошибке памяти при вычислении скалярного произведения.

Есть ли более быстрые способы?

Answer 1

В общем, старайтесь не использовать такие операции, как tile, когда вместо этого вы можете использовать широковещательную рассылку, особенно когда скорость / память являются проблемой.

https://docs.scipy.org/doc/numpy-1.15.0/user/basics.broadcasting.html

С небольшой алгебройВы можете записать это как единое матричное векторное произведение, за которым следует норма.Это может помочь вам избежать временных переменных и сэкономить память.

Вот рабочий пример.Обратите внимание, что в этом примере все трехмерные векторы являются векторами столбцов, поэтому p равно 3x1000 вместо 1000x3.Вам нужно будет транспонировать p, чтобы подключить его к этому примеру.

import numpy as np
# define an example line with unit n
a = np.array([1,2,3])
n = np.array([4,5,6])
norm2n = np.sum(n**2)
n = n/np.sqrt(norm2n)
# get some point data p
p = np.random.randn(3,1000)
# form the projection matrix (see use of None in broadcasting at link above)
P = np.eye(3) - n[:,None]*n[None,:]
# perform the projection using matrix multiplication    
projected = P.dot(a[:,None]-p)
# get the distance
dist = np.sqrt(np.sum(projected**2, axis=0))

Answer 2

Вы можете векторизовать это следующим образом:

temp = np.subtract(a, xyz)  # so we only have to compute this once
dist = np.linalg.norm(np.subtract(temp, np.multiply(np.dot(temp, n)[:, None], n)),
                      axis=-1)
# 220 ms ± 6 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

По сравнению с синхронизацией для вашего первого примера кода выше:

# 30 s ± 1.89 s per loop (mean ± std. dev. of 7 runs, 1 loop each)

Это дает мне тот же результат, что и ваш первый пример кода(проверено с помощью np.array_equal()), и это на пару порядков быстрее.

---

Пояснение

Хитрость заключается в том, чтобы заставить вызов np.multiply() работать правильно, добавивдополнительная ось к результату np.dot(), который я делаю со срезом [:, None] после np.dot().В основном None, используемый в кусочках срезки, является ярлыком для добавления оси, поэтому результат np.dot() для вас должен иметь форму (2521909,), а после скобок с None он будет иметь форму (2521909, 1).Результат np.multiply() (и temp) будет иметь форму (2521909, 3), и мы возьмем норму вдоль последней оси, чтобы получить трехмерное расстояние от линии до каждой из ваших 2521909 точек.