GeoPandas: Temukan titik terdekat di kerangka data lain

20

Saya punya 2 geodataframe:

import geopandas as gpd
from shapely.geometry import Point
gpd1 = gpd.GeoDataFrame([['John',1,Point(1,1)],['Smith',1,Point(2,2)],['Soap',1,Point(0,2)]],columns=['Name','ID','geometry'])
gpd2 = gpd.GeoDataFrame([['Work',Point(0,1.1)],['Shops',Point(2.5,2)],['Home',Point(1,1.1)]],columns=['Place','geometry'])

dan saya ingin mencari nama titik terdekat di gpd2 untuk setiap baris di gpd1:

desired_output = 

    Name  ID     geometry  Nearest
0   John   1  POINT (1 1)     Home
1  Smith   1  POINT (2 2)    Shops
2   Soap   1  POINT (0 2)     Work

Saya sudah mencoba membuatnya berfungsi menggunakan fungsi lambda:

gpd1['Nearest'] = gpd1.apply(lambda row: min_dist(row.geometry,gpd2)['Place'] , axis=1)

dengan

def min_dist(point, gpd2):

    geoseries = some_function()
    return geoseries
RedM
sumber
Metode ini bekerja untuk saya: stackoverflow.com/questions/37402046/… lihat tautan
Johnny Cheesecutter

Jawaban:

16

Anda bisa langsung menggunakan fungsi Shapely Poin terdekat (geometri GeoSeries adalah geometri Shapely):

from shapely.ops import nearest_points
# unary union of the gpd2 geomtries 
pts3 = gpd2.geometry.unary_union
def near(point, pts=pts3):
     # find the nearest point and return the corresponding Place value
     nearest = gpd2.geometry == nearest_points(point, pts)[1]
     return gpd2[nearest].Place.get_values()[0]
gpd1['Nearest'] = gpd1.apply(lambda row: near(row.geometry), axis=1)
gpd1
    Name  ID     geometry  Nearest
0   John   1  POINT (1 1)     Home
1  Smith   1  POINT (2 2)    Shops
2   Soap   1  POINT (0 2)     Work

Penjelasan

for i, row in gpd1.iterrows():
    print nearest_points(row.geometry, pts3)[0], nearest_points(row.geometry, pts3)[1]
 POINT (1 1) POINT (1 1.1)
 POINT (2 2) POINT (2.5 2)
 POINT (0 2) POINT (0 1.1)
gen
sumber
Sesuatu tidak bekerja untuk saya dan saya tidak bisa memahaminya. Fungsi mengembalikan GeoSeries kosong meskipun geometrinya padat. Sebagai contoh: sample_point = gpd2.geometry.unary_union[400] / sample_point in gpd2.geometry Ini mengembalikan True. gpd2.geometry == sample_point Ini keluar semua False.
robroc
Tambahan di atas: gpd2.geometry.geom_equals(sample_point)berfungsi.
robroc
13

Jika Anda memiliki kerangka data yang besar, saya telah menemukan bahwa metode scipyindeks spasial cKDTree .querymengembalikan hasil yang sangat cepat untuk pencarian tetangga terdekat. Karena menggunakan indeks spasial, urutan besarnya lebih cepat daripada perulangan melalui dataframe dan kemudian menemukan minimum semua jarak. Ini juga lebih cepat daripada menggunakan shapely's nearest_pointsdengan RTree (metode indeks spasial tersedia melalui geopanda) karena cKDTree memungkinkan Anda untuk membuat vektor pencarian Anda, sedangkan metode lainnya tidak.

Berikut adalah fungsi pembantu yang akan mengembalikan jarak dan 'Nama' tetangga terdekat di gpd2dari setiap titik di gpd1. Ini mengasumsikan kedua gdf memiliki geometrykolom (poin).

import geopandas as gpd
import numpy as np
import pandas as pd

from scipy.spatial import cKDTree
from shapely.geometry import Point

gpd1 = gpd.GeoDataFrame([['John', 1, Point(1, 1)], ['Smith', 1, Point(2, 2)],
                         ['Soap', 1, Point(0, 2)]],
                        columns=['Name', 'ID', 'geometry'])
gpd2 = gpd.GeoDataFrame([['Work', Point(0, 1.1)], ['Shops', Point(2.5, 2)],
                         ['Home', Point(1, 1.1)]],
                        columns=['Place', 'geometry'])

def ckdnearest(gdA, gdB):
    nA = np.array(list(zip(gdA.geometry.x, gdA.geometry.y)) )
    nB = np.array(list(zip(gdB.geometry.x, gdB.geometry.y)) )
    btree = cKDTree(nB)
    dist, idx = btree.query(nA, k=1)
    gdf = pd.concat(
        [gdA, gdB.loc[idx, gdB.columns != 'geometry'].reset_index(),
         pd.Series(dist, name='dist')], axis=1)
    return gdf

ckdnearest(gpd1, gpd2)

Dan jika Anda ingin menemukan titik terdekat ke LineString, berikut ini contoh kerjanya:

import itertools
from operator import itemgetter

import geopandas as gpd
import numpy as np
import pandas as pd

from scipy.spatial import cKDTree
from shapely.geometry import Point, LineString

gpd1 = gpd.GeoDataFrame([['John', 1, Point(1, 1)],
                         ['Smith', 1, Point(2, 2)],
                         ['Soap', 1, Point(0, 2)]],
                        columns=['Name', 'ID', 'geometry'])
gpd2 = gpd.GeoDataFrame([['Work', LineString([Point(100, 0), Point(100, 1)])],
                         ['Shops', LineString([Point(101, 0), Point(101, 1), Point(102, 3)])],
                         ['Home',  LineString([Point(101, 0), Point(102, 1)])]],
                        columns=['Place', 'geometry'])


def ckdnearest(gdfA, gdfB, gdfB_cols=['Place']):
    A = np.concatenate(
        [np.array(geom.coords) for geom in gdfA.geometry.to_list()])
    B = [np.array(geom.coords) for geom in gdfB.geometry.to_list()]
    B_ix = tuple(itertools.chain.from_iterable(
        [itertools.repeat(i, x) for i, x in enumerate(list(map(len, B)))]))
    B = np.concatenate(B)
    ckd_tree = cKDTree(B)
    dist, idx = ckd_tree.query(A, k=1)
    idx = itemgetter(*idx)(B_ix)
    gdf = pd.concat(
        [gdfA, gdfB.loc[idx, gdfB_cols].reset_index(drop=True),
         pd.Series(dist, name='dist')], axis=1)
    return gdf

c = ckdnearest(gpd1, gpd2)
JHuw
sumber
Apakah mungkin untuk memberikan titik terdekat di telepon juga, menggunakan metode ini? Misalnya untuk mengambil lokasi GPS ke jalan terdekat.
hyperknot
Jawaban ini luar biasa! Namun, kode untuk titik terdekat ke baris menghasilkan bug untuk saya. Tampaknya jarak yang benar dari garis terdekat dikembalikan untuk setiap titik, tetapi baris id yang dikembalikan salah. Saya pikir ini adalah perhitungan idx, tapi saya cukup baru untuk Python, jadi saya tidak bisa membungkus kepala saya dengannya.
Shakedk
1

Menemukannya:

def min_dist(point, gpd2):
    gpd2['Dist'] = gpd2.apply(lambda row:  point.distance(row.geometry),axis=1)
    geoseries = gpd2.iloc[gpd2['Dist'].argmin()]
    return geoseries

Tentu saja ada kritik. Saya bukan penggemar menghitung ulang gpd2 ['Dist'] untuk setiap baris gpd1 ...

RedM
sumber
1

Jawaban oleh Gene tidak berhasil untuk saya. Akhirnya saya menemukan bahwa gpd2.geometry.unary_union menghasilkan geometri yang hanya berisi sekitar 30.000 dari total sekitar 150.000 poin. Untuk orang lain yang mengalami masalah yang sama, inilah cara saya menyelesaikannya:

    from shapely.ops import nearest_points
    from shapely.geometry import MultiPoint

    gpd2_pts_list = gpd2.geometry.tolist()
    gpd2_pts = MultiPoint(gpd2_pts_list)
    def nearest(point, gpd2_pts, gpd2=gpd2, geom_col='geometry', src_col='Place'):
         # find the nearest point
         nearest_point = nearest_points(point, gpd2_pts)[1]
         # return the corresponding value of the src_col of the nearest point
         value = gpd2[gpd2[geom_col] == nearest_point][src_col].get_values()[0]
         return value

    gpd1['Nearest'] = gpd1.apply(lambda x: nearest(x.geometry, gpd2_pts), axis=1)
Inske
sumber
0

Bagi siapa pun yang memiliki kesalahan pengindeksan dengan data mereka sendiri saat menggunakan jawaban yang sangat baik dari @ JHuw , masalah saya adalah bahwa indeks saya tidak selaras. Menyetel ulang indeks gdfA dan gdfB memecahkan masalah saya, mungkin ini dapat membantu Anda juga @ Shakedk .

import itertools
from operator import itemgetter

import geopandas as gpd
import numpy as np
import pandas as pd

from scipy.spatial import cKDTree
from shapely.geometry import Point, LineString

gpd1 = gpd.GeoDataFrame([['John', 1, Point(1, 1)],
                         ['Smith', 1, Point(2, 2)],
                         ['Soap', 1, Point(0, 2)]],
                        columns=['Name', 'ID', 'geometry'])
gpd2 = gpd.GeoDataFrame([['Work', LineString([Point(100, 0), Point(100, 1)])],
                         ['Shops', LineString([Point(101, 0), Point(101, 1), Point(102, 3)])],
                         ['Home',  LineString([Point(101, 0), Point(102, 1)])]],
                        columns=['Place', 'geometry'])


def ckdnearest(gdfA, gdfB, gdfB_cols=['Place']):
    # resetting the index of gdfA and gdfB here.
    gdfA = gdfA.reset_index(drop=True)
    gdfB = gdfB.reset_index(drop=True)
    A = np.concatenate(
        [np.array(geom.coords) for geom in gdfA.geometry.to_list()])
    B = [np.array(geom.coords) for geom in gdfB.geometry.to_list()]
    B_ix = tuple(itertools.chain.from_iterable(
        [itertools.repeat(i, x) for i, x in enumerate(list(map(len, B)))]))
    B = np.concatenate(B)
    ckd_tree = cKDTree(B)
    dist, idx = ckd_tree.query(A, k=1)
    idx = itemgetter(*idx)(B_ix)
    gdf = pd.concat(
        [gdfA, gdfB.loc[idx, gdfB_cols].reset_index(drop=True),
         pd.Series(dist, name='dist')], axis=1)
    return gdf

c = ckdnearest(gpd1, gpd2)
Markus Rosenfelder
sumber