简介山东省2020年热门跑步位置点shapefile数据包主要面向GIS专业人员和跑步数据分析爱好者可用于分析山东省内热门跑步路线的空间分布、城市偏好与线路特征。压缩包共5个文件包含核心的山东省_热门跑步位置.shp几何文件、.dbf属性表、.prj投影定义、.shx空间索引以及.shp.xml元数据整体仅244KB便于快速下载和轻型空间分析。目前已有127人学习数据基于真实位置点整理自带热门程度、所在城市、线路类型和距离等字段可直接用ArcGIS或QGIS打开。利用热度字段可识别最受欢迎线路结合城市字段可绘制跑步热点地图也可按线路类型或距离做分组统计为城市健身设施规划或跑步文化研究提供底图数据。1. 一条跑步定位记录如何变成山东省的热门位置点跑步APP导出的原始数据往往是几十万行“设备号、时间戳、经纬度”的CSV而最终要交付的却是一个名叫“山东省2020热门跑步位置点shapefile”的GIS图层。这个标题真正要解决的不是怎么把点画在地图上而是怎么从离散轨迹里找到人群反复停留或经过的空间位置。处理路径通常是先建立点要素的字段模型和坐标系规范再做边界过滤与漂移清洗然后用密度聚类或网格聚合把散点变成热点最后用shapefile交付并接受核查。这套路线适合两类人一类是给政府体育部门做全民健身数据分析的GIS工程师另一类是拿到APP埋点数据但不懂空间分析的业务数据团队。下面是我在类似任务里会直接采用的做法。2. 用shapefile装跑步点字段设计与坐标系选型2.1 点要素的字段模型时间、速度以外的隐藏维度如果一份跑步点数据只存经纬度和时间后面做热门位置识别时会发现缺了三个关键信息轨迹归属、运动状态、数据来源。没有轨迹编号就无法区分“一个人反复跑”和“一个人跑一次”没有速度字段就无法过滤GPS漂移没有来源字段就无法解释设备之间的精度差异。我一般会先设计一个“干净点位”图层字段如下字段名类型说明tidint轨迹编号同一圈跑步的点共享一个tidlonfloat源坐标经度建议保留WGS84度数值latfloat源坐标纬度tsdatetime北京时间记录时刻不要存字符串speedfloat瞬时速度单位km/hheadingint航向角0-360度用于剔除异常sourcestring数据来源编码例如app_a / app_b这里有一个shapefile的长期痛点属性表基于dBase格式字段名长度限制在10个字符以内且不原生支持datetime的精度。所以字段名要克制source不能改成source_device这样过长的名字。ts必须写成ISO格式的字符串或双精度时间戳否则写入后读出来是损坏的。我通常用双精度存储Unix时间戳然后在QGIS里用表达式转换显示既保证跨平台兼容又避免dBase日期解析的时区问题。2.2 坐标系WGS84、GCJ-02与山东省地方投影热门位置点最怕坐标系混用。跑步设备的GPS芯片输出的是WGS84经纬度但很多APP为了让轨迹贴合在线底图会把它转换成GCJ-02甚至BD-09。如果你拿到的是APP端导出的点第一步就要确认它的偏移量。简单做法是选一个已知固定点比如泉城广场和一份高精度路网对比另一个做法是看轨迹与真实道路的吻合程度GCJ-02的偏移在高清卫星图上非常明显。山东省的位置范围大约在北纬34.3°到38.4°、东经114.8°到122.7°。做密度聚类时我建议不要直接用经纬度算欧氏距离因为经度1度的距离在高纬度和低纬度不一样。常见做法是维护两份shapefile一份是EPSG:4490CGCS2000地理坐标系用于与外部数据交换和底图发布另一份是投影坐标系副本例如EPSG:4523CGCS2000 / 3-degree Gauss-Kruger zone 39用于计算距离、面积和聚类。两份文件用同样的字段通过reproject同步即可。2.3 用GDAL/OGR检查shapefile的坐标边界和字段拿到任何一份shapefile我先用ogrinfo看空间参考、要素数量和坐标范围避免在异常数据上直接跑分析。命令如下ogrinfo -so runner_points.shp runner_points关键是看输出里的Extent和Geometry: Point字段。如果Extent的经度跨度超过3度说明数据里可能混入多省点位如果MinX小于73度则必然存在异常记录。更细的检查用Python的osgeo包来做from osgeo import ogr ds ogr.Open(runner_points.shp, 0) lyr ds.GetLayer() print(f要素数: {lyr.GetFeatureCount()}) print(f空间参考: {lyr.GetSpatialRef().ExportToWkt()[:80]}...) # 打印前两条要素的坐标和关键属性 for i, feat in enumerate(lyr): geom feat.GetGeometryRef() lon geom.GetX() lat geom.GetY() tid feat.GetField(tid) ts feat.GetField(ts) print(tid, lon, lat, ts) if i 1: break这段代码用于确认三件事图层是否真的是点要素、坐标系是否与预期一致、字段读取是否正常。GetX()和GetY()在投影坐标系下同样是合法的只是返回的是投影坐标值。如果GetSpatialRef()返回空说明shapefile缺少.prj文件这时候必须手动指定坐标系否则后面所有转换都会静默出错。对山东省2020跑步数据来说没有.prj文件的数据基本等于不可用。3. 清洗山东省2020跑步点的三个关键步骤3.1 先用山东省边界过滤再处理边界上的点热门位置点明确限定在山东省内所以第一步是用行政边界过滤点。这里不能直接用within因为GPS误差可能让真实位置在边界外几十米比如滨州沿海或微山湖边的点。我一般对边界做正向缓冲50米再用intersects判断from shapely.geometry import Point from shapely.ops import unary_union import geopandas as gpd # 读取山东省行政边界确保它和点位数据同坐标系 shandong gpd.read_file(shandong_boundary.shp).to_crs(EPSG:4490) points gpd.read_file(runner_points.shp).to_crs(EPSG:4490) # 对边界做缓冲吸收边缘误差 boundary shandong.geometry.unary_union.buffer(0.0005) # 空间过滤 mask points.intersects(boundary) filtered points[mask].copy()缓冲区大小“0.0005”度约等于50米。如果底图坐标系是GCJ-02这里的缓冲值含义就变了必须先把点位重投影到WGS84或CGCS2000再判断。过滤后的要素建议另存为shandong_runner_points.shp不要覆盖原始图层因为后面排查问题时经常要对比原始点。3.2 去重与GPS漂移不要忽略静止点和乱跳点跑步数据里重复点非常多等待红绿灯、系鞋带、跑步机上的数据都会产生大量点位堆叠。如果用DBSCAN做密度聚类这些静止点会制造假热点。但直接删除所有距离小于1米的点又会破坏轨迹的连续性所以我的清洗策略是按“轨迹内时距”去重而不是全局去重。常见做法是按轨迹tid排序保留同一轨迹内相邻两点距离大于5米的点同时检查速度与航向的突变import pandas as pd import numpy as np df pd.read_csv(shandong_raw.csv, parse_dates[ts]) df df.sort_values([tid, ts]) # 计算相邻点距离与速度用简化的球面距离 def haversine(lat1, lon1, lat2, lon2): R 6371e3 p1, p2 np.radians(lat1).copy(), np.radians(lat2) dphi np.radians(lat2 - lat1) dlambda np.radians(lon2 - lon1) a np.sin(dphi/2)**2 np.cos(p1)*np.cos(p2)*np.sin(dlambda/2)**2 return 2 * R * np.arcsin(np.sqrt(a)) df[prev_lat] df.groupby(tid)[lat].shift() df[prev_lon] df.groupby(tid)[lon].shift() df df.dropna(subset[prev_lat]) df[dist] df.apply(lambda r: haversine(r.lat, r.lon, r.prev_lat, r.prev_lon), axis1) df[dt] (df[ts] - df[ts].groupby(df[tid]).shift()).dt.total_seconds().fillna(0) df[speed] df[dist] / df[dt].replace(0, np.nan) * 3.6 # 剔除静止点距离5米和速度突变点速度25km/h设备跳变 cleaned df[(df[dist] 5) | (df[dist].isna())] cleaned cleaned[(cleaned[speed].isna()) | (cleaned[speed] 25)]这里的25km/h远高于跑步速度目的是过滤GPS芯片在隧道、高架下产生的瞬移点而不是过滤真实冲刺。注意dt为0时速度计算会出现除零所以先用replace转成NaN。清洗完成后点位数量通常会减少30%到40%但热度分布的形态不会明显改变这正好说明清洗去的是噪声而不是信号。3.3 坐标转点从CSV或原始表批量生成shapefile“坐标转点”在GIS里有两个层次一个是把Excel/CSV里的经纬度加载成临时地图点另一个是生成一个带完整属性结构的shapefile。网上很多教程只教你用QGIS的“添加XY数据”但那是临时图层通常在项目里保存不了字段别名和坐标系定义。我习惯用ogr直接写点文件from osgeo import ogr src_ds ogr.Open(shandong_points.gpkg) # 或打开原有shapefile src_lyr src_ds.GetLayer() driver ogr.GetDriverByName(ESRI Shapefile) if driver.DeleteDataSource(hot_points.shp): print(删除旧文件) ds driver.CreateDataSource(hot_points.shp) lyr ds.CreateLayer(hot_points, src_lyr.GetSpatialRef(), ogr.wkbPoint) # 从源图层复制字段定义 src_defn src_lyr.GetLayerDefn() for i in range(src_defn.GetFieldCount()): fld_defn src_defn.GetFieldDefn(i) lyr.CreateField(fld_defn) feat_defn lyr.GetLayerDefn() for feat in src_lyr: out_feat ogr.Feature(feat_defn) geom feat.GetGeometryRef().Clone() out_feat.SetGeometry(geom) for i in range(src_defn.GetFieldCount()): out_feat.SetField(i, feat.GetField(i)) lyr.CreateFeature(out_feat) out_feat None ds None这段代码把GPKG中的点位原样落成shapefile保留原始属性。lyr.CreateLayer的第三个参数指定几何类型为点而wbkPoint是2D点如果源数据带Z值要用wkbPoint25D。写完务必调用ds None释放文件句柄否则Windows下下一步访问会锁定文件。shapefile生成后会附带.dbf、.prj、.shx等文件压缩交付时要把同一基名的所有后缀都带上只发.shp会丢失坐标系和字段定义。4. 把点聚成热门位置DBSCAN与网格采样的取舍4.1 为什么简单覆盖度不能说明“热门”一个常规误区是直接统计每个1公里格网的点数然后取Top10。这种做法的致命问题是轨迹内相关性同一个跑者在一公里内留下100个点和100个跑者各留下1个点在网格计数上完全一样但语义截然不同。热门跑步位置应该体现“不同个体、多次到访”的空间而不仅仅是“定位记录密集”。所以我在聚合前做两层归一化第一层是轨迹内抽稀保证一条轨迹的局部点数量不影响整体热度第二层是访问身份归一同一用户在同一网格一天最多贡献1次到访记录。去重字段可以是匿名的设备号也可以用轨迹起点生成的轨迹ID。这样处理后密度值才真正对应“访问人次”的近似代理。4.2 DBSCAN聚类的参数选择eps和min_samples的实操值网格统计擅长表达热度却不能区分“滨河步道沿线连续热门”和“公园出口单点热门”。DBSCAN的优势是不需要预设簇的数量只关注点的空间密度。对山东省2020跑步点我建议在投影坐标系EPSG:4523上直接跑from sklearn.cluster import DBSCAN from sklearn.preprocessing import StandardScaler # 输入清洗后的点位投影坐标用东向和北向两列 coords cleaned[[x_proj, y_proj]].values eps_mu 100 # 100米空间范围 min_pts 20 # 至少20个访问点 labels DBSCAN(epseps_mu, min_samplesmin_pts, metriceuclidean, n_jobs-1).fit_predict(coords) cleaned[cluster_id] labels hot_clusters cleaned[cleaned[cluster_id] 0]这里的100米是一个经验起点适配城市公园、滨水步道的尺度。如果在操场这类小范围场景建议降到30米如果是山区大环线可以提高到300到500米。不要直接用经纬度度数作为eps1度经度在不同纬度实际长度差异太大聚类结果会随位置漂移。min_samples设为20的含义是“半径100米内要存在20个不同访问点”而不是20条GPS记录。为了确保是“不同访问”需要在聚类前把同一用户在同一天、同一200米网格内的点合并。聚类完成后把每个簇的点中心作为热门位置点写入输出shapefile并附带n_visits字段。如果某个簇的边界特别狭长说明热门位置可能是一条跑步径而不是场地这时候可以继续对簇内点做线拟合识别跑道方向。4.3 网格聚合把点变成面并按面积均分布点DBSCAN给出的热门点往往参差不齐但很多应用场景需要“均匀布点”比如要布设监测设备、规划补给站希望在每个热门区域覆盖相同数量的候选位。这时我会用网格聚合先把点变成面再做“按面积均分布点”的采样。方案是这样用山东省边界生成1公里方格统计每个格网中去重后的访问人次然后对指定的Top N区域做空间分层抽样。每个区域内按面积比例分配抽样点数而不是按热度比例分配这样能保证地理覆盖均匀不会让三个采样点挤在公园入口import geopandas as gpd import numpy as np # 生成1公里格网假设热点区域已投影 hull hot_clusters.unary_union.envelope minx, miny, maxx, maxy hull.bounds grid_cells [] cell_size 1000 # 米 for ix in np.arange(minx, maxx, cell_size): for iy in np.arange(miny, maxy, cell_size): cell gpd.GeoDataFrame({geometry: [shapely.geometry.box(ix, iy, ix cell_size, iy cell_size)]}, crsEPSG:4523) if cell.intersects(hot_boundary).any(): grid_cells.append(cell) grid gpd.pd.concat(grid_cells) # 格网与访问点连接 grid gpd.sjoin(grid, hot_clusters[[geometry, n_visits]], howleft, predicatecontains) counts grid.dissolve(bygrid_id, aggfunccount).rename(columns{n_visits_right:visit_cnt})这个格网面文件本身就是shapefile中的Polygon图层可以直接用于展示“热门区域”。如果你需要的是点就在每个格网内含要素的范围内随机布点例如每个格网固定放2个候选点再根据实际道路距离调整。这样生成的点位兼具三个特性空间均匀、覆盖热门区域、数量可控。注意grid_id要在生成格网时手动写入否则dissolve后没有聚合键。5. 核查与交付批量出图和二次验证5.1 用QGIS打开shapefile快速核查热点位置交付前先在QGIS里加载山东省底图和热点shapefile。用图层样式里的“按数量分级”以n_visits字段做渐变圆点同时叠加路网数据。重点看两类异常一是热点落在农田或大型立交互通中二是热点成一条直线贴在道路边上这往往是坐标偏移后的结果。遇到第二种情况可以对热点做50米缓冲区与OSM路网线相交计算每个热点缓冲区内道路长度占缓冲区面积的比例比例过小则标红待查。5.2 PyQGIS脚本批量出图把热门位置点变成可读地图人工逐个导出热点地图浪费时间我用QGIS的printLayout机制批量出图保证每张图的比例尺、图例和指北针完全一致。脚本核心如下from qgis.core import QgsProject, QgsLayout, QgsLayoutItemMap, QgsLayoutItemLegend from qgis.PyQt.QtCore import QSizeF project QgsProject.instance() layout QgsLayout(project) layout.initializeDefaults() layout.setPageSize(QSizeF(297, 210)) # A4横向 # 添加地图到布局 map_item QgsLayoutItemMap(layout) map_item.setRect(20, 20, 260, 170) map_item.setFrameEnabled(True) map_item.zoomToExtent(hot_layer.extent().buffered(500)) layout.addLayoutItem(map_item)循环内只需要对每个热点区域执行zoomToExtent然后调用QgsLayoutExporter.exportToPng导出。批量出图最容易忽略的是投影layout里的底图默认跟随项目工程坐标系如果你的工程是WGS84而热点图层是投影坐标缩放范围会用错。我在脚本里强制把画布和布局地图设为相同CRS再计算热点中心点。出图后建议把PNG文件名的前缀设置为热点ID例如hot_42.png方便后续拼接巡检报告。5.3 用轨迹方向验证热门位置点的真实性最后一个验证技巧是看热点处的轨迹方向分布。一个真实跑步热点应该是双向或多方向通过的如果所有轨迹方向完全一致说明这个热点可能是某条单向直行道路的坐标偏移而不是一个场地。提取每个热点中心周边200米范围内的轨迹起点和终点计算航向角差值并画成角度直方图。方向极径越均匀热点越可信。如果发现某个热点方向高度集中我一般会检查两点一是该处是否引用了高架桥下方的GPS漂移点二是是否没有过滤掉骑行数据骑行轨迹往往沿道路单向高速移动方向一致性很高。这时候回到原始数据按speed字段剔除大于12km/h的记录再重新跑一轮。这一步虽然简单但往往能过滤掉20%以上的假热门点也让最终交付的“山东省2020热门跑步位置点shapefile”经得起数据复核。本文还有配套的精品资源点击获取 GPS轨迹如何变成热门跑步位置点？用shapefile与DBSCAN做空间聚类分析
简介山东省2020年热门跑步位置点shapefile数据包主要面向GIS专业人员和跑步数据分析爱好者可用于分析山东省内热门跑步路线的空间分布、城市偏好与线路特征。压缩包共5个文件包含核心的山东省_热门跑步位置.shp几何文件、.dbf属性表、.prj投影定义、.shx空间索引以及.shp.xml元数据整体仅244KB便于快速下载和轻型空间分析。目前已有127人学习数据基于真实位置点整理自带热门程度、所在城市、线路类型和距离等字段可直接用ArcGIS或QGIS打开。利用热度字段可识别最受欢迎线路结合城市字段可绘制跑步热点地图也可按线路类型或距离做分组统计为城市健身设施规划或跑步文化研究提供底图数据。1. 一条跑步定位记录如何变成山东省的热门位置点跑步APP导出的原始数据往往是几十万行“设备号、时间戳、经纬度”的CSV而最终要交付的却是一个名叫“山东省2020热门跑步位置点shapefile”的GIS图层。这个标题真正要解决的不是怎么把点画在地图上而是怎么从离散轨迹里找到人群反复停留或经过的空间位置。处理路径通常是先建立点要素的字段模型和坐标系规范再做边界过滤与漂移清洗然后用密度聚类或网格聚合把散点变成热点最后用shapefile交付并接受核查。这套路线适合两类人一类是给政府体育部门做全民健身数据分析的GIS工程师另一类是拿到APP埋点数据但不懂空间分析的业务数据团队。下面是我在类似任务里会直接采用的做法。2. 用shapefile装跑步点字段设计与坐标系选型2.1 点要素的字段模型时间、速度以外的隐藏维度如果一份跑步点数据只存经纬度和时间后面做热门位置识别时会发现缺了三个关键信息轨迹归属、运动状态、数据来源。没有轨迹编号就无法区分“一个人反复跑”和“一个人跑一次”没有速度字段就无法过滤GPS漂移没有来源字段就无法解释设备之间的精度差异。我一般会先设计一个“干净点位”图层字段如下字段名类型说明tidint轨迹编号同一圈跑步的点共享一个tidlonfloat源坐标经度建议保留WGS84度数值latfloat源坐标纬度tsdatetime北京时间记录时刻不要存字符串speedfloat瞬时速度单位km/hheadingint航向角0-360度用于剔除异常sourcestring数据来源编码例如app_a / app_b这里有一个shapefile的长期痛点属性表基于dBase格式字段名长度限制在10个字符以内且不原生支持datetime的精度。所以字段名要克制source不能改成source_device这样过长的名字。ts必须写成ISO格式的字符串或双精度时间戳否则写入后读出来是损坏的。我通常用双精度存储Unix时间戳然后在QGIS里用表达式转换显示既保证跨平台兼容又避免dBase日期解析的时区问题。2.2 坐标系WGS84、GCJ-02与山东省地方投影热门位置点最怕坐标系混用。跑步设备的GPS芯片输出的是WGS84经纬度但很多APP为了让轨迹贴合在线底图会把它转换成GCJ-02甚至BD-09。如果你拿到的是APP端导出的点第一步就要确认它的偏移量。简单做法是选一个已知固定点比如泉城广场和一份高精度路网对比另一个做法是看轨迹与真实道路的吻合程度GCJ-02的偏移在高清卫星图上非常明显。山东省的位置范围大约在北纬34.3°到38.4°、东经114.8°到122.7°。做密度聚类时我建议不要直接用经纬度算欧氏距离因为经度1度的距离在高纬度和低纬度不一样。常见做法是维护两份shapefile一份是EPSG:4490CGCS2000地理坐标系用于与外部数据交换和底图发布另一份是投影坐标系副本例如EPSG:4523CGCS2000 / 3-degree Gauss-Kruger zone 39用于计算距离、面积和聚类。两份文件用同样的字段通过reproject同步即可。2.3 用GDAL/OGR检查shapefile的坐标边界和字段拿到任何一份shapefile我先用ogrinfo看空间参考、要素数量和坐标范围避免在异常数据上直接跑分析。命令如下ogrinfo -so runner_points.shp runner_points关键是看输出里的Extent和Geometry: Point字段。如果Extent的经度跨度超过3度说明数据里可能混入多省点位如果MinX小于73度则必然存在异常记录。更细的检查用Python的osgeo包来做from osgeo import ogr ds ogr.Open(runner_points.shp, 0) lyr ds.GetLayer() print(f要素数: {lyr.GetFeatureCount()}) print(f空间参考: {lyr.GetSpatialRef().ExportToWkt()[:80]}...) # 打印前两条要素的坐标和关键属性 for i, feat in enumerate(lyr): geom feat.GetGeometryRef() lon geom.GetX() lat geom.GetY() tid feat.GetField(tid) ts feat.GetField(ts) print(tid, lon, lat, ts) if i 1: break这段代码用于确认三件事图层是否真的是点要素、坐标系是否与预期一致、字段读取是否正常。GetX()和GetY()在投影坐标系下同样是合法的只是返回的是投影坐标值。如果GetSpatialRef()返回空说明shapefile缺少.prj文件这时候必须手动指定坐标系否则后面所有转换都会静默出错。对山东省2020跑步数据来说没有.prj文件的数据基本等于不可用。3. 清洗山东省2020跑步点的三个关键步骤3.1 先用山东省边界过滤再处理边界上的点热门位置点明确限定在山东省内所以第一步是用行政边界过滤点。这里不能直接用within因为GPS误差可能让真实位置在边界外几十米比如滨州沿海或微山湖边的点。我一般对边界做正向缓冲50米再用intersects判断from shapely.geometry import Point from shapely.ops import unary_union import geopandas as gpd # 读取山东省行政边界确保它和点位数据同坐标系 shandong gpd.read_file(shandong_boundary.shp).to_crs(EPSG:4490) points gpd.read_file(runner_points.shp).to_crs(EPSG:4490) # 对边界做缓冲吸收边缘误差 boundary shandong.geometry.unary_union.buffer(0.0005) # 空间过滤 mask points.intersects(boundary) filtered points[mask].copy()缓冲区大小“0.0005”度约等于50米。如果底图坐标系是GCJ-02这里的缓冲值含义就变了必须先把点位重投影到WGS84或CGCS2000再判断。过滤后的要素建议另存为shandong_runner_points.shp不要覆盖原始图层因为后面排查问题时经常要对比原始点。3.2 去重与GPS漂移不要忽略静止点和乱跳点跑步数据里重复点非常多等待红绿灯、系鞋带、跑步机上的数据都会产生大量点位堆叠。如果用DBSCAN做密度聚类这些静止点会制造假热点。但直接删除所有距离小于1米的点又会破坏轨迹的连续性所以我的清洗策略是按“轨迹内时距”去重而不是全局去重。常见做法是按轨迹tid排序保留同一轨迹内相邻两点距离大于5米的点同时检查速度与航向的突变import pandas as pd import numpy as np df pd.read_csv(shandong_raw.csv, parse_dates[ts]) df df.sort_values([tid, ts]) # 计算相邻点距离与速度用简化的球面距离 def haversine(lat1, lon1, lat2, lon2): R 6371e3 p1, p2 np.radians(lat1).copy(), np.radians(lat2) dphi np.radians(lat2 - lat1) dlambda np.radians(lon2 - lon1) a np.sin(dphi/2)**2 np.cos(p1)*np.cos(p2)*np.sin(dlambda/2)**2 return 2 * R * np.arcsin(np.sqrt(a)) df[prev_lat] df.groupby(tid)[lat].shift() df[prev_lon] df.groupby(tid)[lon].shift() df df.dropna(subset[prev_lat]) df[dist] df.apply(lambda r: haversine(r.lat, r.lon, r.prev_lat, r.prev_lon), axis1) df[dt] (df[ts] - df[ts].groupby(df[tid]).shift()).dt.total_seconds().fillna(0) df[speed] df[dist] / df[dt].replace(0, np.nan) * 3.6 # 剔除静止点距离5米和速度突变点速度25km/h设备跳变 cleaned df[(df[dist] 5) | (df[dist].isna())] cleaned cleaned[(cleaned[speed].isna()) | (cleaned[speed] 25)]这里的25km/h远高于跑步速度目的是过滤GPS芯片在隧道、高架下产生的瞬移点而不是过滤真实冲刺。注意dt为0时速度计算会出现除零所以先用replace转成NaN。清洗完成后点位数量通常会减少30%到40%但热度分布的形态不会明显改变这正好说明清洗去的是噪声而不是信号。3.3 坐标转点从CSV或原始表批量生成shapefile“坐标转点”在GIS里有两个层次一个是把Excel/CSV里的经纬度加载成临时地图点另一个是生成一个带完整属性结构的shapefile。网上很多教程只教你用QGIS的“添加XY数据”但那是临时图层通常在项目里保存不了字段别名和坐标系定义。我习惯用ogr直接写点文件from osgeo import ogr src_ds ogr.Open(shandong_points.gpkg) # 或打开原有shapefile src_lyr src_ds.GetLayer() driver ogr.GetDriverByName(ESRI Shapefile) if driver.DeleteDataSource(hot_points.shp): print(删除旧文件) ds driver.CreateDataSource(hot_points.shp) lyr ds.CreateLayer(hot_points, src_lyr.GetSpatialRef(), ogr.wkbPoint) # 从源图层复制字段定义 src_defn src_lyr.GetLayerDefn() for i in range(src_defn.GetFieldCount()): fld_defn src_defn.GetFieldDefn(i) lyr.CreateField(fld_defn) feat_defn lyr.GetLayerDefn() for feat in src_lyr: out_feat ogr.Feature(feat_defn) geom feat.GetGeometryRef().Clone() out_feat.SetGeometry(geom) for i in range(src_defn.GetFieldCount()): out_feat.SetField(i, feat.GetField(i)) lyr.CreateFeature(out_feat) out_feat None ds None这段代码把GPKG中的点位原样落成shapefile保留原始属性。lyr.CreateLayer的第三个参数指定几何类型为点而wbkPoint是2D点如果源数据带Z值要用wkbPoint25D。写完务必调用ds None释放文件句柄否则Windows下下一步访问会锁定文件。shapefile生成后会附带.dbf、.prj、.shx等文件压缩交付时要把同一基名的所有后缀都带上只发.shp会丢失坐标系和字段定义。4. 把点聚成热门位置DBSCAN与网格采样的取舍4.1 为什么简单覆盖度不能说明“热门”一个常规误区是直接统计每个1公里格网的点数然后取Top10。这种做法的致命问题是轨迹内相关性同一个跑者在一公里内留下100个点和100个跑者各留下1个点在网格计数上完全一样但语义截然不同。热门跑步位置应该体现“不同个体、多次到访”的空间而不仅仅是“定位记录密集”。所以我在聚合前做两层归一化第一层是轨迹内抽稀保证一条轨迹的局部点数量不影响整体热度第二层是访问身份归一同一用户在同一网格一天最多贡献1次到访记录。去重字段可以是匿名的设备号也可以用轨迹起点生成的轨迹ID。这样处理后密度值才真正对应“访问人次”的近似代理。4.2 DBSCAN聚类的参数选择eps和min_samples的实操值网格统计擅长表达热度却不能区分“滨河步道沿线连续热门”和“公园出口单点热门”。DBSCAN的优势是不需要预设簇的数量只关注点的空间密度。对山东省2020跑步点我建议在投影坐标系EPSG:4523上直接跑from sklearn.cluster import DBSCAN from sklearn.preprocessing import StandardScaler # 输入清洗后的点位投影坐标用东向和北向两列 coords cleaned[[x_proj, y_proj]].values eps_mu 100 # 100米空间范围 min_pts 20 # 至少20个访问点 labels DBSCAN(epseps_mu, min_samplesmin_pts, metriceuclidean, n_jobs-1).fit_predict(coords) cleaned[cluster_id] labels hot_clusters cleaned[cleaned[cluster_id] 0]这里的100米是一个经验起点适配城市公园、滨水步道的尺度。如果在操场这类小范围场景建议降到30米如果是山区大环线可以提高到300到500米。不要直接用经纬度度数作为eps1度经度在不同纬度实际长度差异太大聚类结果会随位置漂移。min_samples设为20的含义是“半径100米内要存在20个不同访问点”而不是20条GPS记录。为了确保是“不同访问”需要在聚类前把同一用户在同一天、同一200米网格内的点合并。聚类完成后把每个簇的点中心作为热门位置点写入输出shapefile并附带n_visits字段。如果某个簇的边界特别狭长说明热门位置可能是一条跑步径而不是场地这时候可以继续对簇内点做线拟合识别跑道方向。4.3 网格聚合把点变成面并按面积均分布点DBSCAN给出的热门点往往参差不齐但很多应用场景需要“均匀布点”比如要布设监测设备、规划补给站希望在每个热门区域覆盖相同数量的候选位。这时我会用网格聚合先把点变成面再做“按面积均分布点”的采样。方案是这样用山东省边界生成1公里方格统计每个格网中去重后的访问人次然后对指定的Top N区域做空间分层抽样。每个区域内按面积比例分配抽样点数而不是按热度比例分配这样能保证地理覆盖均匀不会让三个采样点挤在公园入口import geopandas as gpd import numpy as np # 生成1公里格网假设热点区域已投影 hull hot_clusters.unary_union.envelope minx, miny, maxx, maxy hull.bounds grid_cells [] cell_size 1000 # 米 for ix in np.arange(minx, maxx, cell_size): for iy in np.arange(miny, maxy, cell_size): cell gpd.GeoDataFrame({geometry: [shapely.geometry.box(ix, iy, ix cell_size, iy cell_size)]}, crsEPSG:4523) if cell.intersects(hot_boundary).any(): grid_cells.append(cell) grid gpd.pd.concat(grid_cells) # 格网与访问点连接 grid gpd.sjoin(grid, hot_clusters[[geometry, n_visits]], howleft, predicatecontains) counts grid.dissolve(bygrid_id, aggfunccount).rename(columns{n_visits_right:visit_cnt})这个格网面文件本身就是shapefile中的Polygon图层可以直接用于展示“热门区域”。如果你需要的是点就在每个格网内含要素的范围内随机布点例如每个格网固定放2个候选点再根据实际道路距离调整。这样生成的点位兼具三个特性空间均匀、覆盖热门区域、数量可控。注意grid_id要在生成格网时手动写入否则dissolve后没有聚合键。5. 核查与交付批量出图和二次验证5.1 用QGIS打开shapefile快速核查热点位置交付前先在QGIS里加载山东省底图和热点shapefile。用图层样式里的“按数量分级”以n_visits字段做渐变圆点同时叠加路网数据。重点看两类异常一是热点落在农田或大型立交互通中二是热点成一条直线贴在道路边上这往往是坐标偏移后的结果。遇到第二种情况可以对热点做50米缓冲区与OSM路网线相交计算每个热点缓冲区内道路长度占缓冲区面积的比例比例过小则标红待查。5.2 PyQGIS脚本批量出图把热门位置点变成可读地图人工逐个导出热点地图浪费时间我用QGIS的printLayout机制批量出图保证每张图的比例尺、图例和指北针完全一致。脚本核心如下from qgis.core import QgsProject, QgsLayout, QgsLayoutItemMap, QgsLayoutItemLegend from qgis.PyQt.QtCore import QSizeF project QgsProject.instance() layout QgsLayout(project) layout.initializeDefaults() layout.setPageSize(QSizeF(297, 210)) # A4横向 # 添加地图到布局 map_item QgsLayoutItemMap(layout) map_item.setRect(20, 20, 260, 170) map_item.setFrameEnabled(True) map_item.zoomToExtent(hot_layer.extent().buffered(500)) layout.addLayoutItem(map_item)循环内只需要对每个热点区域执行zoomToExtent然后调用QgsLayoutExporter.exportToPng导出。批量出图最容易忽略的是投影layout里的底图默认跟随项目工程坐标系如果你的工程是WGS84而热点图层是投影坐标缩放范围会用错。我在脚本里强制把画布和布局地图设为相同CRS再计算热点中心点。出图后建议把PNG文件名的前缀设置为热点ID例如hot_42.png方便后续拼接巡检报告。5.3 用轨迹方向验证热门位置点的真实性最后一个验证技巧是看热点处的轨迹方向分布。一个真实跑步热点应该是双向或多方向通过的如果所有轨迹方向完全一致说明这个热点可能是某条单向直行道路的坐标偏移而不是一个场地。提取每个热点中心周边200米范围内的轨迹起点和终点计算航向角差值并画成角度直方图。方向极径越均匀热点越可信。如果发现某个热点方向高度集中我一般会检查两点一是该处是否引用了高架桥下方的GPS漂移点二是是否没有过滤掉骑行数据骑行轨迹往往沿道路单向高速移动方向一致性很高。这时候回到原始数据按speed字段剔除大于12km/h的记录再重新跑一轮。这一步虽然简单但往往能过滤掉20%以上的假热门点也让最终交付的“山东省2020热门跑步位置点shapefile”经得起数据复核。本文还有配套的精品资源点击获取