For one of my projects, I needed a shapefile of Federal Reserve Districts for counties across the continental United States. This wasn’t easily available online, so I took the county file I had, went through each Fed Reserve District website (there are 12 in total), and just tagged counties according to the district they resided within. (I did this in GeoDa using brushing and linking.) Several districts cut across states in unexpected ways, so our data needs to be at the county level, despite districts expanding across states. This shapefile will be shared once I clean up the data a bit more; I just needed information matching the 1930 census, so Wyoming and a handful of counties will need to be added in for a more contemporary version.
Meanwhile, I wanted to share some easy code for visualizing your data in categorical maps. You can do it in QGIS or GeoDa, but with more people trying out map viz and spatial analytics in R and python, I wanted to share the “easy” way of doing it. I was surprised at how complicated some of the online tutorials are for this, so here are quick tutorials to get to your goal in 5 lines or less — assuming your data is a shapefile already. (And keep in mind, your “shapefile” will actually be a group of files: .shx, .shp, .dbf and .prj file).
Visualizing the districts in R, using rgdal
library(rgdal) # Import library setwd("~/Data/Bank_study/Counties") #Set working directory fedRes<-readOGR(".","FedReserveDistricts") #Read the shapefile names(fedRes) # Check the column names. spplot(fedRes,"FEDDIST") #Make a categorical map of US fed districts
This takes a while in R because of the number of counties involved. It shouldn’t, and there’s likely a better library to use, but here’s this for now. The other libraries I tried out were more complicated to get this far, and I prefer “simple” for simple tasks.
Visualizing the districts in python, using pysal
import pysal as ps import numpy as np from pysal.contrib.viz import mapping as maps values = np.array(ps.open(('~/Counties/FedReserveDistricts.dbf')).by_col('FEDDIST')) maps.plot_choropleth(()'~/Counties/FedReserveDistricts.shp'), values, 'unique_values')
This was pretty fast, but the figure is a bit distorted as matplotlib doesn’t nicely account for projections. You can change the figure manually with styling, of course.