Old 'Hood, New 'Hood
Published: 2021-08-24 17:27:51

---
title: Old 'Hood, New 'Hood
published: true
date: 2021-08-11 00:00:00 UTC
tags: 
canonical_url: https://datawookie.dev/blog/2021/08/old-hood-new-hood/
---

I recently moved from suburban South Africa to rural England. I’m figuring out my new environment. Making some maps seemed to be a good way to get familiar with the surroundings.

In the process I wanted to figure out two things:

- how to get maps with a consistent aspect ratio at different latitudes; and
- how to overlay a partially transparent map layer.

To make things more interesting I’ll create maps of both my old and new locations.

![Old 'Hood, New 'Hood](https://datawookie.dev/blog/2021/08/old-hood-new-hood/featured.png)<figcaption>
Image adapted from the cover of “Old Hat New Hat” by Dr Seuss.
</figcaption>

Let’s pull in a few libraries.

```
library(geosphere)
library(ggmap)
```

## Distance and Latitude

Let’s figure out the aspect ratio issue first. This arises because lines of longitude converge at the poles. As a result, a map with the same latitudinal and longitudinal extent, say 1° by 1°, looks very different at 60° latitude than it does at the equator.

We can illustrate this by estimating the distances using the [Haversine formula](https://en.wikipedia.org/wiki/Haversine_formula) (which assumes a spherical Earth).

```
# 1° of longitude on the Equator.
distHaversine(c(0, 0), c(1, 0))

[1] 111319.5

# 1° of longitude at 60° latitude.
distHaversine(c(0, 60), c(1, 60))

[1] 55659.22
```

Twice as long at the equator! A plot illustrates this effect rather nicely.

![](https://datawookie.dev/blog/2021/08/old-hood-new-hood/index_files/figure-html/unnamed-chunk-3-1.png)

Why is this important? Well, it means that 1° by 1° maps at the equator and at 60° latitude will represent very different areas and have different [aspect ratios](https://en.wikipedia.org/wiki/Aspect_ratio_(image)). 💡 For presentation purposes, this effect can be mitigated by choosing an appropriate map projection.

## Mercator & Pseudo-Mercator

![](https://datawookie.dev/blog/2021/08/old-hood-new-hood/globe.png)

The [Pseudo-Mercator](https://en.wikipedia.org/wiki/Web_Mercator_projection) (or Web Mercator) projection is a variation on the Mercator projection which is commonly used for maps on the internet. The maps served by Google Maps and OpenStreetMap, for example, use the Pseudo-Mercator projection. The Pseudo-Mercator projection treats the Earth as a sphere, while the Mercator projection takes into account the fact that the Earth is actually an oblate, ellipsoidal or “squashed” sphere. It’s not _very_ squashed, but it’s squashed enough to make a difference.

Both projections produce maps which look familiar: lines of latitude and longitude are straight, North is at the top, area gets bigger at higher latitudes. The differences are fairly subtle and really only relate to how angles on the surface of the Earth differ from angles on the map. For many applications this is not terribly important.

Why is the Pseudo-Mercator projection used for maps on the internet? It seems that this is ostensibly due to the different computational requirements: the maths for a spherical Earth is simpler (and quicker) than the maths for an ellipsoidal Earth.

## Aspect Ratio

I wanted a simple way to calculate the aspect ratio of a 1° by 1° piece of the Earth’s surface as a function of latitude.

<!-- This can also be implemented using ggplot2:::dist_central_angle(). -->
```
map_aspect <- function(lat, delta = 0.5) {
  if (lat + delta > 90) {
    NA
  } else {
    lon <- distHaversine(c(-delta, lat), c(delta, lat))
    lat <- distHaversine(c(0, lat - delta), c(0, lat + delta))
    lat / lon
  }
}
```

How does that vary with latitude?

![](https://datawookie.dev/blog/2021/08/old-hood-new-hood/index_files/figure-html/unnamed-chunk-5-1.png)

The scale on the vertical axis is logarithmic. The aspect ratio does really get large closer to the poles, where the 1° by 1° region is very narrow indeed.

How am I planning on using this? I want to plot maps for two regions at different latitudes. But I also want the maps to have the same aspect ratio. I need to specify a latitude-longitude bounding box for each map. I start by choosing the height I want the map to be (difference in latitude between top and bottom edges), then I use the aspect ratio to calculate the appropriate width (difference in longitude between left and right edges). Let’s give it a try.

First a map of Durban, South Africa.

![](https://datawookie.dev/blog/2021/08/old-hood-new-hood/index_files/figure-html/unnamed-chunk-8-1.png)

Then a map of West Berkshire, England.

![](https://datawookie.dev/blog/2021/08/old-hood-new-hood/index_files/figure-html/unnamed-chunk-9-1.png)

Both of these maps use the _terrain-background_ tiles from [Stamen Design](http://maps.stamen.com/).

That feels a little underwhelming? Was all this messing around with aspect ratios worth it? I think it was. I now have two maps with the _same shape_ and which represent the same distance along the vertical axis. And that’s pretty cool. To me, anyway.

## Transparent Map Layers

I wanted to be able to stack multiple map layers with varying levels of transparency. Although `ggmap()` has a `darken` argument, this doesn’t do the job. I did some research and did not find a good solution. So I rolled my own.

```
#' Add transparency to a ggmap object
#'
#' Add an alpha (transparency) layer to a \code{ggmap} object. This makes it
#' possible to superimpose multiple map layers and have the lower layers shine
#' through.
#'
#' @param map A \code{ggmap} object.
#' @param alpha The required transparency (a number between 0 and 1).
#'
#' @return A \code{ggmap} object.
#' @export
#'
#' @examples
#' zoom <- 12
#' bbox <- c(-0.3275, 51.407222, 0.0725, 51.607222)
#'
#' base <- get_map(bbox, zoom, maptype = "terrain-lines")
#' overlay <- get_map(bbox, zoom, maptype = "watercolor")
#'
#' # Plot the base map.
#' ggmap(base)
#' # Plot the overlay map (100% opacity).
#' ggmap(overlay)
#'
#' # Plot the base map with the overlay superimposed at 25% opacity.
#' ggmap(base) +
#' inset_ggmap(
#' set_map_alpha(overlay, 0.25)
#' )
set_map_alpha <- function(map, alpha) {
  if(class(map)[1] != "ggmap") {
    stop("map must be a ggmap object", call. = FALSE)
  }
  if(class(alpha) != "numeric" || alpha < 0 || alpha > 1) {
    stop("alpha must be a number between 0 and 1", call. = FALSE)
  }
  
  # Record the attributes & dimensions of the map object.
  map_attributes <- attributes(map)
  map_dimensions <- dim(map)
  
  # Add an alpha channel.
  map <- grDevices::adjustcolor(map, alpha)
  
  # Add back the dimensions (convert vector to matrix).
  dim(map) <- map_dimensions
  # Add back the attributes.
  attributes(map) <- map_attributes
  
  # Convert from matrix to map.
  class(map) <- c("ggmap", "raster")
  
  map
}
```

The `set_map_alpha()` function above accepts a map (returned by `get_map()`) and an `alpha` (between 0 and 1), then creates an alpha channel on the map.

The approach is nothing magical:

1. accept a `ggmap` object;
2. record the attributes and (pixel) dimensions;
3. use the builtin `adjustcolor()` function to add an alpha layer (the underlying data structure of the map is essentially a matrix of RGB strings);
4. the previous step results in a vector, so add back the dimensions, turning it into a matrix;
5. add back the attributes; and
6. convert to a `ggmap` object by assigning to `class()`.

## The Maps

Let’s try this out on those maps. I’ve created two versions of the Durban map, `durban_terrain_lines` and `durban_terrain_background`, which use the _terrain-lines_ and _terrain-background_ tiles from Stamen Design. For the latter I’ve used `set_map_alpha()` to add an alpha channel. To create an overlaid map I use `ggmap()` to create the base layer and then superimpose another later using `inset_ggmap()`.

```
ggmap(durban_terrain_lines) +
  inset_ggmap(durban_terrain_background)
```

How does that look? This is the Durban map with a _terrain-lines_ base and a _terrain-background_ overlay at 50% opacity.

![](https://datawookie.dev/blog/2021/08/old-hood-new-hood/index_files/figure-html/unnamed-chunk-12-1.png)

And now the West Berkshire map with a _terrain-lines_ base and a _watercolor_ overlay at 30% opacity.

![](https://datawookie.dev/blog/2021/08/old-hood-new-hood/index_files/figure-html/unnamed-chunk-13-1.png)

I really like the detail one can see on the _terrain-lines_ map with a hint of colour from the overlaying, transparent layers.

It’s apparent that my new environment is substantially different to where I was before. My old ’hood was great, but my new ’hood is better.

## Epilogue

Finally, just going to throw this in here because I think it looks cool: a map of London with multiple layers.

![](https://datawookie.dev/blog/2021/08/old-hood-new-hood/index_files/figure-html/unnamed-chunk-14-1.png)