【问题标题】:r - Convert output from sf::st_within to vectorr - 将 sf::st_within 的输出转换为向量
【发布时间】:2018-08-23 23:56:06
【问题描述】:

我正在尝试使用 R 中的 sf 包来查看 sf 对象是否在另一个具有 st_within 函数的 sf 对象中。我的问题是这个函数的输出,它是稀疏几何二进制谓词-sgbp,我需要一个向量作为输出,以便之后可以使用dplyr 包进行过滤。这是一个简化的例子:

# object 1: I will test if it is inside object 2
df <- data.frame(lon = c(2.5, 3, 3.5), lat = c(2.5, 3, 3.5), var = 1) %>% 
st_as_sf(coords = c("lon", "lat"), dim = "XY") %>% st_set_crs(4326) %>%
  summarise(var = sum(var), do_union = F) %>% st_cast("LINESTRING")

# object 2: I will test if it contains object 1
box <- data.frame(lon = c(2, 4, 4, 2, 2), lat = c(2, 2, 4, 4,2), var = 1) %>%
  st_as_sf(coords = c("lon", "lat"), dim = "XY") %>% st_set_crs(4326) %>% 
  summarise(var = sum(var), do_union = F) %>% st_cast("POLYGON")

# test 1
df$indicator <- st_within(df$geometry, box$geometry) # gives geometric binary predicate on pairs of sf sets which cannot be used 
df <- df %>% filter(indicator == 1)

这给出了错误:列indicator 必须是一维原子向量或列表。

我尝试在下面解决这个问题:

# test 2
df$indicator <- st_within(df$geometry, box$geometry, sparse = F) %>% 
  diag() # gives matrix that I convert with diag() into vector
df <- df %>% filter(indicator == FALSE)

这行得通,它会删除包含 TRUE 值的行,但由于我的真实数据包含许多观察结果,因此制作矩阵的过程对于我的计算来说非常缓慢。有没有办法让st_within 的输出成为一个字符向量,或者有没有办法将sgbp 转换为与dplyr 兼容的字符向量,而无需制作矩阵?

【问题讨论】:

    标签: r dataframe dplyr gis sf


    【解决方案1】:

    以下是从稀疏几何二元谓词中获取逻辑向量的方法:

    df$indicator <- st_within(df, box) %>% lengths > 0
    

    或不创建新变量的子集:

    df <- df[st_within(df, box) %>% lengths > 0,]
    

    很遗憾,我无法在您的大型数据集上进行测试,但如果它比矩阵方法更快,请告诉我。

    【讨论】:

      【解决方案2】:

      is_within 的结果实际上是一个列表列,因此您可以使用您的 通过“取消列出”它来摆脱这种情况。这样的事情会起作用:

      library(dplyr)
      library(sf)
      
      # object 1: I will test if it is inside object 2 - to make this more interesting
      # I added a second not-contained line
      df <- data.frame(lon = c(2.5, 3, 3.5), lat = c(2.5, 3, 3.5), var = 1) %>% 
        st_as_sf(coords = c("lon", "lat"), dim = "XY") %>% st_set_crs(4326) %>%
        summarise(var = sum(var), do_union = F) %>% st_cast("LINESTRING")
      
      df2 <- data.frame(lon = c(4.5, 5, 6), lat = c(4.5, 5, 6), var = 2) %>% 
        st_as_sf(coords = c("lon", "lat"), dim = "XY") %>% st_set_crs(4326) %>%
        summarise(var = sum(var), do_union = F) %>% st_cast("LINESTRING")
      df3 <- rbind(df, df2)
      
      # object 2: I will test if it contains object 1
      box <- data.frame(lon = c(2, 4, 4, 2, 2), lat = c(2, 2, 4, 4,2), var = 1) %>%
        st_as_sf(coords = c("lon", "lat"), dim = "XY") %>% st_set_crs(4326) %>% 
        summarise(var = sum(var), do_union = F) %>% st_cast("POLYGON")
      
      plot(df3) 
      plot(st_geometry(box), add = TRUE)
      

      # see if the lines are within the box and build a data frame with results
      is_within <- st_within(df3$geometry, box$geometry) %>% 
        lapply(FUN = function(x) data.frame(ind = length(x))) %>% 
        bind_rows()
      
      # add the "indicator" to df3
      df3 <- dplyr::mutate(df3, indicator = is_within$ind) 
      df3
      #> Simple feature collection with 2 features and 2 fields
      #> geometry type:  LINESTRING
      #> dimension:      XY
      #> bbox:           xmin: 2.5 ymin: 2.5 xmax: 6 ymax: 6
      #> epsg (SRID):    4326
      #> proj4string:    +proj=longlat +datum=WGS84 +no_defs
      #>   var indicator                       geometry
      #> 1   3         1 LINESTRING (2.5 2.5, 3 3, 3...
      #> 2   6         0 LINESTRING (4.5 4.5, 5 5, 6 6)
      

      HTH

      reprex package (v0.2.0) 于 2018 年 3 月 15 日创建。

      【讨论】:

        【解决方案3】:

        不要直接使用st_within 函数,而是尝试使用spatial join。 查看以下示例 st_joins 的工作原理

        library(sf)
        library(tidyverse)
        
        lines <-
        data.frame(id=gl(3,2), x=c(-3,2,6,11,7,10), y=c(-1,6,-5,-9,10,5)) %>%
          st_as_sf(coords=c("x","y"), remove=F) %>% 
          group_by(id) %>% 
          summarise() %>%
          st_cast("LINESTRING")
        
        yta10 <-
            st_point(c(0, 0)) %>%
            st_buffer(dist = 10) %>%
            st_sfc() %>%
            st_sf(yta = "10m")
        

        使用左连接保留所有线,但您可以看到哪些线位于多边形内

        lines %>% st_join(yta10, left=TRUE)
        

        内部联接(left = FALSE)只保留内部联接

        lines %>% st_join(yta10, left=FALSE)
        

        后者也可以通过

        lines[yta10,]
        

        【讨论】:

        • st_join 只是合并线和框,但是它不会检查一个几何图形是否在另一个几何图形内,这就是我要问的。我需要一个逻辑或字符向量作为st_within 或可能来自另一个函数的输出
        • @Hans Gardfjell st_join 只检查 y 中的哪些几何图形完全在(包含在)x 内,还是检查重叠/相交的几何图形? R 中 he 'sf' 包的文档相当稀少。
        猜你喜欢
        • 2015-03-26
        • 2023-03-17
        • 2020-11-14
        • 1970-01-01
        • 1970-01-01
        • 2018-01-12
        • 1970-01-01
        • 1970-01-01
        • 2020-08-25
        相关资源
        最近更新 更多