3

Sorry if the title of this post confuses you, I could not come up with a more precise one. I am trying to overlay Voronoi polygons, generated using Houston crime data, on Houston city map at a certain zoom level. When I plot the polygons by themselves and restrict the visible area to the map's bounding box usingcoord_cartesian(xlim =, ylim =),the plot looks fine (see 1st plot below). But when they are plotted over the city map, polygons near the box edges are clipped and distorted (2nd plot), which as I understand from this cheatsheet should not happen. What am I doing wrong?

# load Houston crime data
suppressMessages(library(ggmap))
data(crime)
set.seed(42)
crime <- crime[sample(1:nrow(crime), 2000), ] # truncated for fast run

# convert to SpatialPointsDataFrame
suppressMessages(library(sp))
coords <- SpatialPoints(crime[, c("lon", "lat")])
crime_spdf <- SpatialPointsDataFrame(coords, crime)
proj4string(crime_spdf) <- CRS("+proj=longlat +ellps=WGS84")

# create Voronoi polygons
suppressMessages(library(spatstat))
suppressMessages(library(maptools))
vor_pp <- as(dirichlet(as.ppp(crime_spdf)), "SpatialPolygons")
proj4string(vor_pp) <- CRS("+proj=longlat +ellps=WGS84")

# get Houston map
houston_map <- get_map(location = geocode("Houston"),
                       zoom = 16,
                       maptype = "satellite")

# get map bounding box
xlim <- bb2bbox(attr(houston_map, "bb"))[c(1, 3)]
ylim <- bb2bbox(attr(houston_map, "bb"))[c(2, 4)]

# create Voronoi polygon plot
vor_df <- fortify(vor_pp)
plt1 <- ggplot(data = vor_df) +
  geom_polygon(aes(x = long, y = lat, group = group),
               color = "black",
               fill = NA) +
  coord_cartesian(xlim = xlim, ylim = ylim) +
  ggtitle("Voronoi Polygon Plot") +
  theme(aspect.ratio = 1,
        axis.title = element_blank(),
        axis.text = element_blank(),
        axis.ticks = element_blank())

# create map + polygon plot
plt2 <- ggmap(houston_map) +
  geom_polygon(data = vor_df,
               aes(x = long, y = lat, group = group),
               color = "gray80",
               fill = "red",
               alpha = 0.2) +
  coord_cartesian(xlim = xlim, ylim = ylim) +
  ggtitle("Map + Polygon Plot") +
  theme(aspect.ratio = 1,
        axis.title = element_blank(),
        axis.text = element_blank(),
        axis.ticks = element_blank())

# plot side-by-side
suppressMessages(library(grid))
pushViewport(viewport(layout = grid.layout(1,2)))
print(plt1, vp=viewport(layout.pos.col = 1, layout.pos.row = 1))
print(plt2, vp=viewport(layout.pos.col = 2, layout.pos.row = 1))

enter image description here Thanks.

Claus Wilke
  • 16,992
  • 7
  • 53
  • 104
Manojit
  • 611
  • 1
  • 8
  • 18

1 Answers1

3

The issue can be solved by adding the arguments extent = "normal" and maprange = FALSE to ggmap:

ggmap(houston_map, extent = "normal", maprange = FALSE) +
  geom_polygon(data = vor_df,
               aes(x = long, y = lat, group = group),
               color = "gray80",
               fill = "red",
               alpha = 0.2) +
  coord_cartesian(xlim = xlim, ylim = ylim) +
  ggtitle("Map + Polygon Plot") +
  theme(aspect.ratio = 1,
        axis.title = element_blank(),
        axis.text = element_blank(),
        axis.ticks = element_blank())

enter image description here

maprange: logical for use with base_layer; should the map define the x and y limits?

extent: how much of the plot should the map take up? "normal", "device", or "panel" (default)

adding expand = 0 to coord_cartesian results in:

enter image description here

missuse
  • 19,056
  • 3
  • 25
  • 47
  • Thanks! So, am I correct in understanding that the`extent`of polygons remains default "panel", while that of the map is shrunk to "normal"? – Manojit Feb 11 '18 at 18:25
  • 1
    @Manojit they key option is `maprange` where you specify the map should not define the plot limits, since when it does it overrides `coord_cartesian`. When the map defines the limits it does so in a `scale_x_continuous`/`scale_y_continuous` way which removes points outside of the graph resulting in cropped polygons. I am unsure why `extent` affects this but only `normal` works as intended. Check edit about `expand`. – missuse Feb 11 '18 at 18:38
  • Perfect!`expand = 0`did the trick of removing the unneeded strip around the box - thanks again. – Manojit Feb 11 '18 at 18:48