{"id":3563,"date":"2015-07-26T11:33:20","date_gmt":"2015-07-26T16:33:20","guid":{"rendered":"http:\/\/rud.is\/b\/?p=3563"},"modified":"2018-03-10T07:54:15","modified_gmt":"2018-03-10T12:54:15","slug":"making-staticinteractive-voronoi-map-layers-in-ggplotleaflet","status":"publish","type":"post","link":"https:\/\/rud.is\/b\/2015\/07\/26\/making-staticinteractive-voronoi-map-layers-in-ggplotleaflet\/","title":{"rendered":"Making Static\/Interactive Voronoi Map Layers In ggplot\/leaflet"},"content":{"rendered":"

Despite having shown various ways to overcome D3 cartographic envy<\/a>, there are always more examples that can cause the green monster to rear it’s ugly head.<\/p>\n

Take the Voronoi Arc Map<\/a> example.<\/p>\n


\n\"Voronoi_Arc_Map\"<\/a>
\n<\/center><\/p>\n

For those in need of a primer, a Voronoi tesslation\/diagram is:<\/p>\n

\n …a partitioning of a plane into regions based on distance to points in a specific subset of the plane. That set of points (called seeds, sites, or generators) is specified beforehand, and for each seed there is a corresponding region consisting of all points closer to that seed than to any other.<\/em> Wikipedia<\/a><\/sup>\n<\/p><\/blockquote>\n

We can overlay a Voronoi tessalation on top of a map in R as well thanks to the deldir<\/code> package (which has been around since the “S” days!). Let’s get (most of) the package requirements cruft out of the way, first:<\/p>\n

library(sp)\nlibrary(rgdal)\nlibrary(deldir)\nlibrary(dplyr)\nlibrary(ggplot2)\nlibrary(ggthemes)<\/pre>\n
Now we’ll [ab]use the data from the Arc Map example:<\/p>\n
flights <- read.csv(\"http:\/\/bl.ocks.org\/mbostock\/raw\/7608400\/flights.csv\", stringsAsFactors=FALSE)\nairports <- read.csv(\"http:\/\/bl.ocks.org\/mbostock\/raw\/7608400\/airports.csv\", stringsAsFactors=FALSE)<\/pre>\n
Since the D3 example cheats and only uses the continental US (CONUS) we'll do the same and we'll also filter out only those airports mentioned in the flights<\/code> data and get the total # of incoming\/outgoing flights for each airport:<\/p>\n
conus <- state.abb[!(state.abb %in% c(\"AK\", \"HI\"))]\nairports <- filter(airports,\n                   state %in% conus,\n                   iata %in% union(flights$origin, flights$destination))\norig <- select(count(flights, origin), iata=origin, n1=n)\ndest <- select(count(flights, destination), iata=destination, n2=n)\nairports <- left_join(airports,\n                      select(mutate(left_join(orig, dest),\n                                    tot=n1+n2),\n                             iata, tot)) %>% \n            filter(!is.na(tot))<\/pre>\n
Since we're going to initially plot polygons in ggplot (and, eventually, in leaflet), we'll need to work with Spatial<\/code> objects, so let's make those airport lat\/lon pairs into a SpatialPointsDataFrame<\/code>:<\/p>\n
vor_pts <- SpatialPointsDataFrame(cbind(airports$longitude,\n                                        airports$latitude),\n                                  airports, match.ID=TRUE)<\/pre>\n
The deldir<\/code> function returns a pretty complex object. Thankfully, the authors of the package realized that one might just want the polygons from the computation and pre-made a function: tile.list<\/code> for computing\/extracting them. Those polygons aren't, however, closed and we really want to keep the airport data associatd with them, so we need to close the polygons and associate the data. Since we're likely going to repeat this task, let's make it a (very badly named) function:<\/p>\n
SPointsDF_to_voronoi_SPolysDF <- function(sp) {\n\n  # tile.list extracts the polygon data from the deldir computation\n  vor_desc <- tile.list(deldir(sp@coords[,1], sp@coords[,2]))\n\n  lapply(1:(length(vor_desc)), function(i) {\n\n    # tile.list gets us the points for the polygons but we\n    # still have to close them, hence the need for the rbind\n    tmp <- cbind(vor_desc[[i]]$x, vor_desc[[i]]$y)\n    tmp <- rbind(tmp, tmp[1,])\n\n    # now we can make the Polygon(s)\n    Polygons(list(Polygon(tmp)), ID=i)\n\n  }) -> vor_polygons\n\n  # hopefully the caller passed in good metadata!\n  sp_dat <- sp@data\n\n  # this way the IDs _should_ match up w\/the data & voronoi polys\n  rownames(sp_dat) <- sapply(slot(SpatialPolygons(vor_polygons),\n                                  'polygons'),\n                             slot, 'ID')\n\n  SpatialPolygonsDataFrame(SpatialPolygons(vor_polygons),\n                           data=sp_dat)\n\n}<\/pre>\n
Before we can make the plots, we need to put the Spatial<\/code> objects into the proper form for ggplot2 (and get the U.S. state map):<\/p>\n
vor <- SPointsDF_to_voronoi_SPolysDF(vor_pts)\n\nvor_df <- fortify(vor)\n\nstates <- map_data(\"state\")<\/pre>\n
Now<\/em> we can have some fun. Let's try to mimic the D3 example map as closely as possible. We'll lay down the CONUS map, add a points layer for the the airports, sizing & styling them just like the D3 example. Note that we order the points so that the smallest ones appear on top (so we can still see them).<\/p>\n
We'll then lay down our newly created Voronoi layer. We'll also use the same projection (Albers) that the D3 examples uses:<\/p>\n
gg <- ggplot()\n# base map\ngg <- gg + geom_map(data=states, map=states,\n                    aes(x=long, y=lat, map_id=region),\n                    color=\"white\", fill=\"#cccccc\", size=0.5)\n# airports layer\ngg <- gg + geom_point(data=arrange(airports, desc(tot)),\n                      aes(x=longitude, y=latitude, size=sqrt(tot)),\n                      shape=21, color=\"white\", fill=\"steelblue\")\n# voronoi layer\ngg <- gg + geom_map(data=vor_df, map=vor_df,\n                    aes(x=long, y=lat, map_id=id),\n                    color=\"#a5a5a5\", fill=\"#FFFFFF00\", size=0.25)\ngg <- gg + scale_size(range=c(2, 9))\ngg <- gg + coord_map(\"albers\", lat0=30, lat1=40)\ngg <- gg + theme_map()\ngg <- gg + theme(legend.position=\"none\")\ngg<\/pre>\n
\"ggplot-1\"<\/a><\/p>\n
While that's pretty, it's not exactly useful<\/em>. I'm sure there are times when it's important to show the Voronoi polygons, but they are especially useful when they are used to help with user interface interactions.<\/p>\n
In the case of this map, some airport \"bubbles\" are very small and many overlap, making a \"click\" (or even \"hover\") a potentially painstaking task for someone looking to get more data out of the visualization. The D3 example uses Voronoi polygons to make it super-easy for the user to hover over a map area and get more info about the flights for the closest airport to the mouse pointer.<\/p>\n
We'll use the leaflet htmlwidget to do something similar. Until I can figure out \"hover\" events for R+leaflet, you'll have to live with \"click\".<\/p>\n
First we'll need some additional packages:<\/p>\n
library(leaflet)\nlibrary(rgeos)\nlibrary(htmltools)<\/pre>\n
And, we'll also need a U.S. shapefile (which we simplify since the polygons are pretty detailed and that's not necessary for this vis):<\/p>\n
url <- \"http:\/\/eric.clst.org\/wupl\/Stuff\/gz_2010_us_040_00_500k.json\"\nfil <- \"gz_2010_us_040_00_500k.json\"\n\nif (!file.exists(fil)) download.file(url, fil, cacheOK=TRUE)\n\nstates_m <- readOGR(\"gz_2010_us_040_00_500k.json\", \n                    \"OGRGeoJSON\", verbose=FALSE)\nstates_m <- subset(states_m, \n                   !NAME %in% c(\"Alaska\", \"Hawaii\", \"Puerto Rico\"))\ndat <- states_m@data # gSimplify whacks the data bits\nstates_m <- SpatialPolygonsDataFrame(gSimplify(states_m, 0.05,\n                                               topologyPreserve=TRUE),\n                                     dat, FALSE)<\/pre>\n
The leaflet vis idiom is similar to the ggplot idiom. I'm using a base tile layer since I was too lazy to figure out how to change the leaflet default gray background map color. The map polygons are added, then the circles\/bubbles (note that you work in meters<\/em> with addCircles<\/code> which lets leaflet scale the bubbles as you zoom in\/out). Finally, the Voronoi layer is added. I kept the stroke visible purely for demonstration purposes. You need<\/em> to keep fill=TRUE<\/code> otherwise the Voronoi layer won't get click\/hover events and once I figure out how to trigger popups on hover and use a static popup layer, this will let users hover around the map to get the underlying airport flight information.<\/p>\n
leaflet(width=900, height=650) %>%\n  # base map\n  addProviderTiles(\"Hydda.Base\") %>%\n  addPolygons(data=states_m,\n              stroke=TRUE, color=\"white\", weight=1, opacity=1,\n              fill=TRUE, fillColor=\"#cccccc\", smoothFactor=0.5) %>%\n  # airports layer\n  addCircles(data=arrange(airports, desc(tot)),\n             lng=~longitude, lat=~latitude,\n             radius=~sqrt(tot)*5000, # size is in m for addCircles O_o\n             color=\"white\", weight=1, opacity=1,\n             fillColor=\"steelblue\", fillOpacity=1) %>%\n  # voronoi (click) layer\n  addPolygons(data=vor,\n              stroke=TRUE, color=\"#a5a5a5\", weight=0.25,\n              fill=TRUE, fillOpacity = 0.0,\n              smoothFactor=0.5, \n              popup=sprintf(\"Total In\/Out: %s\",\n                            as.character(vor@data$tot)))<\/pre>\n