如何有效地将 lat/long 对映射到邮政编码（没有 API）
How to efficiently map lat/long pairs to zipcodes (without an API)
我正在将大量 lat/long 坐标对映射到相关的邮政编码。由于调用限制，记录数量太大，无法使用 API 之类的 Google 地图或地理名称。
我有一个查找 table，其中包含邮政编码和每个邮政编码的 lat/long 质心。您可以在此处查找 table：
# zipcode data with lat/lon coordinates
url <- "http://www.boutell.com/zipcodes/zipcode.zip"
fil <- "ziplatlong.zip"
# download an unzip
if (!file.exists(fil)) { download.file(url, fil) }
unzip(fil, exdir="zips")
library(readr)
ziplkp<-read_csv("zips/zipcode.csv")
对于我数据中的每个 lat/long 对，我想通过查找 lat/long 对与我查找中的每个质心之间的最小绝对差，将其与最近的邮政编码质心相匹配table。
将这种 "lookup" 函数按行应用于大量记录的最有效方法是什么？
示例数据：lat/long 坐标列表：
latlongdata <-
structure(list(dropoff_longitude = c(-73.981705, -73.993553,
-73.973305, -73.988823, -73.938484, -74.015503, -73.95472, -73.9571,
-73.971298, -73.99794), dropoff_latitude = c(40.760559, 40.756348,
40.762646, 40.777504, 40.684692, 40.709881, 40.783371, 40.776657,
40.752148, 40.720535)), row.names = c(8807218L, 9760613L, 3175671L,
10878727L, 2025038L, 5345659L, 14474481L, 1650223L, 684883L,
9129975L), class = "data.frame", .Names = c("dropoff_longitude",
"dropoff_latitude"))
print(latlongdata)
dropoff_longitude dropoff_latitude
8807218 -73.98171 40.76056
9760613 -73.99355 40.75635
3175671 -73.97330 40.76265
10878727 -73.98882 40.77750
2025038 -73.93848 40.68469
5345659 -74.01550 40.70988
14474481 -73.95472 40.78337
1650223 -73.95710 40.77666
684883 -73.97130 40.75215
9129975 -73.99794 40.72053
**ZipLooker 函数：查找从输入坐标对到最近的邮政编码质心和 returns 那个 zip
的最小绝对距离
library(dplyr)
ZipLooker<-function(dropoff_longitude,dropoff_latitude){
if(is.na(dropoff_longitude)|is.na(dropoff_latitude)){
z<-NA_character_
} else {
tryCatch({
x<-ziplkp1
x$latdiff=abs(dropoff_latitude-x$Latitude)
x$londiff=abs(dropoff_longitude-x$Longitude)
x$totdiff=x$latdiff+x$londiff
z<-head(top_n(x,1,-totdiff),n=1)$Postal
return(z)
}, error=function(e) NA)
}
}
使用 dplyr 的 rowwsie() 函数应用 Ziplooker 函数
latlongdata %>%
rowwise() %>%
mutate(zipcode=ZipLooker(dropoff_longitude,dropoff_latitude)
)
我使用这些方法将经纬度转换为包含多边形：
library(maptools)
points.file<-readShapePoints("path.to.pts.shp")
poly.file<-read.ShapePoy("path.to.poly.shp")
points.file %over% poly.file
这里有一个完整的解决方案：
library(sp)
library(maptools)
library(zipcode)
# grab the zip code boundaries
url <- "http://www2.census.gov/geo/tiger/GENZ2014/shp/cb_2014_us_zcta510_500k.zip"
fil <- "ztca.zip"
# don't waste bandwidth
if (!file.exists(fil)) { download.file(url, fil) }
unzip(fil, exdir="ztca")
# read them in (this takes a bit)
ztca <- readShapePoly("ztca/cb_2014_us_zcta510_500k.shp", verbose=TRUE)
# extract NY
ny <- ztca[as.character(ztca$ZCTA5CE10) %in% as.character(zipcode[zipcode$state=="NY",]$zip),]
# your points
latlongdata <-
structure(list(dropoff_longitude = c(-73.981705, -73.993553,
-73.973305, -73.988823, -73.938484, -74.015503, -73.95472, -73.9571,
-73.971298, -73.99794), dropoff_latitude = c(40.760559, 40.756348,
40.762646, 40.777504, 40.684692, 40.709881, 40.783371, 40.776657,
40.752148, 40.720535)), row.names = c(8807218L, 9760613L, 3175671L,
10878727L, 2025038L, 5345659L, 14474481L, 1650223L, 684883L,
9129975L), class = "data.frame", .Names = c("dropoff_longitude",
"dropoff_latitude"))
# make them all super spatial-like (must be in lon,lat format)
pts <- SpatialPoints(as.matrix(latlongdata[,1:2]))
# figure out where they are (this can take a bit)
dat <- pts %over% ny
# merge your data back in (there are many ways to do this)
dat$lon <- latlongdata$dropoff_longitude
dat$lat <- latlongdata$dropoff_latitude
rownames(dat) <- rownames(latlongdata)
# boom
dat
## ZCTA5CE10 AFFGEOID10 GEOID10 ALAND10 AWATER10 lon lat
## 8807218 10019 8600000US10019 10019 1780742 0 -73.98171 40.76056
## 9760613 10018 8600000US10018 10018 836253 0 -73.99355 40.75635
## 3175671 10022 8600000US10022 10022 1107169 0 -73.97330 40.76265
## 10878727 10069 8600000US10069 10069 249044 0 -73.98882 40.77750
## 2025038 11221 8600000US11221 11221 3582803 0 -73.93848 40.68469
## 5345659 10280 8600000US10280 10280 300652 38759 -74.01550 40.70988
## 14474481 10128 8600000US10128 10128 1206195 0 -73.95472 40.78337
## 1650223 10028 8600000US10028 10028 811363 0 -73.95710 40.77666
## 684883 10017 8600000US10017 10017 820953 0 -73.97130 40.75215
## 9129975 10013 8600000US10013 10013 1425085 0 -73.99794 40.72053
比较 dplyr rowwise 选项与 hrbrmstr 的 maptools 解决方案的速度，看起来 dplyr 胜出（至少在较小的数据集上）
# Test lat/long data
latlongdata <-
structure(list(dropoff_longitude = c(-73.981705, -73.993553,
-73.973305, -73.988823, -73.938484, -74.015503, -73.95472, -73.9571,
-73.971298, -73.99794), dropoff_latitude = c(40.760559, 40.756348,
40.762646, 40.777504, 40.684692, 40.709881, 40.783371, 40.776657,
40.752148, 40.720535)), row.names = c(8807218L, 9760613L, 3175671L,
10878727L, 2025038L, 5345659L, 14474481L, 1650223L, 684883L,
9129975L), class = "data.frame", .Names = c("dropoff_longitude",
"dropoff_latitude"))
# zipcode data with lat/lon coordinates
url <- "http://www.boutell.com/zipcodes/zipcode.zip"
fil <- "ziplatlong.zip"
# download an unzip
if (!file.exists(fil)) { download.file(url, fil) }
unzip(fil, exdir="zips")
library(readr)
ziplkp<-read_csv("zips/zipcode.csv")
# Method 1: dplyr + ZipLooker function
ZipLooker<-function(dropoff_longitude,dropoff_latitude){
if(is.na(dropoff_longitude)|is.na(dropoff_latitude)){
z<-NA_character_
} else {
tryCatch({
x<-ziplkp
x$latdiff=abs(dropoff_latitude-x$latitude)
x$londiff=abs(dropoff_longitude-x$longitude)
x$totdiff=x$latdiff+x$londiff
z<-head(top_n(x,1,-totdiff),n=1)$zip
return(z)
}, error=function(e) NA_character_)
}
}
latlongdata %>%
rowwise() %>%
mutate(zipcode=ZipLooker(dropoff_longitude,dropoff_latitude)
)
# Method 2: maptools + sp
library(sp)
library(maptools)
# grab the zip code boundaries
url <- "http://www2.census.gov/geo/tiger/GENZ2014/shp/cb_2014_us_zcta510_500k.zip"
fil <- "ztca.zip"
# don't waste bandwidth
if (!file.exists(fil)) { download.file(url, fil) }
unzip(fil, exdir="ztca")
# read them in (this takes a bit)
ztca <- readShapePoly("ztca/cb_2014_us_zcta510_500k.shp", verbose=TRUE)
# extract NY
ny <- ztca[as.character(ztca$ZCTA5CE10) %in% as.character(ziplkp[ziplkp$state=="NY",]$zip),]
# make them all super spatial-like (must be in lon,lat format)
pts <- SpatialPoints(as.matrix(latlongdata[,1:2]))
# figure out where they are (this can take a bit)
dat <- pts %over% ny
# merge your data back in (there are many ways to do this)
dat$lon <- latlongdata$dropoff_longitude
dat$lat <- latlongdata$dropoff_latitude
rownames(dat) <- rownames(latlongdata)
# comparing the two (only the bulkiest parts)
library(microbenchmark)
microbenchmark(
dat <- pts %over% ny
,
latlongdata %>%
rowwise() %>%
mutate(zipcode=ZipLooker(dropoff_longitude,dropoff_latitude)
)
,times = 10)
输出：
Unit: milliseconds
expr
dat <- pts %over% ny
latlongdata %>% rowwise() %>% mutate(zipcode = ZipLooker(dropoff_longitude, dropoff_latitude))
min lq median uq max neval
275.89494 286.38187 297.9254 421.8727 445.7165 10
95.18166 97.09873 101.8102 108.8677 122.0515 10
我正在将大量 lat/long 坐标对映射到相关的邮政编码。由于调用限制，记录数量太大，无法使用 API 之类的 Google 地图或地理名称。
我有一个查找 table，其中包含邮政编码和每个邮政编码的 lat/long 质心。您可以在此处查找 table：
# zipcode data with lat/lon coordinates
url <- "http://www.boutell.com/zipcodes/zipcode.zip"
fil <- "ziplatlong.zip"
# download an unzip
if (!file.exists(fil)) { download.file(url, fil) }
unzip(fil, exdir="zips")
library(readr)
ziplkp<-read_csv("zips/zipcode.csv")
对于我数据中的每个 lat/long 对，我想通过查找 lat/long 对与我查找中的每个质心之间的最小绝对差，将其与最近的邮政编码质心相匹配table。
将这种 "lookup" 函数按行应用于大量记录的最有效方法是什么？
示例数据：lat/long 坐标列表：
latlongdata <-
structure(list(dropoff_longitude = c(-73.981705, -73.993553,
-73.973305, -73.988823, -73.938484, -74.015503, -73.95472, -73.9571,
-73.971298, -73.99794), dropoff_latitude = c(40.760559, 40.756348,
40.762646, 40.777504, 40.684692, 40.709881, 40.783371, 40.776657,
40.752148, 40.720535)), row.names = c(8807218L, 9760613L, 3175671L,
10878727L, 2025038L, 5345659L, 14474481L, 1650223L, 684883L,
9129975L), class = "data.frame", .Names = c("dropoff_longitude",
"dropoff_latitude"))
print(latlongdata)
dropoff_longitude dropoff_latitude
8807218 -73.98171 40.76056
9760613 -73.99355 40.75635
3175671 -73.97330 40.76265
10878727 -73.98882 40.77750
2025038 -73.93848 40.68469
5345659 -74.01550 40.70988
14474481 -73.95472 40.78337
1650223 -73.95710 40.77666
684883 -73.97130 40.75215
9129975 -73.99794 40.72053
**ZipLooker 函数：查找从输入坐标对到最近的邮政编码质心和 returns 那个 zip
的最小绝对距离library(dplyr)
ZipLooker<-function(dropoff_longitude,dropoff_latitude){
if(is.na(dropoff_longitude)|is.na(dropoff_latitude)){
z<-NA_character_
} else {
tryCatch({
x<-ziplkp1
x$latdiff=abs(dropoff_latitude-x$Latitude)
x$londiff=abs(dropoff_longitude-x$Longitude)
x$totdiff=x$latdiff+x$londiff
z<-head(top_n(x,1,-totdiff),n=1)$Postal
return(z)
}, error=function(e) NA)
}
}
使用 dplyr 的 rowwsie() 函数应用 Ziplooker 函数
latlongdata %>%
rowwise() %>%
mutate(zipcode=ZipLooker(dropoff_longitude,dropoff_latitude)
)
我使用这些方法将经纬度转换为包含多边形：
library(maptools)
points.file<-readShapePoints("path.to.pts.shp")
poly.file<-read.ShapePoy("path.to.poly.shp")
points.file %over% poly.file
这里有一个完整的解决方案：
library(sp)
library(maptools)
library(zipcode)
# grab the zip code boundaries
url <- "http://www2.census.gov/geo/tiger/GENZ2014/shp/cb_2014_us_zcta510_500k.zip"
fil <- "ztca.zip"
# don't waste bandwidth
if (!file.exists(fil)) { download.file(url, fil) }
unzip(fil, exdir="ztca")
# read them in (this takes a bit)
ztca <- readShapePoly("ztca/cb_2014_us_zcta510_500k.shp", verbose=TRUE)
# extract NY
ny <- ztca[as.character(ztca$ZCTA5CE10) %in% as.character(zipcode[zipcode$state=="NY",]$zip),]
# your points
latlongdata <-
structure(list(dropoff_longitude = c(-73.981705, -73.993553,
-73.973305, -73.988823, -73.938484, -74.015503, -73.95472, -73.9571,
-73.971298, -73.99794), dropoff_latitude = c(40.760559, 40.756348,
40.762646, 40.777504, 40.684692, 40.709881, 40.783371, 40.776657,
40.752148, 40.720535)), row.names = c(8807218L, 9760613L, 3175671L,
10878727L, 2025038L, 5345659L, 14474481L, 1650223L, 684883L,
9129975L), class = "data.frame", .Names = c("dropoff_longitude",
"dropoff_latitude"))
# make them all super spatial-like (must be in lon,lat format)
pts <- SpatialPoints(as.matrix(latlongdata[,1:2]))
# figure out where they are (this can take a bit)
dat <- pts %over% ny
# merge your data back in (there are many ways to do this)
dat$lon <- latlongdata$dropoff_longitude
dat$lat <- latlongdata$dropoff_latitude
rownames(dat) <- rownames(latlongdata)
# boom
dat
## ZCTA5CE10 AFFGEOID10 GEOID10 ALAND10 AWATER10 lon lat
## 8807218 10019 8600000US10019 10019 1780742 0 -73.98171 40.76056
## 9760613 10018 8600000US10018 10018 836253 0 -73.99355 40.75635
## 3175671 10022 8600000US10022 10022 1107169 0 -73.97330 40.76265
## 10878727 10069 8600000US10069 10069 249044 0 -73.98882 40.77750
## 2025038 11221 8600000US11221 11221 3582803 0 -73.93848 40.68469
## 5345659 10280 8600000US10280 10280 300652 38759 -74.01550 40.70988
## 14474481 10128 8600000US10128 10128 1206195 0 -73.95472 40.78337
## 1650223 10028 8600000US10028 10028 811363 0 -73.95710 40.77666
## 684883 10017 8600000US10017 10017 820953 0 -73.97130 40.75215
## 9129975 10013 8600000US10013 10013 1425085 0 -73.99794 40.72053
比较 dplyr rowwise 选项与 hrbrmstr 的 maptools 解决方案的速度，看起来 dplyr 胜出（至少在较小的数据集上）
# Test lat/long data
latlongdata <-
structure(list(dropoff_longitude = c(-73.981705, -73.993553,
-73.973305, -73.988823, -73.938484, -74.015503, -73.95472, -73.9571,
-73.971298, -73.99794), dropoff_latitude = c(40.760559, 40.756348,
40.762646, 40.777504, 40.684692, 40.709881, 40.783371, 40.776657,
40.752148, 40.720535)), row.names = c(8807218L, 9760613L, 3175671L,
10878727L, 2025038L, 5345659L, 14474481L, 1650223L, 684883L,
9129975L), class = "data.frame", .Names = c("dropoff_longitude",
"dropoff_latitude"))
# zipcode data with lat/lon coordinates
url <- "http://www.boutell.com/zipcodes/zipcode.zip"
fil <- "ziplatlong.zip"
# download an unzip
if (!file.exists(fil)) { download.file(url, fil) }
unzip(fil, exdir="zips")
library(readr)
ziplkp<-read_csv("zips/zipcode.csv")
# Method 1: dplyr + ZipLooker function
ZipLooker<-function(dropoff_longitude,dropoff_latitude){
if(is.na(dropoff_longitude)|is.na(dropoff_latitude)){
z<-NA_character_
} else {
tryCatch({
x<-ziplkp
x$latdiff=abs(dropoff_latitude-x$latitude)
x$londiff=abs(dropoff_longitude-x$longitude)
x$totdiff=x$latdiff+x$londiff
z<-head(top_n(x,1,-totdiff),n=1)$zip
return(z)
}, error=function(e) NA_character_)
}
}
latlongdata %>%
rowwise() %>%
mutate(zipcode=ZipLooker(dropoff_longitude,dropoff_latitude)
)
# Method 2: maptools + sp
library(sp)
library(maptools)
# grab the zip code boundaries
url <- "http://www2.census.gov/geo/tiger/GENZ2014/shp/cb_2014_us_zcta510_500k.zip"
fil <- "ztca.zip"
# don't waste bandwidth
if (!file.exists(fil)) { download.file(url, fil) }
unzip(fil, exdir="ztca")
# read them in (this takes a bit)
ztca <- readShapePoly("ztca/cb_2014_us_zcta510_500k.shp", verbose=TRUE)
# extract NY
ny <- ztca[as.character(ztca$ZCTA5CE10) %in% as.character(ziplkp[ziplkp$state=="NY",]$zip),]
# make them all super spatial-like (must be in lon,lat format)
pts <- SpatialPoints(as.matrix(latlongdata[,1:2]))
# figure out where they are (this can take a bit)
dat <- pts %over% ny
# merge your data back in (there are many ways to do this)
dat$lon <- latlongdata$dropoff_longitude
dat$lat <- latlongdata$dropoff_latitude
rownames(dat) <- rownames(latlongdata)
# comparing the two (only the bulkiest parts)
library(microbenchmark)
microbenchmark(
dat <- pts %over% ny
,
latlongdata %>%
rowwise() %>%
mutate(zipcode=ZipLooker(dropoff_longitude,dropoff_latitude)
)
,times = 10)
输出：
Unit: milliseconds
expr
dat <- pts %over% ny
latlongdata %>% rowwise() %>% mutate(zipcode = ZipLooker(dropoff_longitude, dropoff_latitude))
min lq median uq max neval
275.89494 286.38187 297.9254 421.8727 445.7165 10
95.18166 97.09873 101.8102 108.8677 122.0515 10