матрица расстояний кривых в питоне

У меня есть набор кривых, определенных как 2D-массивы (количество точек, количество координат). Я рассчитываю для них матрицу расстояний, используя расстояние Хаусдорфа. Мой текущий код выглядит следующим образом. К сожалению, это слишком медленно с 500-600 кривыми, каждая из которых имеет 50-100 3D-точек. Есть ли более быстрый способ для этого?

def distanceBetweenCurves(C1, C2):
    D = scipy.spatial.distance.cdist(C1, C2, 'euclidean')

    #none symmetric Hausdorff distances
    H1 = np.max(np.min(D, axis=1))
    H2 = np.max(np.min(D, axis=0))

    return (H1 + H2) / 2.

def distanceMatrixOfCurves(Curves):
    numC = len(Curves)

    D = np.zeros((numC, numC))
    for i in range(0, numC-1):
        for j in range(i+1, numC):
            D[i, j] = D[j, i] = distanceBetweenCurves(Curves[i], Curves[j])

    return D

python numpy scipy distance curves
person ahmethungari    schedule 03.12.2012    source источник
comment
scipy.spatial.distance.cdist медленная часть или двойная петля внутри distanceMatrixOfCurves? Если эти кривые выпуклые, можно было бы оптимизировать первую возможную медленную часть. Пересекаются ли эти кривые или содержатся внутри других? Я чувствую, что вы могли бы повторно использовать ранее найденные расстояния, чтобы ускорить новые вычисления. Конечно, это просто болтовня, у меня есть аналогичная проблема с минимальными (минимальными (..)) мерами, и мне было трудно обобщить эти соображения, которые я здесь излагаю. Вы пробовали или думали о чем-то помимо кода, который вы разместили?   -  person mmgp    schedule 04.12.2012
comment
Я попытался реализовать евклидово расстояние самостоятельно (вместо использования cdist), ничего особенного не изменилось. Я думаю, что проблема в двойной петле. Кривые (некоторые из них) пересекаются и содержатся внутри других...   -  person ahmethungari    schedule 04.12.2012
comment
@ahmethungari, вы должны профилировать свой код, чтобы быть уверенным (cProfile + runnakerun отлично интерпретировать результаты) каково конкретное узкое место. Я не очень разбираюсь в этих вещах, но вы можете обойтись без выделения большой матрицы расстояний, состоящей из всех пар, которую вычисляет cdist — если вы добавите код, который генерирует небольшой пример данных, это будет проще. чтобы помочь вам.   -  person YXD    schedule 04.12.2012
comment
@MrE, хотя, безусловно, можно избавиться от cdist (просто используйте форму распространения волны), что может потребовать больше памяти, а в худшем случае это не поможет. Это не означает, что это не может помочь в сокращении времени выполнения, но это проблематично.   -  person mmgp    schedule 04.12.2012
comment
Конечно, я не думаю, что есть простое решение. Мое единственное другое предложение состояло бы в том, чтобы хранить точки в пространственной структуре данных, такой как kd-дерево для каждой кривой, что, по крайней мере, ускорит поиск ближайшего соседа к данной точке, но в целом может оказаться медленнее. Интересно посмотреть, что люди предлагают здесь.   -  person YXD    schedule 04.12.2012
comment
1. Вам действительно нужна полная матрица D или можно обойтись только верхней или нижней треугольной матрицей? Этот D[i,j] = D[j,i] =... определенно не годится для локальности данных; 2. Пробовали ли вы использовать сжатие списка или map вместо двойных циклов?   -  person ev-br    schedule 04.12.2012
comment
Это никоим образом не снижает вычислительную сложность. В нынешнем виде я думаю, что ОП просто получит незначительные улучшения реализации, т.е. использует библиотеку/язык/пакет X, потому что он запускает Y (где Y — тот же метод, возможно, с некоторыми незначительными штрихами) раз быстрее! и т. д. Надеюсь, я ошибся.   -  person mmgp    schedule 04.12.2012
comment
@Женя, мне не нужна вся матрица, ты прав. Как вы думаете, понимание списка быстрее, чем циклы? Я могу попробовать это...   -  person ahmethungari    schedule 05.12.2012

Ответы (3)


arrow_upward
6
arrow_downward

Ваш вопрос также может быть связан с этим

Это какая-то трудная проблема. Возможным способом было бы реализовать евклидово расстояние самостоятельно, полностью отказаться от scipy и использовать JIT pypy. компилятор. Но, скорее всего, это не принесет вам много пользы.

Лично я бы порекомендовал вам написать подпрограмму на C.

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

  • Предположим, у вас есть 500 кривых, каждая из которых имеет 75 точек. При подходе грубой силы вы в конечном итоге вычисляете евклидово расстояние 500 * 499 * 75 * 75 = 1 403 437 500 раз. Неудивительно, что этот подход работает вечно.

Я не эксперт в этом, но я знаю, что расстояние Хаусдорфа широко используется при обработке изображений. Я бы посоветовал вам просмотреть литературу по алгоритмам, оптимизированным для скорости. Отправной точкой может быть это или этот документ. Кроме того, в сочетании с расстоянием Хаусдорфа часто упоминается диаграмма Ворони.

Я надеюсь, что эти ссылки могут помочь вам с этой проблемой.

person jojo    schedule 04.12.2012
comment
Первая связанная статья интересна, решает часть проблемы. Теперь другую проблему можно более или менее сформулировать так: для трех множеств ABC, если мы знаем расстояние Хаусдорфа между A и B, а также внутренние детали, полученные в результате такого вычисления, возможно ли найти расстояние Хаусдорфа от A до C с меньшей вычислительной сложностью, чем это было для вычислений для A и B, если есть пути, учитывая некоторую метрику, из A в C, которые пересекают B ? Ну, это было намного дольше, чем я ожидал, когда начал писать. Я надеюсь, что кто-то может это понять. - person mmgp; 04.12.2012
comment
спасибо за все ваши советы. У меня появилась идея: использовать более умный алгоритм. Я надеюсь, что сделаю это после того, как использую текущую версию для бумажного дедлайна :) - person ahmethungari; 05.12.2012
comment
@ahmethungari Мне жаль, что я не могу помочь вам больше, но у этой проблемы нет простого решения. Однако у меня есть одна идея для ускорения вычислений: вместо вычисления для заданной точки в Curve1 euclid.dist. со всеми точками из Curve2 вы можете случайным образом выбрать точку из Curve2, вычислить euclid.dist. (назовем ее r), то для каждой другой точки из Curve2 можно вычислить манхэттонское расстояние (много быстрее вычислять) и проверьте, соответствует ли каждое измерение манхэттонского расстояния. больше, чем r, если так... - person jojo; 05.12.2012
comment
@ahmethungari произвольно выбранная точка ближе, и вам не нужно вычислять euclid.dist. (если вы не уверены, почему: треугольное неравенство). если нет, вычислите euclid.dist. и если он меньше r, обновить r. Вы по-прежнему будете проходить через все точки в Curve2, но часто вам нужно будет просто вычислить манхэттонское расстояние, что быстрее. ...ну, это всего лишь идея, и наверняка есть более изощренные способы ускорить этот процесс. - person jojo; 05.12.2012
comment
@ahmethungari рассмотрите возможность принятия этого ответа, если он вас устраивает. :) - person jojo; 07.08.2014

arrow_upward
3
arrow_downward

Недавно я ответил на аналогичный вопрос здесь: Расстояние Хаусдорфа между трехмерными сетками

Надеюсь, это поможет, я столкнулся с 25 х 25 000 точек в попарном сравнении (всего 25 х 25 х 25 000 точек), и мой код работает от 1 минуты до 3-4 часов (в зависимости от количества точек). Я не вижу много вариантов математически для увеличения скорости.

Альтернативой может быть использование разных языков программирования (C/C++) или перенос этого расчета на GPU (CUDA). Я играю с подходом CUDA прямо сейчас.

Изменить от 12.03.2015:

Я смог ускорить это сравнение, выполнив параллельные вычисления на базе ЦП. Это был самый быстрый путь. Я использовал хороший пример пакета pp (parallel python) и запускал его на трех разных компьютерах и комбинации phython. К сожалению, у меня все время были ошибки памяти с 32-разрядной версией python 2.7, поэтому я установил 64-разрядную версию WinPython 2.7 и несколько экспериментальных 64-разрядных пакетов numpy.

введите здесь описание изображения

Так что для меня это усилие было весьма полезным, и для меня это было не так сложно, как CUDA.... Удачи

person Akos Gulyban    schedule 26.11.2015

arrow_upward
0
arrow_downward

Есть несколько методов, которые вы можете попробовать:

  1. Использование numpy-MKL, в котором вместо numpy используется высокопроизводительная библиотека Math Kernel от Intel;
  2. Использование Bootleneck для функций массива;
  3. Использование Cpython для вычислений.
person foool    schedule 04.12.2012