Как вернуть все пары точек, которые удовлетворяют критерию порога расстояния, используя матрицу расстояний haversine

У меня 39803 различных латинских балла (вот первые 5 баллов)

lat lon
0   27.269987   -82.497004
1   27.537598   -82.508422
2   27.337793   -82.533753
3   27.497719   -82.570261
4   27.512062   -82.722158

Что я рассчитал матрицу расстояния хаверсин для

array([[  0.        ,  29.77832527,   8.36846516, ...,  30.07072193,
      7.60700598,   7.63477669],
   [ 29.77832527,   0.        ,  22.35749757, ...,   2.6836159 ,
     22.17639199,  23.07090099],
   [  8.36846516,  22.35749757,   0.        , ...,  23.07825172,
      3.10333262,   0.75441483],
   ..., 
   [ 30.07072193,   2.6836159 ,  23.07825172, ...,   0.        ,
     22.53766911,  23.75965211],
   [  7.60700598,  22.17639199,   3.10333262, ...,  22.53766911,
      0.        ,   3.03795035],
   [  7.63477669,  23.07090099,   0.75441483, ...,  23.75965211,
      3.03795035,   0.        ]])

таким образом, моя матрица составляет 39 803 на 39 803. Я хотел бы найти способ найти все пары, которые были на расстоянии менее 25 метров. Например, если мы рассмотрим первый массив

[  0.        ,  29.77832527,   8.36846516, ...,  30.07072193,
  7.60700598,   7.63477669]

пара

(lat[0],lon[0])-(lat[2],lon[2])=(27.269987,-82.497004)-(27.337793,-82.533753) = 8.36846516

удовлетворяет этому критерию, но

(lat[0],lon[0])-(lat[1],lon[1])=(27.269987,-82.497004)-(27.537598,-82.508422) = 29.77832527

не. Я хотел бы получить подмножество баллов, которые удовлетворяют этому критерию. Это то, что я до сих пор:

X=df[['lat','lon']].dropna(axis=0) 
coors=np.array(X)

from math import radians, cos, sin, asin, sqrt

from scipy.spatial.distance import pdist, squareform
from sklearn.cluster import DBSCAN

import matplotlib.pyplot as plt
import pandas as pd


def haversine(lonlat1, lonlat2):

    lat1, lon1 = lonlat1
    lat2, lon2 = lonlat2
    lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2])

    # haversine formula 
    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2
    c = 2 * asin(sqrt(a)) 
    r = 6371 # Radius of earth in kilometers. Use 3956 for miles
    return c * r

distance_matrix = squareform(pdist(X, (lambda u,v: haversine(u,v))))

Я был бы очень признателен за любую помощь, как я мог бы потенциально вычислить это

1 ответ

Как насчет использования numpy.argwhere() функционировать? Это можно использовать для поиска индекса строки / столбца записей в матрице, которые удовлетворяют определенному критерию. В вашем случае это может быть что-то вроде:

numpy.argwhere(distances < 25)

Это вернет массив Nx2, дающий индексы строк / столбцов для N точек, которые разделены расстоянием менее 25 метров. Очевидно, что вам может потребоваться дополнительная фильтрация, чтобы избежать дублирования, связанного с вашей симметричной матрицей расстояний, и устранить нулевое расстояние между каждой точкой и самой собой. Возможно, вы могли бы справиться с этим, взяв матрицу расстояний и выполнив следующее преобразование:

distances = distances + 26 * (numpy.tril(numpy.ones_like(distances), k=1))

что заставляет все элементы на диагонали и в верхнем треугольнике матрицы расстояний иметь значение выше вашего 25-метрового порога.

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

indices = [ (i, j) for i in range(Npoints)
                     for j in range(i+1, Npoints)
                     if haversine(points[i], points[j]) < minDistance ]

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

Гибридный подход может включать в себя взятие подмножеств вашего списка точек, формирование матриц расстояний между каждой парой подмножеств и их фильтрацию с использованием приведенного выше numpy.argwhere() метод.

Другие вопросы по тегам