如何利用OSMNX自动标记GPS轨迹中位于步行区/绿地的点位？

我正在处理智能手机APP采集的真实GPS数据，数据集包含每秒一次的经纬度点位信息。目前我使用OSMNX库在地图上绘制轨迹，可通过可视化结合带有“pedestrian area”和“park”标签的OSMNX图层判断GPS点位是否位于步行区或绿地，但希望能通过OSMNX信息与真实数据的比对，自动标记每个GPS点位是否处于这两类区域中。

以下是我当前用于绘制单条轨迹及对应步行区、公园OSM图层的代码：

import networkx as nx
import osmnx as ox
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
import glob
import os
import osrm

#%matplotlib inline
ox.config(log_console=True)
ox.__version__

# 创建步行路网图
G = ox.graph_from_point((41.31367857092018, 2.0233411472684057), dist=1000, network_type='walk')

# 绘制路网并保留绘图轴
fig, ax = ox.plot_graph(G, show=False, close=False)

# 读取GPS数据
df = pd.read_csv('2018-11-05_sgv_0101_PEU.csv')

# 获取公园和步行区的地理数据
place = "Viladecans, Baix Llobregat"   
tags_park = {"leisure": "park"}
tags_pedestrian = {"highway": "pedestrian", "area": True}
gdf_parks = ox.geometries_from_place(place, tags_park)
gdf_ped_areas = ox.geometries_from_place(place, tags_pedestrian)

# 绘制公园（深绿）和步行区（紫色）
gdf_parks.plot(ax=ax, color='darkgreen')
gdf_ped_areas.plot(ax=ax, color='purple')

# 绘制起点和GPS轨迹点
ax.scatter(2.0233411472684057, 41.31367857092018, marker='*', c='yellow', s=200)
ax.scatter(df['longitude'], df['latitude'], c='red', s=1)

plt.show()

实现自动标记的解决方案 #

要实现自动标记，核心是判断每个GPS点是否落在公园或步行区的多边形范围内，可借助geopandas的空间操作功能完成：

1. 准备依赖 #

确保已安装geopandas：

pip install geopandas

2. 添加自动标记逻辑 #

在读取GPS数据后，插入以下代码：

import geopandas as gpd
from shapely.geometry import Point

# 将GPS数据转换为带空间信息的GeoDataFrame
geometry = [Point(lon, lat) for lon, lat in zip(df['longitude'], df['latitude'])]
gdf_gps = gpd.GeoDataFrame(df, geometry=geometry, crs=gdf_parks.crs)  # 匹配OSM图层坐标系

# 方法1：逐个点位判断（适合小数据集）
gdf_gps['in_park'] = gdf_gps.geometry.apply(lambda point: any(gdf_parks.geometry.contains(point)))
gdf_gps['in_pedestrian_area'] = gdf_gps.geometry.apply(lambda point: any(gdf_ped_areas.geometry.contains(point)))

# 方法2：空间连接（适合大数据集，效率更高）
# park_join = gpd.sjoin(gdf_gps, gdf_parks[['geometry']], how='left', predicate='within')
# gdf_gps['in_park'] = ~park_join.index_right.isna()
# ped_join = gpd.sjoin(gdf_gps, gdf_ped_areas[['geometry']], how='left', predicate='within')
# gdf_gps['in_pedestrian_area'] = ~ped_join.index_right.isna()

# 合并标记：是否处于公园或步行区
gdf_gps['in_target_area'] = gdf_gps['in_park'] | gdf_gps['in_pedestrian_area']

3. 可视化标记结果（可选） #

修改原散点绘制代码，根据标记区分点位颜色：

# 绘制标记后的点位
ax.scatter(gdf_gps[gdf_gps['in_target_area']]['longitude'], 
           gdf_gps[gdf_gps['in_target_area']]['latitude'], 
           c='blue', s=2, label='在公园/步行区')
ax.scatter(gdf_gps[~gdf_gps['in_target_area']]['longitude'], 
           gdf_gps[~gdf_gps['in_target_area']]['latitude'], 
           c='red', s=1, label='不在目标区域')
plt.legend()

4. 导出标记后的数据（可选） #

将带标记的GPS数据保存为CSV文件：

gdf_gps.to_csv('gps_with_area_tags.csv', index=False)

内容的提问来源于stack exchange，提问作者Ferran