Reproducible Modifiable Inset Maps

Over the weekend, I starting having a look at generating static maps (rather than interactive web maps) using matplotlib/Basemap, with one eye on the reproducible educational materials production idea, and the ways in which Basemap might be useful to authors creating new materials that are capable of being reused and/or maintained, with modification.

One of the things that struck me was how authors may want to produce different sorts of map. Basemap has importers for several different flavours of map tile that can essentially be treated as map styles, so once you have defined your map, you should be able to tile it in different ways without having to redefine the map.

It’s also worth noting how easy it is to change the projection…

Another thing that struck me was how maps-on-maps (inset maps?) can often help provide a situate the wider geographical context of a region that is being presented in some detail.

There are quite a few hardwired defaults baked in and the generator could be parameterised in all sorts of ways (eg changing plot size, colour themes, etc.)

Also on my to do list for the basic maps is a simple way of adding things like great circle connectors between two points, adding clear named location points, etc etc.

You can find my work in progress/gettingstarted demos on this theme on Azure Notebooks here.

If you’re interested, here’s what I came up with for the inset maps… It’s longer than it needs to be because it incorporates various bits and pieces for rendering default / demo views if no actual regions are specified.

from operator import itemgetter

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1.inset_locator import zoomed_inset_axes
from mpl_toolkits.axes_grid1.inset_locator import mark_inset
import numpy as np

def createInsetMap(base=None, inset=None,
                   zoom=None, loc=1, baseresolution=None,
                   insetresolution=None, connector=True):

    ''' Create an inset map for a particular region. '''

    fig = plt.figure()
    ax = fig.add_subplot(111)

    #Create a basemap

    #Add some default logic to provide some demos
    #If no args, use a World map
    if base is None or base.lower=='world':
        baseargs={'projection':'cyl','lat_0':0, 'lon_0':0}
    #If set UK base, then we can add a more detailed inset...
    elif isinstance(base,str) and base.lower()=='uk':
        baseargs={'llcrnrlon':ukllcrnrlon,'llcrnrlat':ukllcrnrlat,
                  'urcrnrlon':ukurcrnrlon,'urcrnrlat':ukurcrnrlat,
                  'resolution':'l'}
    else:
        #should really check base is a dict
        baseargs=base

    if baseresolution is not None:
        baseargs.update({'resolution':baseresolution})

    map1=Basemap(**baseargs)

    map1.drawmapboundary(fill_color='#f8f8ff')
    map1.drawcoastlines()

    plt.xticks(visible=False)
    plt.yticks(visible=False)

    #Now define the inset map

    #This default is to make for some nice default demos
    #With no explicit settings, inset UK on World map
    if base is None and inset is None:
        insetargs={'llcrnrlon':ukllcrnrlon,'llcrnrlat':ukllcrnrlat,
                  'urcrnrlon':ukurcrnrlon,'urcrnrlat':ukurcrnrlat,
                  'resolution':'l'}
        zoom=10
        loc=3
    #If the base is UK, and no inset is selected, demo an IW inset
    elif (isinstance(base,str) and base.lower()=='uk') and inset is None:
        insetargs={'llcrnrlon':iwllcrnrlon,'llcrnrlat':iwllcrnrlat,
                   'urcrnrlon':iwurcrnrlon,'urcrnrlat':iwurcrnrlat,
                   'resolution':'h'}
        zoom=10
        loc=2
    else:
        #Should really check inset is a dict...
        insetargs=inset

    axins = zoomed_inset_axes(ax, zoom, loc=loc)

    #The following seem to be set automatically?
    #axins.set_xlim(llcrnrlon, llcrnrlat)
    #axins.set_ylim(llcrnrlat, urcrnrlat)

    if insetresolution is not None:
        insetargs.update({'resolution':insetresolution})

    map2 = Basemap(ax=axins,**insetargs)

    #map2.drawmapboundary(fill_color='#7777ff')
    map2.fillcontinents(color='#ddaabb', lake_color='#7777ff', zorder=0)
    map2.drawcoastlines()
    map2.drawcountries()

    #Add connector lines from marked area on large map to zoomed area
    if connector:
        mark_inset(ax, axins, loc1=2, loc2=4, fc="none", ec="0.5")

    plt.show()

So, what other fragments of workflow, or use case, might be handy when creating (reusable / modifiable) maps in an educational context?

Experimenting With R – Point to Point Mapping With Great Circles

I’ve started doodling again… This time, around maps, looking for recipes that make life easier plotting lines to connect points on maps. The most attractive maps seem to use great circles to connect one point with another, these providing the shortest path between two points when you consider the Earth as a sphere.

Here’s one quick experiment (based on the Flowing Data blog post How to map connections with great circles), for an R/Shiny app that allows you to upload a CSV file containing a couple of location columns (at least) and an optional “amount” column, and it’ll then draw lines between the points on each row.

greatcircle map demo

The app requires us to solve several problems, including:

  • how to geocode the locations
  • how to plot the lines as great circles
  • how to upload the CSV file
  • how to select the from and two columns from the CSV file
  • how to optionally select a valid numerical column for setting line thickness

Let’s start with the geocoder. For convenience, I’m going to use the Google geocoder via the geocode() function from the ggmap library.

#Locations are in two columns, *fr* and *to* in the *dummy* dataframe
#If locations are duplicated in from/to columns, dedupe so we don't geocode same location more than once
locs=data.frame(place=unique(c(as.vector(dummy[[fr]]),as.vector(dummy[[to]]))),stringsAsFactors=F)
#Run the geocoder against each location, then transpose and bind the results into a dataframe
cbind(locs, t(sapply(locs$place,geocode, USE.NAMES=F))) 

The locs data is a vector of locations:

                    place
1              London, UK
2            Cambridge,UK
3            Paris,France
4       Sydney, Australia
5           Paris, France
6             New York,US
7 Cape Town, South Africa

The sapply(locs$place,geocode, USE.NAMES=F) function returns data that looks like:

    [,1]       [,2]     [,3]     [,4]      [,5]     [,6]      [,7]     
lon -0.1254872 0.121817 2.352222 151.207   2.352222 -74.00597 18.42406 
lat 51.50852   52.20534 48.85661 -33.86749 48.85661 40.71435  -33.92487

The transpose (t() gives us:

     lon        lat      
[1,] -0.1254872 51.50852 
[2,] 0.121817   52.20534 
[3,] 2.352222   48.85661 
[4,] 151.207    -33.86749
[5,] 2.352222   48.85661 
[6,] -74.00597  40.71435 
[7,] 18.42406   -33.92487

The cbind() binds each location with its lat and lon value:

                    place        lon       lat
1              London, UK -0.1254872  51.50852
2            Cambridge,UK   0.121817  52.20534
3            Paris,France   2.352222  48.85661
4       Sydney, Australia    151.207 -33.86749
5           Paris, France   2.352222  48.85661
6             New York,US  -74.00597  40.71435
7 Cape Town, South Africa   18.42406 -33.92487

Code that provides a minimal example for uploading the data from a CSV file on the desktop to the Shiny app, then creating dynamic drop lists containing column names, can be found here: Simple file geocoder (R/shiny app).

The following snippet may be generally useful for getting a list of column names from a data frame that correspond to numerical columns:

#Get a list of column names for numerical columns in data frame df
nums <- sapply(df, is.numeric)
names(nums[nums])

The code for the full application can be found as a runnable gist in RStudio from here: R/Shiny app – great circle mapping. [In RStudio, install.packages(“shiny”); library(shiny); runGist(9690079). The gist contains a dummy data file if you want to download it to try it out…]

Here’s the code explicitly…

The global.R file loads the necessary packages, installing them if they are missing:

#global.R

##This should detect and install missing packages before loading them - hopefully!
list.of.packages <- c("shiny", "ggmap","maps","geosphere")
new.packages <- list.of.packages[!(list.of.packages %in% installed.packages()[,"Package"])]
if(length(new.packages)) install.packages(new.packages)
lapply(list.of.packages,function(x){library(x,character.only=TRUE)}) 

The ui.R file builds the Shiny app’s user interface. The drop down column selector lists are populated dynamically with the names of the columns in the data file once it is uploaded. An optional Amount column can be selected – the corresponding list only displays the names of numerical columns. (The lists of location columns to be geocoded should really be limited to non-numerical columns.) The action button prevents the geocoding routines firing until the user is ready – select the columns appropriately before geocoding (error messages are not handled very nicely;-)

#ui.R
shinyUI(pageWithSidebar(
  headerPanel("Great Circle Map demo"),
  
  sidebarPanel(
    #Provide a dialogue to upload a file
    fileInput('datafile', 'Choose CSV file',
              accept=c('text/csv', 'text/comma-separated-values,text/plain')),
    #Define some dynamic UI elements - these will be lists containing file column names
    uiOutput("fromCol"),
    uiOutput("toCol"),
    #Do we want to make use of an amount column to tweak line properties?
    uiOutput("amountflag"),
    #If we do, we need more options...
    conditionalPanel(
      condition="input.amountflag==true",
      uiOutput("amountCol")
    ),
    conditionalPanel(
      condition="input.amountflag==true",
      uiOutput("lineSelector")
    ),
    #We don't want the geocoder firing until we're ready...
    actionButton("getgeo", "Get geodata")
    
  ),
  mainPanel(
    tableOutput("filetable"),
    tableOutput("geotable"),
    plotOutput("geoplot")
  )
))

The server.R file contains the server logic for the app. One thing to note is the way we isolate some of the variables in the geocoder reactive function. (Reactive functions fire when one of the external variables they contain changes. To prevent the function firing when a variable it contains changes, we need to isolate it. (See the docs for me; for example, Shiny Lesson 7: Reactive outputs or Isolation: avoiding dependency.)

#server.R

shinyServer(function(input, output) {

  #Handle the file upload
  filedata <- reactive({
    infile <- input$datafile
    if (is.null(infile)) {
      # User has not uploaded a file yet
      return(NULL)
    }
    read.csv(infile$datapath)
  })

  #Populate the list boxes in the UI with column names from the uploaded file  
  output$toCol <- renderUI({
    df <-filedata()
    if (is.null(df)) return(NULL)
    
    items=names(df)
    names(items)=items
    selectInput("to", "To:",items)
  })
  
  output$fromCol <- renderUI({
    df <-filedata()
    if (is.null(df)) return(NULL)
    
    items=names(df)
    names(items)=items
    selectInput("from", "From:",items)
  })
  
  #If we want to make use of an amount column, we need to be able to say so...
  output$amountflag <- renderUI({
    df <-filedata()
    if (is.null(df)) return(NULL)
    
    checkboxInput("amountflag", "Use values?", FALSE)
  })

  output$amountCol <- renderUI({
    df <-filedata()
    if (is.null(df)) return(NULL)
    #Let's only show numeric columns
    nums <- sapply(df, is.numeric)
    items=names(nums[nums])
    names(items)=items
    selectInput("amount", "Amount:",items)
  })
  
  #Allow different line styles to be selected
  output$lineSelector <- renderUI({
    radioButtons("lineselector", "Line type:",
                 c("Uniform" = "uniform",
                   "Thickness proportional" = "thickprop",
                   "Colour proportional" = "colprop"))
  })
  
  #Display the data table - handy for debugging; if the file is large, need to limit the data displayed [TO DO]
  output$filetable <- renderTable({
    filedata()
  })
  
  #The geocoding bit... Isolate variables so we don't keep firing this...
  geodata <- reactive({
    if (input$getgeo == 0) return(NULL)
    df=filedata()
    if (is.null(df)) return(NULL)
    
    isolate({
      dummy=filedata()
      fr=input$from
      to=input$to
      locs=data.frame(place=unique(c(as.vector(dummy[[fr]]),as.vector(dummy[[to]]))),stringsAsFactors=F)      
      cbind(locs, t(sapply(locs$place,geocode, USE.NAMES=F))) 
    })
  })

  #Weave the goecoded data into the data frame we made from the CSV file
  geodata2 <- reactive({
    if (input$getgeo == 0) return(NULL)
    df=filedata()
    if (input$amountflag != 0) {
      maxval=max(df[input$amount],na.rm=T)
      minval=min(df[input$amount],na.rm=T)
      df$b8g43bds=10*df[input$amount]/maxval
    }
    gf=geodata()
    df=merge(df,gf,by.x=input$from,by.y='place')
    merge(df,gf,by.x=input$to,by.y='place')
  })
  
  #Preview the geocoded data
  output$geotable <- renderTable({
    if (input$getgeo == 0) return(NULL)
    geodata2()
  })
  
  #Plot the data on a map...
  output$geoplot<- renderPlot({
    if (input$getgeo == 0) return(map("world"))
    #Method pinched from: http://flowingdata.com/2011/05/11/how-to-map-connections-with-great-circles/
    map("world")
    df=geodata2()
    
    pal <- colorRampPalette(c("blue", "red"))
    colors <- pal(100)
    
    for (j in 1:nrow(df)){
      inter <- gcIntermediate(c(df[j,]$lon.x[[1]], df[j,]$lat.x[[1]]), c(df[j,]$lon.y[[1]], df[j,]$lat.y[[1]]), n=100, addStartEnd=TRUE)

      #We could possibly do more styling based on user preferences?
      if (input$amountflag == 0) lines(inter, col="red", lwd=0.8)
      else {
        if (input$lineselector == 'colprop') {
          colindex <- round( (df[j,]$b8g43bds[[1]]/10) * length(colors) )
          lines(inter, col=colors[colindex], lwd=0.8)
        } else if (input$lineselector == 'thickprop') {
          lines(inter, col="red", lwd=df[j,]$b8g43bds[[1]])
        } else lines(inter, col="red", lwd=0.8)
      } 
    } 
  })

})

So that’s the start of it… this app could be further developed in several ways, for example allowing the user to filter or colour displayed lines according to factor values in a further column (commodity type, for example), or produce a lattice of maps based on facet values in a column.

I also need to figure how to to save maps, and maybe produce zoomable ones. If geocoded points all lay within a blinding box limited to a particular geographical area, scaling the map view to show just that area might be useful.

Other techniques might include using proportional symbols (circles) at line landing points to show the sum of values incoming to that point, or some of values outgoing, or the difference between the two; (maybe use green for incoming outgoing, then size by the absolute difference?)

Further Innovations in Campus Mapping

Almost a year or so ago I posted some quick round-ups of some recent innovations (at the time) in UK HE’s use of online campus maps (Open Data Powered Location Based Services in UK Higher Education, Innovations in Campus Mapping). Retweeting a post from Mike Nolan (Mapping the campus) about Edge Hill’s use of OpenStreetMap as the basis of an intereactive campus map (a nice example of how local benefits can support the common good?) I got a couple of tweeted responses about use of OpenStreetMap elsewhere:

– @danmcquillan mentioned how “my social computing 2nd years added detail & edits to @openstreetmap for @GoldsmithsUoL” [GoldsmithsUoL on OpenStreetMap]

– @julieallinson sent a link to the blog post that announced the University of York’s latest interactive map release, as well as to a behind the scenes post, which describes their use of custom tiles from CloudMade with the Google Map API providing the interactive map infrastructure. Data for location markers is stored in a Google Spreadsheet and then transformed into an appropriate JSON format via YQL which is glued into the map using a dash of JQuery. Magic:-)

I was also interested to see how the York U Estates department was releasing information about their building codes.

A quick trawl turned up a couple of other approaches to campus mapping that I don’t think I’ve mentioned before:

– the University of Warwick interactive map appears to use Bing maps, and as part of the offering provides a crude room level search:

– From the GoGeo team at Edina, a how-to post on creating simple campus maps using Digimap ROAM. Digimap ROAM is a service that allows you to annotate Ordnance Survey maps that can then be “printed off” as PDF documents.

The GoGeo blog post ends with a tease – “In the next post, we will look at some more advanced uses of Digimap data in campus maps” – but I can’t find any evidence of such a follow on post appearing?

So – any other innovations in campus based interactive maps out there (and in particular, web team blog posts about some of the technical details, including links to github repos containing the relevant code, perhaps?!;-)

PS Further OSM based campus maps: University of Cambridge

Innovations in Campus Mapping

A post on the ever interesting Google Maps Mania blog alerted me to a new interactive campus map produced by the Marketing folk at Loughborough University built on top of Google maps (I wonder if this is partly driven by Loughborough’s seemingly close ties with the Google education folk? If so, I wonder if they’re doing anything interesting relating to search in education too?)

Loughborough campus map maps.lboro.ac.uk

The layers are reminiscent of the layers in the Southampton campus map (e.g. Open Data Powered Location Based Services in UK Higher Education).

A simpler approach is taken by Liverpool University, who have bookmarked a few markers onto a map:

Liverpool campus map

Bristol University’s map appears to fall somewhere between the two in terms of complexity:

Bristol U campu map

Have any other UK HEIs rolled out interactive maps built on top of a third party API such as Google’s or maybe OpenStreetmap?

I know that Lincoln and Southampton are also currently battling it out in the “who can do the best 3D models of university buildings in Google Earth” (e.g. Lincoln’s KML file – h/t @alexbilbie;-). Any other UK HEIs with 2.5/3D building model layers? [h/t re: 2.5D to @gothwin] (Why’s this useful? Because if you can identify buildings placed as 3D models in Google Earth, you also have the lat/long data for those building profiles;-) Of course, just having building markers on a Google Map would also give you that geo-location data.

Lincoln U 3d layers/kml in google earth

I wonder whether any universities have Google Streetmap (or equivalent) style navigation? Where campuses are town centre based, and built on public roads, I guess they may be? I also wonder how many universities have started marking up locations as Google Places, or places on other location based services (I notice the Lincoln KML layer links through to a venue/place definition on Foursquare). This sort of detail would then support services based around checkins. For a homebrew example of what this might look like, see the OU’s check-in app: wayOU – mobile location tracking app using linked data.

A crude alternative, for promo purposes at least, is the campus tour, of course. Here’s one from several years ago looking at the OU’s Walton Hall campus:

Any others out there? Or any other map/location based initiatives being carried out by UK HEIs, whether similar to the above mentioned ones, or maybe even something completely different?!;-)

UPDATE: I’m guessing interactive maps are tending towards the norm by now, but as an when I come across them, I;ll try to add links here:
Uni of York interactive campus map (Google map; bootstrapped from an SU app? Makes me think of Student as Producer and the LNCD initiative…;-)