在R中,如何使用Leaflet和SF创建基于某些变量进行颜色编码的空间网格

iugsix8n  于 2023-11-14  发布在  其他
关注(0)|答案(1)|浏览(98)

我正在和美国的一些州合作。我有经纬度数据(包含两个不同的列),我有第三列,我们称之为L1,包含从1到10.我想在州内创建一个正方形网格,然后使用色标填充网格,其中颜色由落入该正方形的所有地理点的L1平均值确定在电网上。
在运行下面的代码后,我收到一条错误消息,说“Error in to_ring.default(x):Don't know how to get polygon data from object of class XY,POINT,sfg”,这让我很困惑,因为我以为我传递的是网格坐标,每个网格都是一个多边形。但显然我传递的是点而不是多边形。请帮助。

# Create a simple feature object from your data
points <- st_as_sf(data, coords = c("Long", "Lat"), crs = 4326)

# Create a grid over the bounding box of your points
bb <- st_bbox(points)
grid <- st_make_grid(st_as_sfc(bb), cellsize = 0.1)
grid <- st_sf(grid, crs = st_crs(points))

# Calculate the average of L1 for each grid cell
grid$L1 <- aggregate(points["L1"], grid, function(x) mean(as.numeric(x), na.rm = TRUE))$L1

grid_sp <- as(grid, "Spatial")

# Create a color palette function
pal <- colorNumeric(palette = "viridis", domain = grid$L1)

# Create a map widget by calling leaflet()
m <- leaflet() %>%
  # Add a Tiles layer to the map
  addTiles() %>%
  # Add state boundaries to the map
  addPolygons(data = us_map, weight = 1, color = "#444444", fill = NA) %>%
  # Add the grid layer to the map
  addPolygons(data = grid, fillColor = ~pal(L1), fillOpacity = 0.8, 
              color = "#BDBDC3", weight = 1) %>%
  # Add a legend
  addLegend(pal = pal, values = ~L1, opacity = 0.8, title = "Average L1")

# Print the map
m

字符串

xzv2uavs

xzv2uavs1#

我认为你可以通过使用terra::rasterize()函数比使用sf::st_make_grid()更容易地完成你所追求的。注意,使用terra::方法,你将创建一个栅格网格,而不是一个矢量多边形网络,所以你对leaflet::addPolygons()的第二次调用应该是leaflet::addRasterImage()sf::st_make_grid()调用需要what="polygons"参数来创建一个多边形。我想这就是你收到的错误的来源。下面是一个使用terra::rasterize()工作流的可重现示例。在该示例中,我首先使用美国俄勒冈州的矢量版本,并创建一个空白栅格以适合其边界框和所需分辨率。然后创建一个随机点网络,并随机为它们分配1到10的L1值。最后使用terra::rasterize()函数将点聚合到栅格像元中,对每个栅格单元内的所有点取L1的平均值。我调整了对leaflet::的调用,以适应addRasterImage()函数的栅格。希望这是你想要的。

library(sf)
library(terra)
library(leaflet)
library(spData)

data("us_states")

##Get a simple feature of the US State of Oregon## 
oregon<-st_transform(us_states[us_states$NAME=="Oregon",], crs=4326)#transformed to WGS84, EPSG:4326

##Creating points within the State##
pts<-st_as_sf(st_sample(oregon, 50000, type="random"))
pts[,"L1"]<-sample(1:10, 50000, replace=TRUE)

##Create blank raster to populate##
OR_rast<-rast(extent=ext(oregon), ncols=500, nrows=200, crs="EPSG:4326")

#Rasterizing points' L1 values using the terra::rasterize() function##
L1_rast<-rasterize(x=pts, y=OR_rast, field="L1", fun="mean")

# Create a color palette function
pal <- colorNumeric(palette = "viridis", domain = 1:10)

# Create a map widget by calling leaflet()
m <- leaflet() %>%
  # Add a Tiles layer to the map
  addTiles() %>%
  # Add state boundaries to the map
  addPolygons(data = oregon, weight = 1, color = "#444444", fill = NA) %>%
  # Add the grid layer to the map
  addRasterImage(x = L1_rast, colors=pal, opacity = 0.8) %>%
  # Add a legend
  addLegend(pal = pal, values = 1:10, opacity = 0.8, title = "Average L1")

# Print the map
m

字符串

相关问题