如何将 lat/long 映射到 R 中 polygon/shapefile 的最近点
How can I map a lat/long to the nearest point of my polygon/shapefile in R
我正在尝试将一组 lats/longs 映射到我的形状文件以获得 CSDUID/census 分区信息。
但是,有一些异常值或 lat/longs 在 shapefile 之外。我想知道是否可以将这些点映射到 shapefile 的最近点,以便我可以为我的分析的后期部分获取划分信息?
可重现的例子
假设我有一个类似 gr
的 shapefile
nc <- st_transform(st_read(system.file("shape/nc.shp", package="sf")), 2264)
gr = st_sf(
label = apply(expand.grid(1:10, LETTERS[10:1])[,2:1], 1, paste0, collapse = " "),
geom = st_make_grid(nc))
# gr
# Simple feature collection with 100 features and 1 field
# geometry type: POLYGON
# dimension: XY
# bbox: xmin: 406262.2 ymin: 48374.87 xmax: 3052887 ymax: 1044158
# epsg (SRID): 2264
# proj4string: +proj=lcc +lat_1=36.16666666666666 +lat_2=34.33333333333334 +lat_0=33.75 +lon_0=-79 +x_0=609601.2192024384 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs
# First 10 features:
# label geom
# 1 J 1 POLYGON ((406262.2 48374.87...
# 2 J 2 POLYGON ((670924.7 48374.87...
# 3 J 3 POLYGON ((935587.2 48374.87...
# 4 J 4 POLYGON ((1200250 48374.87,...
# 5 J 5 POLYGON ((1464912 48374.87,...
# 6 J 6 POLYGON ((1729575 48374.87,...
# 7 J 7 POLYGON ((1994237 48374.87,...
# 8 J 8 POLYGON ((2258900 48374.87,...
# 9 J 9 POLYGON ((2523562 48374.87,...
# 10 J 10 POLYGON ((2788225 48374.87,...
我有一个点在边界框外:
point = st_point(c(106160,28370))
library(ggplot2)
ggplot() + geom_sf(data = gr) + geom_sf(data = point)
如何将点映射到最近的多边形?
您在 sf
程序包中具有强大的空间分析功能,尤其是空间连接功能。假设您的点数据被命名为 point
并且您的多边形被命名为.... polygons
(都是 sf
对象),您可以使用 sf::st_join
进行空间连接。
st_join
接受不同类型的空间映射。我认为您对 st_nearest_feature
合并感兴趣(首先删除落在多边形内的点,例如第一个 st_join
使用 st_within
)
library(sf)
st_join(
points,
polygons,
join = st_nearest_feature
)
在编辑后的示例中,将给出:
point = st_sf(st_sfc(point), crs = st_crs(gr))
st_join(point, gr, join = st_nearest_feature)
请参阅 ?st_join
文档。
我正在尝试将一组 lats/longs 映射到我的形状文件以获得 CSDUID/census 分区信息。
但是,有一些异常值或 lat/longs 在 shapefile 之外。我想知道是否可以将这些点映射到 shapefile 的最近点,以便我可以为我的分析的后期部分获取划分信息?
可重现的例子
假设我有一个类似 gr
nc <- st_transform(st_read(system.file("shape/nc.shp", package="sf")), 2264)
gr = st_sf(
label = apply(expand.grid(1:10, LETTERS[10:1])[,2:1], 1, paste0, collapse = " "),
geom = st_make_grid(nc))
# gr
# Simple feature collection with 100 features and 1 field
# geometry type: POLYGON
# dimension: XY
# bbox: xmin: 406262.2 ymin: 48374.87 xmax: 3052887 ymax: 1044158
# epsg (SRID): 2264
# proj4string: +proj=lcc +lat_1=36.16666666666666 +lat_2=34.33333333333334 +lat_0=33.75 +lon_0=-79 +x_0=609601.2192024384 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=us-ft +no_defs
# First 10 features:
# label geom
# 1 J 1 POLYGON ((406262.2 48374.87...
# 2 J 2 POLYGON ((670924.7 48374.87...
# 3 J 3 POLYGON ((935587.2 48374.87...
# 4 J 4 POLYGON ((1200250 48374.87,...
# 5 J 5 POLYGON ((1464912 48374.87,...
# 6 J 6 POLYGON ((1729575 48374.87,...
# 7 J 7 POLYGON ((1994237 48374.87,...
# 8 J 8 POLYGON ((2258900 48374.87,...
# 9 J 9 POLYGON ((2523562 48374.87,...
# 10 J 10 POLYGON ((2788225 48374.87,...
我有一个点在边界框外:
point = st_point(c(106160,28370))
library(ggplot2)
ggplot() + geom_sf(data = gr) + geom_sf(data = point)
如何将点映射到最近的多边形?
您在 sf
程序包中具有强大的空间分析功能,尤其是空间连接功能。假设您的点数据被命名为 point
并且您的多边形被命名为.... polygons
(都是 sf
对象),您可以使用 sf::st_join
进行空间连接。
st_join
接受不同类型的空间映射。我认为您对 st_nearest_feature
合并感兴趣(首先删除落在多边形内的点,例如第一个 st_join
使用 st_within
)
library(sf)
st_join(
points,
polygons,
join = st_nearest_feature
)
在编辑后的示例中,将给出:
point = st_sf(st_sfc(point), crs = st_crs(gr))
st_join(point, gr, join = st_nearest_feature)
请参阅 ?st_join
文档。