Created
March 8, 2014 12:27
-
-
Save chrishanretty/9429801 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
library(rgeos) | |
library(sp) ## the best package for handling shapefiles | |
library(spdep) ## for the function poly2nb | |
library(maptools) ## for the function readShapePoly | |
library(RColorBrewer) ## for pretty colour palettes | |
### Read in the shapefile | |
x <- readShapePoly("data/westminster_const_region.shp") | |
## Add the Press Association number | |
name2pa <- read.csv("data/name2pa.csv",header=T) | |
x$refno <- name2pa$refno[match(as.character(x$NAME),name2pa$ShapeFileName)] | |
# Create spatial polygons object | |
gpclibPermit() | |
sp.const <- unionSpatialPolygons(x, x$refno) ## may take a while | |
# Lose some detail | |
sp.const <- gSimplify(sp.const,100) | |
# Find out how many neighbours each constituency has | |
neighbouring.districts <- poly2nb(sp.const) | |
## Let's turn this into a data frame we can use | |
nb.df <- data.frame(refno = names(sp.const), | |
numNeighbours = card(neighbouring.districts)) | |
## Add some other info | |
nb.df <- merge(nb.df,name2pa,all.x=T,all.y=F) | |
## Check islands as a sanity check | |
nb.df[which(nb.df$numNeighbours==0),] | |
# Create spatial polygon data frame | |
# Use function SpatialPolygonsDataFrame from package "sp" (Pebesma 2005). | |
# This matches external data in with geodata from sp.const: | |
# Rownames must match | |
rownames(nb.df) <- nb.df$refno | |
area.spdf <- SpatialPolygonsDataFrame(sp.const, nb.df) | |
# Now plot this | |
manual.colours <- cut(area.spdf$numNeighbours, | |
breaks = quantile(area.spdf$numNeighbours,probs=seq(0,1,length.out=7))) | |
manual.colours.numeric <- as.numeric(manual.colours) | |
my.palette <- brewer.pal(7,name="PuBuGn") | |
manual.colours.hex <- my.palette[manual.colours.numeric] | |
plot(area.spdf, col = manual.colours.hex) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment