From d3a1ccf036bfc4c0d260b82de7e63d514fa8ba75 Mon Sep 17 00:00:00 2001
From: Jeremy Kidwell <bots@jeremykidwell.info>
Date: Fri, 25 Jan 2019 16:55:16 +0000
Subject: [PATCH] working on graphs, split into 2

---
 mapping_draft.Rmd       | 489 ++++++----------------------------------
 mapping_draft_part2.Rmd | 271 ++++++++++++++++++++++
 2 files changed, 346 insertions(+), 414 deletions(-)
 create mode 100644 mapping_draft_part2.Rmd

diff --git a/mapping_draft.Rmd b/mapping_draft.Rmd
index 7b5151d..012368b 100644
--- a/mapping_draft.Rmd
+++ b/mapping_draft.Rmd
@@ -7,7 +7,7 @@ affiliation: University of Birmingham
 institute: University of Birmingham
 e-mail: “[j.kidwell@bham.ac.uk](mailto:j.kidwell@bham.ac.uk)”
 date: “`r Sys.Date()`”
-bibliography: everything.bib
+bibliography: /Users/jeremy/Dropbox/bibtex/everything.bib
 linkcolor: black
 fontsize: 12pt
 ---
@@ -18,21 +18,32 @@ knitr::opts_chunk$set(fig.path='figures/', warning=FALSE, echo=FALSE, message=FA
 ```
 
 ```{r load_packages, include=FALSE}
+## Default repo
+setwd("/Users/jeremy/gits/mapping_environmental_action")
+
+local({r <- getOption("repos")
+       r["CRAN"] <- "https://cloud.r-project.org" 
+       options(repos=r)
+})
 # TODO: need to alter this to use new sf data class as sp is deprecated
 require(sf) # new simplefeature data class, supercedes sp in many ways
 require(rgdal)
 require(GISTools) # deprecated by sf()
+require(sp) # needed for proj4string, deprecated by sf()
 require(ggplot2)
 require(broom) # required for tidying SPDF to data.frame for ggplot2
 require(rgeos)
-require(sp) # needed for proj4string, deprecated by sf()
 require(ggmap)
 require(maptools)
 require(RCurl)
+require(tibble) # using for grouped bar plot
+require(tidyr)  # using for grouped bar plot
 require(dplyr)
+library(reshape2)  # using for grouped bar plot
 require(pander)
 require(scales)
 require(sqldf) # using sqldf to filter while loading very large data sets
+require(plotly) # allows for export of plots to dynamic web pages
 
 # Set up local workspace:
 if (dir.exists("data") == FALSE) {
@@ -73,7 +84,7 @@ coordinates(ecs) <- c("X", "Y")
 proj4string(ecs) = CRS(bng)
 ```
 
-There are `length(ecs)` eco-congregations in Scotland. By some measurements, particularly in terms of individual sites and possibly also with regards to volunteers, this makes Eco-Congregation Scotland one of the largest environmental third-sector groups in Scotland.[^159141043] 
+There are `r length(ecs)` eco-congregations in Scotland. By some measurements, particularly in terms of individual sites and possibly also with regards to volunteers, this makes Eco-Congregation Scotland one of the largest environmental third-sector groups in Scotland.[^159141043] 
 
 In seeking to conduct GIS and statistical analysis of ECS, it is important to note that there some ways in which these sites are statistically opaque. Our research conducted through interviews at a sampling of sites and analysis of a variety of documents suggests that there is a high level of diversity both in terms of the number of those participating in environmental action and the types of action underway at specific sites. Work at a particular site can also ebb and flow over the course of time. Of course, as research into other forms of activism and secular environmental NGOs has shown, this is no different from any other third sector volunteer group. Variability is a regular feature of groups involved in activism and/or environmental concern. 
 
@@ -103,9 +114,10 @@ ecs$registration <- as.Date(ecs$registration, "%Y-%m-%d")
 # ecs$award2 <- as.Date(ecs$award2)
 # ecs$award3 <- as.Date(ecs$award3)
 # ecs$award4 <- as.Date(ecs$award4)
+ecs_complete_cases <- ecs[complete.cases(ecs$year_begun),]
 ```
 
-The programme launched officially in 2001 at Dunblane Cathedral in Stirling and by 2005 the project had `length(ecs_complete_cases[ecs_complete_cases$year_begun < 2006, ])` congregations registered to be a part of the programme and 25 which had completed the curriculum successfully and received an Eco-Congregation award. By 2011, the number of registrations had tripled to `length(ecs_complete_cases[ecs_complete_cases$year_begun < 2012, ])` and the number of awarded congregations had quadrupled to `sum(ecs$award1 < "2012-01-01", na.rm=TRUE)`. This process of taking registrations and using a tiered award or recognition scheme is common to many voluntary organisations. The ECS curriculum was developed in part by consulting the Eco-Congregation England and Wales materials which had been released just a year earlier in 1999, though it has been subsequently revised, particularly with a major redesign in 2010. In the USA, a number of similar groups take a similar approach including Earth Ministry (earthministry.org) and Green Faith (greenfaith.org). 
+The programme launched officially in 2001 at Dunblane Cathedral in Stirling and by 2005 the project had `r length(ecs_complete_cases[ecs_complete_cases$year_begun < 2006, ])` congregations registered to be a part of the programme and 25 which had completed the curriculum successfully and received an Eco-Congregation award. By 2011, the number of registrations had tripled to `r length(ecs_complete_cases[ecs_complete_cases$year_begun < 2012, ])` and the number of awarded congregations had quadrupled to `r sum(ecs$award1 < "01/01/2012", na.rm=TRUE)`. This process of taking registrations and using a tiered award or recognition scheme is common to many voluntary organisations. The ECS curriculum was developed in part by consulting the Eco-Congregation England and Wales materials which had been released just a year earlier in 1999, though it has been subsequently revised, particularly with a major redesign in 2010. In the USA, a number of similar groups take a similar approach including Earth Ministry (earthministry.org) and Green Faith (greenfaith.org). 
 
 In the case of Eco-Congregation Scotland, congregations are invited to begin by "registering" their interest in the programme by completing a basic one-sided form. The next step requires the completion of an award application, which includes a facilitated curriculum called a "church check-up" and after an application is submitted, the site is visited and assessed by third-party volunteer assessors. Sites are invited to complete additional applications for further awards which are incremental (as is the application process). Transition communities, at least in the period reflected on their map, go through a similar process (though this does not involve the use of a supplied curriculum) by which they are marked first as "interested," become "active" and then gain "official" status.[^1554162]
 
@@ -114,23 +126,28 @@ In the case of Eco-Congregation Scotland, congregations are invited to begin by
 # read in polygon for Scottish admin boundaries
 # TODO: upload bundle of admin data to new zenodo repository and alter below to use new URLs
 # TODO: need to remove readOGR below once st_read is confirmed to be working as sf
+
+if (file.exists("data/scotland_ca_2010.shp") == FALSE) {
 download.file("https://borders.ukdataservice.ac.uk/ukborders/easy_download/prebuilt/shape/Scotland_ca_2010.zip", 
               destfile = "data/Scotland_ca_2010.zip")
 unzip("data/Scotland_ca_2010.zip", exdir = "data")
+}
 admin_lev1 <- readOGR("./data", "scotland_ca_2010") 
 admin_lev1_sf <- st_read("data/scotland_ca_2010.shp")
 
 # read in polygon for intermediate admin boundary layers
+if (file.exists("data/scotland_parlcon_2011.shp") == FALSE) {
 download.file("http://census.edina.ac.uk/ukborders/easy_download/prebuilt/shape/Scotland_parlcon_2011.zip", 
               destfile = "data/Scotland_parlcon_2011.zip")
 unzip("data/Scotland_parlcon_2011.zip", exdir = "data")
+}
 admin_lev2 <- readOGR("./data", "scotland_parlcon_2011")
 
-# TODO note - need to optimise here, pre-filter sqldf? shapefile on "geo_code"" or "label" begins with "S"
-download.file("http://census.edina.ac.uk/ukborders/easy_download/prebuilt/shape/infuse_ward_lyr_2011.zip", 
-              destfile = "data/infuse_ward_lyr_2011.zip")
-unzip("data/infuse_ward_lyr_2011.zip", exdir = "data")
-admin_lev3 <- readOGR("./data", "infuse_ward_lyr_2011")
+## TODO note - need to optimise here, pre-filter sqldf? shapefile on "geo_code"" or "label" begins with "S"
+# download.file("http://census.edina.ac.uk/ukborders/easy_download/prebuilt/shape/infuse_ward_lyr_2011.zip", 
+#              destfile = "data/infuse_ward_lyr_2011.zip")
+# unzip("data/infuse_ward_lyr_2011.zip", exdir = "data")
+# admin_lev3 <- readOGR("./data", "infuse_ward_lyr_2011")
 
 # read in Transition Towns data and turn it into a SpatialPointsDataFrame
 transition <- read.csv(text=getURL("https://zenodo.org/record/165519/files/SCCAN_1.4.csv"))
@@ -171,6 +188,8 @@ proj4string(permaculture) <- proj4string(admin_lev1)
 
 # This code will generate a table of frequencies for each spatialpointsdataframe in admin
 # calculate count of ECS for fields in admin and provide percentages
+# JK Note: need to convert from poly.counts, which uses sp() to st_covers() which uses sf() - cf. https://stackoverflow.com/questions/45314094/equivalent-of-poly-counts-to-count-lat-long-pairs-falling-inside-of-polygons-w#45337050
+
 admin_lev1$ecs_count <- poly.counts(ecs,admin_lev1)
 admin_lev1$ecs_percent<- prop.table(admin_lev1$ecs_count)
 # calculate count of places of worship in PointX db for fields in admin and provide percentages
@@ -204,34 +223,34 @@ admin_lev2$dtas_percent<- prop.table(admin_lev2$dtas_count)
 admin_lev2$permaculture_count <- poly.counts(permaculture,admin_lev2)
 admin_lev2$permaculture_percent<- prop.table(admin_lev2$permaculture_count)
 
-# run totals for intermediate boundaries level 3
-# This code will generate a table of frequencies for each spatialpointsdataframe in admin_lev3
-# calculate count of ECS for fields in admin_lev3 and provide percentages
-admin_lev3$ecs_count <- poly.counts(ecs,admin_lev3)
-admin_lev3$ecs_percent<- prop.table(admin_lev3$ecs_count)
-# calculate count of places of worship in PointX db for fields in admin_lev3 and provide percentages
-admin_lev3$pow_count <- poly.counts(pow_pointX,admin_lev3)
-admin_lev3$pow_percent<- prop.table(admin_lev3$pow_count)
-# calculate count of Transition for fields in admin_lev3 and provide percentages
-admin_lev3$transition_count <- poly.counts(transition,admin_lev3)
-admin_lev3$transition_percent<- prop.table(admin_lev3$transition_count)
-# calculate count of dtas for fields in admin_lev3 and provide percentages
-admin_lev3$dtas_count <- poly.counts(dtas,admin_lev3)
-admin_lev3$dtas_percent<- prop.table(admin_lev3$dtas_count)
-# calculate count of permaculture for fields in admin_lev3 and provide percentages
-admin_lev3$permaculture_count <- poly.counts(permaculture,admin_lev3)
-admin_lev3$permaculture_percent<- prop.table(admin_lev3$permaculture_count)
+## run totals for intermediate boundaries level 3
+## This code will generate a table of frequencies for each spatialpointsdataframe in admin_lev3
+## calculate count of ECS for fields in admin_lev3 and provide percentages
+#admin_lev3$ecs_count <- poly.counts(ecs,admin_lev3)
+#admin_lev3$ecs_percent<- prop.table(admin_lev3$ecs_count)
+## calculate count of places of worship in PointX db for fields in admin_lev3 and provide percentages
+#admin_lev3$pow_count <- poly.counts(pow_pointX,admin_lev3)
+#admin_lev3$pow_percent<- prop.table(admin_lev3$pow_count)
+## calculate count of Transition for fields in admin_lev3 and provide percentages
+#admin_lev3$transition_count <- poly.counts(transition,admin_lev3)
+#admin_lev3$transition_percent<- prop.table(admin_lev3$transition_count)
+## calculate count of dtas for fields in admin_lev3 and provide percentages
+#admin_lev3$dtas_count <- poly.counts(dtas,admin_lev3)
+#admin_lev3$dtas_percent<- prop.table(admin_lev3$dtas_count)
+## calculate count of permaculture for fields in admin_lev3 and provide percentages
+#admin_lev3$permaculture_count <- poly.counts(permaculture,admin_lev3)
+#admin_lev3$permaculture_percent<- prop.table(admin_lev3$permaculture_count)
 
 # Import csv with population data for each level of administrative subdivision and join to spatialdataframe
-
 # Add working directory for parish level data (not sure why I've done this)
 
-if (dir.exists("parishes") == FALSE) {
-  dir.create("parishes")
-}
-download.file("http://census.edina.ac.uk/ukborders/easy_download/prebuilt/shape/England_cp_1991.zip", destfile = "parishes/parishes-1991.zip")
-unzip("parishes/parishes-1991.zip", exdir = "parishes")
-parishes <- rgdal::readOGR(dsn = "parishes", "england_cp_1991")
+# JK Note: not using parish data below - have temporarily removed the following.
+# if (dir.exists("parishes") == FALSE) {
+#  dir.create("parishes")
+# }
+# download.file("http://census.edina.ac.uk/ukborders/easy_download/prebuilt/shape/England_cp_1991.zip", destfile = "parishes/parishes-1991.zip")
+# unzip("parishes/parishes-1991.zip", exdir = "parishes")
+# parishes <- rgdal::readOGR(dsn = "parishes", "england_cp_1991")
 
 # Load population statistics for normalising data by population on admin_lev1
 # TODO - consider adapting to use ONS mid-year statistics: https://www.ons.gov.uk/peoplepopulationandcommunity/populationandmigration/populationestimates/datasets/populationestimatesforukenglandandwalesscotlandandnorthernireland
@@ -249,11 +268,11 @@ admin_joined_lev2 <- merge(x=admin_lev2, y=admin_lev2_pop, by.x = "code", by.y =
 # admin_joined_lev3$Postcode.Population <- as.numeric(gsub("[,-]","",admin_joined_lev3$Postcode.Population))
 ```
 
-Perhaps the first important question to ask of these groups is, where are they? I calculated the spread of eco-congregations and transition groups across each of the 32 council areas in Scotland. Every council area in Scotland has at least one eco-congregation or transition group). The most are located in `as.character(admin_lev1$NAME_2[which.max(admin_lev1$ecs_count)])`, with `max(admin_lev1$ecs_count)`, whereas the mean among all the 32 council areas is `mean(admin_lev1$ecs_count)`, with a median of `median(admin_lev1$ecs_count)`, standard deviation of `sd(admin_lev1$ecs_count)`, and interquartile range of `IQR(admin_lev1$ecs_count)`. The following choropleth maps show the relative concentration of eco-congregations (indicated by yellow to red). 
+Perhaps the first important question to ask of these groups is, where are they? I calculated the spread of eco-congregations and transition groups across each of the 32 council areas in Scotland. Every council area in Scotland has at least one eco-congregation or transition group). The most are located in `r as.character(admin_lev1$NAME_2[which.max(admin_lev1$ecs_count)])`, with `r max(admin_lev1$ecs_count)`, whereas the mean among all the 32 council areas is `r mean(admin_lev1$ecs_count)`, with a median of `r median(admin_lev1$ecs_count)`, standard deviation of `r sd(admin_lev1$ecs_count)`, and interquartile range of `r IQR(admin_lev1$ecs_count)`. The following choropleth maps show the relative concentration of eco-congregations (indicated by yellow to red). 
 
 (*TODO: need to implement*) Though there are too few eco-congregations and transition groups for a numerically significant representation in any of the intermediate geographies, mapping the concentration of sites by agricultural parishes allows for a more granular visual and I include this for comparison sake. Note, for the sake of a more accurate visual communication, we have also marked out areas of Scotland that are uninhabited with hash marks on the map of agricultural parishes. (*TODO: this will be done in the final draft, once I get my image masking fixed!*).[^15571030]
 
-```{r create_admin_ecs_choropleth}
+```{r admin_ecs_choropleth}
 # TODO: 
 # 1. Need to rewrite create_admin_ecs_choropleth to use ggplot
 # 2. Need to augment plots to display in 2x2 plots (side-by-side comparison 
@@ -261,119 +280,17 @@ Perhaps the first important question to ask of these groups is, where are they?
 # Row 1 plot using polygons from admin_lev1 and row 2 plot using ploygons from admin_lev2
 # 3. Need to clip choropleth polygons to buildings shapefile
 
-savePlot <- function(myPlot) {
-        pdf(file="figures/admin_choropleth_ecs.pdf")
-        print(myPlot)
-        dev.off()
-}
+plot(admin_lev1)
 
 myplot <- ggplot() + geom_sf(data = admin_lev1_sf) + geom_point(data=as.data.frame(ecs), aes(x=X, y=Y))
 #  + geom_point(data=as.data.frame(ecs), aes(x=X, y=Y))
 
-savePlot(myplot)
-
-
-
-
-var01 <- admin_lev1$ecs_count
-bins <- unique(quantile(var01, seq(0,1,length.out=30)))
-admin_lev1$binId01 <- findInterval(var01, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin, col=colSet01[admin_lev1$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-# use the following to add points to the map (with transparency)
-# points(ecs, pch='.', col=rgb(0,0,0,alpha=0.15))
-# use the following to add region names to the map
-# regionnames <- admin_lev1$NAME_2
-# title(main="Figure 1", sub="Eco-Congregation Scotland\nconcentrations by administrative region", cex.main=0.75, cex.sub=0.5)
-# save to file
-# dev.copy(png,'figures/admin_choropleth_1.png')
-dev.off()
-
-# run choropleth for lev2
-pdf(file="figures/admin_lev2_choropleth_ecs.pdf", width=4, height=4)
-var01 <- admin_lev2$ecs_count
-bins <- unique(quantile(var01, seq(0,1,length.out=30)))
-admin_lev2$binId01 <- findInterval(var01, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin_lev2, col=colSet01[admin_lev2$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-dev.off()
-
-# run choropleth1 for lev3
-
-pdf(file="figures/admin_lev3_choropleth_ecs.pdf", width=4, height=4)
-var01 <- admin_lev3$ecs_count
-bins <- unique(quantile(var01, seq(0,1,length.out=30)))
-admin_lev3$binId01 <- findInterval(var01, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin_lev3, col=colSet01[admin_lev3$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-dev.off()
-
-# New section for population normalised choropleth maps ------------------
-
-# Run plot for level 1
-pdf(file="figures/admin_choropleth_ecs_popnorm.pdf", width=4, height=4)
-var02 <- admin_joined$ecs_count / admin_joined$X2011_pop
-bins <- unique(quantile(var02, seq(0,1,length.out=30)))
-admin_joined$binId01 <- findInterval(var02, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin_joined, col=colSet01[admin_joined$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-dev.off()
-
-# Run plot for level 2
-pdf(file="figures/admin_lev2_choropleth_ecs_popnorm.pdf", width=4, height=4)
-var02 <- admin_joined_lev2$ecs_count / admin_joined_lev2$Postcode.Population
-bins <- unique(quantile(var02, seq(0,1,length.out=30)))
-admin_joined_lev2$binId01 <- findInterval(var02, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin_joined_lev2, col=colSet01[admin_joined_lev2$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-dev.off()
-
-# Run plot for level 3
-pdf(file="figures/admin_lev3_choropleth_ecs_popnorm.pdf", width=4, height=4)
-var02 <- admin_joined_lev3$ecs_count / admin_joined_lev3$Postcode.Population
-bins <- unique(quantile(var02, seq(0,1,length.out=30)))
-admin_joined_lev3$binId01 <- findInterval(var02, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin_joined_lev3, col=colSet01[admin_joined_lev3$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-dev.off()
-
-# New section for overall places of worship count normalised choropleth maps
-
-pdf(file="figures/admin_choropleth_ecs_pownorm.pdf", width=4, height=4)
-var03 <- admin_lev1$ecs_count / admin_lev1$pow_count
-bins <- unique(quantile(var03, seq(0,1,length.out=30)))
-admin_lev1$binId01 <- findInterval(var03, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin, col=colSet01[admin_lev1$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-dev.off()
-
-pdf(file="figures/admin_lev2_choropleth_ecs_pownorm.pdf", width=4, height=4)
-var03 <- admin_lev2$ecs_count / admin_lev2$pow_count
-bins <- unique(quantile(var03, seq(0,1,length.out=30)))
-admin_lev2$binId01 <- findInterval(var03, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin_lev2, col=colSet01[admin_lev2$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-dev.off()
-
-pdf(file="figures/admin_lev3_choropleth_ecs_pownorm.pdf", width=4, height=4)
-var03 <- admin_lev3$ecs_count / admin_lev3$pow_count
-bins <- unique(quantile(var03, seq(0,1,length.out=30), 'na.rm'=TRUE))
-admin_lev3$binId01 <- findInterval(var03, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(admin_lev3, col=colSet01[admin_lev3$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-dev.off()
+ggsave("figures/admin_choropleth_ecs.pdf")
 ```
 
-Given the way population and places of worship are unevenly distributed across Scotland it is important to represent data in terms of relative distribution. For this study, we attempted to "normalise" our data in two different ways, (1) as shown by Figure 2 above, by taking population figures from the 2011 census (see data sheet in Appendix A) and (2) by adjusting relative to the number of places of worship in each council region.[^15914204] The latter of these two can yield particularly unexpected results. Thus, of the `length(pow_pointX)` "places of worship" in Scotland, the highest concentration is actually the `as.character(admin_lev1$NAME_2[which.max(admin_lev1$pow_count)])` region, with `max(admin_lev1$pow_count)`, second is `max( admin_lev1$pow_count[admin_lev1$pow_count!=max(admin_lev1$pow_count)] )` (`as.character(admin_lev1$NAME_2[which.max( admin_lev1$pow_count[admin_lev1$pow_count!=max(admin_lev1$pow_count)])] )`). Rank of Council Areas by population and number of places of worship is also included in Appendix A.
+![][Figure2]
+
+Given the way population and places of worship are unevenly distributed across Scotland it is important to represent data in terms of relative distribution. For this study, we attempted to "normalise" our data in two different ways, (1) as shown by Figure 2 above, by taking population figures from the 2011 census (see data sheet in Appendix A) and (2) by adjusting relative to the number of places of worship in each council region.[^15914204] The latter of these two can yield particularly unexpected results. Thus, of the `r length(pow_pointX)` "places of worship" in Scotland, the highest concentration is actually the `r as.character(admin_lev1$NAME_2[which.max(admin_lev1$pow_count)])` region, with `r max(admin_lev1$pow_count)`, second is `r max( admin_lev1$pow_count[admin_lev1$pow_count!=max(admin_lev1$pow_count)] )` (`r as.character(admin_lev1$NAME_2[which.max( admin_lev1$pow_count[admin_lev1$pow_count!=max(admin_lev1$pow_count)])] )`). Rank of Council Areas by population and number of places of worship is also included in Appendix A.
 
 ```{r create_admin_proportions}
 # Calculate factors by which ECS representation exceeds rep by population and total pow counts
@@ -384,17 +301,28 @@ admin_lev1$ecs_dtas_factor <- ((admin_lev1$ecs_percent - admin_lev1$dtas_percent
 
 ```
 
-We can use this data to normalise our figures regarding Eco-Congregation Scotland communities and this draws the presence in Edinburgh of ECS communities into even sharper relief, as Edinburgh, though ranked second in terms of population and fifth in terms of places of worship, ranks first for the presence of all ECS congregations and awarded ECS congregations. However, taking population as the basis for normalisation first, we find that Edinburgh is far from the most prominent outlier. In trying to communicate this difference for a lay-audience, we have chosen to list this difference as a multiplier (i.e. there are 2.x times as many congregations as their share of population and an average figure of congregations might allow for) as this conveys the difference in a straight-forward way. Outliers where the disparity between their relative share of the total ECS footprint and their relative share of population is different by a positive ratio of more than double include the Orkney Islands (3.7 times more eco-congregations than their expected average share based on population), Argyll and Bute (`admin_lev1[CODE=S12000023]$ecs_pop_factor)4.2x), Stirling (2.76x), and Perthshire and Kinross (2.18x). Interestingly, there are no outliers whose relative share of the total footprint of ECS is double or more in the negative direction (see Appendix A chart for full numbers). 
+We can use this data to normalise our figures regarding Eco-Congregation Scotland communities and this draws the presence in Edinburgh of ECS communities into even sharper relief, as Edinburgh, though ranked second in terms of population and fifth in terms of places of worship, ranks first for the presence of all ECS congregations and awarded ECS congregations. However, taking population as the basis for normalisation first, we find that Edinburgh is far from the most prominent outlier. In trying to communicate this difference for a lay-audience, we have chosen to list this difference as a multiplier (i.e. there are 2.x times as many congregations as their share of population and an average figure of congregations might allow for) as this conveys the difference in a straight-forward way. Outliers where the disparity between their relative share of the total ECS footprint and their relative share of population is different by a positive ratio of more than double include the Orkney Islands (3.7 times more eco-congregations than their expected average share based on population), Argyll and Bute (`r admin_lev1[CODE=S12000023]$ecs_pop_factor)4.2x), Stirling (2.76x), and Perthshire and Kinross (2.18x). Interestingly, there are no outliers whose relative share of the total footprint of ECS is double or more in the negative direction (see Appendix A chart for full numbers). 
 
-Turning to the total of `length(pow_pointX)` "places of worship" in Scotland, we find a slightly different picture of the relative concentration of Eco-Congregations in Scotland. In this case, the outliers are 
+Turning to the total of `r length(pow_pointX)` "places of worship" in Scotland, we find a slightly different picture of the relative concentration of Eco-Congregations in Scotland. In this case, the outliers are 
 
 Whereas our initial measurements indicated a prominent lead for Edinburgh, by normalising our data in this way we can highlight the stronger-than-expected presence of several others that might otherwise escape notice because they lie in a region with significantly lower population or numerically less places of worship. Taking the PointX data on "places of worship" in Scotland, we find a less dramatic picture, but also a slightly different one. The positive outliers include East Renfrewshire (3.4x) Edinburgh (2.9x), Stirling (2.2), West Lothian (1.9x) and Aberdeen (1.5x). Again, negative outliers are far less dramatic, with only Midlothian possessing a ratio of more than 100% negative difference from the number of "places of worship" at 1.5x *fewer*.
 
 ```{r create_admin_barplot}
 # comvert admin back to dataframe for analysis
 admin.df<-data.frame(admin_lev1)
-# Note to self: still doesn't work! Need to fix this bar plot
-qplot(admin_lev1$ecs_count, admin_lev1$transition_count data=admin, geom = "bar")
+# Goal here is to generate a grouped bar plot; https://www.r-graph-gallery.com/48-grouped-barplot-with-ggplot2/
+# Need to flatten admin_lev1 based on all the count columns and generate using ggplot
+
+admin_lev1_flattened <- melt(admin_lev1, id.vars='name')
+
+as_tibble(admin_lev1) %>%
+  rownames_to_column(var = "name") %>%
+  gather(ecs_count, transition_count, -name)
+
+ggplot(admin_lev1, aes(fill=condition, y=value, x=specie)) + 
+    geom_bar(position="dodge", stat="identity")
+ 
+
 ```
 
 We can compare the representation in these various regions against our comparison groups to see how other community-based organisations cluster in Scottish administrative districts. Here there are some significant contrasts. Scottish Community Development trusts are most intensely concentrated in the Highlands and Argyll & Bute. But, this is consistent with all the other categories, Eco-Congregations, Places of Worship, and dtass are all over-represented in this area, varying only by the degree. Edinburgh is different, here we find that Eco-Congregations and Transition projects are over-represented, while dtass are under-represented. Finally, the highlands are another strong contrast, here we find a very strong over-representation by transition towns and dtass while the representation of Eco-Congregations is relatively close to the population share for that area. The two areas of greatest contrast for Eco-Congregations from the other groups are unsurprising, Edinburgh is the location of the ECS offices, while Stirling is the area in which ECS first began (see Appendix B for full data). 
@@ -405,259 +333,6 @@ We can compare the representation in these various regions against our compariso
 
 [Figure 4, normalised by PointX data; also including agricultural parishes etc. as above]
 
-# Christian Denominations #
-
-Eco-Congregation Scotland describes itself as an "ecumenical movement helping local groups of Christians link environmental issues to their faith, reduce their environmental impact and engage with their local community." There are several ties to the Church of Scotland, as the denomination provides office space to Eco-Congregation Scotland in the Church of Scotland complex at 121 George Street in Edinburgh and provides funding for one full-time member of staff. In spite of this, ECS has, from the start, attempted to emphasise its ecumenical aspirations and this is reflected in a wide variety of ways. The name "eco-congregation" is meant to be tradition neutral (in interviews, staff noted how they have sought to avoid names such as "eco-kirk" which would be the more obvious Presbyterian title, or "eco-community" or "eco-church" which might indicate allegiance towards another). Further, the group has a environmental chaplain on their staff whose position is funded by the United Reformed Church, and other members of staff are funded by the Scottish government, and as such, carry no formal affiliation with a religious institution. This diversity and ecumenicism is reflected in a membership which is, though dominated by the Church of Scotland, nevertheless, made up of a range of Christian traditions. 
-
-Though these are not numerically significant, it is important to note that some member congregations describe themselves as ecumenical communities, and others are hybrids reflecting the merging of two traditions. As this ecumenical/hybrid designation involves a small number of the overall total, for the sake of this research, these have been combined into a category called "ecumenical." Further, as research conducted by Church of Scotland statistician Fiona Tweedie has shown, in many Scottish communities with only one church, members of this church will specify their denominational affiliation in a variety of ways (Roman Catholic, Quaker, Methodist, etc.) even though the church and its minister are formally affiliated with the Church of Scotland.[^159142242] So, we should be careful not to assume that the various denominational affiliations of eco-congregations are indicative in an absolute way.
-
-A wide variety of historians and sociologists of religion have noted the regional significance of different Christian denominations in Scotland so we sought to assess the relative distribution and concentration of eco-congregations by denomination. Finding comparative statistics is a complex task, made more complicated by several factors. First, most demographic data on religious belonging in Scotland comes in the form of the 2011 census and as such is far more atomised than this data-set which identifies groups at the level of "congregations" rather than individuals. Equating these two is also complex, as participation by members of congregations can be measured in a variety of ways, there are often a small number of active participants in each eco-congregation group, but may also be a large scale, but passive, support by the wider community. 
-
-So why provide this kind of data (i.e. at the level of individual churches) when more granular data (i.e. at the level of individuals persons) is available in the form of the census and related parallel publications such as the 2008 Scottish Environmental Attitudes survey? We believe that mapping places of worship provides a useful intermediate level of analysis and may complement our more atomised understanding of EA which has been assessed at the level of individual persons to date. Because representation within some administrative areas of Scotland, can lead to a small number of data points, we have kept analysis to a National level and have not provided more specific administrative-area level calculations.
-
-```{r create_ecs_denomination_table}
-# TODO: Need to find a prettier way to do this:
-table(ecs@data[["denomination"]])
-# TODO: Add dataframe with overall church numbers for the UK from "ScottishChurches" dataset by JK
-```
-
-![][Figure5]
-
-As one might expect, there is a strong representation of the Church of Scotland, almost 74% of eco-congregations, with this number remaining the same when we only count awarded sites. We can confirm, on the basis of this analysis that ECS has a disproportional representation by Church of Scotland churches. At the 2002 church census count, it only represented 40.20% of Scottish churches (1666 of 4144 total churches). Similarly, on the 2011 Scottish census, only 32.44% of persons claimed to be members of the Church of Scotland. We can adjust this representation to 60%, if one excludes the 2,445,204 persons (46% of the total on the census) who reported either "no religion" or adherence to a religious tradition not currently represented among the eco-congregation sites. There is a slight over-representation by the United Reformed church, though this seems considerably more dramatic when one takes into account the fact that this is a trebling or more of their overall share of Scottish churches. The URC makes up only sightly more than 1% of church buildings in Scotland and a tiny 0.04% of respondents to the 2011 census. The Scottish Episcopal church hovers right around a proportional representation within ECS. More concerning are the significant underrepresentation by Roman Catholic churches, Baptists, the Free Church of Scotland, and other independent churches. 
-
-While Roman Catholic churches make up just over 10% of the church buildings in Scotland, less than 5% of churches registered as eco-congregations are RC. Even more dramatic is the quartering of baptist churches, and the non-existent representation among the significant group of independent churches and small denominations. These make up nearly 25% of all Scottish churches (over a thousand) and yet only 4 have registered as eco-congregations. We provide several tentative advisories in response to these under-representations in the final section of this paper.
-
-# Eco-Congregations, Urban, Rural and Remote
-
-```{r ur8fold}
-# read in relevant polygons for UR8fold scale
-download.file("http://sedsh127.sedsh.gov.uk/Atom_data/ScotGov/ZippedShapefiles/SG_UrbanRural_2016.zip", 
-              destfile = "data/SG_UrbanRural_2016.zip")
-unzip("data/SG_UrbanRural_2016.zip", exdir = "data")
-urbanrural <- readOGR("./data", "SG_UrbanRural_2016") 
-
-# TODO: worth considering uploading data to zenodo for long-term reproducibility as ScotGov shuffles this stuff around periodically breaking URLs
-
-# This code will generate a table of frequencies for each spatialpointsdataframe in urbanrural
-# calculate count of ECS for fields in urbanrural
-urbanrural$ecs_count <- poly.counts(ecs,urbanrural)
-urbanrural$ecs_percent<- prop.table(urbanrural$ecs_count)
-# calculate count of places of worship in PointX db for fields in urbanrural and provide percentages
-urbanrural$pow_count <- poly.counts(pow_pointX,urbanrural)
-urbanrural$pow_percent<- prop.table(urbanrural$pow_count)
-# calculate count of Transition for fields in urbanrural
-urbanrural$transition_count <- poly.counts(transition,urbanrural)
-urbanrural$transition_percent<- prop.table(urbanrural$transition_count)
-# calculate count of dtas for fields in urbanrural
-urbanrural$dtas_count <- poly.counts(dtas,urbanrural)
-urbanrural$dtas_percent<- prop.table(urbanrural$dtas_count)
-# calculate count of permaculture for fields in urbanrural
-urbanrural$permaculture_count <- poly.counts(permaculture,urbanrural)
-urbanrural$permaculture_percent<- prop.table(urbanrural$permaculture_count)
-
-# Create dataframe for analysis
-urbanrural_df<-data.frame(urbanrural)
-```
-
-Rather than bifurcate congregations into an urban/rural dichotomy, for this study we used the Scottish Government's eight-point remoteness scale to categorise eco-congregations along a spectrum of highly populated to remote areas. This 8-fold scale (calculated biennially) offers a more nuanced measurement that combines measurements of remoteness and population along the following lines:
-
-1. Large Urban Areas - Settlements of over 125,000 people.
-2. Other Urban Areas - Settlements of 10,000 to 125,000 people. 
-3. Accessible Small Towns - Settlements of between 3,000 and 10,000 people, and within a 30 minute drive time of a Settlement of 10,000 or more. 
-4. Remote Small Towns - Settlements of between 3,000 and 10,000 people, and with a drive time between 30 and 60 minutes to a Settlement of 10,000 or more.
-5. Very Remote Small Towns - Settlements of between 3,000 and 10,000 people, and with a drive time of over 60 minutes to a Settlement of 10,000 or more.
-6. Accessible Rural Areas - Areas with a population of less than 3,000 people, and within a drive time of 30 minutes to a Settlement of 10,000 or more.
-7. Remote Rural Areas - Areas with a population of less than 3,000 people, and with a drive time of between 30 and 60 minutes to a Settlement of 10,000 or more.
-8. Very Remote Rural Areas - Areas with a population of less than 3,000 people, and with a drive time of over 60 minutes to a Settlement of 10,000 or more.
-
-The key question which this analysis seeks to answer is whether ECS, or the other groups surveyed, are more concentrated in Urban or Rural areas, so as is the case below with our analysis of deprivation, we are concerned with the outer conditions, i.e. the urban areas (items 1-2) and remote areas (items 7-8).
-
-Of all the groups surveyed in this study, Eco-Congregation Scotland is the most heavily concentrated in large urban areas (33.53%), exceeding by almost 50% the rate for all places of worship (22.96% in large urban areas). Transition is a much more modest 20% and development trusts a bit lower at 15%. It is interesting to note that the rate of ECS concentration in these large urban areas matches the level of overall population distribution (34.5%). On the other end of the scale, Eco-Congregation Scotland is the least concentrated in remote rural areas (with 3.93% on level 7 and 5.44% on level 8 on the urban-rural scale), though again, they correlate roughly to the general population distribution (3.2% and 2.9% respectively). Places of worship outpace both the population of Scotland and the footprint of Eco-Congregation Scotland, with 14.98% in very remote rural areas, but this is exceeded by transition at 16.47% and both by Scottish community development trusts at 32.14%. So while Eco-Congregation Scotland correlates roughly with Scottish population distribution across the urban-rural scale, it has a considerably more urban profile than either of the other two groups surveyed.
-
-```{r create_urbanrural_ecs_chart_choropleth}
-# Create bar chart
-# TODO: use ggplot here to generate stacked bar chart showing representation across 8 urban/rural bands from urbanrural_df
-# TODO: generate map based on urbanrural polygons, shade based on 8-fold categories, draw transparent dots onto polygons
-
-pdf(file="figures/urbanrural_choropleth_ecs.pdf", width=4, height=4)
-bins <- unique(quantile(urbanrural$ecs_count, seq(0,1,length.out=30)))
-urbanrural$binId01 <- findInterval(urbanrural$ecs_count, bins)
-colSet01 <- rev(heat.colors(length(bins)))
-plot(urbanrural, col=colSet01[urbanrural$binId01], border="grey", lwd=0.25)
-par(mar=c(5,3,2,2)+0.1)
-# use the following to add points to the map (with transparency)
-# points(ecs, pch='.', col=rgb(0,0,0,alpha=0.15))
-title(main="Figure 5", sub="Eco-Congregation Scotland\nconcentrations in Urban Rural 8-fold classifications", cex.main=0.75, cex.sub=0.5)
-# save to file
-# dev.copy(png,'figures/urbanrural_choropleth_1.png')
-dev.off()
-```
-
-![][Figure6 map of ECS concentration by UR Areas Urban - Rural via Brown to Green]
-
-
-# Wealth, Employment, and Literacy
-
-```{r simd}
-# read in relevant polygons, Scottish Index of Multiple deprivation
-download.file("http://simd.scot/2016/data/simd2016_withgeog.zip", 
-              destfile = "data/simd2016_withgeog.zip")
-unzip("data/simd2016_withgeog.zip", exdir = "data", junkpaths = TRUE)
-simd_shapes <- readOGR("./data", "sc_dz_11") 
-simd_indicators <- read.csv("./data/simd2016_withinds.csv")
-simd <- merge(x=simd_shapes, y=simd_indicators, by="Data_Zone")
-  
-# commenting out old dataset in light of reproducible version above
-# simd <- readOGR("data", "simd_04-12_all_data")
-
-# STAGE 1, augment each dataset with relevant (geolocated) columns from SIMD
-# examine which ecs fall within each SIMD classification
-cbind(ecs@data, over(ecs, simd))
-# assign combined table with SIMD columns to attribute table slot of ecs table
-ecs@data=cbind(ecs@data,over(ecs,simd))
-# examine where pointX falls within each SIMD classification
-cbind(pow_pointX@data, over(pow_pointX, simd))
-# assign combined table with SIMD columns to attribute table slot of ecs table
-pow_pointX@data=cbind(pow_pointX@data,over(pow_pointX,simd))
-
-# examine which transition fall within each SIMD classifications
-cbind(transition@data, over(transition, simd))
-# assign combined table with SIMD columns to attribute table slot of transition table
-transition@data=cbind(transition@data,over(transition,simd))
-# examine which permaculture fall within each SIMD classifications
-cbind(permaculture@data, over(permaculture, simd))
-# assign combined table with SIMD columns to attribute table slot of permaculture table
-permaculture@data=cbind(permaculture@data,over(permaculture,simd))
-# examine which dtas fall within each SIMD classifications
-cbind(dtas@data, over(dtas, simd))
-# assign combined table with SIMD columns to attribute table slot of dtas table
-dtas@data=cbind(dtas@data,over(dtas,simd))
-
-# STAGE 2, extract NULL cells from each data set to prevent errors in stage 3
-
-# convert back to data frame for null cell extraction
-ecs<-data.frame(ecs)
-# split out null and normal cells
-ecs_clean<-ecs[complete.cases(ecs),]
-ecs_null<-ecs[!complete.cases(ecs),]
-# convert back to spatialpointdataframe
-coordinates(ecs) <- c("X", "Y")
-proj4string(ecs) <- proj4string(admin_lev1)
-
-# convert back to data frame for null cell extraction
-transition<-data.frame(transition)
-# split out null and normal cells
-transition_clean<-transition[complete.cases(transition),]
-transition_null<-transition[!complete.cases(transition),]
-# convert back to spatialpointdataframe
-coordinates(transition) <- c("X", "Y")
-proj4string(transition) <- proj4string(admin_lev1)
-
-# convert back to data frame for null cell extraction
-permaculture<-data.frame(permaculture)
-# split out null and normal cells
-permaculture_clean<-permaculture[complete.cases(permaculture),]
-permaculture_null<-permaculture[!complete.cases(permaculture),]
-# convert back to spatialpointdataframe
-coordinates(permaculture) <- c("X", "Y")
-proj4string(permaculture) <- proj4string(admin_lev1)
-
-# convert back to data frame for null cell extraction
-dtas<-data.frame(dtas)
-# split out null and normal cells
-dtas_clean<-dtas[complete.cases(dtas),]
-dtas_null<-dtas[!complete.cases(dtas),]
-# convert back to spatialpointdataframe
-coordinates(dtas) <- c("X", "Y")
-proj4string(dtas) <- proj4string(admin_lev1)
-
-# STAGE 3a, calculate sums based on SIMD12R columns and generate new integer sets with quintile count data
-simd_rownames = c("Quintile 1","Quintile 2","Quintile 3","Quintile 4","Quintile 5")
-simdr12_ecs = c((sum(ecs_clean$SIMDR12<1301)), (sum(ecs_clean$SIMDR12 > 1300 & ecs_clean$SIMDR12 < 2602)), (sum(ecs_clean$SIMDR12 > 2601 & ecs_clean$SIMDR12 < 3903)), (sum(ecs_clean$SIMDR12 > 3902 & ecs_clean$SIMDR12 < 5204)), (sum(ecs_clean$SIMDR12 > 5203 & ecs_clean$SIMDR12 < 6505)))
-# names(simdr12_ecs) <- simd_rownames
-simdr12_transition = c((sum(transition_clean$SIMDR12<1301)), (sum(transition_clean$SIMDR12 > 1300 & transition_clean$SIMDR12 < 2602)), (sum(transition_clean$SIMDR12 > 2601 & transition_clean$SIMDR12 < 3903)), (sum(transition_clean$SIMDR12 > 3902 & transition_clean$SIMDR12 < 5204)), (sum(transition_clean$SIMDR12 > 5203 & transition_clean$SIMDR12 < 6505)))
-# names(simdr12_transition) <- simd_rownames
-simdr12_permaculture = c((sum(permaculture_clean$SIMDR12<1301)), (sum(permaculture_clean$SIMDR12 > 1300 & permaculture_clean$SIMDR12 < 2602)), (sum(permaculture_clean$SIMDR12 > 2601 & permaculture_clean$SIMDR12 < 3903)), (sum(permaculture_clean$SIMDR12 > 3902 & permaculture_clean$SIMDR12 < 5204)), (sum(permaculture_clean$SIMDR12 > 5203 & permaculture_clean$SIMDR12 < 6505)))
-# names(simdr12_permaculture) <- simd_rownames
-simdr12_dtas = c((sum(dtas_clean$SIMDR12<1301)), (sum(dtas_clean$SIMDR12 > 1300 & dtas_clean$SIMDR12 < 2602)), (sum(dtas_clean$SIMDR12 > 2601 & dtas_clean$SIMDR12 < 3903)), (sum(dtas_clean$SIMDR12 > 3902 & dtas_clean$SIMDR12 < 5204)), (sum(dtas_clean$SIMDR12 > 5203 & dtas_clean$SIMDR12 < 6505)))
-# names(simdr12_dtas) <- simd_rownames
-
-# STAGE 3b, calculate sums based on INCR12 columns and generate new integer sets with quintile count data
-incr12_ecs = c((sum(ecs_clean$INCR12<1301)), (sum(ecs_clean$INCR12 > 1300 & ecs_clean$INCR12 < 2602)), (sum(ecs_clean$INCR12 > 2601 & ecs_clean$INCR12 < 3903)), (sum(ecs_clean$INCR12 > 3902 & ecs_clean$INCR12 < 5204)), (sum(ecs_clean$INCR12 > 5203 & ecs_clean$INCR12 < 6505)))
-incr12_transition = c((sum(transition_clean$INCR12<1301)), (sum(transition_clean$INCR12 > 1300 & transition_clean$INCR12 < 2602)), (sum(transition_clean$INCR12 > 2601 & transition_clean$INCR12 < 3903)), (sum(transition_clean$INCR12 > 3902 & transition_clean$INCR12 < 5204)), (sum(transition_clean$INCR12 > 5203 & transition_clean$INCR12 < 6505)))
-incr12_permaculture = c((sum(permaculture_clean$INCR12<1301)), (sum(permaculture_clean$INCR12 > 1300 & permaculture_clean$INCR12 < 2602)), (sum(permaculture_clean$INCR12 > 2601 & permaculture_clean$INCR12 < 3903)), (sum(permaculture_clean$INCR12 > 3902 & permaculture_clean$INCR12 < 5204)), (sum(permaculture_clean$INCR12 > 5203 & permaculture_clean$INCR12 < 6505)))
-incr12_dtas = c((sum(dtas_clean$INCR12<1301)), (sum(dtas_clean$INCR12 > 1300 & dtas_clean$INCR12 < 2602)), (sum(dtas_clean$INCR12 > 2601 & dtas_clean$INCR12 < 3903)), (sum(dtas_clean$INCR12 > 3902 & dtas_clean$INCR12 < 5204)), (sum(dtas_clean$INCR12 > 5203 & dtas_clean$INCR12 < 6505)))
-
-# STAGE 3c, calculate sums based on EMPR12 columns and generate new integer sets with quintile count data
-empr12_ecs = c((sum(ecs_clean$EMPR12<1301)), (sum(ecs_clean$EMPR12 > 1300 & ecs_clean$EMPR12 < 2602)), (sum(ecs_clean$EMPR12 > 2601 & ecs_clean$EMPR12 < 3903)), (sum(ecs_clean$EMPR12 > 3902 & ecs_clean$EMPR12 < 5204)), (sum(ecs_clean$EMPR12 > 5203 & ecs_clean$EMPR12 < 6505)))
-empr12_transition = c((sum(transition_clean$EMPR12<1301)), (sum(transition_clean$EMPR12 > 1300 & transition_clean$EMPR12 < 2602)), (sum(transition_clean$EMPR12 > 2601 & transition_clean$EMPR12 < 3903)), (sum(transition_clean$EMPR12 > 3902 & transition_clean$EMPR12 < 5204)), (sum(transition_clean$EMPR12 > 5203 & transition_clean$EMPR12 < 6505)))
-empr12_permaculture = c((sum(permaculture_clean$EMPR12<1301)), (sum(permaculture_clean$EMPR12 > 1300 & permaculture_clean$EMPR12 < 2602)), (sum(permaculture_clean$EMPR12 > 2601 & permaculture_clean$EMPR12 < 3903)), (sum(permaculture_clean$EMPR12 > 3902 & permaculture_clean$EMPR12 < 5204)), (sum(permaculture_clean$EMPR12 > 5203 & permaculture_clean$EMPR12 < 6505)))
-empr12_dtas = c((sum(dtas_clean$EMPR12<1301)), (sum(dtas_clean$EMPR12 > 1300 & dtas_clean$EMPR12 < 2602)), (sum(dtas_clean$EMPR12 > 2601 & dtas_clean$EMPR12 < 3903)), (sum(dtas_clean$EMPR12 > 3902 & dtas_clean$EMPR12 < 5204)), (sum(dtas_clean$EMPR12 > 5203 & dtas_clean$EMPR12 < 6505)))
-
-# STAGE 3d, calculate sums based on HER12 columns and generate new integer sets with quintile count data
-her12_ecs = c((sum(ecs_clean$HER12<1301)), (sum(ecs_clean$HER12 > 1300 & ecs_clean$HER12 < 2602)), (sum(ecs_clean$HER12 > 2601 & ecs_clean$HER12 < 3903)), (sum(ecs_clean$HER12 > 3902 & ecs_clean$HER12 < 5204)), (sum(ecs_clean$HER12 > 5203 & ecs_clean$HER12 < 6505)))
-her12_transition = c((sum(transition_clean$HER12<1301)), (sum(transition_clean$HER12 > 1300 & transition_clean$HER12 < 2602)), (sum(transition_clean$HER12 > 2601 & transition_clean$HER12 < 3903)), (sum(transition_clean$HER12 > 3902 & transition_clean$HER12 < 5204)), (sum(transition_clean$HER12 > 5203 & transition_clean$HER12 < 6505)))
-her12_permaculture = c((sum(permaculture_clean$HER12<1301)), (sum(permaculture_clean$HER12 > 1300 & permaculture_clean$HER12 < 2602)), (sum(permaculture_clean$HER12 > 2601 & permaculture_clean$HER12 < 3903)), (sum(permaculture_clean$HER12 > 3902 & permaculture_clean$HER12 < 5204)), (sum(permaculture_clean$HER12 > 5203 & permaculture_clean$HER12 < 6505)))
-her12_dtas = c((sum(dtas_clean$HER12<1301)), (sum(dtas_clean$HER12 > 1300 & dtas_clean$HER12 < 2602)), (sum(dtas_clean$HER12 > 2601 & dtas_clean$HER12 < 3903)), (sum(dtas_clean$HER12 > 3902 & dtas_clean$HER12 < 5204)), (sum(dtas_clean$HER12 > 5203 & dtas_clean$HER12 < 6505)))
-
-# STAGE 3e, calculate sums based on ESTR12 columns and generate new integer sets with quintile count data
-estr12_ecs = c((sum(ecs_clean$ESTR12<1301)), (sum(ecs_clean$ESTR12 > 1300 & ecs_clean$ESTR12 < 2602)), (sum(ecs_clean$ESTR12 > 2601 & ecs_clean$ESTR12 < 3903)), (sum(ecs_clean$ESTR12 > 3902 & ecs_clean$ESTR12 < 5204)), (sum(ecs_clean$ESTR12 > 5203 & ecs_clean$ESTR12 < 6505)))
-estr12_transition = c((sum(transition_clean$ESTR12<1301)), (sum(transition_clean$ESTR12 > 1300 & transition_clean$ESTR12 < 2602)), (sum(transition_clean$ESTR12 > 2601 & transition_clean$ESTR12 < 3903)), (sum(transition_clean$ESTR12 > 3902 & transition_clean$ESTR12 < 5204)), (sum(transition_clean$ESTR12 > 5203 & transition_clean$ESTR12 < 6505)))
-estr12_permaculture = c((sum(permaculture_clean$ESTR12<1301)), (sum(permaculture_clean$ESTR12 > 1300 & permaculture_clean$ESTR12 < 2602)), (sum(permaculture_clean$ESTR12 > 2601 & permaculture_clean$ESTR12 < 3903)), (sum(permaculture_clean$ESTR12 > 3902 & permaculture_clean$ESTR12 < 5204)), (sum(permaculture_clean$ESTR12 > 5203 & permaculture_clean$ESTR12 < 6505)))
-estr12_dtas = c((sum(dtas_clean$ESTR12<1301)), (sum(dtas_clean$ESTR12 > 1300 & dtas_clean$ESTR12 < 2602)), (sum(dtas_clean$ESTR12 > 2601 & dtas_clean$ESTR12 < 3903)), (sum(dtas_clean$ESTR12 > 3902 & dtas_clean$ESTR12 < 5204)), (sum(dtas_clean$ESTR12 > 5203 & dtas_clean$ESTR12 < 6505)))
-
-# STAGE 4a - calculate percentages
-simdr12_ecs_percent<- prop.table(simdr12_ecs)
-simdr12_transition_percent<- prop.table(simdr12_transition)
-simdr12_permaculture_percent<- prop.table(simdr12_permaculture)
-simdr12_dtas_percent<- prop.table(simdr12_dtas)
-incr12_ecs_percent<- prop.table(incr12_ecs)
-incr12_transition_percent<- prop.table(incr12_transition)
-incr12_permaculture_percent<- prop.table(incr12_permaculture)
-incr12_dtas_percent<- prop.table(incr12_dtas)
-empr12_ecs_percent<- prop.table(empr12_ecs)
-empr12_transition_percent<- prop.table(empr12_transition)
-empr12_permaculture_percent<- prop.table(empr12_permaculture)
-empr12_dtas_percent<- prop.table(empr12_dtas)
-her12_ecs_percent<- prop.table(her12_ecs)
-her12_transition_percent<- prop.table(her12_transition)
-her12_permaculture_percent<- prop.table(her12_permaculture)
-her12_dtas_percent<- prop.table(her12_dtas)
-estr12_ecs_percent<- prop.table(estr12_ecs)
-estr12_transition_percent<- prop.table(estr12_transition)
-estr12_permaculture_percent<- prop.table(estr12_permaculture)
-estr12_dtas_percent<- prop.table(estr12_dtas)
-
-# STAGE 4b, generate data frame using integer sets
-simd = data.frame(simdr12_ecs, simdr12_ecs_percent, incr12_ecs, incr12_ecs_percent, empr12_ecs, empr12_ecs_percent, her12_ecs, her12_ecs_percent, estr12_ecs, estr12_ecs_percent, simdr12_transition, simdr12_transition_percent, incr12_transition, incr12_transition_percent, empr12_transition, empr12_transition_percent, her12_transition, her12_transition_percent, estr12_transition, estr12_transition_percent, simdr12_permaculture, simdr12_permaculture_percent, incr12_permaculture, incr12_permaculture_percent, empr12_permaculture, empr12_permaculture_percent, her12_permaculture, her12_permaculture_percent, estr12_permaculture, estr12_permaculture_percent, simdr12_dtas, simdr12_dtas_percent, incr12_dtas, incr12_dtas_percent, empr12_dtas, empr12_dtas_percent, her12_dtas, her12_dtas_percent, estr12_dtas, estr12_dtas_percent)
-write.csv(simd, "derivedData/simd.csv", row.names=FALSE)
-
-
-simd_percents_only = data.frame(simd_rownames, simdr12_ecs_percent, incr12_ecs_percent, empr12_ecs_percent, her12_ecs_percent, estr12_ecs_percent, simdr12_transition_percent, incr12_transition_percent, empr12_transition_percent, her12_transition_percent, estr12_transition_percent, simdr12_permaculture_percent, incr12_permaculture_percent, empr12_permaculture_percent, her12_permaculture_percent, estr12_permaculture_percent, simdr12_dtas_percent, incr12_dtas_percent, empr12_dtas_percent, her12_dtas_percent, estr12_dtas_percent)
-write.csv(simd_percents_only, "derivedData/simd_percents_only.csv", row.names=FALSE)
-
-
-# STAGE 5, generate cool charts
-
-# clustered bar charts
-# convert to long format
-library(reshape2)
-simd_percents_only_long <- melt(simd_percents_only, id.vars = "simd_rownames", 
-              measure.vars = grep("^12", names(simd_percents_only), value = TRUE))
-
-qplot(data=simd_percents_only_long , geom="bar", fill=(factor(simd_rownames)))
-
-# jitterplot option, from Teutonico 2015, p. 63
-
-```
-
-Another crucial point of assessment relates to the relation of Eco-Congregation communities to the Scottish Index of Multiple Deprivation. This instrument aggregates a large variety of factors which can lead to deprivation including crime rates, employment levels, access to services (implicating remoteness), and literacy. By assessing ECS, Transition, and dtas against the deprivation scale, we can assess whether eco-congregations fall within particular demographics and also whether the fully aggregated SIMD measurement provides a useful point of comparison for our purposes. The SIMD essentially divides Scotland into 6407 geographic zones and then ranks them based on their relative deprivation. This data set can be split into any number of groups, but for our purposes we have settled on Quintiles, splitting the SIMD data set at every 1302 entries. We then measured where each transition group, ECS, and dtas fell within these zones and calculated how they fell into these five quintiles, from more to least deprived.
-
-The first, and most compelling finding is that, in general Eco-Congregation Scotland and Transition Scotland are both roughly the same and match the level of population distribution in the lowest quintile of the general SIMD measurement. 8% of transition groups and eco-congregation groups which have received awards and 9% of the population are located within this quintile. However, taken in relation to the distribution of places of worship in the lowest quintile, we find that eco-congregations are located at half the rate that places of worship are (15%) and dtass match this much more closely at 14%. Turning towards the top quintile, this pattern also holds, here both transition groups (21%) and eco-congregations (21% and 29% of awarded congregations) depart from the population distribution in this upper quintile (which is 10%). Again, general places of worship (at 11%) and DTASs (at 5%) take the opposite direction. We can say decisively that in communities which have been identified as good candidates for intervention to reduce deprivation, ECS and Transition are less likely, and they are over-represented at the areas which fall into the least deprived quintile.
-
-We can find divergence between transition communities and eco-congregation when we split out SIMD domains. In the lowest quartile, measuring exclusively for the income domain, ECS is more represented (11%) - roughly the same as DTAS (12%), and transition is less (6%) represented. In general (as shown on the chart in Appendix D), these trends hold when representation of our groups are measured within other non-remoteness domains of the SIMD. Our basic conclusion is that transition towns are least likely to operate within the lowest quartile of SIMD and DTASs are most likely, with ECS somewhere in the middle. Given the general disparity against the presence of places of worship, it seems fair to suggest that this might be an area for improvement, perhaps even worth developing a special programme which might target areas in SIMD quartile 1 for eco-congregation outreach. This might be considered particularly in light of the starkest underrepresentation of ECS and transition within the SIMD domain of education, skills, and training.
-
-
 # Appendix A
 
 ```{r pander_admin_table}
@@ -665,7 +340,7 @@ admin.shortened <- admin_lev1[,c(7,12:22)]
 # Output CSV files for various levels of admin
 write.csv(admin, "derivedData/admin.csv", row.names=FALSE)
 write.csv(admin_lev2, "derivedData/admin_lev2.csv", row.names=FALSE)
-write.csv(admin_lev3, "derivedData/admin_lev3.csv", row.names=FALSE)
+# write.csv(admin_lev3, "derivedData/admin_lev3.csv", row.names=FALSE)
 write.csv(admin.shortened, "derivedData/admin.csv", row.names=FALSE)
 write.csv(ecs, "derivedData/ecs.csv", row.names=FALSE)
 write.csv(transition, "derivedData/transition.csv", row.names=FALSE)
@@ -679,20 +354,6 @@ panderOptions("digits", 2)
 pander(admin.shortened)
 ```
 
-# Appendix B
-
-(JK note to self: same as above, but augmented with multipliers by which categories are different from one another)
-
-# Appendix C - Data by Urban / Rural Classification
-
-```{r pander_urbanrural_table}
-urbanrural.shortened <- urbanrural[,c(2,6,9:18)]
-write.csv(urbanrural, "derivedData/urbanrural.csv", row.names=FALSE)
-write.csv(urbanrural.shortened, "derivedData/urbanrural.csv", row.names=FALSE)
-urbanrural.shortened<-data.frame(urbanrural.shortened)
-panderOptions("digits", 2)
-pander(urbanrural.shortened)
-```
 
 
 # Citations
diff --git a/mapping_draft_part2.Rmd b/mapping_draft_part2.Rmd
new file mode 100644
index 0000000..d70385b
--- /dev/null
+++ b/mapping_draft_part2.Rmd
@@ -0,0 +1,271 @@
+
+# Christian Denominations #
+
+Eco-Congregation Scotland describes itself as an "ecumenical movement helping local groups of Christians link environmental issues to their faith, reduce their environmental impact and engage with their local community." There are several ties to the Church of Scotland, as the denomination provides office space to Eco-Congregation Scotland in the Church of Scotland complex at 121 George Street in Edinburgh and provides funding for one full-time member of staff. In spite of this, ECS has, from the start, attempted to emphasise its ecumenical aspirations and this is reflected in a wide variety of ways. The name "eco-congregation" is meant to be tradition neutral (in interviews, staff noted how they have sought to avoid names such as "eco-kirk" which would be the more obvious Presbyterian title, or "eco-community" or "eco-church" which might indicate allegiance towards another). Further, the group has a environmental chaplain on their staff whose position is funded by the United Reformed Church, and other members of staff are funded by the Scottish government, and as such, carry no formal affiliation with a religious institution. This diversity and ecumenicism is reflected in a membership which is, though dominated by the Church of Scotland, nevertheless, made up of a range of Christian traditions. 
+
+Though these are not numerically significant, it is important to note that some member congregations describe themselves as ecumenical communities, and others are hybrids reflecting the merging of two traditions. As this ecumenical/hybrid designation involves a small number of the overall total, for the sake of this research, these have been combined into a category called "ecumenical." Further, as research conducted by Church of Scotland statistician Fiona Tweedie has shown, in many Scottish communities with only one church, members of this church will specify their denominational affiliation in a variety of ways (Roman Catholic, Quaker, Methodist, etc.) even though the church and its minister are formally affiliated with the Church of Scotland.[^159142242] So, we should be careful not to assume that the various denominational affiliations of eco-congregations are indicative in an absolute way.
+
+A wide variety of historians and sociologists of religion have noted the regional significance of different Christian denominations in Scotland so we sought to assess the relative distribution and concentration of eco-congregations by denomination. Finding comparative statistics is a complex task, made more complicated by several factors. First, most demographic data on religious belonging in Scotland comes in the form of the 2011 census and as such is far more atomised than this data-set which identifies groups at the level of "congregations" rather than individuals. Equating these two is also complex, as participation by members of congregations can be measured in a variety of ways, there are often a small number of active participants in each eco-congregation group, but may also be a large scale, but passive, support by the wider community. 
+
+So why provide this kind of data (i.e. at the level of individual churches) when more granular data (i.e. at the level of individuals persons) is available in the form of the census and related parallel publications such as the 2008 Scottish Environmental Attitudes survey? We believe that mapping places of worship provides a useful intermediate level of analysis and may complement our more atomised understanding of EA which has been assessed at the level of individual persons to date. Because representation within some administrative areas of Scotland, can lead to a small number of data points, we have kept analysis to a National level and have not provided more specific administrative-area level calculations.
+
+```{r create_ecs_denomination_table}
+# TODO: Need to find a prettier way to do this:
+table(ecs@data[["denomination"]])
+# TODO: Add dataframe with overall church numbers for the UK from "ScottishChurches" dataset by JK
+```
+
+![][Figure5]
+
+As one might expect, there is a strong representation of the Church of Scotland, almost 74% of eco-congregations, with this number remaining the same when we only count awarded sites. We can confirm, on the basis of this analysis that ECS has a disproportional representation by Church of Scotland churches. At the 2002 church census count, it only represented 40.20% of Scottish churches (1666 of 4144 total churches). Similarly, on the 2011 Scottish census, only 32.44% of persons claimed to be members of the Church of Scotland. We can adjust this representation to 60%, if one excludes the 2,445,204 persons (46% of the total on the census) who reported either "no religion" or adherence to a religious tradition not currently represented among the eco-congregation sites. There is a slight over-representation by the United Reformed church, though this seems considerably more dramatic when one takes into account the fact that this is a trebling or more of their overall share of Scottish churches. The URC makes up only sightly more than 1% of church buildings in Scotland and a tiny 0.04% of respondents to the 2011 census. The Scottish Episcopal church hovers right around a proportional representation within ECS. More concerning are the significant underrepresentation by Roman Catholic churches, Baptists, the Free Church of Scotland, and other independent churches. 
+
+While Roman Catholic churches make up just over 10% of the church buildings in Scotland, less than 5% of churches registered as eco-congregations are RC. Even more dramatic is the quartering of baptist churches, and the non-existent representation among the significant group of independent churches and small denominations. These make up nearly 25% of all Scottish churches (over a thousand) and yet only 4 have registered as eco-congregations. We provide several tentative advisories in response to these under-representations in the final section of this paper.
+
+# Eco-Congregations, Urban, Rural and Remote
+
+```{r ur8fold}
+# read in relevant polygons for UR8fold scale
+download.file("http://sedsh127.sedsh.gov.uk/Atom_data/ScotGov/ZippedShapefiles/SG_UrbanRural_2016.zip", 
+              destfile = "data/SG_UrbanRural_2016.zip")
+unzip("data/SG_UrbanRural_2016.zip", exdir = "data")
+urbanrural <- readOGR("./data", "SG_UrbanRural_2016") 
+
+# TODO: worth considering uploading data to zenodo for long-term reproducibility as ScotGov shuffles this stuff around periodically breaking URLs
+
+# This code will generate a table of frequencies for each spatialpointsdataframe in urbanrural
+# calculate count of ECS for fields in urbanrural
+urbanrural$ecs_count <- poly.counts(ecs,urbanrural)
+urbanrural$ecs_percent<- prop.table(urbanrural$ecs_count)
+# calculate count of places of worship in PointX db for fields in urbanrural and provide percentages
+urbanrural$pow_count <- poly.counts(pow_pointX,urbanrural)
+urbanrural$pow_percent<- prop.table(urbanrural$pow_count)
+# calculate count of Transition for fields in urbanrural
+urbanrural$transition_count <- poly.counts(transition,urbanrural)
+urbanrural$transition_percent<- prop.table(urbanrural$transition_count)
+# calculate count of dtas for fields in urbanrural
+urbanrural$dtas_count <- poly.counts(dtas,urbanrural)
+urbanrural$dtas_percent<- prop.table(urbanrural$dtas_count)
+# calculate count of permaculture for fields in urbanrural
+urbanrural$permaculture_count <- poly.counts(permaculture,urbanrural)
+urbanrural$permaculture_percent<- prop.table(urbanrural$permaculture_count)
+
+# Create dataframe for analysis
+urbanrural_df<-data.frame(urbanrural)
+```
+
+Rather than bifurcate congregations into an urban/rural dichotomy, for this study we used the Scottish Government's eight-point remoteness scale to categorise eco-congregations along a spectrum of highly populated to remote areas. This 8-fold scale (calculated biennially) offers a more nuanced measurement that combines measurements of remoteness and population along the following lines:
+
+1. Large Urban Areas - Settlements of over 125,000 people.
+2. Other Urban Areas - Settlements of 10,000 to 125,000 people. 
+3. Accessible Small Towns - Settlements of between 3,000 and 10,000 people, and within a 30 minute drive time of a Settlement of 10,000 or more. 
+4. Remote Small Towns - Settlements of between 3,000 and 10,000 people, and with a drive time between 30 and 60 minutes to a Settlement of 10,000 or more.
+5. Very Remote Small Towns - Settlements of between 3,000 and 10,000 people, and with a drive time of over 60 minutes to a Settlement of 10,000 or more.
+6. Accessible Rural Areas - Areas with a population of less than 3,000 people, and within a drive time of 30 minutes to a Settlement of 10,000 or more.
+7. Remote Rural Areas - Areas with a population of less than 3,000 people, and with a drive time of between 30 and 60 minutes to a Settlement of 10,000 or more.
+8. Very Remote Rural Areas - Areas with a population of less than 3,000 people, and with a drive time of over 60 minutes to a Settlement of 10,000 or more.
+
+The key question which this analysis seeks to answer is whether ECS, or the other groups surveyed, are more concentrated in Urban or Rural areas, so as is the case below with our analysis of deprivation, we are concerned with the outer conditions, i.e. the urban areas (items 1-2) and remote areas (items 7-8).
+
+Of all the groups surveyed in this study, Eco-Congregation Scotland is the most heavily concentrated in large urban areas (33.53%), exceeding by almost 50% the rate for all places of worship (22.96% in large urban areas). Transition is a much more modest 20% and development trusts a bit lower at 15%. It is interesting to note that the rate of ECS concentration in these large urban areas matches the level of overall population distribution (34.5%). On the other end of the scale, Eco-Congregation Scotland is the least concentrated in remote rural areas (with 3.93% on level 7 and 5.44% on level 8 on the urban-rural scale), though again, they correlate roughly to the general population distribution (3.2% and 2.9% respectively). Places of worship outpace both the population of Scotland and the footprint of Eco-Congregation Scotland, with 14.98% in very remote rural areas, but this is exceeded by transition at 16.47% and both by Scottish community development trusts at 32.14%. So while Eco-Congregation Scotland correlates roughly with Scottish population distribution across the urban-rural scale, it has a considerably more urban profile than either of the other two groups surveyed.
+
+```{r create_urbanrural_ecs_chart_choropleth}
+# Create bar chart
+# TODO: use ggplot here to generate stacked bar chart showing representation across 8 urban/rural bands from urbanrural_df
+# TODO: generate map based on urbanrural polygons, shade based on 8-fold categories, draw transparent dots onto polygons
+
+pdf(file="figures/urbanrural_choropleth_ecs.pdf", width=4, height=4)
+bins <- unique(quantile(urbanrural$ecs_count, seq(0,1,length.out=30)))
+urbanrural$binId01 <- findInterval(urbanrural$ecs_count, bins)
+colSet01 <- rev(heat.colors(length(bins)))
+plot(urbanrural, col=colSet01[urbanrural$binId01], border="grey", lwd=0.25)
+par(mar=c(5,3,2,2)+0.1)
+# use the following to add points to the map (with transparency)
+# points(ecs, pch='.', col=rgb(0,0,0,alpha=0.15))
+title(main="Figure 5", sub="Eco-Congregation Scotland\nconcentrations in Urban Rural 8-fold classifications", cex.main=0.75, cex.sub=0.5)
+# save to file
+# dev.copy(png,'figures/urbanrural_choropleth_1.png')
+dev.off()
+```
+
+![][Figure6 map of ECS concentration by UR Areas Urban - Rural via Brown to Green]
+
+
+# Wealth, Employment, and Literacy
+
+```{r simd}
+# read in relevant polygons, Scottish Index of Multiple deprivation
+download.file("http://simd.scot/2016/data/simd2016_withgeog.zip", 
+              destfile = "data/simd2016_withgeog.zip")
+unzip("data/simd2016_withgeog.zip", exdir = "data", junkpaths = TRUE)
+simd_shapes <- readOGR("./data", "sc_dz_11") 
+simd_indicators <- read.csv("./data/simd2016_withinds.csv")
+simd <- merge(x=simd_shapes, y=simd_indicators, by="Data_Zone")
+  
+# commenting out old dataset in light of reproducible version above
+# simd <- readOGR("data", "simd_04-12_all_data")
+
+# STAGE 1, augment each dataset with relevant (geolocated) columns from SIMD
+# examine which ecs fall within each SIMD classification
+cbind(ecs@data, over(ecs, simd))
+# assign combined table with SIMD columns to attribute table slot of ecs table
+ecs@data=cbind(ecs@data,over(ecs,simd))
+# examine where pointX falls within each SIMD classification
+cbind(pow_pointX@data, over(pow_pointX, simd))
+# assign combined table with SIMD columns to attribute table slot of ecs table
+pow_pointX@data=cbind(pow_pointX@data,over(pow_pointX,simd))
+
+# examine which transition fall within each SIMD classifications
+cbind(transition@data, over(transition, simd))
+# assign combined table with SIMD columns to attribute table slot of transition table
+transition@data=cbind(transition@data,over(transition,simd))
+# examine which permaculture fall within each SIMD classifications
+cbind(permaculture@data, over(permaculture, simd))
+# assign combined table with SIMD columns to attribute table slot of permaculture table
+permaculture@data=cbind(permaculture@data,over(permaculture,simd))
+# examine which dtas fall within each SIMD classifications
+cbind(dtas@data, over(dtas, simd))
+# assign combined table with SIMD columns to attribute table slot of dtas table
+dtas@data=cbind(dtas@data,over(dtas,simd))
+
+# STAGE 2, extract NULL cells from each data set to prevent errors in stage 3
+
+# convert back to data frame for null cell extraction
+ecs<-data.frame(ecs)
+# split out null and normal cells
+ecs_clean<-ecs[complete.cases(ecs),]
+ecs_null<-ecs[!complete.cases(ecs),]
+# convert back to spatialpointdataframe
+coordinates(ecs) <- c("X", "Y")
+proj4string(ecs) <- proj4string(admin_lev1)
+
+# convert back to data frame for null cell extraction
+transition<-data.frame(transition)
+# split out null and normal cells
+transition_clean<-transition[complete.cases(transition),]
+transition_null<-transition[!complete.cases(transition),]
+# convert back to spatialpointdataframe
+coordinates(transition) <- c("X", "Y")
+proj4string(transition) <- proj4string(admin_lev1)
+
+# convert back to data frame for null cell extraction
+permaculture<-data.frame(permaculture)
+# split out null and normal cells
+permaculture_clean<-permaculture[complete.cases(permaculture),]
+permaculture_null<-permaculture[!complete.cases(permaculture),]
+# convert back to spatialpointdataframe
+coordinates(permaculture) <- c("X", "Y")
+proj4string(permaculture) <- proj4string(admin_lev1)
+
+# convert back to data frame for null cell extraction
+dtas<-data.frame(dtas)
+# split out null and normal cells
+dtas_clean<-dtas[complete.cases(dtas),]
+dtas_null<-dtas[!complete.cases(dtas),]
+# convert back to spatialpointdataframe
+coordinates(dtas) <- c("X", "Y")
+proj4string(dtas) <- proj4string(admin_lev1)
+
+# STAGE 3a, calculate sums based on SIMD12R columns and generate new integer sets with quintile count data
+simd_rownames = c("Quintile 1","Quintile 2","Quintile 3","Quintile 4","Quintile 5")
+simdr12_ecs = c((sum(ecs_clean$SIMDR12<1301)), (sum(ecs_clean$SIMDR12 > 1300 & ecs_clean$SIMDR12 < 2602)), (sum(ecs_clean$SIMDR12 > 2601 & ecs_clean$SIMDR12 < 3903)), (sum(ecs_clean$SIMDR12 > 3902 & ecs_clean$SIMDR12 < 5204)), (sum(ecs_clean$SIMDR12 > 5203 & ecs_clean$SIMDR12 < 6505)))
+# names(simdr12_ecs) <- simd_rownames
+simdr12_transition = c((sum(transition_clean$SIMDR12<1301)), (sum(transition_clean$SIMDR12 > 1300 & transition_clean$SIMDR12 < 2602)), (sum(transition_clean$SIMDR12 > 2601 & transition_clean$SIMDR12 < 3903)), (sum(transition_clean$SIMDR12 > 3902 & transition_clean$SIMDR12 < 5204)), (sum(transition_clean$SIMDR12 > 5203 & transition_clean$SIMDR12 < 6505)))
+# names(simdr12_transition) <- simd_rownames
+simdr12_permaculture = c((sum(permaculture_clean$SIMDR12<1301)), (sum(permaculture_clean$SIMDR12 > 1300 & permaculture_clean$SIMDR12 < 2602)), (sum(permaculture_clean$SIMDR12 > 2601 & permaculture_clean$SIMDR12 < 3903)), (sum(permaculture_clean$SIMDR12 > 3902 & permaculture_clean$SIMDR12 < 5204)), (sum(permaculture_clean$SIMDR12 > 5203 & permaculture_clean$SIMDR12 < 6505)))
+# names(simdr12_permaculture) <- simd_rownames
+simdr12_dtas = c((sum(dtas_clean$SIMDR12<1301)), (sum(dtas_clean$SIMDR12 > 1300 & dtas_clean$SIMDR12 < 2602)), (sum(dtas_clean$SIMDR12 > 2601 & dtas_clean$SIMDR12 < 3903)), (sum(dtas_clean$SIMDR12 > 3902 & dtas_clean$SIMDR12 < 5204)), (sum(dtas_clean$SIMDR12 > 5203 & dtas_clean$SIMDR12 < 6505)))
+# names(simdr12_dtas) <- simd_rownames
+
+# STAGE 3b, calculate sums based on INCR12 columns and generate new integer sets with quintile count data
+incr12_ecs = c((sum(ecs_clean$INCR12<1301)), (sum(ecs_clean$INCR12 > 1300 & ecs_clean$INCR12 < 2602)), (sum(ecs_clean$INCR12 > 2601 & ecs_clean$INCR12 < 3903)), (sum(ecs_clean$INCR12 > 3902 & ecs_clean$INCR12 < 5204)), (sum(ecs_clean$INCR12 > 5203 & ecs_clean$INCR12 < 6505)))
+incr12_transition = c((sum(transition_clean$INCR12<1301)), (sum(transition_clean$INCR12 > 1300 & transition_clean$INCR12 < 2602)), (sum(transition_clean$INCR12 > 2601 & transition_clean$INCR12 < 3903)), (sum(transition_clean$INCR12 > 3902 & transition_clean$INCR12 < 5204)), (sum(transition_clean$INCR12 > 5203 & transition_clean$INCR12 < 6505)))
+incr12_permaculture = c((sum(permaculture_clean$INCR12<1301)), (sum(permaculture_clean$INCR12 > 1300 & permaculture_clean$INCR12 < 2602)), (sum(permaculture_clean$INCR12 > 2601 & permaculture_clean$INCR12 < 3903)), (sum(permaculture_clean$INCR12 > 3902 & permaculture_clean$INCR12 < 5204)), (sum(permaculture_clean$INCR12 > 5203 & permaculture_clean$INCR12 < 6505)))
+incr12_dtas = c((sum(dtas_clean$INCR12<1301)), (sum(dtas_clean$INCR12 > 1300 & dtas_clean$INCR12 < 2602)), (sum(dtas_clean$INCR12 > 2601 & dtas_clean$INCR12 < 3903)), (sum(dtas_clean$INCR12 > 3902 & dtas_clean$INCR12 < 5204)), (sum(dtas_clean$INCR12 > 5203 & dtas_clean$INCR12 < 6505)))
+
+# STAGE 3c, calculate sums based on EMPR12 columns and generate new integer sets with quintile count data
+empr12_ecs = c((sum(ecs_clean$EMPR12<1301)), (sum(ecs_clean$EMPR12 > 1300 & ecs_clean$EMPR12 < 2602)), (sum(ecs_clean$EMPR12 > 2601 & ecs_clean$EMPR12 < 3903)), (sum(ecs_clean$EMPR12 > 3902 & ecs_clean$EMPR12 < 5204)), (sum(ecs_clean$EMPR12 > 5203 & ecs_clean$EMPR12 < 6505)))
+empr12_transition = c((sum(transition_clean$EMPR12<1301)), (sum(transition_clean$EMPR12 > 1300 & transition_clean$EMPR12 < 2602)), (sum(transition_clean$EMPR12 > 2601 & transition_clean$EMPR12 < 3903)), (sum(transition_clean$EMPR12 > 3902 & transition_clean$EMPR12 < 5204)), (sum(transition_clean$EMPR12 > 5203 & transition_clean$EMPR12 < 6505)))
+empr12_permaculture = c((sum(permaculture_clean$EMPR12<1301)), (sum(permaculture_clean$EMPR12 > 1300 & permaculture_clean$EMPR12 < 2602)), (sum(permaculture_clean$EMPR12 > 2601 & permaculture_clean$EMPR12 < 3903)), (sum(permaculture_clean$EMPR12 > 3902 & permaculture_clean$EMPR12 < 5204)), (sum(permaculture_clean$EMPR12 > 5203 & permaculture_clean$EMPR12 < 6505)))
+empr12_dtas = c((sum(dtas_clean$EMPR12<1301)), (sum(dtas_clean$EMPR12 > 1300 & dtas_clean$EMPR12 < 2602)), (sum(dtas_clean$EMPR12 > 2601 & dtas_clean$EMPR12 < 3903)), (sum(dtas_clean$EMPR12 > 3902 & dtas_clean$EMPR12 < 5204)), (sum(dtas_clean$EMPR12 > 5203 & dtas_clean$EMPR12 < 6505)))
+
+# STAGE 3d, calculate sums based on HER12 columns and generate new integer sets with quintile count data
+her12_ecs = c((sum(ecs_clean$HER12<1301)), (sum(ecs_clean$HER12 > 1300 & ecs_clean$HER12 < 2602)), (sum(ecs_clean$HER12 > 2601 & ecs_clean$HER12 < 3903)), (sum(ecs_clean$HER12 > 3902 & ecs_clean$HER12 < 5204)), (sum(ecs_clean$HER12 > 5203 & ecs_clean$HER12 < 6505)))
+her12_transition = c((sum(transition_clean$HER12<1301)), (sum(transition_clean$HER12 > 1300 & transition_clean$HER12 < 2602)), (sum(transition_clean$HER12 > 2601 & transition_clean$HER12 < 3903)), (sum(transition_clean$HER12 > 3902 & transition_clean$HER12 < 5204)), (sum(transition_clean$HER12 > 5203 & transition_clean$HER12 < 6505)))
+her12_permaculture = c((sum(permaculture_clean$HER12<1301)), (sum(permaculture_clean$HER12 > 1300 & permaculture_clean$HER12 < 2602)), (sum(permaculture_clean$HER12 > 2601 & permaculture_clean$HER12 < 3903)), (sum(permaculture_clean$HER12 > 3902 & permaculture_clean$HER12 < 5204)), (sum(permaculture_clean$HER12 > 5203 & permaculture_clean$HER12 < 6505)))
+her12_dtas = c((sum(dtas_clean$HER12<1301)), (sum(dtas_clean$HER12 > 1300 & dtas_clean$HER12 < 2602)), (sum(dtas_clean$HER12 > 2601 & dtas_clean$HER12 < 3903)), (sum(dtas_clean$HER12 > 3902 & dtas_clean$HER12 < 5204)), (sum(dtas_clean$HER12 > 5203 & dtas_clean$HER12 < 6505)))
+
+# STAGE 3e, calculate sums based on ESTR12 columns and generate new integer sets with quintile count data
+estr12_ecs = c((sum(ecs_clean$ESTR12<1301)), (sum(ecs_clean$ESTR12 > 1300 & ecs_clean$ESTR12 < 2602)), (sum(ecs_clean$ESTR12 > 2601 & ecs_clean$ESTR12 < 3903)), (sum(ecs_clean$ESTR12 > 3902 & ecs_clean$ESTR12 < 5204)), (sum(ecs_clean$ESTR12 > 5203 & ecs_clean$ESTR12 < 6505)))
+estr12_transition = c((sum(transition_clean$ESTR12<1301)), (sum(transition_clean$ESTR12 > 1300 & transition_clean$ESTR12 < 2602)), (sum(transition_clean$ESTR12 > 2601 & transition_clean$ESTR12 < 3903)), (sum(transition_clean$ESTR12 > 3902 & transition_clean$ESTR12 < 5204)), (sum(transition_clean$ESTR12 > 5203 & transition_clean$ESTR12 < 6505)))
+estr12_permaculture = c((sum(permaculture_clean$ESTR12<1301)), (sum(permaculture_clean$ESTR12 > 1300 & permaculture_clean$ESTR12 < 2602)), (sum(permaculture_clean$ESTR12 > 2601 & permaculture_clean$ESTR12 < 3903)), (sum(permaculture_clean$ESTR12 > 3902 & permaculture_clean$ESTR12 < 5204)), (sum(permaculture_clean$ESTR12 > 5203 & permaculture_clean$ESTR12 < 6505)))
+estr12_dtas = c((sum(dtas_clean$ESTR12<1301)), (sum(dtas_clean$ESTR12 > 1300 & dtas_clean$ESTR12 < 2602)), (sum(dtas_clean$ESTR12 > 2601 & dtas_clean$ESTR12 < 3903)), (sum(dtas_clean$ESTR12 > 3902 & dtas_clean$ESTR12 < 5204)), (sum(dtas_clean$ESTR12 > 5203 & dtas_clean$ESTR12 < 6505)))
+
+# STAGE 4a - calculate percentages
+simdr12_ecs_percent<- prop.table(simdr12_ecs)
+simdr12_transition_percent<- prop.table(simdr12_transition)
+simdr12_permaculture_percent<- prop.table(simdr12_permaculture)
+simdr12_dtas_percent<- prop.table(simdr12_dtas)
+incr12_ecs_percent<- prop.table(incr12_ecs)
+incr12_transition_percent<- prop.table(incr12_transition)
+incr12_permaculture_percent<- prop.table(incr12_permaculture)
+incr12_dtas_percent<- prop.table(incr12_dtas)
+empr12_ecs_percent<- prop.table(empr12_ecs)
+empr12_transition_percent<- prop.table(empr12_transition)
+empr12_permaculture_percent<- prop.table(empr12_permaculture)
+empr12_dtas_percent<- prop.table(empr12_dtas)
+her12_ecs_percent<- prop.table(her12_ecs)
+her12_transition_percent<- prop.table(her12_transition)
+her12_permaculture_percent<- prop.table(her12_permaculture)
+her12_dtas_percent<- prop.table(her12_dtas)
+estr12_ecs_percent<- prop.table(estr12_ecs)
+estr12_transition_percent<- prop.table(estr12_transition)
+estr12_permaculture_percent<- prop.table(estr12_permaculture)
+estr12_dtas_percent<- prop.table(estr12_dtas)
+
+# STAGE 4b, generate data frame using integer sets
+simd = data.frame(simdr12_ecs, simdr12_ecs_percent, incr12_ecs, incr12_ecs_percent, empr12_ecs, empr12_ecs_percent, her12_ecs, her12_ecs_percent, estr12_ecs, estr12_ecs_percent, simdr12_transition, simdr12_transition_percent, incr12_transition, incr12_transition_percent, empr12_transition, empr12_transition_percent, her12_transition, her12_transition_percent, estr12_transition, estr12_transition_percent, simdr12_permaculture, simdr12_permaculture_percent, incr12_permaculture, incr12_permaculture_percent, empr12_permaculture, empr12_permaculture_percent, her12_permaculture, her12_permaculture_percent, estr12_permaculture, estr12_permaculture_percent, simdr12_dtas, simdr12_dtas_percent, incr12_dtas, incr12_dtas_percent, empr12_dtas, empr12_dtas_percent, her12_dtas, her12_dtas_percent, estr12_dtas, estr12_dtas_percent)
+write.csv(simd, "derivedData/simd.csv", row.names=FALSE)
+
+
+simd_percents_only = data.frame(simd_rownames, simdr12_ecs_percent, incr12_ecs_percent, empr12_ecs_percent, her12_ecs_percent, estr12_ecs_percent, simdr12_transition_percent, incr12_transition_percent, empr12_transition_percent, her12_transition_percent, estr12_transition_percent, simdr12_permaculture_percent, incr12_permaculture_percent, empr12_permaculture_percent, her12_permaculture_percent, estr12_permaculture_percent, simdr12_dtas_percent, incr12_dtas_percent, empr12_dtas_percent, her12_dtas_percent, estr12_dtas_percent)
+write.csv(simd_percents_only, "derivedData/simd_percents_only.csv", row.names=FALSE)
+
+
+# STAGE 5, generate cool charts
+
+# clustered bar charts
+# convert to long format
+library(reshape2)
+simd_percents_only_long <- melt(simd_percents_only, id.vars = "simd_rownames", 
+              measure.vars = grep("^12", names(simd_percents_only), value = TRUE))
+
+qplot(data=simd_percents_only_long , geom="bar", fill=(factor(simd_rownames)))
+
+# jitterplot option, from Teutonico 2015, p. 63
+
+```
+
+Another crucial point of assessment relates to the relation of Eco-Congregation communities to the Scottish Index of Multiple Deprivation. This instrument aggregates a large variety of factors which can lead to deprivation including crime rates, employment levels, access to services (implicating remoteness), and literacy. By assessing ECS, Transition, and dtas against the deprivation scale, we can assess whether eco-congregations fall within particular demographics and also whether the fully aggregated SIMD measurement provides a useful point of comparison for our purposes. The SIMD essentially divides Scotland into 6407 geographic zones and then ranks them based on their relative deprivation. This data set can be split into any number of groups, but for our purposes we have settled on Quintiles, splitting the SIMD data set at every 1302 entries. We then measured where each transition group, ECS, and dtas fell within these zones and calculated how they fell into these five quintiles, from more to least deprived.
+
+The first, and most compelling finding is that, in general Eco-Congregation Scotland and Transition Scotland are both roughly the same and match the level of population distribution in the lowest quintile of the general SIMD measurement. 8% of transition groups and eco-congregation groups which have received awards and 9% of the population are located within this quintile. However, taken in relation to the distribution of places of worship in the lowest quintile, we find that eco-congregations are located at half the rate that places of worship are (15%) and dtass match this much more closely at 14%. Turning towards the top quintile, this pattern also holds, here both transition groups (21%) and eco-congregations (21% and 29% of awarded congregations) depart from the population distribution in this upper quintile (which is 10%). Again, general places of worship (at 11%) and DTASs (at 5%) take the opposite direction. We can say decisively that in communities which have been identified as good candidates for intervention to reduce deprivation, ECS and Transition are less likely, and they are over-represented at the areas which fall into the least deprived quintile.
+
+We can find divergence between transition communities and eco-congregation when we split out SIMD domains. In the lowest quartile, measuring exclusively for the income domain, ECS is more represented (11%) - roughly the same as DTAS (12%), and transition is less (6%) represented. In general (as shown on the chart in Appendix D), these trends hold when representation of our groups are measured within other non-remoteness domains of the SIMD. Our basic conclusion is that transition towns are least likely to operate within the lowest quartile of SIMD and DTASs are most likely, with ECS somewhere in the middle. Given the general disparity against the presence of places of worship, it seems fair to suggest that this might be an area for improvement, perhaps even worth developing a special programme which might target areas in SIMD quartile 1 for eco-congregation outreach. This might be considered particularly in light of the starkest underrepresentation of ECS and transition within the SIMD domain of education, skills, and training.
+
+
+
+
+
+
+# Appendix B
+
+(JK note to self: same as above, but augmented with multipliers by which categories are different from one another)
+
+# Appendix C - Data by Urban / Rural Classification
+
+```{r pander_urbanrural_table}
+urbanrural.shortened <- urbanrural[,c(2,6,9:18)]
+write.csv(urbanrural, "derivedData/urbanrural.csv", row.names=FALSE)
+write.csv(urbanrural.shortened, "derivedData/urbanrural.csv", row.names=FALSE)
+urbanrural.shortened<-data.frame(urbanrural.shortened)
+panderOptions("digits", 2)
+pander(urbanrural.shortened)
\ No newline at end of file