使用 ggplot2 按邻域绘制收入图答案

【问题标题】：Mapping Income by Neighborhood with ggplot2使用 ggplot2 按邻域绘制收入图
【发布时间】：2015-10-09 05:27:20
【问题描述】：

我有俄克拉荷马城不同部分的平均收入数据，我正在尝试创建该城市的等值线。这是数据。

library(acs)
library(ggplot2)
library(ggmap)
library(UScensus2010)
library(spdep)
library(RColorBrewer)

# Use your census API key (You'll need to get your own API key)
#http://api.census.gov/data/key_signup.html
api.key.install(key="c369cd6ed053a84332caa62301eb8afe98bed825")

# Load in Shape File (You'll need to download this file from the census)
#ftp://ftp2.census.gov/geo/tiger/TIGER2013/TRACT/tl_2013_40_tract.zip
geodat<-readShapePoly((""), proj4string=CRS('+proj=longlat +datum=NAD83'))

geodat<-geodat[geodat$COUNTYFP==109,]

# Tract Level Data
mdat<-geodat@data

# American Community Survey Data: Median HH Income for OK Census Tracts
ok.counties=geo.make(state="OK", county="Oklahoma", tract="*")
ok.income<-acs.fetch(geography=ok.counties, table.number="B19013", endyear=2013)

# Merge Data Sets 
geo_dat<-geography(ok.income)
var_dat<-as.data.frame(estimate(ok.income))
acs_data<-cbind(geo_dat,var_dat)
mdat2<-merge(mdat,acs_data,by.x="TRACTCE", by.y= "tract")

#cut data
mdat2$cut <- cut(mdat2$B19013_001, breaks = c(5000, 200000, by = 10000))
map <- mdat2[, 11:12]

#create map
 okc <- ggplot(mdat2, aes(x = INTPTLON, y = INTPTLAT, fill = B19013_001)) +
                  geom_map(aes(fill=cut), map=map)

我是 ggplot2 的新手，似乎无法让它正常工作。我可以使用 spplot 创建一个等值线，但我想用 ggplot 制作一个，这样我就可以将它与其他地图分层。

编辑：忘记包含我正在尝试创建一个以中等收入作为填充的等值线。 mdat2 数据框中的收入中位数为 B19013_001。 INTPTLON 矢量是城市不同部分的经度坐标，INTPTLAT 是纬度坐标。 Here 是我在 sppolot 中创建的一个绘图示例。我正在尝试在 ggplot 中制作类似的东西

以下是 mdat2 的前 10 个观察结果：

  TRACTCE STATEFP COUNTYFP       GEOID NAME.x          NAMELSAD MTFCC FUNCSTAT   ALAND AWATER    INTPTLAT
1  100100      40      109 40109100100   1001 Census Tract 1001 G5020        S 2141969      0 +35.5207581
2  100200      40      109 40109100200   1002 Census Tract 1002 G5020        S 3349288      0 +35.5039838
3  100300      40      109 40109100300   1003 Census Tract 1003 G5020        S 2098254      0 +35.5123370
4  100400      40      109 40109100400   1004 Census Tract 1004 G5020        S 2593440      0 +35.5004315
5  100500      40      109 40109100500   1005 Census Tract 1005 G5020        S 2734867      0 +35.5006179
6  100600      40      109 40109100600   1006 Census Tract 1006 G5020        S  551754      0 +35.5040878
      INTPTLON                                       NAME.y state county B19013_001
1 -097.5315537 Census Tract 1001, Oklahoma County, Oklahoma    40    109      33440
2 -097.5546474 Census Tract 1002, Oklahoma County, Oklahoma    40    109      29338
3 -097.5259348 Census Tract 1003, Oklahoma County, Oklahoma    40    109      85592
4 -097.4855750 Census Tract 1004, Oklahoma County, Oklahoma    40    109      21084
5 -097.5037840 Census Tract 1005, Oklahoma County, Oklahoma    40    109      23411
6 -097.5172200 Census Tract 1006, Oklahoma County, Oklahoma    40    109      62857

【问题讨论】：

到底是什么问题？并请包括您的一些数据（可能是模拟的）。很少有人会获得密码/API 密钥来回答您的问题。一点：不要在 aes-calls 中使用data$column。这是不必要的，尤其是当一切都来自 mdat2 时。
@Heroka 我添加了一张使用 spplot 创建的等值线图像。我正在尝试在 ggplot 中重新创建类似的东西。我怎样才能为这样的事情添加示例数据？如果这样会更容易，我可以包含我的 API 密钥。 mdat2数据框中的重要向量是INTPTLON、INTPTLAT和B19013_001，分别代表经度、纬度和收入中位数。
那么，什么是错误/挂断？
看起来你的 INTPLON / INTPLAT 是点——你需要映射完整的多边形几何图形，而不仅仅是点（使用正确的数据结构，ggplot 可以在没有显式 x / y 参数的情况下自动绘制）。网上有几个关于如何准备空间数据的教程，比如goo.gl/PkeHR6--see，如果你能让geom_map工作的话。

标签： r ggplot2 choropleth

【解决方案1】：

这是一个适用于 ggplot 的解决方案：

library(acs)
library(ggplot2)
library(ggmap)
library(UScensus2010)
library(spdep)
library(RColorBrewer)


library(dplyr)
library(scales)

# Use your census API key (You'll need to get your own API key)
#http://api.census.gov/data/key_signup.html
api.key.install(key="c369cd6ed053a84332caa62301eb8afe98bed825")

# Load in Shape File (You'll need to download this file from the census)
#ftp://ftp2.census.gov/geo/tiger/TIGER2013/TRACT/tl_2013_40_tract.zip

## load, subset shapefile
geodat<-readShapePoly("/home/aksel/Downloads/ohio//tl_2013_40_tract.shp", proj4string=CRS('+proj=longlat +datum=NAD83'))
geodat<-geodat[geodat$COUNTYFP==109,]

## fortify for ggplot digestion
geodat.f<-fortify(geodat,region="GEOID")

# American Community Survey Data: Median HH Income for OK Census Tracts
ok.counties=geo.make(state="OK", county="Oklahoma", tract="*")
ok.income<-acs.fetch(geography=ok.counties, table.number="B19013", endyear=2013)


# Merge Data Sets 
geo_dat<-geography(ok.income)
var_dat<-as.data.frame(estimate(ok.income))
acs_data<-cbind(geo_dat,var_dat)
acs_data$id<- paste("40109", acs_data$tract, sep = "")

## from dplyr
mapdata<-left_join(geodat.f,acs_data)

ggplot() +
  geom_polygon(data = mapdata, aes(x = long, y = lat, group = group,
                                    fill = B19013_001), color = "black", size = 0.5)+
  scale_fill_distiller(palette = "Reds", labels = comma,
                       breaks = pretty_breaks(n = 10), values = c(1,0)) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_nothing(legend = TRUE) +
ggtitle('Map of 40109')

产生：

【讨论】：

供您参考，我正在开发一个包，该包使用单个函数 geom_boundary( )、github.com/GL-Li/ggtiger 在 ggmap 上绘制人口普查边界。希望对你有用。