Plot over OpenStreetMap with ggplot2

For those who work with geo-referenced data, it is useful to represent a certain magnitude, such as bias error or correlation of a data group, for a set of localisations or coordinates. This is relatively easy to do and there are many software for it. One way to do this is by using R and combine the OpenStreetMap (OSM) layers with the power of ggplot2, which will give us a world of possibilities for doing something really cool.

First we need to load the libraries ggplot2 and OpenStreetMap. The latter depends on rJava, so be sure that you have Java devel libraries installed (java-openjdk-devel).

library(OpenStreetMap)
library(ggplot2)

In this example I will use OpenStreetMap packages, there is another similar package for Google Maps, but as you probably know, Google limits the maximum tile area to be downloaded for a given level of zoom, for large areas you should use lower zoom levels otherwise you will get just a part of the area requested. The function openmap provides the OpenStreetMap tiles given a level of zoom (if NULL, it is determined automatically) or a minimum number of tiles. We provide the maximum and minimum Longitude and Latitude defining the upper-left and lower-right corners of your area.

LAT1 =  30 ; LAT2 = 50
LON1 = -10 ; LON2 = 10

map <- openmap(c(LAT2,LON1), c(LAT1,LON2), zoom = NULL,
               type = c("osm", "stamen-toner", "stamen-terrain","stamen-watercolor", "esri","esri-topo")[6],
               mergeTiles = TRUE)

Check the parameter type for selecting different types of maps layers.

As Google maps, OSM uses Mercator projection. In this example, for convenience, I prefer to work in standard WGS84 longlat coordinates. The line below re-project the original OSM map.

## OSM CRS :: "+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs"

map.latlon <- openproj(map, projection = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")

Now lets generate synthetic data for plotting over OSM. In this case I have defined a data.frame including random LON, LAT coordinates, a point ID and the BIAS as variable to be represented in the map. For real data, use read.table or any other read-from-file procedures to import your data to R.

bias.dtf <-  data.frame(LON  = runif(10, LON1, LON2),
                        LAT  = runif(10, LAT1, LAT2),
                        BIAS = runif(10, 5, 10),
                        ID   = paste("Pt.",sprintf("%02d",1:10)))

No it is time to use ggplot2 to overlap our data into the map. The map object we get from OSM is not directly in the format for ggplot2, we should then apply autoplot function. Afterwards, we can use all the power of ggplot2 to include points, labels, paths, etc. In this case I have plotted the coordinates points and the associated bias error as colors.

mytheme <- theme(plot.title = element_text(face = "bold",size = rel(1.2), hjust = 0.5),
                 panel.background = element_rect(colour = NA),
                 plot.background = element_rect(colour = NA),
                 axis.title = element_text(face = "bold",size = rel(1)),
                 axis.title.y = element_text(angle=90,vjust =2),
                 axis.title.x = element_text(vjust = -0.2))

OSMap <- autoplot(map.latlon)  +
    labs(title = "Plot over OpenStreetMap", subtitle = "Bias [%]",x = "Longitude", y="Latitude")+
    geom_label(data=bias.dtf, aes(x=LON ,y=LAT, label=ID, fontface=7), hjust=1, vjust=0, size=4) +
    geom_point(data=bias.dtf, aes(x=LON, y=LAT, fill=BIAS), size=4, shape=21) +
    scale_fill_gradientn(colours = rainbow(10)) +
    mytheme

This is a very simple plot (I know!), but this is good starting point for doing much more complex graphics by including density maps, pie or bar charts on each point, etc.

You can also use the package plotly for generating interactive graphics adding just a few of code lines.