Lab 9: Mapping with Leaflet

Author

Lindsay Poirier

Introduction

In this lab, we will build a map that visualizes toxic emissions in Louisiana’s Cancer Alley - an 85-mile stretch of land from Baton Rouge to New Orleans, where some of the most serious health disparities in the U.S. are reported. In doing so, we will gain practice in producing point maps in Leaflet.

Learning Goals

  • Transform point data to an appropriate CRS
  • Map point data in Leaflet
  • Create palettes for points on maps
  • Add legends and labels to a map

Review of Key Terms

Tip

You may wish to reference this reference guide when completing this lab.

Coordinate Reference System (CRS)

a system used to locate geographic points on different spatial projections

Toxic Release Inventory

The EPA’s Emergency Planning and Community Right to Know Act (EPCRA) of 1986 established the Toxic Release Inventory as a mechanism to monitor and inform the public of toxic emissions released in their communities. Every year, certain U.S. industrial facilities are required to report to the EPA the amounts of certain chemical on-site and off-site releases in pounds. Facilities required to report include those that employ more than 10 individuals, release more than a certain threshold of a TRI-regulated chemical, and are classified by a specified set of Standard Industrial Codes, including mining, utilities, manufacturing, publishing, and hazardous waste.

Today, we are going to examine TRI data from 2023 in Louisiana’s Cancer Alley. Cancer Alley is an 85-mile stretch of land along the Mississippi River where more than 200 fossil fuel and petrochemical operate. It is an area of the country facing some of the highest cancer risks in the United States, and it is inhabited predominantly by Black communities, raising significant concerns about environmental injustices. Next week, we are going to look at the demographics of the community’s in cancer alley and consider its racial history. This week, we are going to focus on chemical emissions in this area.

Setting Up Your Environment

  1. Install the following packages in your Environment:
  • install.packages("leaflet")
  • install.packages("sf")
  1. Run the code below to load today’s data frames into your environment.
library(tidyverse)
library(leaflet)
library(RColorBrewer)
library(sf)

tri_la_2023 <- read_csv("https://raw.githubusercontent.com/SDS-192-Intro/public-website-fall-25/refs/heads/main/data/2023_la.csv", name_repair = make.names) |>
  filter(X50..UNIT.OF.MEASURE == "Pounds") |>
  arrange(desc(X107..TOTAL.RELEASES)) |>
  mutate(chemical_total_text = paste(X37..CHEMICAL, X107..TOTAL.RELEASES, sep = ": ")) |>
  group_by(X2..TRIFD) |>
  summarize(X4..FACILITY.NAME = first(X4..FACILITY.NAME),
            X12..LATITUDE = first(X12..LATITUDE),
            X13..LONGITUDE = first(X13..LONGITUDE),
            X23..INDUSTRY.SECTOR = first(X23..INDUSTRY.SECTOR),
            X107..TOTAL.RELEASES = sum(X107..TOTAL.RELEASES),
            X..CHEMICAL.TEXT = paste(chemical_total_text, collapse = " <br> "))

Mapping in Leaflet

We are going to work with the leaflet() package to design a map of the toxic emissions released in Louisiana in 2023. The map will color and size the points by the total emissions. At any point during these exercises, you can reference the leaflet documentation to help you build out these maps.

Initialize your Map

There are three steps to initializing a map in leaflet:

  1. Call the leaflet() function to create a map widget.
  2. Call the setView(lat = 0.00, lng = 0.00, zoom = 0). This function determines where the map will initially focus. We provide a set of coordinates that will be the map’s center point when we load it, and a number (from 1 to 16) to indicate the level at which to zoom in on the map.
  3. Call addProviderTiles() to load up provider tiles that constitute a base map. A number of different map providers have provider tiles that we can reference here. A few examples of the arguments we can supply to this function include include:
  • providers$OpenStreetMap
  • providers$Stamen.Toner
  • providers$CartoDB.Positron
  • providers$Esri.NatGeoWorldMap

Run the code below to initialize a map.

map1 <- leaflet(width = "100%") %>%
  setView(lat = 0.00, lng = 0.00, zoom = 0) %>%
  addProviderTiles(providers$Stamen.Toner)

map1
Question

Adjust the code below to center the map on Louisiana. You’ll need to look up the coordinates for LA, keeping in mind that South and West coordinates will be negative. Adjust the zoom to keep the whole state in view (setting the zoom level between 1 and 16). When you are happy with the View, switch out the provider tiles to a base map that won’t distract from the data points we will layer on top of this map. Keep in mind our discussions regarding Visualization Aesthetics here.

la_map <- leaflet(width = "100%") %>%
  setView(lat = 0.00, lng = 0.00, zoom = 0) %>%
  addProviderTiles(providers$Stamen.Toner)

# Uncomment below to view the map!
# la_map

Set Geometry and Transform Data CRS

We will be looking to convert our data into an object with coordinate-based geometry so that we can map it. We can use the st_as_sf() function to do this. st_as_sf() takes two arguments:

  1. The names of the columns in our data frame containing geographic coodinates in the format coords = c("<longitude column names>", "<latitude column name>")
  2. The coordinate reference system that the dataset currently uses in the format crs = <crs number>.

After adding this geometry column with st_as_sf(), we need to make sure that the CRS for our coordinates is consistent with the CRS of the basemap we will be placing the points on. The coordinates in tri_la_2023 data are in NAD 83 (EPSG:4269) and our basemap is in WGS 84 (EPSG:4326).

Question

In the code block below, create a new geometry column with the given coordinates, using the function st_as_sf(). The column for longitude will go in the first coordinate position, and latitude will go in the second. Be sure to set the data’s current CRS (4269) in that function. Then transform the CRS to 4326 using st_transform().

# Uncomment and fill in the blanks below. The column for longitude will go in the first coordinate position, and latitude will go in the second. 

#tri_la_2023 <- tri_la_2023 %>%
#  st_as_sf(coords = c("____", "____"), crs = ____) %>%
#  st_transform(____)

Add Layers

Once we have our initial map, we can add layers to the map that display different forms of geospatial data. There are a number of different functions in leaflet that we can use to add layers. For instance, we can add markers to the map at a certain geo-coordinates using the addMarkers() function. We can also add polygons and rectangles to the map using the addRectangles() function or the addPolygons() function. Today we are going to work exclusively with the addCircleMarkers() function. This allows us to add a circle at the latitude and longitude for each row in our dataset. It also allows us to adjust the circle’s color and size according to values in our dataset.

Question

In the code block below, pipe addCircleMarkers() onto the basemap and list data = tri_la_2023 as an argument in that function. It should look like my map below.

# Uncomment and add the function

#la_map %>%

Map isn’t so legible/beautiful at this point, right?

Styling the Map

Question

Copy the map you just created into the code chunk below, and pipe the addCircleMarkers() function onto your code. Check out the help pages for the addCircleMarkers() function, and add some arguments to help with the map’s legibility. At the very least, you should adjust the radius, weight, color, fillColor, and fillOpacity, and understand how each of these arguments will change the style of the map. For now you can set the color to “black” and the fillColor to “white”. See if you can get your map to match mine below.

# Copy previous map here! 

Sizing the Circles

Question

Copy and paste the last map that you created below. We are going to size each circle by the total chemical releases at that facility. Total releases is stored in the column X107..TOTAL.RELEASES in our dataset, so we could size the circles by setting the radius to ~X107..TOTAL.RELEASES. However, this would lead to some massive circles on our map as certain facilities have released hundreds of thousands of pounds in emissions, and the value we supply to radius will determine the pixels of the circle on our map. To deal with this, we are going to divide the total releases by 1000000. Set the radius to ~X107..TOTAL.RELEASES/1000000 below.

# Copy map colored by bins here!

Creating Color Palettes

Now that our map is looking more legible, let’s color the circles by the type of industry for each facility. Remembering back to our lesson on Understanding Datasets and Visualization Aesthetics, we should keep in mind that X23..INDUSTRY.SECTOR is a categorical variable, and therefore, we will create a qualitative color palette to map it.

There is one function in leaflet for creating a categorical color palette:

  1. colorFactor(): Creates a palette by assigning a color to each unique value

This function takes the following arguments:

  • palette: the colors we wish to map to the data. We will use preset palettes from the RColorBrewer. You can call display.brewer.all() to see the list of palettes or reference here: http://applied-r.com/rcolorbrewer-palettes/
  • domain: the values we wish to apply the palette to. Here we reference the column from the data frame we wish to color the points by using the accessor $.
Question

Create a categorical palettes below using colorFactor(), setting the palette to “Set2”, and the domain to tri_la_2023$X23..INDUSTRY.SECTOR.

#Uncomment below to create the palette

#pal_fact <- 

Now that we have created a palette, we can map the colors to the points in our map.

Question

Copy and paste the last map you created into the code chunk below. We are going to adjust the fillColor by setting it to the variable from the dataset that we wish to color (X23.INDUSTRY.SECTOR), colored by one of the palettes that we created. We can do this by setting the fill color argument equal to ~pal_fact(X23.INDUSTRY.SECTOR).

We also need to add a legend to explain the colors. Add a pipe to the end of the addCircleMarkers() and then add the function addLegend(). Consult the help pages for the addLegend() function to determine how to add a legend for the meaning of the colors represented on the map. At the very least, you will need an argument for data, pal, and values.

Your map should look like mine below.

# Copy the map here!

Now there are 19 unique industry sectors in this dataset, and that may be too many unique values to be able to readily recognize differentiation in color. Let’s change up how we approach coloring this map. Instead, let’s color the circles by the 2023 total releases at each facility on the map. Keep in mind that X107..TOTAL.RELEASES is a numeric variable, and therefore, we will create a sequential color palette to map it.

There are three functions in leaflet for creating a sequential color palette:

  1. colorNumeric(): Creates a palette by assigning numbers to different colors on a spectrum
  2. colorBin(): Creates a palette by grouping numbers into a specified number of equally-spaced intervals (e.g. 0-10, >10-20, >20-30)
  3. colorQuantile(): Creates a palette by grouping numbers into a specified number of equally-sized quantiles

Each of these functions takes the arguments palette and domain, just as we did when calling colorFactor(). In addition, the colorBin() function takes the argument bins, which indicates the number of color intervals to be created. The colorQuantile() functions takes the argument n, which indicates the number of quantiles to map data into.

Question

Create three palettes below (one using each function colorNumeric(), colorBin(), and colorQuantile()), setting the palette to “YlOrRd”, and the domain to tri_la_2023$X107..TOTAL.RELEASES. For colorBin(), set the bins to 6, and for colorQuantile(), set n to 4.

pal_num <- colorNumeric(palette="YlOrRd", 
                        domain = tri_la_2023$X107..TOTAL.RELEASES)
#Uncomment below to create the other two palettes.

#pal_bin <- 

#pal_quant <- 

Question

Copy and paste the last map you created into the code chunk below, three times. For each, we are going to adjust the fillColor by setting it to the variable from the dataset that we wish to color (X107..TOTAL.RELEASES), colored by one of the palettes that we created. We can do this by setting the fill color argument equal to:

  • ~pal_num(X107..TOTAL.RELEASES): for coloring the points according the pal_num() palette we created on the X107..TOTAL.RELEASES variable
  • ~pal_bin(X107..TOTAL.RELEASES): for coloring the points according the pal_bin() palette we created on the X107..TOTAL.RELEASES variable
  • ~pal_quant(X107..TOTAL.RELEASES): for coloring the points according the pal_quant palette we created on the X107..TOTAL.RELEASES() variable

Adjust the legends for each of the maps.

Your maps should look like mine below.

# Copy map first time here!

# Copy map second time here!

# Copy map third time here!

Question

Below explain the following: Why do the colors appear differently on each map? Which map best represents the distribution of values in X107..TOTAL.RELEASES?

Popups

Leaflet popups display information about a point we click on it. They are configured via Hypertext Markup Language (HTML). HTML is the markup language used to create webpages. It is akin to RMarkdown; for both HTML and RMarkdown we use certain characters to signal how we want our text formatted and styled.

There are a few HTML tags that you need to know in order to style a popup:

  • <br> creates a newline
  • <b></b> makes the enclosing text bold
  • <i></i> makes the enclosing text italicized

So let’s say I had my name stored in the variable name and my institution stored in the variable institution. If I wanted to create a popup with my name in bold and my institution beneath it italics, I would use the following HTML tags: paste("<b>", name, "</b><br><i>", institution, "</i>"). Here I’ve pasted together my name (wrapped in the bold tags), a break, and my instituion (wrapped in the italic tags).

Question

Now, we are going to add popups on the map. Copy and paste the code from the previous step into the code chunk below. Then do the following:

  1. Add the popup argument to the addCircleMarkers() function. Using the paste() function, create a popup that includes the facility’s name (in bold), the facility’s industry sector (in italics), and the chemicals released by the facility. (Note that the last the data point is stored in X..CHEMICAL.TEXT). Note that you will need to place a ~ in front of the paste() function to communicate that there are variables in your argument.

  2. Add the popupOptions argument to the addCircleMarkers() function. Set the popupOptions to popupOptions(maxHeight = 200, closeOnClick = TRUE)). This will ensure that really tall popups will scroll vertically rather than taking up an entire page.

Your popups should look like mine in the map below.

# Copy previous map here!

Ethical Considerations

Notably, while the Emergency Planning and Community Right-to-Know Act mandates reporting of emissions, it does not mandate monitoring of emissions. While other environmental regulations do set certain monitoring standards for specific TRI chemicals and pollution activities, for all other chemicals and activities, facilities are required to report based on a “reasonable estimate” of releases and other waste management quantities. In the early 1990s the National Wildlife Federation uncovered that certain facilities showing tremendous reductions in emissions from one year to the next had not actually cut down on pollution, but instead developed different ways of estimating their emissions. The public doesn’t see how emissions are calculated; they only see the final reported emissions. Given this history of certain facilities reporting “phantom reductions”, what kinds of legal, practice-based, or social interventions might ensure that this data remains reliable and accountable to the communities that use it to assess risks? What might make these interventions challenging to implement; what would need to change in order for these interventions to be realized? Share your ideas on our sds-192-discussions Slack channel.