2011 Open Atlas Project: Output Area Key Statistics

1. Download and consolidate 2011 Output Area Level Key Statistics

The first stage is to create a list of the tables that you want to download - these are constructed using the Nomis naming conventions for the download files - e.g. ks102ew_2011_oa for Key Statistics table 102 England and Wales at Output Area level.


# Create a list of tables
table_list <- c("KS101EW", "KS102EW", "KS103EW", "KS104EW", "KS105EW", "KS106EW", 
    "KS107EW", "KS201EW", "KS202EW", "KS204EW", "KS205EW", "KS206EW", "KS207WA", 
    "KS208WA", "KS209EW", "KS301EW", "KS401EW", "KS402EW", "KS403EW", "KS404EW", 
    "KS405EW", "KS501EW", "KS601EW", "KS602EW", "KS603EW", "KS604EW", "KS605EW", 
    "KS606EW", "KS607EW", "KS608EW", "KS609EW", "KS610EW", "KS611EW", "KS612EW", 
    "KS613EW")

table_list <- tolower(paste(table_list, "_2011_oa", sep = ""))

Download and unzip the files for each of the OA tables.

# Set a download folder for the census CSV files setwd('/Volumes/Macintosh
# HD 2/Dropbox/Projects/2011_Census/2011_Data')
setwd("/Volumes/Macintosh HD 2/Dropbox/Projects/2011_Census/2011_Data")

# Download Files
for (n in 1:length(table_list)) {
    file <- as.character(table_list[n])
    temp <- tempfile(fileext = ".zip")
    download.file(paste("http://www.nomisweb.co.uk/output/census/2011/", file, 
        ".zip", sep = ""), temp)
    unzip(temp, junkpaths = TRUE)
    unlink(temp)
    csv_file <- paste(toupper(sub("_2011_oa", "", file)), "DATA.CSV", sep = "")
    assign(file, read.csv(csv_file))
}

Next, compile a single object containing the table variable names and their corresponding descriptions.

Variable_Desc <- NULL
for (n in 1:length(table_list)) {
    file <- as.character(table_list[n])
    csv_file <- paste(toupper(sub("_2011_oa", "", file)), "DESC0.CSV", sep = "")
    temp <- read.csv(csv_file)
    Variable_Desc <- rbind(Variable_Desc, temp)
    remove(temp)
}

Finally, we will create a single file containing all of the 2011 key statistics using the tables that cover both England and Wales.

# Remove Welsh only tables from list
table_list <- table_list[regexpr("wa", table_list) < 0]

# Create lookup
Key_Statistics_2011 <- as.data.frame(ks101ew_2011_oa[, 1])
colnames(Key_Statistics_2011) <- "GeographyCode"

# Merge Loop
for (n in 1:length(table_list)) {
    Key_Statistics_2011 <- data.frame(Key_Statistics_2011, get(as.character(table_list[n]))[match(get(as.character(table_list[n]))[, 
        "GeographyCode"], get(as.character(table_list[n]))[, "GeographyCode"]), 
        ])
    Key_Statistics_2011$GeographyCode.1 <- NULL
}

# Export CSV
setwd("/Volumes/Macintosh HD 2/Dropbox/Projects/2011_Census")
write.csv(Key_Statistics_2011, "Census_2011OA_Key_Statistics.csv")

2. Create local authority district map books

In this next section we will run through the method of creating an OA level map of the 2011 Census data for each variable in the Key Statistics tables.

# Load Libraries
library(rgdal)
library(RColorBrewer)
library("sp")
library(GISTools)
library(classInt)
library(maptools)


# Get 2011 OA
temp <- tempfile(fileext = ".zip")
download.file("http://www.ons.gov.uk/ons/external-links/other-ns-online/census-geography/2011-oa-boundary/oa-2011-ew-bgc-shp.html", 
    temp)
unzip(temp)
OA <- readOGR(".", "OA_2011_EW_BGC")
unlink(temp)

# Get 2011 Merged Wards
temp <- tempfile(fileext = ".zip")
download.file("http://www.ons.gov.uk/ons/external-links/other-ns-online/census-geography/2011-cmw-boundary/cmwd-2011-ew-bgc-shp.html", 
    temp)
unzip(temp)
WARD <- readOGR(".", "CMWD_2011_EW_BGC")
unlink(temp)

# Get 2011 OA higher geog lookup
temp <- tempfile(fileext = ".zip")
download.file("http://www.ons.gov.uk/ons/external-links/other-ns-online/census-geography/exact-fit-2011-lookup/oa11-lsoa11-msoa11-lad11-ew-lu-csv.html", 
    temp)
unzip(temp)
oa_higher_lookup <- read.csv("OA11_LSOA11_MSOA11_LAD11_EW_LUv2.csv")
unlink(temp)

# Create key stats table List
a <- c("KS101EW", "KS102EW", "KS103EW", "KS104EW", "KS105EW", "KS106EW", "KS107EW", 
    "KS201EW", "KS202EW", "KS204EW", "KS205EW", "KS206EW", "KS207WA", "KS208WA", 
    "KS209EW", "KS301EW", "KS401EW", "KS402EW", "KS403EW", "KS404EW", "KS405EW", 
    "KS501EW", "KS601EW", "KS602EW", "KS603EW", "KS604EW", "KS605EW", "KS606EW", 
    "KS607EW", "KS608EW", "KS609EW", "KS610EW", "KS611EW", "KS612EW", "KS613EW")
b <- c("Usual resident population", "Age structure", "Marital and civil partnership status", 
    "Living arrangements", "Household composition", "Adults not in employment and dependent children and persons with long-term health problem or disability for all households", 
    "Lone parent households with dependent children", "Ethnic group", "National identity", 
    "Country of birth", "Passports held", "Household language", "Welsh language skills", 
    "Welsh language profile", "Religion", "Health and provision of unpaid care", 
    "Dwellings, household spaces and accommodation type", "Tenure", "Rooms, bedrooms and central heating", 
    "Car or van availability", "Communal establishment residents", "Qualifications and students", 
    "Economic activity", "Economic activity - Males", "Economic activity - Females", 
    "Hours worked", "Industry", "Industry - Males", "Industry - Females", "Occupation", 
    "Occupation - Males", "Occupation - Females", "NS-SeC", "NS-SeC - Males", 
    "NS-SeC - Females")


Key_stat_tables <- data.frame(a, b)
colnames(Key_stat_tables) <- c("KS2011_Tab_Code", "KS2011_Tab_Name")
# remove welsh only tables
Key_stat_tables <- Key_stat_tables[regexpr("WA", Key_stat_tables$KS2011_Tab_Code) < 
    0, ]

The next stage is to merge the the 2011 Key Statistics data with the Output Area spatial polygons data frame.

# Join key statistics data onto 2011 OA
OA@data = data.frame(OA@data, Key_Statistics_2011[match(OA@data[, "OA11CD"], 
    Key_Statistics_2011[, "GeographyCode"]), ])

# Join higher geog onto OA
OA@data = data.frame(OA@data, oa_higher_lookup[match(OA@data[, "OA11CD"], oa_higher_lookup[, 
    "OA11CD"]), ])

Ok, we are now ready to do some mapping. This code creates a map book for each local authority district. The maps feature the percentage values from each Key Statistic table; and choropleth maps are presented for all those local authority districts where there were more than 8 output areas with values. This restriction was made for aesthetic restrictions related to the manipulation of the percentage scores into the coloured bins used to shade the choropleth.

# Set working directory to output maps
setwd("/Volumes/Macintosh HD 2/Dropbox/Projects/2011_Census/2011_Output_Area_Maps_LA")

# Map Colours
my_colours <- brewer.pal(7, "Blues")

# Create LAD list
LAD_LIST <- as.data.frame(table(OA@data$LAD11CD, OA@data$LAD11NM))
LAD_LIST <- LAD_LIST[LAD_LIST$Freq > 0, c("Var1", "Var2")]
colnames(LAD_LIST) <- c("LAD11CD", "LAD11NM")
# LAD_LIST <- LAD_LIST[c(sample(1:nrow(LAD_LIST),
# 1),sample(1:nrow(LAD_LIST), 1),sample(1:nrow(LAD_LIST), 1)),1:2]

# Create Map list
map_list <- as.character(Variable_Desc[Variable_Desc$ColumnVariableMeasurementUnit != 
    "Count", 1])
# Remove Welsh only variables from map list
map_list <- map_list[regexpr("WA", map_list) < 0]

# map_list <-
# c('KS405EW0015','KS405EW0020','KS405EW0025','KS608EW0012','KS105EW0024')


# Loop through each LAD
for (i in 17:nrow(LAD_LIST)) {
    LAD <- as.character(LAD_LIST[i, 1])
    LAD_NAME <- as.character(LAD_LIST[i, 2])

    pdf(paste(LAD, ".pdf", sep = ""), width = 10, height = 10, title = paste(LAD_NAME), 
        paper = "a4")

    WARD_to_Map <- WARD[WARD@data$LAD11CD == LAD, ]

    for (x in 1:length(map_list)) {

        # Get a name for the variable
        MAP_NAME <- Variable_Desc[Variable_Desc$ColumnVariableCode == paste(as.character(map_list[x])), 
            ]
        OA_to_Map <- OA[OA@data$LAD11CD == LAD, ]
        tmp <- as.numeric(OA_to_Map@data[, paste(as.character(map_list[x]))])
        OA_to_Map@data$mapvar <- tmp

        # Get table name, and add to MAP_NAMES
        MAP_NAME$KS2011_Tab_Code <- substring(MAP_NAME$ColumnVariableCode, 1, 
            7)
        MAP_NAME = data.frame(MAP_NAME, Key_stat_tables[match(MAP_NAME[, "KS2011_Tab_Code"], 
            Key_stat_tables[, "KS2011_Tab_Code"]), ])



        # Mapping block set breaks for 6 categories and then add the values to a
        # single list - style = fisher is jenks a test is made to check that all
        # the values aren't zero
        if ((table(OA_to_Map@data$mapvar)[1] - length(OA_to_Map@data$mapvar)) %in% 
            c(-1, -2, -3, -4, -5, -6, -7, -8, 0)) {


        } else {

            breaks <- classIntervals(OA_to_Map@data$mapvar, n = 7, style = "fisher", 
                unique = TRUE)
            breaks <- breaks$brks
            plot(OA_to_Map, col = my_colours[findInterval(OA_to_Map@data$mapvar, 
                breaks, all.inside = TRUE)], axes = FALSE, border = NA)

            plot(WARD_to_Map, border = "#707070", add = TRUE)

            # Add on text labels for the wards
            pointLabel(coordinates(WARD_to_Map)[, 1], coordinates(WARD_to_Map)[, 
                2], labels = WARD_to_Map@data$CMWD11NM, cex = 0.5)

            # Add the legend
            legend("bottomleft", legend = leglabs(round(breaks, digits = 2), 
                between = " to "), fill = my_colours, bty = "n", cex = 0.7)

            mtext("Contains National Statistics data © Crown copyright and database right 2013; Contains Ordnance Survey data © Crown copyright and database right 2013", 
                side = 1, line = 2, cex = 0.5)
            mtext("Map created by Alex Singleton (http://www.alex-singleton.com)", 
                side = 1, line = 1, cex = 0.5)
            title(LAD_NAME)
            mtext(paste("Table: ", as.character(MAP_NAME[, 7]), " \n ", as.character(MAP_NAME[, 
                1]), " (", as.character(MAP_NAME[, 4]), ")", sep = ""), side = 3, 
                cex = 0.7)
            remove(list = c("MAP_NAME", "OA_to_Map", "tmp", "breaks"))

        }

    }

    dev.off()
    remove(list = c("WARD_to_Map", "LAD", "LAD_NAME"))

}