Geopandas to_crs转换后记录缺失问题求助

我有一个包含POLYGON和MULTIPOLYGON几何数据的Geopandas GeoDataFrame，尝试将其从其他坐标参考系统（CRS）转换为EPSG:4326。由于该GeoDataFrame约有20万条记录，我采取了以下步骤：

• 将完整GeoDataFrame拆分为200个各约1000条记录的小GeoDataFrame
• 执行small_gdf.to_crs('epsg:4326', inplace=True)进行转换
• 导出每个小GeoDataFrame：small_gdf.to_file(f'small_gdf_{filecounter}.shp')

转换过程耗时约2天，但使用pd.concat合并所有小GeoDataFrame后，结果仅保留了原数据约60%的记录。请问是否会因为to_crs转换失败导致记录丢失？
同时，我打算为每个小GeoDataFrame新增一列，重新执行to_crs操作以追踪转换过程中丢失的记录。

原代码示例 #

import geopandas as gpd
gdf = gpd.read_file('bigShapefilePath.shp')
n_records = len(gdf)

# 创建分块的起止索引元组
chunksize = 1000
i=0
list_start_end_idx_tuples = []
for start in range(i, n_records, chunksize):
    end = start+999
    if end > n_records:
        end = n_records-1

    start_end_idx_tuple = (start, end)
    list_start_end_idx_tuples.append(start_end_idx_tuple)

# 分块转换CRS并导出
parts_folderpath = <parts_folderpath>
file_counter=1
for each_start_end in list_start_end_idx_tuples:
    start, end = each_start_end
    small_gdf = gdf.iloc[start:end+1]
    small_gdf['WITHIN_PART_IDX'] = range(len(small_gdf))
    small_gdf.to_crs('epsg:4326', inplace=True)
    small_gdf.to_file(f'{parts_folderpath}/small_gdf_part{file_counter}.shp')

    file_counter+=1


# 查找分块文件
full_folderpath = <full_folderpath>
i=0
list_smallgdf_filename = []
list_smallgdf_filenamenext = []

for dir, subdir, filenames in os.walk(parts_folderpath):
    for filenamenext in filenames:
        if ('.shp' in filenamenext) and ('.xml' not in filenamenext):
            filename = filenamenext.split('.')[0]
            i+=1
            list_smallgdf_filename.append(filename)
            list_smallgdf_filenamenext.append(filenamenext)


# 合并分块文件
i=0
for filenamenext in list_smallgdf_filenamenext:
    small_gdf = gpd.read_file(f'{parts_folderpath}/{filenamenext}')
    small_filename = list_smallgdf_filename[i]
    part_num = small_filename.split('_')[-1]
    small_gdf['PART_NUM'] = int(part_num)
    
    if i<1:
        concat_gdf = small_gdf
    else:
        concat_gdf = pd.concat([concat_gdf, small_gdf])
    i+=1

concat_gdf.to_file(f'{full_folderpath}/concat_gdf.shp')

问题解答 #

1. to_crs转换是否会导致记录丢失？ #

是的，to_crs转换失败确实可能引发记录丢失，核心原因包括：

• 几何拓扑错误：原数据中的POLYGON/MULTIPOLYGON存在自相交、无效环等问题时，坐标转换会失败，生成空几何对象。而Shapefile格式默认会过滤掉空几何记录，导致导出后这部分记录直接丢失。
• CRS参数不匹配：如果原CRS的定义不准确（比如缺少投影参数），转换过程会出现异常，同样可能产生空几何。
• 分块IO异常：磁盘空间不足、文件写入中断等IO错误会导致部分分块文件不完整，合并时自然缺失对应记录。

2. 追踪丢失记录的改进方案 #

核心思路

保留原数据的唯一标识（比如原始索引），在转换前后对比索引定位丢失记录；同时添加几何有效性检查，标记转换失败的记录而非直接丢弃。

优化后的代码

import geopandas as gpd
import pandas as pd
import os

# 读取原始数据
gdf = gpd.read_file('bigShapefilePath.shp')
n_records = len(gdf)
chunksize = 1000
parts_folderpath = "<parts_folderpath>"
full_folderpath = "<full_folderpath>"

# 简化分块索引生成
chunk_ranges = [(i, min(i + chunksize, n_records)) for i in range(0, n_records, chunksize)]

# 分块转换并导出，保留原始索引
file_counter = 1
for start, end in chunk_ranges:
    small_gdf = gdf.iloc[start:end].copy()
    # 保留原始索引，用于后续对比
    small_gdf["ORIGINAL_INDEX"] = small_gdf.index
    # 转换CRS（避免inplace=True引发的潜在内存问题）
    small_gdf = small_gdf.to_crs("epsg:4326")
    # 标记几何有效性，排查转换失败的记录
    small_gdf["GEOM_VALID"] = small_gdf.geometry.is_valid
    # 导出所有记录（包括无效几何，避免Shapefile自动过滤）
    small_gdf.to_file(f"{parts_folderpath}/small_gdf_part{file_counter}.shp")
    file_counter += 1

# 合并分块文件（按文件名排序，确保顺序正确）
concat_gdf = gpd.GeoDataFrame()
for filename in sorted(os.listdir(parts_folderpath)):
    if filename.endswith(".shp") and not filename.endswith(".xml.shp"):
        small_gdf = gpd.read_file(f"{parts_folderpath}/{filename}")
        concat_gdf = pd.concat([concat_gdf, small_gdf], ignore_index=True)

# 对比原始索引，找出丢失的记录
original_indexes = set(gdf.index)
processed_indexes = set(concat_gdf["ORIGINAL_INDEX"])
missing_indexes = original_indexes - processed_indexes

print(f"丢失记录数：{len(missing_indexes)}")
print(f"丢失记录索引示例：{list(missing_indexes)[:10]}")

# 导出合并后的数据
concat_gdf.to_file(f"{full_folderpath}/concat_gdf.shp")

额外优化建议 #

• 预处理修复几何错误：转换前先检查并修复拓扑错误，比如用gdf.geometry = gdf.geometry.buffer(0)修复自相交多边形，减少转换失败的概率。
• 替换Shapefile为GeoParquet：Shapefile有单文件大小限制（最大2GB），且对空几何支持差，GeoParquet更适合大数据量存储，转换和读写速度更快。
• 添加日志记录：分块处理时添加日志（比如用logging模块），记录每个分块的处理状态、记录数，方便排查异常。

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