Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
R Choropleth with shape file
###
### Barry Rowlingson, Lancaster University
###
## needed for shapefiles:
require(rgdal)
## needed for colour mapping - not on CRAN:
## http://r-forge.r-project.org/projects/colourscheme/
## try:
## install.packages("colourschemes", repos="http://R-Forge.R-project.org")
require(colourschemes)
## pretty colours:
require(RColorBrewer)
## read the unemployment data:
unem=read.table("http://datasets.flowingdata.com/unemployment09.csv",sep=",",
as.is=TRUE,
colClass=c(rep("character",8),"numeric"),
col.names=c("code","f1","f2","name","yr","x1","x2","x3","percentage")
)
## compute the FIPS code:
unem$fips=paste(unem$f1,unem$f2,sep="")
## read the county data:
## source: http://www.census.gov/geo/www/cob/co2000.html
## this is in lat-long, with Alaska, Hawaii, and Puerto Rico
## in their proper geographic place. So I'm just going to
## plot the lower 48
county=readOGR("Maps/","co99_d00")
county$fips = paste(county$STATE,county$COUNTY,sep="")
## work out how the fips codes match up
m = match(county$fips,unem$fips)
## add the unemployment data to the county map:
county$unem = unem$percentage[m]
## fix a couple of missing data:
county$unem[is.na(county$unem)]=0
## use the purple-red palette:
colours = brewer.pal(6,"PuRd")
## make a colour scheme:
sd = data.frame(col=colours,values=c(1,3,5,7,9,11))
sc = nearestScheme(sd)
## set the plot region to the lower 48:
## the aspect ratio may not be right. Meh. This is a
## quick exercise! If I was doing this for real I'd find
## a better shapefile!
plot(c(-129,-61),c(21,53),type="n",axes=FALSE,xlab="",ylab="")
## add the counties
plot(county,col=sc(county$unem),add=TRUE,border="white",lwd=0.2)
@reubano

This comment has been minimized.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.