Объединение картограммы, сделанной в ggplot и ggmap

Создал картограмму с помощью ggplot2. Вот код ggplot

okc <- 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')

Вот пример данных из mapdata:

       long      lat order  hole piece         group          id
1 -97.54285 35.51951     1 FALSE     1 40109100100.1 40109100100
2 -97.54282 35.51954     2 FALSE     1 40109100100.1 40109100100
3 -97.54280 35.51963     3 FALSE     1 40109100100.1 40109100100
4 -97.54276 35.51976     4 FALSE     1 40109100100.1 40109100100
5 -97.54270 35.51993     5 FALSE     1 40109100100.1 40109100100
6 -97.54266 35.52016     6 FALSE     1 40109100100.1 40109100100
                                          NAME state county  tract B19013_001
1 Census Tract 1001, Oklahoma County, Oklahoma    40    109 100100      33440
2 Census Tract 1001, Oklahoma County, Oklahoma    40    109 100100      33440
3 Census Tract 1001, Oklahoma County, Oklahoma    40    109 100100      33440
4 Census Tract 1001, Oklahoma County, Oklahoma    40    109 100100      33440
5 Census Tract 1001, Oklahoma County, Oklahoma    40    109 100100      33440
6 Census Tract 1001, Oklahoma County, Oklahoma    40    109 100100      33440

Он произвел этот сюжет.

Карта 40109

Я также создал карту дорог с помощью ggmap. Вот код:

map <- get_map(location = c(lon = mean(mapdata$lon), lat = mean(mapdata$lat))
               , zoom = 10
               , maptype = "roadmap"
               , color = "bw")
p <- ggmap(map) +
  scale_x_continuous(limits = c(min(mapdata$lon), max(mapdata$lon)), expand = c(0, 0)) +
  scale_y_continuous(limits = c(min(mapdata$lat), max(mapdata$lat)), expand = c(0, 0))
p

И вот карта, которую он производит.

карта улиц

Когда я пытаюсь объединить их, я получаю сообщение об ошибке. Вот код, который я использую, чтобы объединить их и ошибку:

okc <- okc + p

Error in p + o : non-numeric argument to binary operator
In addition: Warning message:
Incompatible methods ("+.gg", "Ops.data.frame") for "+"

Я не уверен, почему я получаю эту ошибку. Это потому, что карты не масштабируются одинаково? Я не мог понять, как еще масштабировать ggmap, кроме как с помощью очень неточной функции масштабирования. Если у кого-нибудь есть идеи, как наложить хороплет поверх ggmap, я буду очень благодарен.

Вот остальная часть кода для воссоздания хороплет ggplot.

    library(acs)
    library(ggplot2)
    library(ggmap)
    library(UScensus2010)
    library(RColorBrewer)
    library(dplyr)
    library(scales)

    #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("insert shapefile here", 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)

    okc <- 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 OKC')

r ggplot2 choropleth ggmap
person Mike Tibb    schedule 16.10.2015    source источник
comment
С предоставленными данными я действительно многого не могу сделать. Поскольку вы хотите рисовать многоугольники поверх карты, вам нужно сделать что-то вроде этого: ggmap(map) + geom_polygon(data = mapdata, aes(x = long, y = lat, group = group, fill = B19013_001), color = "black", size = 0.5) Я думаю, вы хотите предоставить базовый слой, который является картой, используя ggmap(), затем вы рисуете полигоны поверх него.   -  person jazzurro    schedule 16.10.2015
comment
@jazzurro Да, это именно то, что я пытаюсь сделать. Код, который вы предоставили, является хорошим началом, но он укладывает картограмму поверх карты дорог. Есть ли способ изменить непрозрачность картограммы, чтобы вы действительно могли видеть, что находится под ней?   -  person Mike Tibb    schedule 18.10.2015
comment
Понимаю. В этом случае вы хотите использовать alpha в geom_polygon(). Ты можешь сделать; ggmap(map) + geom_polygon(data = mapdata, aes(x = long, y = lat, group = group, fill = B19013_001), color = "black", size = 0.5, alpha = 0.5). Поиграйте со значением альфа (от 0 до 1) и посмотрите, какое значение дает вам правильное изображение.   -  person jazzurro    schedule 19.10.2015
comment
Кстати, вы можете предоставить шейп-файл, который вы используете. Просто сказать "insert shapefile here" не позволит пользователям SO помочь вам.   -  person jazzurro    schedule 19.10.2015
comment
@jazzurro, как я могу добавить это напрямую? Я вставил ссылку на шейп-файл выше, но не должен ли кто-то, пытающийся запустить код, загрузить его и сохранить в своей собственной системе?   -  person Mike Tibb    schedule 19.10.2015
comment
Если все материалы есть, постараюсь посмотреть и помочь. Если вы уже нашли свой путь, вам не нужно. :)   -  person jazzurro    schedule 20.10.2015

Ответы (1)


arrow_upward
3
arrow_downward

На самом деле это намного лучше сделано в Leaflet. Он выглядит лучше эстетически, а также с точки зрения кода, он более интуитивно понятен.

library(leaflet)
library(rgdal)
library(RColorBrewer)

pal <- colorNumeric("OrRd", domain = new$pct_minority_popn)

leaflet(mapdata) %>%
 addTiles %>%
 addPolygons(stroke=T, fillOpacity=.5, smoothFactor=.5, color=~pal(B19013_001)) %>%
 addLegend("bottomright", pal=pal, values=~B19013_001, title="Legend Title", opacity=.8)

Вы можете изменить нижнюю карту, заменив команду addTiles чем-то вроде addProviderTiles("CartoDB.Positron"). Вы можете увидеть остальные параметры и дополнительную информацию в буклете по адресу: https://rstudio.github.io/leaflet/basemaps.html

person Abdalah El-Barrad    schedule 08.11.2015