Crime Mapping In R with Leaflet

Introduction

The Leaflet R package is a powerful tool that makes it easy to build interactive maps using the open-source JavaScript library, Leaflet. This tutorial demonstrates how to use leaflet to create beautiful, interactive crime maps.

Getting Started

To begin, we’ll install leaflet: install.packages("leaflet"). We’ll also install leaflet.extras for some additional map modification tools: install.packages("leaflet.extras). Then we’ll load the packages as well as a handful of other useful packages that we’ll need to carry out the following analyses.

library(leaflet)
library(leaflet.extras)
library(ggplot2)
library(tidyverse)
library(reshape2)

Data Preparation

We will be mapping crime data from Columbus, Ohio. The panel data set, which we call crime_data contains geo-coded crime incident data from 2013 to 2019. yearly crime rates in each geographic unit were calculated by dividing raw counts of crimes within each block group by the block group’s total population and multiplying by 1000.

Re-shaping the Data

Our goal is to analyze how the violent crime rate changed from 2013 to 2019 within block groups. We’ll begin by filtering our data so only information from 2013 and 2019 and re-shaping the data to wide form. To calculate the change in the violent crime rate from 2013 to 2019, we simply subtract the 2013 rate from the 2019 rate.

crime_data_wide <- 
  #Filters out 2013 and 2019 and retains the variables needed for the analysis
  crime_data %>% 
  filter(year==2013|year==2019) %>% 
  dplyr::select(c("block_group","year","viol_index_crime_rate_rep")) %>% 
  mutate(block_group = as.character(block_group)) %>% 
  rename(viol_crime_rate = viol_index_crime_rate_rep) %>% 
  #Casts data to wide form
  pivot_wider(., names_from = year, values_from = viol_crime_rate, names_prefix = "viol_crime_rate_") %>% 
  #Calculate change in violent crime from 2013 to 2019
  rowwise() %>% 
  mutate(viol_crime_change = viol_crime_rate_2019 - viol_crime_rate_2013) %>% 
  #Remove rows with any missing values
  na.omit()

After re-shaping the data we are left with the following data structure:

block_group	viol_crime_rate_2013	viol_crime_rate_2019	viol_crime_change
390490001101	0.000000	0.0000000	0.000000
390490001102	1.114827	0.0000000	-1.114827
390490001103	2.597403	0.0000000	-2.597403
390490001104	4.566210	1.3947001	-3.171510
390490001201	0.000000	1.2019231	1.201923
390490001202	0.000000	0.0000000	0.000000
390490001203	1.602564	0.0000000	-1.602564
390490002101	2.378121	0.9930487	-1.385073
390490002102	0.000000	0.0000000	0.000000
390490002103	0.000000	0.0000000	0.000000

Converting Data to a Shapefile

We’ll now need to convert the data to a shapefile. To do so, we’ll first need to get GIS information for block groups in Columbus, Ohio. To achieve this, we’ll use the block_groups() function from the tigris package. If your data already contains spatial polygon information, this step can be skipped. In this case, however, we’ll gather block groups from Franklin County, Ohio and merge the geometries with our crime rate data.

#install.packages("tigris)
library(tigris)

fc_block_groups <- 
  block_groups("Ohio", "Franklin",year=2023) %>% 
  dplyr::select(c("GEOID")) %>% 
  rename(block_group = "GEOID")

#Merge shapefiie with crime rate data
crime_data_geo <- left_join(crime_data_wide ,fc_block_groups,by="block_group")

After merging the block group polygons onto the crime rate data, we can convert the data to an sfsimple feature object. We can do this by using the st_as_sf() function from the sf package. Sf objects are a special class of spatial data that contain geographic features, like polygons and points.

#install.packages("sf)
library(sf)

#create shapefile
crime_data_sf <- st_as_sf(crime_data_geo)
#Set coordinate reference system
crime_data_sf <- st_transform(crime_data_sf,"+proj=longlat +datum=WGS84")

print(class(crime_data_sf))

[1] "sf"         "rowwise_df" "tbl_df"     "tbl"        "data.frame"

Mapping with Leaflet

Now that our data is ready we can begin mapping violent crime using leaflet.

Creating Base Map:

The first step to creating a leaflet map is creating the base map layer.

crime_map <-
  leaflet() %>% 
  addTiles() 

crime_map

Adding Block Groups

Next, we’ll add the block groups contained in crime_data_sf to the map. We’ll do so with the addPolygons command. Note that new layers are added to the map using the tidyverse piping operator (%>%).

We’ll set the outline color of the polygons to white; the fill color to grey; the fill opacity to 0.5; the opacity of the outline to 1; and the weight of the outline to 1. These parameters can be adjusted as needed.

crime_map %>% 
  addPolygons(data = crime_data_sf,
              color = "white",
              fillColor = "grey",
              opacity = 1,
              weight = 1,
              fillOpacity = 0.5)

Coloring Map According to Crime Rate

While adding the polygons to the map is a good start, what we really want to accomplish is mapping violent crime. To achieve this, we need to add color to the map corresponding to different levels of violent crime. Adjusting the color according to levels of violent crime can be accomplished by replacing the fillColor argument with a color palette that corresponds to the violent crime rate and changes in violent crime.

Naturally, then, the next step is to create a color palette that captures changes in violent crime in a visually appealing way. Creating a color palette for the violent crime is relatively simple. In our data, the violent crime rate for 2019 is distributed as follows:

crime_data_sf %>% 
  ggplot(.,aes(x=viol_crime_rate_2019)) + 
  geom_histogram(binwidth = 0.5) +
  theme_minimal() +
  xlab("Violent Crime Rate in 2019")

Let’s also look at the summary statistics

summary(crime_data_sf$viol_crime_rate_2019)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  0.000   0.000   1.110   4.104   5.310 129.167 

The distribution informs the cut points needed to create a color palette using colorBin. More information on R color palettes can be found here. In our case, we’ll be using the “YlorRd” (Yellow-Red) palette to map violent crime rates from 2019. We’ll apply the domain to crime_data_sf$viol_crime_rate_2019 and set the bins based to the distribution.

#| warning: false
#| message: false
#| include: true

viol_rate_2019_palette <- 
  colorBin("YlOrRd", 
           domain = crime_data_sf$viol_crime_rate_2019,
           bins = c(0,0.25,0.5,0.75,1,1.25,1.5,2,3,5,10,20))

After we’ve created the color palette, we can set the fillColor argument accordingly. Doing so applies the color palette that we created to the values of violent crimes associated with each layer of block group polygons. We’ll also add a corresponding legend using addLegend.

crime_map %>% 
  addPolygons(data = crime_data_sf,
              color = "white",
              fillColor = viol_rate_2019_palette(crime_data_sf$viol_crime_rate_2019),
              opacity = 1,
              weight = 1,
              fillOpacity = 0.5) %>% 
  addLegend(pal = viol_rate_2019_palette, 
            values = crime_data_sf$viol_crime_rate_2019, 
            opacity = 0.5, 
            title = "Violent Crime Rate 2019",
            position="topleft")

Coloring Map According to Change in Violent Crime

While we’ve shown how to map the violent crime rate by color, mapping changes in the violent crime rate is somewhat more complicated. This is because change values can take on negative and positive values. From a visual perspective, it is beneficial to color negative values differently than positive values. A useful tool for selecting a color palette that provides such a contrast is ColorBrewer 2.

The code below creates a palette for values of change in violent crime. The palette applies a blue color scale to negative values (reflecting a decline in violent crime) and a red color palette to positive values (reflecting increases in violent crime). First, we’ll explore the distribution of the change variable to determine appropriate cut points.

crime_data_sf %>% 
  ggplot(.,aes(x=viol_crime_change)) + 
  geom_histogram(binwidth = 0.5) +
  theme_minimal() +
  xlab("Change in Violent Crime Rate")

We’ll also take a look at the summary statistics:

summary(crime_data_sf$viol_crime_change)

    Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
-48.3044  -1.5819   0.0000  -1.1274   0.2743  72.9486 

It looks like there are some problematic outliers. We’ll remove them and check the summary statistic again before creating the color palette.

crime_data_sf <-
  crime_data_sf %>% 
  filter(viol_crime_change<=25&viol_crime_change>=-25)

summary(crime_data_sf$viol_crime_change)

    Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
-24.8807  -1.4503   0.0000  -0.8884   0.2810  18.8959 

Based on the distribution, we’ll now create a color palette for the change in violent crime. We’ll get a suitable “diverging” color palette from Color Brewer 2. We’ll select the 11-class “RdBu” palette and apply this palette to violent_crime_change. We need to ensure that the number of bins specified with the bins argument in colorBin matches the number of colors in the Color Brewer palette. Ideally, zero will take on the middle neutral color. Finally, we’ll use rev() to reverse the order of the color vector so that red represents increases in violent crime.

change_palette <- 
  colorBin(palette = rev(c('#67001f','#b2182b','#d6604d','#f4a582','#fddbc7','#f7f7f7','#d1e5f0','#92c5de','#4393c3','#2166ac','#053061')),
           domain = crime_data_sf$viol_crime_change,
           bins =  c(-25,-15,-10,-6,-3,-1,0,1,3,6,10,15,25))

Now that we’ve specified the palette, we can map the change in violent crime on our leaflet map.

crime_map %>% 
  addPolygons(data = crime_data_sf,
              color = "white",
              fillColor = change_palette(crime_data_sf$viol_crime_change),
              opacity = 1,
              weight = 1,
              fillOpacity = 0.8,) %>% 
  addLegend(pal = change_palette, 
            values = crime_data_sf$viol_crime_change, 
            opacity = 0.5, 
            title = "Change in Violent Crime",
            position="topleft")

To further improve the map, we can add pop-up labels that show the change in violent crime in a given block group. To do so, we first create a label using HTML syntax and then apply the label with label and modify the label’s position using labelOptions. We also add highlightOptions which highlight the block group over which the courser is hovering.

#Create label:
change_label <- sprintf(
  "Change in Violent Crime Rate: <strong>%g</strong><br/>",
  round(crime_data_sf$viol_crime_change,2)
) %>% lapply(htmltools::HTML)

crime_map %>% 
  addPolygons(data = crime_data_sf,
              color = "white",
              fillColor = change_palette(crime_data_sf$viol_crime_change),
              opacity = 1,
              weight = 1,
              fillOpacity = 0.8,
              highlightOptions = highlightOptions(
              weight = 5,
              color = "#666",
              fillOpacity = 0.7,
              bringToFront = FALSE),
              label = change_label,
              labelOptions = labelOptions(
              style = list("font-weight" = "normal", padding = "3px 8px"),
              textsize = "15px",
              direction = "auto")) %>% 
  addLegend(pal = change_palette, 
            values = crime_data_sf$viol_crime_change, 
            opacity = 0.5, 
            title = "Change in Violent Crime",
            position="topleft")