How Color Great Lakes Blue from GADM maps

Question

I am using getData("GADM") to obtain map data for the US and Canada. Hudson's Bay and the Atlantic are transparent so I can add a blue background that shows through. However, all the Great Lakes have disappeared. Can I get them to be transparent as well to let the blue show through?

The data from GADM is detailed and the plotting of the maps is very slow in RStudio on my Mac. I have added code found here to speed up the plotting.

library(tidyverse)
library(raster)

# Convert the polygons into data frames so we can make lines
poly2df <- function(poly) {
  # Convert the polygons into data frames so we can make lines
  # Number of regions
  n_regions <- length(poly@polygons)

  # Get the coords into a data frame
  poly_df <- c()
  for(i in 1:n_regions) {
    # Number of polygons for first region
    n_poly <- length(poly@polygons[[i]]@Polygons)
    print(paste("There are",n_poly,"polygons"))
    # Create progress bar
    pb <- txtProgressBar(min = 0, max = n_poly, style = 3)
    for(j in 1:n_poly) {
      poly_df <- rbind(poly_df, NA, 
                       poly@polygons[[i]]@Polygons[[j]]@coords)
      # Update progress bar
      setTxtProgressBar(pb, j)
    }
    close(pb)
    print(paste("Finished region",i,"of",n_regions))
  }
  poly_df <- data.frame(poly_df)
  names(poly_df) <- c('lon','lat')
  return(poly_df)
}

# Get the main polygons, will determine by area.
getSmallPolys <- function(poly, minarea = 0.01) {
  # Get the areas
  areas <- lapply(poly@polygons, 
                  function(x) sapply(x@Polygons, function(y) y@area))

  # Quick summary of the areas
  print(quantile(unlist(areas)))

  # Which are the big polygons?
  bigpolys <- lapply(areas, function(x) which(x > minarea))
  length(unlist(bigpolys))

  # Get only the big polygons and extract them
  for(i in 1:length(bigpolys)){
    if(length(bigpolys[[i]]) >= 1 && bigpolys[[i]] >= 1){
      poly@polygons[[i]]@Polygons <- poly@polygons[[i]]@Polygons[bigpolys[[i]]]
      poly@polygons[[i]]@plotOrder <- 1:length(poly@polygons[[i]]@Polygons)
    }
  }
  return(poly)
}

provinces <- c("Ontario", "Manitoba", "Québec")

# Get ALL provincial boundaries as spatial polygons data frame
canada <- getData("GADM", country = "CAN", level = 1)
usa <- getData("GADM", country = "USA", level = 1)

# Narrow down to just the provinces of interest
ca.provinces <- canada[canada$NAME_1 %in% provinces,]

ca.provinces_small <- getSmallPolys(ca.provinces) %>% 
  gSimplify(tol = 0.01, topologyPreserve = TRUE) %>% 
  poly2df()

usa_small <- getSmallPolys(usa) %>% 
  gSimplify(tol = 0.01, topologyPreserve = TRUE) %>% 
  poly2df()


map_theme <- theme_void() +
  theme(
    panel.background = element_rect(fill = "steelblue")
  )


ggplot() +

  # Crop area
  coord_fixed(xlim = c(-95, -65),  ylim = c(41, 55), ratio = 1.6) +

  # Canada
  geom_polygon(data = ca.provinces_small,  aes(lon, lat), 
               fill = "dark grey", size = 0.3, color = "black") +
  geom_polygon(data = usa_small,  aes(lon, lat), 
               fill = "dark grey", size = 0.3, color = "black") + 
  map_theme

Answer 1

Using the lakes shapefile referenced in the other answer, you can subtract the lakes from the boundaries to get the land region for each country. Working with sf objects helps greatly too:

library(raster)
canada <- getData("GADM", country = "CAN", level = 1)
usa <- getData("GADM", country = "USA", level = 1)
provinces <- c("Ontario", "Manitoba", "Québec")
ca.provinces <- canada[canada$NAME_1 %in% provinces,]

library(sp)
library(sf)

Now get the lakes as an sf object, and convert the provinces and usa to sf objects:

lakes = st_read("./Great_Lakes.shp")

ca.provinces = st_as_sf(ca.provinces)
usa = st_as_sf(usa)

Subtract the lakes from the admin regions:

ca.land = st_difference(ca.provinces,lakes)
usa.land = st_difference(usa,lakes)

Now we use ggplot2 with geom_sf and your theme, and coord_sf to set the limits:

library(ggplot2)

ggplot() + geom_sf(data=usa.land) + geom_sf(data=ca.land) + coord_sf(xlim = c(-95, -65),  ylim = c(41, 55)) + map_theme

Note how the background shows through, which you don't get by simply splatting a lake polygon shape over the top with the right colour.

If things are still slow then use st_simplify to reduce the complexity of polygons. But this isn't slow on my fairly old hardware.

Answer 2

One solution would be download a Great Lakes shapefile (e.g. the shapes here) and add them as a last layer on your map, the same color as the background.

Edit: You can get one good Great Lakes shapefile here rather than shapes for each individual link as above.