使用显式二元谓词子集化sfc对象时的报错问题求助

问题背景 #

运行以下R代码时，最后两个子集化操作报错：

library(sf)
#> Linking to GEOS 3.11.2, GDAL 3.7.2, PROJ 9.3.0; sf_use_s2() is TRUE

p1 <- st_sfc(st_point(c(0, 0)))
p2 <- st_sfc(st_point(c(1, 1)))

p1[p2]
#> Geometry set for 0 features 
#> Bounding box:  xmin: NA ymin: NA xmax: NA ymax: NA
#> CRS:           NA
p1[p2, op = st_intersects]
#> Error in `[.default`(p1, p2, op = st_intersects): incorrect number of dimensions
p1[p2, , op = st_intersects]
#> Error in `[.default`(p1, p2, , op = st_intersects): incorrect number of dimensions

查看[.sfc的定义可知，其默认操作就是st_intersects：

sf:::`[.sfc`
#> function (x, i, j, ..., op = st_intersects) 
#> {
#>     if (!missing(i) && (inherits(i, "sf") || inherits(i, "sfc") || 
#>         inherits(i, "sfg"))) 
#>         i = lengths(op(x, i, ...)) != 0
#>     st_sfc(NextMethod(), crs = st_crs(x), precision = st_precision(x), 
#>         dim = if (length(x)) 
#>             class(x[[1]])[1]
#>         else "XY")
#> }
#> <bytecode: 0x00000255ad682df8>
#> <environment: namespace:sf>

报错原因 #

报错核心是参数传递的位置混淆：

• [.sfc的参数顺序为(x, i, j, ..., op = st_intersects)，其中j是针对二维结构的索引参数，但sfc是一维几何集合，不需要j参数。
• 当写p1[p2, op = st_intersects]时，R会将op = st_intersects解析为位置参数j的取值，而非命名参数op。此时j被赋予函数值，默认的[.default无法处理这种维度错误，因此抛出异常。
• 即使写p1[p2, , op = st_intersects]，旧版本sf的[.sfc逻辑仍会触发默认索引方法的维度检查，导致报错。

正确用法 #

方法1：提前生成逻辑索引 #

绕过参数歧义问题，手动生成子集化所需的逻辑向量：

# 用指定空间关系生成逻辑索引
idx <- lengths(st_intersects(p1, p2)) != 0
p1[idx]

如果要使用其他空间操作（比如st_disjoint），只需替换对应函数：

idx <- lengths(st_disjoint(p1, p2)) != 0
p1[idx]

方法2：直接传递op参数（修复后版本可用） #

sf包已对此问题进行调整修复，修复后的版本可直接通过命名参数传递op：

p1[p2, op = st_intersects]

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