如何矢量化比较来自两个不同数据帧的匀称对象的函数?

问题描述

我有一个 Pandas 数据框和一个 geopandas 数据框。在 Pandas 数据框中,我有一列 Points 包含 shapely.geometry Point 对象。 geopandas 框架中的 geometry 列包含 polygon 对象。我想做的是在 Pandas 框架中取一个 Point 并测试它是否是 geopandas 框架中的 within any polygon 对象.

在熊猫框架的新列中,我想要以下内容。如果 Point 在给定的 polygon 内(即 within 调用返回 True),我希望 Point 行中的新列的值是 geopandas 框架中 polygon 行中不同列的值。

我有一个解决这个问题的有效方法,但它没有被矢量化。是否可以对其进行矢量化?

示例:

import geopandas as gpd
import pandas as pd
from shapely.geometry import Point,polygon

# Create random frame,geometries are supposed to be mutually exclusive
gdf = gpd.GeoDataFrame({'A': [1,2],'geometry': [polygon([(10,5),(5,6)]),polygon([(1,2),(2,5))]})

# Create random pandas
df = pd.DataFrame({'Foo': ['bar','Bar'],'Points': [Point(4,Point(1,2)]})

# My non-vectorized solution
df['new'] = ''
for i in df.index:
    for j in gdf.index:
        if df.at[i,'Points'].within(gdf.at[j,'geometry']):
            df.at[i,'new'] = gdf.at[j,'A'] 

这很好用,因此当点位于多边形内时,df['new'] 将包含列 gdf['A'] 中的任何内容。我希望有一种方法可以让我对这个操作进行矢量化。

解决方法

您可以计算 PointsPolygon 的所有点之间的欧几里得距离。而且,只要距离等于 0,就会为您提供一个交点。我的方法如下。请注意,我将从数据框中获取所有点和多边形点的部分留给了您。可能,像 pandas.Series.toList 这样的函数应该提供这一点。

import numpy as np
from scipy.spatial.distance import cdist

polygon = [[10,5],[5,6],[1,2],[2,5]]
points = [[4,2]]

# return distances between all the items of the two arrays
distances = cdist(polygon,points) 

print(distances)
[[6.         9.48683298]
 [1.41421356 5.65685425]
 [4.24264069 0.        ]
 [2.         3.16227766]]

我们现在要做的就是获取数组中 0 的索引。如您所见,我们的交点在第 3 行第 2 列,即多边形的第 3 项或点的第 2 项。


for i,dist in enumerate(distances.flatten()):
    if dist==0:
        intersect_index = np.unravel_index(i,shape=distances.shape)
        intersect_point = polygon[intersect_index[0]]
        print(intersect_point)
[1,2]

这应该会为您提供您正在寻找的矢量化形式。

,

我找到了适合我的目的的解决方案。不是最优雅的,但仍然比循环快得多。

def within_vectorized(array,point):
# Create array of False and True values 
    _array = np.array([point.within(p) for p in array])
# When the first element of np.where tuple is not empty
    if np.where(_array)[0].size != 0:
        return np.where(_array)[0][0]
    else:
        return -1

# Create dummy value row geopandas frame
# This will have an empty Polygon object in the geometry column and NaN's everywhere else
dummy_values = np.empty((1,gdf.shape[1]))
dummy_values[:] = np.nan
dummy_values = dummy_values.tolist()[0]
dummy_values[-1] = Polygon()
gdf.loc[-1] = dummy_values

# Use loc where index is retrieved by calling vectorized function
df['A'] = gdf.loc[df['Point'].apply(lambda x: within_vectorized(gdf['geometry'],x)),'A'].to_list()