library(sf)
library(spData)
## pointsDF: A data.frame whose first column contains longitudes and
## whose second column contains latitudes.
##
## states: An sf MULTIPOLYGON object with 50 states plus DC.
##
## name_col: Name of a column in `states` that supplies the states'
## names.
lonlat_to_state <- function(pointsDF,
states = spData::us_states,
name_col = "NAME") {
## Convert points data.frame to an sf POINTS object
pts <- st_as_sf(pointsDF, coords = 1:2, crs = 4326)
## Transform spatial data to some planar coordinate system
## (e.g. Web Mercator) as required for geometric operations
states <- st_transform(states, crs = 3857)
pts <- st_transform(pts, crs = 3857)
## Find names of state (if any) intersected by each point
state_names <- states[[name_col]]
ii <- as.integer(st_intersects(pts, states))
state_names[ii]
}
## Test the function with points in Wisconsin, Oregon, and France
testPoints <- data.frame(x = c(-90, -120, 0), y = c(44, 44, 44))
lonlat_to_state(testPoints)
## [1] "Wisconsin" "Oregon" NA
library(sp)
library(maps)
library(maptools)
# The single argument to this function, pointsDF, is a data.frame in which:
# - column 1 contains the longitude in degrees (negative in the US)
# - column 2 contains the latitude in degrees
lonlat_to_state_sp <- function(pointsDF) {
# Prepare SpatialPolygons object with one SpatialPolygon
# per state (plus DC, minus HI & AK)
states <- map('state', fill=TRUE, col="transparent", plot=FALSE)
IDs <- sapply(strsplit(states$names, ":"), function(x) x[1])
states_sp <- map2SpatialPolygons(states, IDs=IDs,
proj4string=CRS("+proj=longlat +datum=WGS84"))
# Convert pointsDF to a SpatialPoints object
pointsSP <- SpatialPoints(pointsDF,
proj4string=CRS("+proj=longlat +datum=WGS84"))
# Use 'over' to get _indices_ of the Polygons object containing each point
indices <- over(pointsSP, states_sp)
# Return the state names of the Polygons object containing each point
stateNames <- sapply(states_sp@polygons, function(x) x@ID)
stateNames[indices]
}
# Test the function using points in Wisconsin and Oregon.
testPoints <- data.frame(x = c(-90, -120), y = c(44, 44))
lonlat_to_state_sp(testPoints)
[1] "wisconsin" "oregon" # IT WORKS
library(sp)
library(rgdal)
#lat and long
Lat <- 57.25
Lon <- -9.41
#make a data frame
coords <- as.data.frame(cbind(Lon,Lat))
#and into Spatial
points <- SpatialPoints(coords)
#SpatialPolygonDataFrame - I'm using a shapefile of UK counties
counties <- readOGR(".", "uk_counties")
#assume same proj as shapefile!
proj4string(points) <- proj4string(counties)
#get county polygon point is in
result <- as.character(over(points, counties)$County_Name)
# library(remotes)
# install_github("JVAdams/jvamisc") # if you are installing this for the first time you will need to load the remotes package
library(jvamisc)
library(stringr)
#> Warning: package 'stringr' was built under R version 4.2.3
# Example Data
data <- data.frame(
longitude = c(-74.28000,-80.62036,-77.43923),
latitude = c(40.99194,33.82849,37.54588))
# Use function latlong2 from library(jvamisc) to convert lat and long points to state
data$state <- latlong2(data, to = 'state')
# Use function str_to_title form library(stringr) to make the first letter of each state uppercase
data$state <- str_to_title(data$state)
# Convert state name to state abbreviation
data$state_abb <- state.abb[match(data$state, state.name)]
data
#> longitude latitude state state_abb
#> 1 -74.28000 40.99194 New Jersey NJ
#> 2 -80.62036 33.82849 South Carolina SC
#> 3 -77.43923 37.54588 Virginia VA
library(maps)
library(sf)
## Get the states map, turn into sf object
US <- st_as_sf(map("state", plot = FALSE, fill = TRUE))
## Test the function using points in Wisconsin and Oregon
testPoints <- data.frame(x = c(-90, -120), y = c(44, 44))
# Make it a spatial dataframe, using the same coordinate system as the US spatial dataframe
testPoints <- st_as_sf(testPoints, coords = c("x", "y"), crs = st_crs(US))
#.. and perform a spatial join!
st_join(testPoints, US)
ID geometry
1 wisconsin POINT (-90 44)
2 oregon POINT (-120 44)
6条答案
按热度按时间6tdlim6h1#
这里有两个选项,一个使用sf,一个使用sp包函数。sf是用于分析空间数据的更现代的(并且在2020年,推荐)包,但如果它仍然有用,我将留下我最初的2012年答案,展示如何使用sp相关函数来完成此操作。
方法一(使用sf):
如果需要更高分辨率的州边界,请使用
sf::st_read()
或其他方法将自己的矢量数据作为sf
对象读取。一个不错的选择是安装rnaturalearth包,并使用它从rnaturalearthhires加载状态向量层。然后使用我们刚刚定义的lonlat_to_state()
函数,如下所示:要获得非常准确的结果,您可以从this page下载包含GADM维护的美国行政边界的地理包。然后,加载状态边界数据并像这样使用它们:
方法二(使用sp):
下面是一个函数,它在较低的48个状态中获取lat-long的 Dataframe ,并为每个点返回它所在的状态。
该函数的大部分内容只是准备
sp
包中的over()
函数所需的SpatialPoints
和SpatialPolygons
对象,该函数执行计算点和多边形的“交集”的真实的繁重工作:fnx2tebb2#
你可以用几行R来实现。
o3imoua43#
参见sp包中的
?over
。您需要将州边界设置为 SpatialPolygonsDataFrame。
mhd8tkvw4#
示例数据(多边形和点)
使用光栅::提取
yvt65v4c5#
这里是一个快速和简单的方法来转换纬度和经度到美国。州和美国州缩写
创建于2023-06-20使用reprex v2.0.2
lpwwtiir6#
使用
sf
非常简单: