在R语言中将纬度经度坐标转换为州代码

Question

在R语言中将纬度经度坐标转换为州代码

rlatitude-longitudegoogle-geocoder

23

有没有一种快速的方法在 R 中将纬度和经度坐标转换成州代码？我一直在使用 zipcode 包作为查找表，但是当我查询大量的纬度/经度值时它太慢了。

如果不能在 R 中实现，是否可以使用 google 地理编码或任何其他类型的快速查询服务来完成？

谢谢！

- Michael

4

请参见我在这里的答案，使用 ggmap::revgeocode：https://stackoverflow.com/questions/46150851/how-to-get-california-county-location-from-latitude-and-longitude-information/46151310#46151310。 - moodymudskipper

6个回答

13

你可以用几行R代码完成它。

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)

- tingenek

3

请参考sp包中的?over。

您需要将州边界作为SpatialPolygonsDataFrame。

- mdsumner

0

这是一个快速简便的方法，可以将纬度和经度转换为美国州名和美国州名缩写。

# 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

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

- tassones

0

示例数据（多边形和点）

library(raster)
pols <- shapefile(system.file("external/lux.shp", package="raster"))
xy <- coordinates(p)

使用 raster::extract。

extract(p, xy)

#   point.ID poly.ID ID_1       NAME_1 ID_2           NAME_2 AREA
#1         1       1    1     Diekirch    1         Clervaux  312
#2         2       2    1     Diekirch    2         Diekirch  218
#3         3       3    1     Diekirch    3          Redange  259
#4         4       4    1     Diekirch    4          Vianden   76
#5         5       5    1     Diekirch    5            Wiltz  263
#6         6       6    2 Grevenmacher    6       Echternach  188
#7         7       7    2 Grevenmacher    7           Remich  129
#8         8       8    2 Grevenmacher   12     Grevenmacher  210
#9         9       9    3   Luxembourg    8         Capellen  185
#10       10      10    3   Luxembourg    9 Esch-sur-Alzette  251
#11       11      11    3   Luxembourg   10       Luxembourg  237
#12       12      12    3   Luxembourg   11           Mersch  233

- Robert Hijmans

-1

使用 sf 非常简单：

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)

- HAVB

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Josh O'Brien · Accepted Answer

这里提供了两个选项，一个使用sf包函数，另一个使用sp包函数。在2020年，推荐使用更现代的sf包来分析空间数据，但以防有用，我将保留我的原始2012年的回答，展示如何使用sp相关函数来完成此操作。

方法1 (使用sf)：

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

如果您需要更高分辨率的状态边界，请使用sf::st_read()或其他方式将自己的矢量数据读入为sf对象。一个不错的选择是安装rnaturalearth包，并从rnaturalearthhires加载状态矢量层。然后按照这里所示使用刚刚定义的lonlat_to_state()函数：

library(rnaturalearth)
us_states_ne <- ne_states(country = "United States of America",
                          returnclass = "sf")
lonlat_to_state(testPoints, states = us_states_ne, name_col = "name")
## [1] "Wisconsin" "Oregon"    NA

对于非常精确的结果，您可以从此页面下载包含由GADM维护的美国行政边界的地理数据包。然后，加载州边界数据并像这样使用它们：

USA_gadm <- st_read(dsn = "gadm36_USA.gpkg", layer = "gadm36_USA_1")
lonlat_to_state(testPoints, states = USA_gadm, name_col = "NAME_1")
## [1] "Wisconsin" "Oregon"    NA

方法二（使用sp）：

这里有一个函数，它接受一个包含在美国下48个州范围内的经纬度数据框，并为每个点返回所在的州名。

大部分函数只是准备需要由sp软件包中的over()函数使用的SpatialPoints和SpatialPolygons对象，该函数执行真正繁重的工作，即计算点和多边形之间的“交集”：

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