Ottimizzazione della distanza GeoPandas

Question

Ho creato la seguente funzione che misura la distanza da un punto in un GDF a tutti i punti in un altro GDF e restituisce una tabella con la distanza più breve per ogni punto. Funziona bene per un punto, ma ho trascurato il fatto che ho una tabella di 4000 punti e quindi ci vogliono 10 minuti. L'ho eseguito in PostGIS e posso ridurlo a meno di un secondo. C'è un modo per farlo in Python che potrebbe corrispondere alla velocità di PostGIS?

def get_distance_to(gdf_in, aoi_df, aoi):
   dist_df_list = list()
   for row in range(len(gdf_in)):
      single_row = gdf_in.iloc[row]
      distances = aoi_df.geometry.distance(single_row.geometry)
      dist_list = distances.to_list()
      closest_aoi = min(dist_list)
      single_row["dist_to_"+aoi] = closest_aoi
      df = single_row.to_frame().T
      dist_df_list.append(df)
completed_distances = pd.concat(dist_df_list, ignore_index=True, sort=False)
return completed_distances

le mie tabelle di input hanno un aspetto simile a questo

e la tabella di output ha questo aspetto

Answer 1

Per qualsiasi operazione spaziale di questo tipo, dovresti sempre provare a utilizzare l'indice spaziale. Se sei interessato solo alla distanza minima, quanto segue dovrebbe darti un'opzione relativamente performante.

import geopandas as gpd
from shapely.geometry import Point
import pandas as pd
import random

gdf = gpd.GeoDataFrame(geometry=[Point(random.randint(0, 1000), random.randint(0, 1000)) for _ in range(1000)])
gdf2 = gpd.GeoDataFrame(geometry=[Point(random.randint(0, 1000), random.randint(0, 1000)) for _ in range(1000)])

def get_nearest_distance(left, right, initial_buffer):
    """get distance from left to right"""
    buffered = left.buffer(initial_buffer)

    distances = []
    for i in range(len(buffered)):
        geom = buffered.geometry.iloc[i]
        query = right.sindex.query(geom)
        while query.size == 0:
            query = right.sindex.query(geom.buffer(b))
            b += initial_buffer
        distances.append(right.iloc[query].distance(left.geometry.iloc[i]).min())

    return pd.Series(distances, index=left.index)

gdf['distance_to_x'] = get_nearest_distance(gdf, gdf2, 50)

Da 1000 a 1000 punti, è meno di un secondo, rispetto a circa un minuto impiegato dal codice del gene @.

Per renderlo efficiente, dovresti indovinare la initial_bufferdistanza in cui pensi saranno solo pochi punti. Se non è presente nessuno, espande il buffer finché non ne colpisce alcuni.

In generale, se vuoi le migliori prestazioni da GeoPandas, dovresti usare l'ultima versione (questo codice richiede 0.8) e pygeos di dipendenza opzionale (https://geopandas.readthedocs.io/en/latest/getting_started/install.html#using-the-optional-pygeos-dependency), che può accelerare il codice sopra in ordine di grandezza.

Answer 2

Iterare su righe in un DataFrame (Geo) in (Geo) Pandas è molto lento, vedi Approccio ottimale per l'iterazione su un DataFrame, ad esempio

L'iterazione in Pandas è un anti-pattern ed è qualcosa che dovresti fare solo quando hai esaurito ogni altra opzione. ( Come iterare su righe in un DataFrame in Pandas )

Puoi provare a usare (Geo)DataFrame.apply()e in forma: punto più vicino come in GeoPandas: trova il punto più vicino in un altro dataframe senza foriterazione (vedi commento)

import geopandas as gpd
from shapely.geometry import Point
from shapely.ops import nearest_points
gpd1 = gpd.read_file("point1.shp") # red points
gpd2 = gpd.read_file("point2.shp") # blue points
pts3 = gpd2.geometry.unary_union
 def near(point, pts=pts3):
     # find the nearest point and return the corresponding value
     nearest = gpd2.geometry == nearest_points(point, pts)[1]
     return gpd2[nearest].id.values[0],gpd2[nearest].geometry.values[0]
 gpd1['Nearest'] = gpd1.apply(lambda row: near(row.geometry)[0], axis=1)
 gpd1['geom2'] = gpd1.apply(lambda row: near(row.geometry)[1], axis=1)
 print(gpd1)
         id         geometry         Nearest        geom2
 0   1   POINT (-0.99013 0.48096)      3     POINT (-0.77574 0.64739)
 1   2   POINT (-1.00987 0.08039)      4     POINT (-0.73060 0.10860)
 2   3   POINT (-0.71932 -0.13117)     5     POINT (-0.57827 -0.08039)
 3   4   POINT (-0.90268 -0.28914)     5     POINT (-0.57827 -0.08039)

Calcola la distanza

gpd1['distance'] = gpd1.apply(lambda row: row.geometry.distance(row.geom2), axis=1)
gpd1.drop('geom2', axis=1, inplace=True)
print(gpd1)
    id         geometry             Nearest     distance
 0   1   POINT (-0.99013 0.48096)      3        0.271406
 1   2   POINT (-1.00987 0.08039)      4        0.280688
 2   3   POINT (-0.71932 -0.13117)     5        0.149905
 3   4   POINT (-0.90268 -0.28914)     5        0.385759