Python中合并GeoDataFrame多边形以满足人口阈值需求

1. 成熟算法 #

• 贪心区域聚合算法：此类场景的主流方案。核心逻辑是每次定位人口最少且未达阈值的区域，将其与人口最多的相邻区域合并，重复操作直到所有区域满足阈值。该策略能最大化保留高人口区域的独立性，保证结果粒度尽可能细。
• 自底向上层次聚类：将单个区域作为初始簇，计算相邻簇的合并成本（如人口增量、几何复杂度变化），优先合并成本最低的相邻簇，直到所有簇的人口达标。适合对合并“代价”有特定要求的场景。

2. Python可用工具 #

• GeoPandas：核心空间数据处理库，提供邻接查询（sjoin）、多边形合并（dissolve/unary_union）等基础能力，是实现的核心依赖。
• libpysal：专业空间分析库，可快速计算Rook/Queen邻接矩阵（共享边/顶点即为相邻），省去手动判断邻接的繁琐，大幅提升效率。
• pandas：配合GeoPandas完成人口数据的聚合、筛选等操作。

3. 最佳实践与注意事项 #

• 预构建空间索引：对GeoDataFrame的geometry列构建RTree索引（gdf.sindex），查询相邻区域时用sindex.query()替代全表遍历，显著提升速度。
• 增量更新邻接关系：每次合并后，无需重新计算全量邻接，仅更新被合并区域及其邻居的邻接信息，减少计算开销。
• 优先处理最小人口区域：贪心策略下，先合并人口最少的区域，能最大程度保留高人口区域的独立性，避免过度合并。
• 修复拓扑问题：合并后可能出现无效几何，用shapely.make_valid()或geometry.buffer(0)修复，保证结果的空间有效性。
• 用向量化操作替代循环：尽量使用GeoPandas内置的向量化方法（如dissolve、groupby），避免Python循环，提升处理效率。

代码实现示例 #

import geopandas as gpd
import libpysal
from shapely.geometry import Polygon

# 构造示例GeoDataFrame（替换为你的实际数据）
sample_data = {
    'population': [1000, 1500, 800, 2000, 900],
    'geometry': [
        Polygon([(0,0), (0,1), (1,1), (1,0)]),
        Polygon([(1,0), (1,1), (2,1), (2,0)]),
        Polygon([(0,1), (0,2), (1,2), (1,1)]),
        Polygon([(2,0), (2,1), (3,1), (3,0)]),
        Polygon([(1,1), (1,2), (2,2), (2,1)])
    ]
}
gdf = gpd.GeoDataFrame(sample_data, crs="EPSG:4326")

POP_THRESHOLD = 2500  # 人口阈值

def merge_by_pop_threshold(gdf, threshold):
    merged = gdf.copy().reset_index(drop=True)
    merged['cluster_id'] = range(len(merged))
    
    while True:
        # 计算每个簇的总人口和合并后的几何
        cluster_summary = merged.groupby('cluster_id').agg(
            total_pop=('population', 'sum'),
            geometry=('geometry', 'unary_union')
        ).reset_index()
        
        # 筛选未达阈值的簇
        under_threshold = cluster_summary[cluster_summary['total_pop'] < threshold]
        if under_threshold.empty:
            break
        
        # 选取人口最少的簇作为合并起点
        target_cluster = under_threshold.sort_values('total_pop').iloc[0]
        target_id = target_cluster['cluster_id']
        
        # 计算当前所有簇的邻接关系（Rook规则：共享边为相邻）
        w = libpysal.weights.Rook.from_dataframe(merged)
        neighbors = w.neighbors.get(target_id, [])
        
        if not neighbors:
            # 孤立区域无法合并，直接保留
            break
        
        # 选择人口最多的邻居进行合并（可根据需求调整合并策略）
        neighbor_pops = merged[merged['cluster_id'].isin(neighbors)].groupby('cluster_id')['population'].sum()
        best_neighbor_id = neighbor_pops.sort_values(ascending=False).index[0]
        
        # 合并簇ID
        merged.loc[merged['cluster_id'] == best_neighbor_id, 'cluster_id'] = target_id
        
        # 重新聚合簇数据
        merged = merged.groupby('cluster_id').agg(
            population=('population', 'sum'),
            geometry=('geometry', 'unary_union')
        ).reset_index()
    
    # 修复无效几何
    merged['geometry'] = merged['geometry'].make_valid()
    return merged

# 执行合并
result_gdf = merge_by_pop_threshold(gdf, POP_THRESHOLD)
print(result_gdf)

额外说明 #

• 代码中使用Rook邻接规则，若需要将共享顶点的区域视为相邻，替换为libpysal.weights.Queen.from_dataframe(merged)。
• 合并策略可灵活调整：比如选择合并后总人口最接近阈值的邻居，而非人口最多的，只需修改best_neighbor_id的筛选逻辑。
• 针对大规模数据集，可优化邻接关系的计算逻辑，仅针对活跃簇（未达阈值的簇）计算邻接，进一步提升效率。

内容的提问来源于stack exchange，提问作者J.Tmr