地质灾害点数据全解析:25万+灾害点Shapefile清洗与风险制图实战

发布时间:2026/9/9 15:22:10
地质灾害点数据全解析:25万+灾害点Shapefile清洗与风险制图实战 做全国地质灾害风险底图的时候最头疼的其实不是画图而是找不到一套能直接用的地质灾害点数据。我前后折腾了大半年最后整理出全国25万地质灾害点空间分布数据1949—2019覆盖7大类灾害类型统一转成Shapefile矢量格式GIS软件直接打开就能叠图、查询、做密度分析。这篇文章就把这套数据背后的门道讲透包括7大类灾害分别是什么、Shapefile格式到底有什么讲究、拿到点数据后应该怎么清洗校验以及从25万个点里提炼结论的实操思路。不管你是搞GIS分析、地质灾害防治、城乡规划还是做科研底图这套方法都值得参考。1. 这套数据到底有什么25万灾害点的家底拆解先说数据本体。25万不是随便凑出来的数量级它意味着文件加载进ArcGIS或者QGIS之后一次属性表全选就会有明显卡顿也意味着统计口径稍微变一下某类灾害数量可能差出几千条。所以在用数据之前要先把它的“底细”摸清楚。1.1 七大类灾害类型逐一拆解标题里说的“7大类灾害类型”对应的是我国地质灾害防治工作中最常见的几类滑坡、崩塌、泥石流、地面塌陷、地裂缝、地面沉降以及不稳定斜坡。这7类虽然都叫地质灾害但形成机制完全不同空间分布规律差异也很大混在一起分析等于白做。我这里画个表方便你理解每一类的本质差异灾害类型主要成因空间分布典型特征常见伴生情况滑坡斜坡岩土体在重力、降雨、地震等作用下整体滑动多分布在山区中缓坡、公路铁路边坡、库区两岸常与降雨同期发生可能堵江、毁路崩塌陡坡岩土体在自重作用下突然崩落高陡边坡、断层带、采矿削坡区域块体运动速度快预警时间短泥石流山区沟谷中饱含水沙石碎屑的流动体沟谷、冲沟、植被破坏严重的山区暴雨后集中暴发具有群发性地面塌陷地下空洞或松散土层失稳导致地表塌落岩溶区、矿山采空区、地下隧道沿线可能形成天坑影响建筑安全地裂缝地壳活动或地下水、矿产开采引发地表开裂活动断裂带、地下水超采区、矿区延伸长穿过农田、公路和居民区地面沉降长期过量开采地下水或地下资源导致地面标高降低冲积平原、沿海城市、水源地周边缓变过程累积影响排水和建筑不稳定斜坡已出现变形迹象但尚未失稳的斜坡体滑坡、崩塌隐患点周边高风险隐患点需要持续监测这7类灾害在数据属性表里一般会有专门的“灾害类型”字段但不同来源的字段命名五花八门有的是汉字“滑坡”“崩塌”有的是拼音首字母“HP”“BT”还有的是代码“01、02、03”。如果你拿到手的数据里类型字段不规范第一件事就是先统一编码否则后面所有筛选统计都会错。1.2 1949—2019这个时间跨度意味着什么1949—2019年恰好是70年这个时间窗口不是随便切的。国内成体系的地质灾害记录、调查和数据库建设大致就是从这个时间节点逐步积累起来的。早年间记录比较稀疏属性信息往往只有地点和灾害类型到了近20年记录逐渐完整甚至能包含规模等级、威胁人数、直接经济损失等要素。所以在分析时要注意一个关键问题早期的点少不代表当年灾害少而是记录缺失严重。如果你直接拿1949年的数据点和2019年的数据点做时间序列对比会得出“灾害数量逐年暴增”的结论这是不严谨的。实际上2010年以后记录更完整加之大范围排查和群测群防体系建立大量小型地质灾害被纳入记录数据量增长明显与调查力度有关未必完全等于实际灾情增长。从数据字段看除了灾害类型常见的核心字段至少包括灾害点编号、位置名称、所在省市区县、经度、纬度、发生年份、规模等级、威胁人口、威胁财产、灾害成因、数据来源、备注说明。其中经纬度字段是空间定位的命根子建议在拿到数据后第一时间检查是否有缺值、零值或越界值因为这类问题在高维数据集中太常见了。2. Shapefile矢量格式的实战细节别急着双击加载Shapefile是GIS界最普及的矢量格式之一但它和常见的Excel、CSV不一样它不是一个文件而是一组文件。很多人第一次使用Shapefile时只看到一个.shp后缀文件就把它单独复制走了结果打开提示缺失数据。理解这套文件结构是使用这套地质灾害点数据的第一步。2.1 一个Shapefile是一组兄弟文件少一个都不行Shapefile格式至少包括三个必要文件.shp几何信息文件存储点、线、面的坐标位置。.shx形状索引文件用来快速定位几何对象的位置相当于目录索引。.dbf属性表文件存储每个要素的非空间属性比如灾害类型、发生年份、地点名称等。此外还有几个可选但很重要的辅助文件.prj投影信息文件记录了坐标系统、投影方式。没有这个文件软件只能靠猜容易导致位置偏移。.cpg字符编码声明文件用来告诉软件.dbf属性表是什么编码。没有或编码错误中文属性会乱码成“口口口”或“鍦熸柟”。.sbn/.sbx空间索引文件大数据量下能提升查询速度。在分发或传输这套全国25万地质灾害点数据时一定要把整个文件夹打包成zip而不是只拖出某一个文件。我实测过只传.shp文件给对方对方用ArcGIS打开就会提示“无法打开”或者图形出现但属性表全空。正确做法是压缩时选中整个文件夹或者把同名的所有后缀文件一起选中再打包发送。如果你的数据缺少.prj文件也不必太慌可以通过经纬度范围判断大致坐标系。中国境内的点数据经度范围大约东经73°—135°纬度范围北纬18°—53°。如果数值在这个范围内通常就是WGS84或CGCS2000经纬度坐标如果经度是8位左右的大数比如400000米量级那很可能是经过投影后的平面坐标这类情况需要额外反算。2.2 坐标系不查清楚就叠图偏出去几十公里很正常全国范围的点数据最稳妥的坐标系是WGS84地理坐标系或CGCS2000地理坐标系也就是用经纬度直接表示位置。这种坐标系的优点是全国范围都能统一使用不会因为分带问题出现断裂缺点是长度和面积计算会有畸变做距离、面积量算时不准确。我见过不少新手直接把WGS84经纬度的点数据和地方政府的2000投影坐标边界叠加结果点全跑到地图外。原因就是数据结构不匹配。处理思路分两步第一查看数据的.prj文件确认当前坐标系。在QGIS中可以直接右键图层属性查看在ArcGIS中则是图层属性里的“源”选项卡。第二如果需要做投影转换建议先用Python的pyproj或geopandas统一转到WGS84地理坐标之后再做局部区域的投影分析。举个例子如果要做喀斯特区域的落水洞点分析需要转成高斯-克吕格三度带投影才能精确计算距离。可以用下面的思路import geopandas as gpd gdf gpd.read_file(地质灾害点.shp, encodingutf-8) # 确保当前是WGS84经纬度 gdf gdf.to_crs(epsg4326) # 保存一份WGS84版本便于全球配图 gdf.to_file(地质灾害点_wgs84.shp, encodingutf-8) # 如果需要投影以EPSG:4547CGCS2000 / 3-degree Gauss-Kruger zone 37为例 gdf_proj gdf.to_crs(epsg4547) gdf_proj.to_file(地质灾害点_4547.shp, encodingutf-8)这里特别提一句不是坐标数字大就精度高。很多数据标注为“GCS_WGS_1984”但实际导出时经纬度被保留成整数比如107.123456精度看起来还行但如果字段是浮点型但保留3位小数点位误差可能已经超过100米。对于区域级风险分析100米的误差可接受但如果你要针对某个村庄的房前屋后做隐患点排查这个精度就不够了必须回到原始调查记录去核对。2.3 属性表的快速读取用GDAL一行命令看全貌生产环境中不一定装了ArcGIS但GDAL是开源GIS绕不开的工具。想快速预览Shapefile属性不用打开软件命令行一条命令就能搞定ogrinfo -so -al 地质灾害点.shp这条命令会输出要素数量、图层几何类型、属性字段列表和坐标范围。想看具体属性记录把-so去掉即可ogrinfo -al 地质灾害点.shp | head -50如果想把属性表导成CSV方便在Python或Excel里处理用ogr2ogrogr2ogr -f CSV 地质灾害点.csv 地质灾害点.shp -lco GEOMETRYAS_XY加上GEOMETRYAS_XY导出的CSV会附带X和Y坐标字段直接用pandas读取也不会有乱码问题。正常情况下.dbf默认编码可能是GBK或UTF-8在ogr2ogr转换时可以用-lco ENCODINGUTF-8来统一编码减少后续文件乱码的麻烦。3. 数据处理的第一步清洗、校验与空间化质控拿到25万个点我建议不要急着做任何华丽的可视化先做清洗。地质灾害数据来源复杂历史记录、群众上报、遥感解译、野外调查都可能有同一个灾害点被重复记录、经纬度填反、时间格式混乱等情况非常普遍。3.1 重复点、无效点、类型缺失这些坑必须在分析前填平重复点是最常见的问题。比如某地滑坡在历史台账里算一条在遥感解译里又算一条两套数据合并后就出现坐标几乎完全一致但属性不同的重叠点。处理重复点的逻辑不应该是简单的“坐标相同就删”而要考虑业务含义。我建议至少做两层去重第一层完全重复灾害类型、经纬度、发生年份、地点名称完全相同直接保留一条。第二层疑似重复坐标差值小于100米且灾害类型相同、年份相同但名称或来源字段不一致。这类数据需要通过人工抽查或空间邻接分析确认不建议直接批量删除。无效点主要分三类经纬度为0或空值、经纬度超出中国国界范围、经纬度中经度和纬度填反。经纬度填反的问题很隐蔽全国点数据如果省字段存在可以通过粗定位来判断比如某点记录的省份是四川但坐标落在新疆那大概率是经纬度写反了。处理方法是先按“经度在73—135纬度在18—53”做初筛再用省界缓冲圈进一步校验。类型缺失的情况也不少见。历史早期记录里可能只写了“山体滑坡”或“泥石流”后来统一归纳时字段值会有细微差异比如“泥石流”被记成“泥流”、“泥石流沟”。这时候需要写一个字段映射字典把所有相似名称归到7大类型中。import pandas as pd df pd.read_csv(地质灾害点.csv) type_map { 滑坡: 滑坡, 山体滑坡: 滑坡, 崩塌: 崩塌, 危岩体: 崩塌, 泥石流: 泥石流, 泥流: 泥石流, 地面塌陷: 地面塌陷, 岩溶塌陷: 地面塌陷, 采空塌陷: 地面塌陷, 地裂缝: 地裂缝, 地面沉降: 地面沉降, 沉降: 地面沉降, 不稳定斜坡: 不稳定斜坡, 滑坡隐患: 不稳定斜坡, } df[灾害类型] df[原始类型].map(type_map) # 检查映射遗漏 print(df[df[灾害类型].isna()][原始类型].unique())这段代码的重点是最后一行单独列出没有映射成功的原始类型然后人工判断归类。这一步多做几分钟后面统计就少很多麻烦。3.2 时间字段标准化1949—2019的时序分析怎么做时间字段是这套数据里的“第二灵魂”。但灾害点的时间记录格式堪称混乱有的只有年份有的写“1983年6月”有的精确到日甚至还有“1980年代”“民国三十八年”这种描述。要做时间维度的趋势分析必须先把时间字段转换成标准格式。我的做法是拆成三个字段年份、月份、日期均保留为整数。只有年份的记录月份和日期设为空值如果是连续时段或模糊描述则标记为特殊编码比如“1980-1985”取起始年份1985但增加一个时间精度字段说明该记录的可靠程度。这样后期如果只想用高精度记录可以一键过滤。import pandas as pd df[发生时间] df[发生时间].astype(str).str.strip() def parse_year(t): # 提取4位年份 import re match re.search(r(19|20)\d{2}, t) return int(match.group()) if match else None df[年份] df[发生时间].apply(parse_year)需要注意的是1949—2019这个跨度横跨了多次行政区划调整尤其是地级市和县级的边界变化很大。做“某县70年灾害数量变化”时如果直接用现在的区划边界去切历史灾害点会得出很奇怪的结论比如老县城的数据被算到邻县头上。稳妥的办法是先不切区划用点数据的“所在县”字段做汇总如果必须用边界至少要在说明中标注“按现行行政区划统计”。3.3 叠加行政区划与易发区让点数据真正落到管理单元上清洗完的点数据下一步通常就是叠加分析。最朴素也最常用的操作是把灾害点和行政区划边界做空间连接得到每个县区的灾害点总数和第1类数量统计。这样做的目的是让决策者能按辖区去认领风险。在QGIS里操作很简单右键点图层 → 属性 → 连接选择行政区划面图层空间关系选“包含”。或者用geopandas做spatial joinimport geopandas as gpd disaster gpd.read_file(地质灾害点_wgs84.shp, encodingutf-8) county gpd.read_file(全国县级行政区划.shp, encodingutf-8) # 需要先确保两个图层坐标系一致 county county.to_crs(disaster.crs) # 空间连接 joined gpd.sjoin(disaster, county, howleft, predicatewithin) # 按区县和灾害类型统计数量 result joined.groupby([县级名称, 灾害类型]).size().unstack(fill_value0) result.to_excel(县级灾害点统计.xlsx)我特别提醒一个小细节predicatewithin表示点必须完全落在面内部如果要兼顾边界情况可用intersects。但全国数据里不少灾害点坐落在县界边缘因为坐标精度不高可能刚好落在界外导致归属到隔壁县。这种情况下街道地址或“所在县”字段的优先级往往比空间位置更高建议以属性字段为主空间连接结果仅作参考。4. 25万点怎么画才能说清问题从密度图到风险判读数据清洗好之后进入分析阶段。25万个点直接画在图上比例尺小一点就会糊成一团黑色。这时候就需要“降维表达”。4.1 核密度分析快速找出“灾害窝”核密度分析是地质灾害点分析最常用的手段之一。它可以基于点事件的分布生成一个连续的表面表示单位面积内的灾害点密集程度。简单理解就是想象每个灾害点向四周扩散“影响力”越密集的区域颜色越深。在QGIS中使用“热力图Kernel Density Estimation”插件核密度半径是个关键参数。对全国尺度建议设置为50—100公里对省尺度设置为10—30公里对县域尺度设置为1—5公里。半径太大会过度平滑看不出区域差异半径太小则容易过拟合看上去全是孤立的亮斑。核密度的输出是栅格文件制作专题图时建议再做一次分位分类而不是用等间距分类。因为灾害点分布极不均匀西部山区和中东部平原的数据密度差距很大等间距分类会导致绝大多数区域颜色一样只有少数密集区是深色看图效果很差。4.2 分类型分时段制图25万点不能画成一个色如果只做一张“全国灾害点分布图”然后打上一堆黑点这张图的信息量就太小了。更好的做法是分类型、分时段分别制图。我的常用做法是做成“多面板”的系列图比如第一张滑坡灾害点密度图第二张崩塌灾害点密度图第三张泥石流灾害点密度图第四张1949—1999年灾害点分布第五张2000—2019年灾害点分布。这样每一张图都在回答一个具体问题“滑坡集中在哪”“近20年的灾害点比前50年多在哪里” 分类型制图时不同灾害类型的色系尽量用对比度高一点的比如滑坡用红、崩塌用橙、泥石流用蓝避免用相近色系否则图例很难分辨。对于全量25万个点如果还想看单点位置建议在出图时只控制到大比例尺局部区域。全国一张图上强行显示所有点既没有可读性也会撑爆符号渲染。一个可行的替代方案是用六边形网格聚合——把全国划分成0.5° x 0.5°的六边形统计每个网格内的灾害点数量再按网格填色这样既能保留空间趋势又能规避点重叠。4.3 风险区划初探灾害点叠加基础地理与人口经济数据灾害点分布不直接等于风险还应该叠加“受体”也就是人口、道路、建筑、耕地的暴露度。常见做法是构建一个简单的打分模型历史灾害点密度0—100分按核密度分位赋分降雨因子多年平均降雨量或暴雨日数0—100分地形因子坡度或地形起伏度0—100分人口或承灾体密度0—100分。将4项权重相加得到风险指数。实际操作时我用的是QGIS栅格计算器把所有因子统一重采样到1公里分辨率然后按权重相加。注意这只是一个基于历史频度和静态因子的极简风险模型不能直接用于实际防灾应急预案但用来做区域优先排序、确定排查重点是很快捷的底图工具。这里还要提醒一下灾害点数据代表“曾经发生过”或“已发现隐患”不代表“未来高风险必然在这些点上”。有不少区域虽然历史灾害点少但地形、降雨、工程活动条件非常脆弱风险同样很高。所以在做风险区划时一定要把灾害点数据与孕灾环境因子结合而不是只依赖点密度。5. 我踩过的几个坑和给后来者的建议最后分享几个实操中遇到的真实问题这些坑在文档里基本不会写但一旦踩进去往往要花大量时间填坑。第一个坑是整库传输的脆弱性。Shapefile三个主文件必须一起走但很多协作平台会默认只解析.shp导致对方收到文件夹后其他文件被过滤掉。后来我学乖了无论发什么矢量数据一律打包成zip再传输。如果对方只需要属性我会顺手导出一份CSV双保险。第二个坑是坐标系标注错误。某次我拿到一个数据源.prj文件写着GCS_WGS_1984但用卫星底图一叠点整体偏离约600米。检查后发现实际是通过火星坐标系偏移过的坐标只是.prj没更新。这种问题只有拿高精度影像底图抽查才能发现。建议在任何正式分析前先在在线卫星图上抽查20个点肉眼看看点是否落在合理的坡体、沟谷或村庄附近。第三个坑是年份字段里的“看不见的异常”。比如Excel处理日期后自动变成数字序列或者“1949-00-00”这种非法日期在转成Shapefile属性时变成空值。更隐蔽的是“2019/12/32”这种完全无效的日期被强制读成了2019年的某些数值。清洗时一定要先把时间字段转成字符串再写正则提取年份避免Excel自动转换带来的二次污染。第四个坑是数据量带来的软件卡顿。25万点在全国范围加一个符号化连平缩放都可能卡成幻灯片。我现在的习惯是保留全量数据作为底库用于分析时按需要切片而不是一次性全量加载到地图工程里。比如做省级分析就只导出该省的点做县级再进一步筛。这样既保数据精度又不卡软件。第五个坑是字段语义理解错误。不同批次的调查资料中“威胁人口”可能包含威胁常住人口、威胁流动人口甚至有的是威胁房屋间接估算人口。直接把这些字段混淆在一起做统计会导致某些市县的威胁人数异常翻倍。建议使用这个数据时先花20分钟读字段定义和源数据说明尤其是“规模等级”“威胁人口”这类与排险决策直接相关的字段不要想当然。这套25万地质灾害点数据配合Shapefile矢量格式非常适合做全国或大区域尺度的灾害分布规律研究、易发区划分和风险底图绘制。但归根结底数据只是工具真正决定分析质量的是清洗流程、坐标系管理和对灾害类型背后地学机制的理解。希望这篇内容能让你避开我走过的弯路拿到数据后少一点踩坑多一点有效产出。