R中数据框的空间子集
Spatial subset of data frame in R
我有一个使用 plot()
绘制的大数据框。然后我使用了:
library(splancs)
polygon_xy = getpoly(quiet=FALSE)
并在 select 我感兴趣的区域的地块上画点。这生成了我绘制的多边形的 x,y 坐标。
我想提取位于多边形内的数据,或将我的 df 子集化以仅包含位于多边形内的点。有什么建议吗?
当您使用 plot
和 getpoly
识别坐标时,点数据需要采用特定格式(即具有 x 和 y 的矩阵)。
library(splancs)
library(tidyverse)
library(sf)
set.seed(543)
xy <-
cbind(x = runif(n = 25, min = -118, max = -117),
y = runif(n = 25, min = 40, max = 42))
plot(xy)
# Draw a polygon for study area.
poly <- getpoly()
# Convert to sf objects.
polysf <- st_as_sf(as.data.frame(poly), coords = c("V1", "V2"), crs = 4326) %>%
dplyr::summarise() %>%
st_cast("POLYGON") %>%
st_convex_hull()
xysf <- st_as_sf(as.data.frame(xy), coords = c("x", "y"), crs = 4326)
# Do an intersection to keep only points inside the drawn polygon.
xy_intersect <- st_intersection(polysf, xysf)
输出
Simple feature collection with 9 features and 0 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: -117.7913 ymin: 40.82405 xmax: -117.4264 ymax: 41.7448
Geodetic CRS: WGS 84
geometry
1 POINT (-117.4264 41.18712)
2 POINT (-117.5756 41.7448)
3 POINT (-117.7913 40.82405)
4 POINT (-117.7032 41.15077)
5 POINT (-117.5634 41.23936)
6 POINT (-117.7441 40.84163)
7 POINT (-117.692 41.27514)
8 POINT (-117.6864 40.98462)
9 POINT (-117.5759 40.88477)
用 mapview::mapview(xy_intersect)
从 library(mapview)
绘制
但是,如果您想从原始数据框中提取行,那么这里还有另一种技巧,用于提取落在绘制的多边形内的点(例如,当多边形坐标看起来像 0.003456 时)。
library(splancs)
library(tidyverse)
set.seed(543)
xy <-
cbind(x = runif(n = 25, min = -118, max = -117),
y = runif(n = 25, min = 40, max = 42))
plot(xy)
# Draw a polygon for study area.
poly <- getpoly()
# Plot the results.
plot(xy)
polygon(poly)
# This will return a logical vector for points in the polygon
io <- inout(xy, poly)
points(xy[io,], pch = 16, col = "blue")
# Then, can use the index from io to extract the points that
# are inside the polygon from the original set of points.
extract_points <- as.data.frame(xy)[which(io == TRUE),]
extract_points
输出
x y
2 -117.4506 41.17794
3 -117.4829 40.71030
8 -117.4679 40.71702
19 -117.3354 40.53687
21 -117.5219 40.47077
22 -117.4876 40.18188
25 -117.2015 40.86243
我有一个使用 plot()
绘制的大数据框。然后我使用了:
library(splancs)
polygon_xy = getpoly(quiet=FALSE)
并在 select 我感兴趣的区域的地块上画点。这生成了我绘制的多边形的 x,y 坐标。
我想提取位于多边形内的数据,或将我的 df 子集化以仅包含位于多边形内的点。有什么建议吗?
当您使用 plot
和 getpoly
识别坐标时,点数据需要采用特定格式(即具有 x 和 y 的矩阵)。
library(splancs)
library(tidyverse)
library(sf)
set.seed(543)
xy <-
cbind(x = runif(n = 25, min = -118, max = -117),
y = runif(n = 25, min = 40, max = 42))
plot(xy)
# Draw a polygon for study area.
poly <- getpoly()
# Convert to sf objects.
polysf <- st_as_sf(as.data.frame(poly), coords = c("V1", "V2"), crs = 4326) %>%
dplyr::summarise() %>%
st_cast("POLYGON") %>%
st_convex_hull()
xysf <- st_as_sf(as.data.frame(xy), coords = c("x", "y"), crs = 4326)
# Do an intersection to keep only points inside the drawn polygon.
xy_intersect <- st_intersection(polysf, xysf)
输出
Simple feature collection with 9 features and 0 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: -117.7913 ymin: 40.82405 xmax: -117.4264 ymax: 41.7448
Geodetic CRS: WGS 84
geometry
1 POINT (-117.4264 41.18712)
2 POINT (-117.5756 41.7448)
3 POINT (-117.7913 40.82405)
4 POINT (-117.7032 41.15077)
5 POINT (-117.5634 41.23936)
6 POINT (-117.7441 40.84163)
7 POINT (-117.692 41.27514)
8 POINT (-117.6864 40.98462)
9 POINT (-117.5759 40.88477)
用 mapview::mapview(xy_intersect)
从 library(mapview)
但是,如果您想从原始数据框中提取行,那么这里还有另一种技巧,用于提取落在绘制的多边形内的点(例如,当多边形坐标看起来像 0.003456 时)。
library(splancs)
library(tidyverse)
set.seed(543)
xy <-
cbind(x = runif(n = 25, min = -118, max = -117),
y = runif(n = 25, min = 40, max = 42))
plot(xy)
# Draw a polygon for study area.
poly <- getpoly()
# Plot the results.
plot(xy)
polygon(poly)
# This will return a logical vector for points in the polygon
io <- inout(xy, poly)
points(xy[io,], pch = 16, col = "blue")
# Then, can use the index from io to extract the points that
# are inside the polygon from the original set of points.
extract_points <- as.data.frame(xy)[which(io == TRUE),]
extract_points
输出
x y
2 -117.4506 41.17794
3 -117.4829 40.71030
8 -117.4679 40.71702
19 -117.3354 40.53687
21 -117.5219 40.47077
22 -117.4876 40.18188
25 -117.2015 40.86243