转换几何列并保留其他列

Question

我已将 shapefile (.shp) 导入 R 中以定义一组多边形。 shapefile 采用分组数据框的形式，其中包含一个

geometry

列，其中包含 MULTIPOLYGON XY 坐标列表。我已向此数据框添加了一些额外的数据列：

ID   Value   geometry
A    5       MULTIPOLYGON (((519270 1639...
B    10      MULTIPOLYGON (((519553 1642...
C    25      MULTIPOLYGON (((519000 1666...

这对于使用

ggplot

和

geom_sf

生成静态分区统计图效果很好。但是，要使用

leaflet

生成动态地图，我需要将定义多边形的 XY 坐标转换为纬度/经度坐标。

我尝试使用

st_transform

来执行此操作，虽然它确实生成了一组多边形，但我无法为这些多边形着色，因为没有与它们关联的值：

inter_map2 <- leaflet(df) %>%
  addProviderTiles(providers$CartoDB.Positron) %>%
  addPolygons(data = st_transform(df$geometry, 4326), stroke = TRUE, weight = 0.5, opacity = 1, fillOpacity = 0.05, highlightOptions = highlightOptions(color = "red", weight = 2, bringToFront = TRUE))

我也尝试过使用

extract

但留下了结果的纬度和长列与

NA

值：

df_latlong <- df %>%
  extract(geometry, c('lat', 'lon'), '\\((.*), (.*)\\)', convert = TRUE)

有没有办法做到这一点，同时还保留与多边形关联的

ID

和

Value

列？

Answer 1

sf

对象不仅仅是一个常规的 tibble / 数据框，其坐标位于名为

geometry

的列中，属性和几何图形应该一起使用，在大多数情况下没有理由强制拆分它们。

extract()

不起作用，因为几何图形并不是真正的字符串，打印它只是将几何图形输出为 WKT（众所周知的文本）。通过调用

st_transform(df$geometry, 4326)

，您确实只剩下转换后的几何列表，属性列被删除；它有其用途，但更常见的是变换 sf 对象而不仅仅是几何列，即执行

st_transform(df, 4326)

。

data

中的

leaflet()

参数设置默认数据对象，在您的情况下它是

df

。在以下

addPolygons()

调用中，您不需要覆盖它，除非您想添加更多层，只需使用公式符号 (~) 和引用

df

列。像这样的东西：

library(leaflet)
library(sf)

# sf exampl shapefile, tranform to projected CRS
nc <- st_read(system.file("shape/nc.shp", package="sf")) %>% st_transform(6542)
nc[,1:2]
#> Simple feature collection with 100 features and 2 fields
#> Geometry type: MULTIPOLYGON
#> Dimension:     XY
#> Bounding box:  xmin: 123829.8 ymin: 14740.06 xmax: 930518.6 ymax: 318255.5
#> Projected CRS: NAD83(2011) / North Carolina
#> First 10 features:
#>     AREA PERIMETER                       geometry
#> 1  0.114     1.442 MULTIPOLYGON (((387344.7 27...
#> 2  0.061     1.231 MULTIPOLYGON (((408601.4 29...
#> 3  0.143     1.630 MULTIPOLYGON (((478715.7 27...
#> 4  0.070     2.968 MULTIPOLYGON (((878193.4 28...
#> 5  0.153     2.206 MULTIPOLYGON (((769834.9 27...
#> 6  0.097     1.670 MULTIPOLYGON (((812327.7 27...
#> 7  0.062     1.547 MULTIPOLYGON (((878193.4 28...
#> 8  0.091     1.284 MULTIPOLYGON (((828444.5 29...
#> 9  0.118     1.421 MULTIPOLYGON (((671746.3 27...
#> 10 0.124     1.428 MULTIPOLYGON (((517435.1 27...

# pallete for leaflet from nc$AREA values
pal <- colorNumeric(
  palette = "Greens",
  domain = nc$AREA
)

# pass transformed sf object directly to leaflet and in 
# addPolygons() reference its columns though formula (~) interface:
nc %>% 
  st_transform("WGS84") %>% 
  leaflet() %>%
  addTiles() %>% 
  addPolygons(color = ~pal(AREA), stroke = FALSE, fillOpacity = 1)

^{创建于 2023-06-14，使用 reprex v2.0.2}