13 Day 11: Elevation

13.1 Technologies/Techniques

• Working with R Simple Features `{sf}`
• Dealing with multiple-sourced shapefiles
• Using `{ggplot2}` `geom_sf()` to heavily stylize a map

13.2 Data Source: Maine Elevation Contours

We love Acadia National Park47 and I used today’s challenge to finally chart out all of Mount Desert Island48 (MDI). There are links in the code-blocks to the various data sources where I obtained basemaps and elevation data. Where possible, these elements have been carved up and minified to make it easier to work with them.

``````library(sf)
library(hrbrthemes)
library(tidyverse)``````

The first bit retrieves all of the towns that make up MDI since I could not find a single MDI shapefile (I really should make one after this exercise). We do this work primarily to get the outline of the island so we can use that to extract the elevation data and lake data for that region from two other shapefiles. This is very straightforward in the `{sf}` ecosystem (it “feels” just like working with plain ol’ tidy data frames).

``````if (!all(file.exists(here::here("data", c("mdi.rds", "mdi-border.rds", "mdi-inland-water.rds"))))) {

# http://hub.arcgis.com/datasets/c7fdfaabcffa480c820e8be4a19bb298/data?layer=3

# get all the towns and focus on the ones we care about
filter(TOWN %in% c("Bar Harbor", "Mount Desert", "Tremont", "Southwest Harbor", "Cranberry Isles")) %>%
select(OBJECTID) %>%
filter( # this removes island and meta-borders we don't want
OBJECTID < 3000,
!OBJECTID %in% c(2412, 2519, 2480, 2520, 2530, 2533, 2523, 2266, 2661, 2694, 2255, 2287, 2631, 2674, 2643, 2698)
) %>%
st_union() -> mdi_border # finally, smush them all together

# the actual elevation contours
st_crop(xmin = -68.4719, ymin = 44.2108, xmax = -68.1542, ymax = 44.4666) %>% # reduce scope for faster intersection
st_transform(26919) %>%
st_intersection(mdi_border) -> mdi # only what's in the smushed outline above

# lakes and ponds

# "http://easterndivision.s3.amazonaws.com/Freshwater/LakeData.zip"

filter(STATE_ABBR == "ME") %>% # reduce scope for faster intersection
st_transform(26919) %>%
st_intersection(mdi_border) %>% # only what's in the smushed outline above
filter(!is.na(GNIS_NAME)) -> inland_water # finally, get rid of unnamed water bodies

# save ^^ work out so we don't have to do it again
saveRDS(mdi, here::here("data/mdi.rds"))
saveRDS(mdi_border, here::here("data/mdi-border.rds"))
saveRDS(inland_water, here::here("data/mdi-inland-water.rds"))

}

# read in our hard work

``plot(mdi_border)``

Here’s the elevation contours:

``plot(mdi["ContourEle"], lwd=0.1)``

Here are bodies of water on the island:

``plot(inland_water["GNIS_NAME"])``

We’re going to make our own `{sf}` points geometry object so we can add an annotation layer for some of the mountains on the final map. One key item to note is ensuring we use the same coordinate reference system (CRS) for this `outer_lab` object as we have for the other shapefiles we’ve mangled.

``````tibble(
lat = c(44.3380, 44.3526, 44.3351, 44.3426, 44.3329),
lng = c(-68.3119, -68.2251, -68.2438, -68.2728, -68.2667),
size = c(4.5, 3, 3, 3, 3),
angle = c(-80, 0, 0, 0, 0),
lab = c("Somes Sound", "Cadillac\nMountain", "Pemetic\nMountain", "Sargent\nMountain", "Penobscot\nMountain")
) %>%
st_as_sf(coords = c("lng", "lat"), crs = 4326) %>%
st_transform(26919) -> outer_lab``````

13.3 Drawing the Map

Don’t let this final map fool you. I spent quite a bit of time tweaking aesthetics to get it to this point and I think it still could use a bit of work.

We use the `inland_water` object twice, once to plot the polygons and the second time to label certain lakes/ponds.

``````ggplot() +
geom_sf(data = mdi_border, color = "#999999", fill = NA, size = 2, alpha = 1/20) + # water light outline
geom_sf(data = mdi_border, color = "black", fill = "#f6e3b9", size = 0.5) + # maine border
geom_sf(data = mdi, aes(color = ContourEle), size = 0.1) + # the contours
geom_sf(data = inland_water, color = "black", fill = "#3a4f5b", size = 0.3) + # the lakes and ponds
geom_sf_text( # label inland water we care about
data = inland_water %>%
filter(
GNIS_NAME %in% c("Jordan Pond", "Eagle Lake", "Echo Lake", "Seal Cove Pond") |
(GNIS_NAME == "Long Pond" & ACRES > 100)
),
aes(label = GNIS_NAME),
family = font_es, color = "white",
size = c(3.25, 3.25, 2.5, 2.5, 2.5),  # these align with the ordered factor names
angle = c(87.5, -75, -75, -80, -87.5)
) +
geom_sf_text(
data = outer_lab,
aes(label = lab, size = I(size)+0.125, angle = I(angle)),
family = font_es, color = c("black", "white", "white", "white", "white"), lineheight = 0.875
) +
geom_sf_text(
data = outer_lab,
aes(label = lab, size = I(size), angle = I(angle)),
family = font_es, color = c("white", "black", "black", "black", "black"), lineheight = 0.875
) +
scale_color_viridis_c(
name = "Elevation (ft)", option = "magma", label = scales::comma
) +
coord_sf(datum = NA) +
labs(
x = NULL, y = NULL,
title = "Mount Desert Island / Acadia National Park Elevation Contours",
caption = "Data source: ArcGIS/<maine.gov> • #30DayMapChallenge"
) +
theme_ipsum_es(grid="", plot_title_size = 24) +
theme(plot.background = element_rect(fill = "#3a4f5b", color = "#3a4f5b")) +
theme(panel.background = element_rect(fill = "#3a4f5b", color = "#3a4f5b")) +
theme(legend.title = element_text(color = "white")) +
theme(legend.text = element_text(color = "white")) +
theme(plot.title = element_text(hjust = 0.5, color = "white")) +
theme(plot.caption = element_text(color = "white")) +
theme(legend.position = c(0.9, 0.275))``````

13.4 In Review

We pulled in shapefiles from multiple sources and did a bit of slicing and dicing on and with them. We also used `{ggplot2}` to build a highly stylized final map (vs some of the plain ones we’ve built thus far).

13.5 Try This At Home

Pick one of your favorite hiking, biking, fishing, (etc.) spots and start spelunking for shapefiles you can use to create a similar map. “Googling” for `xyz elevation contours` would be a good start for finding data you can use. Many states/towns are doing a great job providing GIS data these days, so your search should be fruitful and fairly brief.