如何解决R sf包连接空间数据时出现的球面几何失效问题

报错原因 #

sf 1.0及以上版本默认对WGS84（EPSG:4326）坐标系启用s2球面几何引擎，对矢量数据的几何有效性要求远高于之前版本的平面几何计算逻辑。你遇到的报错是因为所用的shapefile多边形存在球面几何无效问题（如自交、多边形环方向不符合s2规范等），之前版本sf未启用s2因此可以正常运行。

解决方案 #

方案1：临时关闭s2引擎（快速适配老逻辑） #

直接退回到旧版sf的平面几何计算逻辑，无需修改原始数据，小范围区域计算误差可忽略，代码如下：

library(sf)
library(readr)
# 关闭s2球面几何引擎
sf_use_s2(FALSE)

shape <- read_sf("data/Provinces_v1_2017.shp")
data <- read_csv("data/data.csv")
pts <- st_as_sf(data, coords = c("dec_lon", "dec_lat"), crs = 4326)

# 执行空间连接
result <- st_join(pts, shape)

方案2：修复矢量几何（适配s2引擎，适合长期使用） #

如果需要保留s2球面计算的精度，可对shapefile的几何进行修复，调整多边形环方向符合s2规范，代码如下：

library(sf)
library(readr)
library(dplyr)

shape <- read_sf("data/Provinces_v1_2017.shp")
# 修复无效几何+调整环方向适配s2要求
shape_fixed <- shape %>%
  st_make_valid() %>%
  st_force_polygon_cw()

data <- read_csv("data/data.csv")
pts <- st_as_sf(data, coords = c("dec_lon", "dec_lat"), crs = 4326)

# 执行空间连接
result <- st_join(pts, shape_fixed)

补充说明 #

• 若修复后仍有报错，可先检查坐标点的经纬度是否在合理范围内：经度取值范围为-180_{180，纬度取值范围为-90}90，排除异常值干扰。
• 大范围跨半球的空间计算场景建议使用方案2，避免平面计算带来的精度误差。

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