Необходимо вычислить расстояния между каждой из 60 000 координат

Для этого вам на самом деле не нужно использовать функцию расстояния geopandas. По сути, единственное, что вам нужно, это scipy. Начните с помещения ваших координат в массив:

import numpy as np
import pandas as pd
import geojson
from scipy.spatial import distance

with open("/Users/katcn/Desktop/Spring 2024/Research/PleaseWork.geojson") as f:
    gj = geojson.load(f)

coords = np.array([feature["geometry"]["coordinates"] for feature in gj['features']])
dist_matrix = distance.cdist(coords, coords, 'euclidean')
dist_df = pd.DataFrame(dist_matrix)
print(dist_df)

Поскольку вы не предоставили данные, я создал образец набора данных с точками континентальной части США (здесь 60 000 точек от KKey West, Флорида до канадской границы и от Западного побережья до Восточного побережья):

import numpy as np
import pandas as pd
from scipy.spatial import distance
import time

np.random.seed(42)


latitudes = np.random.uniform(low=25.0, high=49.0, size=60000)  
longitudes = np.random.uniform(low=-125.0, high=-66.0, size=60000) 
coords = np.column_stack((latitudes, longitudes))

start_time = time.time()

dist_matrix = distance.cdist(coords, coords, 'euclidean')

end_time = time.time()
elapsed_time = end_time - start_time

print("Distance matrix computed in {:.2f} seconds".format(elapsed_time))
print("Shape of the distance matrix:", dist_matrix.shape)

Какие взять

Distance matrix computed in 185.64 seconds
Shape of the distance matrix: (60000, 60000)

Итак, примерно 3 минуты.

ОБНОВЛЕНИЕ: ИСПОЛЬЗОВАНИЕ ГЕОДЕЗИЧЕСКОГО РАССТОЯНИЯ

Это работает, но требует времени:

import numpy as np
import pandas as pd
from geopy.distance import geodesic
import time

np.random.seed(42)

latitudes = np.random.uniform(low=25.0, high=49.0, size=60000)
longitudes = np.random.uniform(low=-125.0, high=-66.0, size=60000)

coords = list(zip(latitudes, longitudes))

def geodesic_distance_matrix(coords):
    n = len(coords)
    dist_matrix = np.zeros((n, n))
    for i in range(n):
        for j in range(i+1, n): 
            dist_matrix[i, j] = geodesic(coords[i], coords[j]).meters
            dist_matrix[j, i] = dist_matrix[i, j]
    return dist_matrix

start_time = time.time()

dist_matrix = geodesic_distance_matrix(coords)

end_time = time.time()
elapsed_time = end_time - start_time

print("Distance matrix computed in {:.2f} seconds".format(elapsed_time))
print("Shape of the distance matrix:", dist_matrix.shape)

Другой подход — использовать haversine следующим образом:

import numpy as np
import pandas as pd
import time

def haversine(lat1, lon1, lat2, lon2):
    lon1, lat1, lon2, lat2 = map(np.radians, [lon1, lat1, lon2, lat2])

    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = np.sin(dlat/2)**2 + np.cos(lat1) * np.cos(lat2) * np.sin(dlon/2)**2
    c = 2 * np.arcsin(np.sqrt(a)) 
    r = 6371  
    return c * r

np.random.seed(42)
latitudes = np.random.uniform(low=25.0, high=49.0, size=60000)
longitudes = np.random.uniform(low=-125.0, high=-66.0, size=60000)
coords = np.column_stack((latitudes, longitudes))

start_time = time.time()

n = coords.shape[0]
dist_matrix = np.zeros((n, n))
for i in range(n):
    dist_matrix[i, :] = haversine(coords[i, 0], coords[i, 1], coords[:, 0], coords[:, 1])

end_time = time.time()
elapsed_time = end_time - start_time

print("Distance matrix computed in {:.2f} seconds".format(elapsed_time))
print("Shape of the distance matrix:", dist_matrix.shape)

который возвращает

Distance matrix computed in 463.75 seconds
Shape of the distance matrix: (60000, 60000)

то есть более чем в два раза больше времени евклидова расстояния.

array([[   0.        , 1683.78776508, 1717.22879531, ..., 1706.77761313,
        3857.52001651, 1579.36733204],
       [1683.78776508,    0.        , 2008.58039812, ..., 2813.5398641 ,
        4379.29512353, 2340.45973263],
       [1717.22879531, 2008.58039812,    0.        , ..., 1151.61134582,
        2399.40193517,  560.14519709],
       ...,
       [1706.77761313, 2813.5398641 , 1151.61134582, ...,    0.        ,
        2244.42005933,  593.632848  ],
       [3857.52001651, 4379.29512353, 2399.40193517, ..., 2244.42005933,
           0.        , 2291.13872501],
       [1579.36733204, 2340.45973263,  560.14519709, ...,  593.632848  ,
        2291.13872501,    0.        ]])

здесь размерность составляет километры. Если вам нужны мили, вам придется поменять

 r = 6371 

к

r = 3956

Этот ответ предназначен для решения следующей проблемы, а именно сохранения полученных результатов в CSV-файл. Учитывая объем данных, pd.to_csv() — не лучший вариант.

Здесь я опираюсь на предыдущее решение, вычисляя матрицу расстояний. Обратите внимание, что матрица на самом деле разрежена, поскольку расстояние между x и y такое же, как расстояние между y и x (поэтому матрица треугольная.

Я создаю фрейм данных, добавляя точки (lat, long). Затем, используя polar, я преобразую фрейм данных pandas в полярный фрейм данных и, наконец, сохраняю его:

import numpy as np
import pandas as pd
import time

def haversine(lat1, lon1, lat2, lon2):
    lon1, lat1, lon2, lat2 = map(np.radians, [lon1, lat1, lon2, lat2])
    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = np.sin(dlat/2)**2 + np.cos(lat1) * np.cos(lat2) * np.sin(dlon/2)**2
    c = 2 * np.arcsin(np.sqrt(a)) 
    r = 6371  
    return c * r

np.random.seed(42)
latitudes = np.random.uniform(low=25.0, high=49.0, size=10000)
longitudes = np.random.uniform(low=-125.0, high=-66.0, size=10000)
coords = np.column_stack((latitudes, longitudes))

start_time = time.time()

n = coords.shape[0]
dist_matrix = np.zeros((n, n))
for i in range(n):
    dist_matrix[i, :] = haversine(coords[i, 0], coords[i, 1], coords[:, 0], coords[:, 1])

end_time = time.time()
elapsed_time = end_time - start_time
print("Distance matrix computed in {:.2f} seconds".format(elapsed_time))
print("Shape of the distance matrix:", dist_matrix.shape)

point_labels = [f"Point {i}" for i in range(n)]
dist_df = pd.DataFrame(dist_matrix, index=point_labels, columns=point_labels)


distance = dist_df.loc['Point 0', 'Point 1']
print(f"Distance from Point 0 to Point 1: {distance} km")

lat_lon_index = pd.MultiIndex.from_arrays([latitudes, longitudes], names=['Latitude', 'Longitude'])
dist_df.index = lat_lon_index
dist_df.columns = lat_lon_index

lat_from, lon_from = latitudes[0], longitudes[0]  
lat_to, lon_to = latitudes[1], longitudes[1]
distance = dist_df.loc[(lat_from, lon_from), (lat_to, lon_to)]
print(f"Distance from ({lat_from}, {lon_from}) to ({lat_to}, {lon_to}): {distance} km")

Теперь dist_df выглядит так

Чтобы сохранить его в CSV, просто сделайте следующее:

dist_df_reset = dist_df.reset_index()

polars_df = pl.from_pandas(dist_df_reset)

polars_df.write_csv('polars_distance_matrix.csv')

Файл будет:

ВАЖНОЕ ПРИМЕЧАНИЕ: Это будет намного быстрее, но вы не можете надеяться, что это займет очень мало времени с таким количеством точек.