使用R的sf包识别多边形边界线：如何排除内部接触线？

问题背景 #

我正在用R的sf包处理新泽西州（NJ）的空间数据，目标是找出构成各普查区（多边形）边界的道路（线要素）——NJ很多普查区的边界完全由道路构成。

我尝试用st_touches()做空间连接，这个函数本应识别仅共享边界点、内部无交集的几何对象。结果确实找到了构成边界的道路（示例中红色的4条），但同时误抓了普查区内部、仅端点接触边界的道路（示例中蓝色的3条）。

核心疑问解答 #

1. 多边形与线的空间连接是否把多边形视为闭合线？ #

不是。sf中多边形是带内部区域的面要素，空间关系判定会明确区分内部、边界、外部三个部分。只要道路的任何部分（端点、线段）接触到多边形的边界，无论道路主体在内部还是外部，只要满足对应空间函数的规则，就会被标记。

2. 用st_touches时，内部的线本不应被识别？ #

理论上是的——st_touches()的官方规则是：两个几何的内部无交集，但边界有交集。如果内部道路的主体在普查区内部，那道路的内部与普查区的内部是相交的，st_touches()应该返回FALSE。

但你遇到的误识别，大概率是以下两种情况：

• 数据精度问题：内部道路的端点刚好落在普查区边界上，sf的浮点数精度判定把“端点接触边界”判定为满足边界交集，同时误判道路内部与普查区内部无交集；
• 道路数据的实际情况：部分内部道路的一小段刚好与普查区边界重合（而非仅端点接触），触发了st_touches()的判定。

解决方案 #

方法1：基于边界线的交集长度过滤 #

直接提取普查区的边界线，通过计算道路与边界线的交集长度，过滤掉仅端点接触（交集长度为0）的道路：

library(sf)
library(dplyr)

# 假设已有：
# tracts：NJ普查区多边形要素（带唯一标识列，比如GEOID）
# roads：道路线要素

# 1. 提取普查区的边界线
tract_bounds <- st_boundary(tracts)

# 2. 先通过st_touches初步关联道路与普查区
road_tract_pre <- st_join(roads, tracts, join = st_touches)

# 3. 过滤：仅保留与对应普查区边界有实际线段重合的道路
road_tract_final <- road_tract_pre %>%
  rowwise() %>%
  # 匹配当前道路对应的普查区边界
  mutate(match_bound = tract_bounds[tract_bounds$GEOID == GEOID, ]) %>%
  # 计算道路与边界线的交集长度
  mutate(intersect_len = st_length(st_intersection(geometry, match_bound))) %>%
  # 过滤掉仅端点接触的情况（交集长度为0）
  filter(intersect_len > units::set_units(0, "meter")) %>%
  ungroup() %>%
  # 清理临时列
  select(-match_bound, -intersect_len)

方法2：用DE-9IM自定义空间关系 #

用st_relate()指定严格的DE-9IM模式，直接筛选出仅与普查区边界有线段重合、完全不在普查区内部的道路：

# DE-9IM模式"FF*F1****"的含义：
# - 道路内部与普查区内部无交集（第一个F）
# - 道路与普查区边界的交集维度为1（即有线段重合，第五个字符1）
boundary_roads <- roads[st_relate(roads, tracts, pattern = "FF*F1****"), ]

# 关联对应的普查区信息
boundary_roads_with_tract <- st_join(boundary_roads, tracts, join = st_intersects)

验证建议 #

可以用st_intersection()查看误识别道路与普查区的交集类型：如果是POINT，说明是仅端点接触；如果是LINESTRING，说明道路确实有一段在普查区边界上。

内容的提问来源于stack exchange，提问作者Aaron R.