7

Есть 733 географических точек F (массив 2*733)
Есть 22 000 географических точек S (массив 2*22 000)

Как быстро посчитать расстояние от каждой точки из F до каждой точки S.

Этот код работает очень долго, около семи минут. Почему?

import scipy.spatial
D = scipy.spatial.distance.cdist(F,S,lambda u, v: haversine(u, v))

А если вычислять евклидову метрику, то быстро.

import scipy.spatial
D2 = scipy.spatial.distance.cdist(F,S,'euclidean')

Вычислять расстояние функцией тоже долго:

import scipy.spatial
D = scipy.spatial.distance.cdist(F,S,lambda u,v:geopy.distance.vincenty(u,v).km) 

Но дело не в сложности функции вычисления расстояния, т.к даже если я буду считать сумму разности координат вычисления займут много времени:

import scipy.spatial
D = scipy.spatial.distance.cdist(F,S,lambda  u,v:(u-v).sum())
5
  • Вы каким haversine пользуетесь? 16 ноя 2017 в 10:54
  • судя по докам - коду нужно просчитать 1.6 миллиона пар координат. Поэтому 7 минут вполне может быть ожидаемым временем. Но если это нужно считать постоянно, может сделать один раз и закешировать?
    – KoVadim
    16 ноя 2017 в 11:04
  • 2
    С scipy не работал, но, насколько я знаю, у него под капотом скомпилированное с++, да ещё заоптимизированное на работу с большими массивами данных. Очевидно, что питоновская лямбда, вызываямая отдельно для каждой пары будет ГОРАЗДО медленнее.
    – Xander
    16 ноя 2017 в 11:42
  • @MaxU haversine 0.4.5 : github.com/mapado/haversine 16 ноя 2017 в 12:26

1 ответ 1

8

Вот рабочий пример.

Данные для примера: http://samplecsvs.s3.amazonaws.com/SacramentocrimeJanuary2006.csv

import numpy as np
import pandas as pd
from sklearn.neighbors import DistanceMetric

# http://scikit-learn.org/stable/modules/generated/sklearn.neighbors.DistanceMetric.html
dist = DistanceMetric.get_metric('haversine')

url = 'http://samplecsvs.s3.amazonaws.com/SacramentocrimeJanuary2006.csv'
df = pd.read_csv(url, usecols=['latitude','longitude'])

первые 5 записей:

In [213]: df.head()
Out[213]:
    latitude   longitude
0  38.550420 -121.391416
1  38.473501 -121.490186
2  38.657846 -121.462101
3  38.506774 -121.426951
4  38.637448 -121.384613

увеличим DataFrame в три раза, чтобы получить примерно такой же объем данных как у вас в вопросе:

df = pd.concat([df] * 3, ignore_index=True).sample(frac=1)
# split our data set into two
F, S = np.split(df, [int(len(df)*0.033)])

получилось два DataFrame'а:

print('shape of [F]: {}'.format(F.shape))
print('shape of [S]: {}'.format(S.shape))

выдает:

shape of [F]: (750, 2)
shape of [S]: (22002, 2)

считаем расстояние:

earth_radius = 6371

D = dist.pairwise(np.radians(F), np.radians(S)) * earth_radius

PS haversine считается по формуле:

2 * arcsin(sqrt(sin^2(0.5*dx) * cos(x1) * cos(x2) * sin^2(0.5*dy)))

и возвращает результат в радианах, который надо умножить на радиус Земли (в км.), чтобы получить расстояние в километрах

Результат:

In [208]: D
Out[208]:
array([[  3.86258192,  12.14112631,  16.1702106 , ...,  11.2399247 ,  13.71076021,   2.59575523],
       [  3.19765341,  15.72580688,  20.35963524, ...,  15.38247654,  17.85975202,   4.60552864],
       [  6.15251463,  19.94251384,  19.77488568, ...,  17.0637875 ,  19.15736473,   6.55354701],
       ...,
       [  7.85715424,   9.52674314,  11.13189495, ...,   6.81544901,   9.16104377,   5.58452969],
       [ 14.02227805,   5.08370143,   6.47846773, ...,   0.57188025,   3.04569744,  11.81577747],
       [ 17.84804533,   3.93791815,   6.91865107, ...,   3.8743883 ,   2.94376778,  15.72668196]])

In [209]: D.shape
Out[209]: (750, 22002)

замер скорости - около 6.5 секунд (для 16.5 миллионов пар точек) на моем стареньком ноуте:

In [210]: %timeit dist.pairwise(np.radians(F), np.radians(S)) * earth_radius
6.45 s ± 156 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

Проверка формулы:

In [204]: lyon = np.array([[45.7597, 4.8422]])

In [205]: paris = np.array([[48.8567, 2.3508]])

In [206]: dist.pairwise(np.radians(lyon), np.radians(paris))  * earth_radius
Out[206]: array([[ 392.21671781]])
1

Ваш ответ

By clicking “Отправить ответ”, you agree to our terms of service and acknowledge you have read our privacy policy.

Всё ещё ищете ответ? Посмотрите другие вопросы с метками или задайте свой вопрос.