R将邮政编码或纬度/经度转换为县答案

【问题标题】：R convert zipcode or lat/long to countyR将邮政编码或纬度/经度转换为县
【发布时间】：2012-10-30 06:02:31
【问题描述】：

我有一个位置列表，其中包含每个位置的城市、州、邮编、纬度和经度。

我单独有一份县级经济指标清单。我玩过zipcode 包、ggmap 包和其他几个免费的地理编码网站，包括US Gazeteer 文件，但似乎找不到匹配这两个部分的方法。

目前是否有任何软件包或其他来源可以做到这一点？

【问题讨论】：

这有帮助吗？ gis.stackexchange.com/questions/28035/…
对于纬度/经度数据，您可能会修改我提供的代码here，以回答有关如何将纬度/经度坐标转换为州代码的问题。
@JoshO'Brien：看起来它可以工作。看起来他们有一张美国所有县的地图。我得把它编码起来看看。感谢您的建议。
@JoshO'Brien：好主意。我刚刚更改了“州”->“县”，效果很好。非常感谢。
哦，你真是太好了，但不用了，谢谢。我很高兴得知该方法与 maps 县数据库和州数据库一样适用。虽然非常适合您发布您使用的内容，然后在几天内接受它。我当然会赞成一个真正的答案，它可能对未来的搜索者有帮助。干杯。

标签： r geolocation geocoding

【解决方案1】：

我最终使用了上面提到的JoshO'Brien 的建议并找到了here。

我拿了他的代码并将state更改为county，如下所示：

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

latlong2county <- function(pointsDF) {
    # Prepare SpatialPolygons object with one SpatialPolygon
    # per county
    counties <- map('county', fill=TRUE, col="transparent", plot=FALSE)
    IDs <- sapply(strsplit(counties$names, ":"), function(x) x[1])
    counties_sp <- map2SpatialPolygons(counties, 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, counties_sp)

    # Return the county names of the Polygons object containing each point
    countyNames <- sapply(counties_sp@polygons, function(x) x@ID)
    countyNames[indices]
}

# Test the function using points in Wisconsin and Oregon.
testPoints <- data.frame(x = c(-90, -120), y = c(44, 44))

latlong2county(testPoints)
[1] "wisconsin,juneau" "oregon,crook" # IT WORKS

【讨论】：

您好，如何修改这段代码以将 lat-long 转换为 FIPS 代码？谢谢。
这段代码给了我一个“未知的椭圆参数”错误，但在删除 SpatialPoints 和 map2SpatialPolygons 中的 proj4string 时对我有用。也许更新包？谢谢你写这篇文章！

【解决方案2】：

将邮政编码与县匹配是困难的。（某些邮政编码跨越多个县，有时跨越一个州。例如 30165）

我不知道有任何特定的 R 包可以为您匹配这些。

不过，您可以从密苏里州人口普查数据中心获得一张不错的表格。
您可以使用以下内容进行数据提取：http://bit.ly/S63LNU

示例输出可能如下所示：

    state,zcta5,ZIPName,County,County2
    01,30165,"Rome, GA",Cherokee AL,
    01,31905,"Fort Benning, GA",Russell AL,
    01,35004,"Moody, AL",St. Clair AL,
    01,35005,"Adamsville, AL",Jefferson AL,
    01,35006,"Adger, AL",Jefferson AL,Walker AL
    ...

注意县2。元数据解释可以在here找到。

    county 
    The county in which the ZCTA is all or mostly contained. Over 90% of ZCTAs fall entirely within a single county.

    county2 
    The "secondary" county for the ZCTA, i.e. the county which has the 2nd largest intersection with it. Over 90% of the time this value will be blank.

另请参阅 ANSI 县代码 http://www.census.gov/geo/www/ansi/ansi.html

【讨论】：

【解决方案3】：

我认为“非人口普查”包很有帮助。

对应的是我用来匹配邮政编码和县的

### code for get county based on zipcode

library(noncensus)
data(zip_codes)
data(counties)

state_fips  = as.numeric(as.character(counties$state_fips))
county_fips = as.numeric(as.character(counties$county_fips))    
counties$fips = state_fips*1000+county_fips    
zip_codes$fips =  as.numeric(as.character(zip_codes$fips))

# test
temp = subset(zip_codes, zip == "30329")    
subset(counties, fips == temp$fips)

【讨论】：

你的前三行，counties$fips的构造，可以用counties$fips <- interaction(counties$state_fips, counties$county_fips)在一行中获得。
看起来非人口普查已从 2020 年 1 月 16 日从 CRAN 中删除

【解决方案4】：

一个简单的选项是使用ggmap 中的geocode() 函数，以及output="more" 或output="all 选项。

这可以采用灵活的输入，例如地址或纬度/经度，并以列表的形式返回地址、城市、县、州、国家、邮政编码等。

require("ggmap")
address <- geocode("Yankee Stadium", output="more")

str(address)
$ lon                        : num -73.9
$ lat                        : num 40.8
$ type                       : Factor w/ 1 level "stadium": 1
$ loctype                    : Factor w/ 1 level "approximate": 1
$ address                    : Factor w/ 1 level "yankee stadium, 1 east 161st street, bronx, ny 10451, usa": 1
$ north                      : num 40.8
$ south                      : num 40.8
$ east                       : num -73.9
$ west                       : num -73.9
$ postal_code                : chr "10451"
$ country                    : chr "united states"
$ administrative_area_level_2: chr "bronx"
$ administrative_area_level_1: chr "ny"
$ locality                   : chr "new york"
$ street                     : chr "east 161st street"
$ streetNo                   : num 1
$ point_of_interest          : chr "yankee stadium"
$ query                      : chr "Yankee Stadium"

另一种解决方案是使用人口普查 shapefile，以及问题中相同的 over() 命令。我在使用 maptools 底图时遇到了一个问题：因为它使用 WGS84 基准，在北美，距离海岸几英里内的点被错误地映射，大约 5% 的数据集不匹配。

试试这个，使用 sp 包和 Census TIGERLine 形状文件

counties <- readShapeSpatial("maps/tl_2013_us_county.shp", proj4string=CRS("+proj=longlat +datum=NAD83"))

# Convert pointsDF to a SpatialPoints object 
pointsSP <- SpatialPoints(pointsDF, proj4string=CRS("+proj=longlat +datum=NAD83"))

countynames <- over(pointsSP, counties)
countynames <- countynames$NAMELSAD

【讨论】：