Begin by preparing this working environment

Load relevant libraries

# For working between R and a postgres database for spatial analysis
library(RPostgreSQL)
## Loading required package: DBI
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(RODBC)
library(rgdal)
## Loading required package: sp
## rgdal: version: 1.2-7, (SVN revision 660)
##  Geospatial Data Abstraction Library extensions to R successfully loaded
##  Loaded GDAL runtime: GDAL 2.0.1, released 2015/09/15
##  Path to GDAL shared files: C:/Users/Blake/Documents/R/win-library/3.3/rgdal/gdal
##  Loaded PROJ.4 runtime: Rel. 4.9.2, 08 September 2015, [PJ_VERSION: 492]
##  Path to PROJ.4 shared files: C:/Users/Blake/Documents/R/win-library/3.3/rgdal/proj
##  Linking to sp version: 1.2-4
# Producing interactive in R
library(leaflet)
library(RColorBrewer)

Establish a connection to a local postgres database

drv.psql <- dbDriver("PostgreSQL")
con <- dbConnect(drv = drv.psql, host="localhost", port="5432", dbname="spatialDB", user="postgres", pass="password")

The next section of code was used for reading in previously scraped Hungry Jack’s store locations, and loading it into postgres. The lat lons had to be converted into geometry data types for quick visualisation in QGIS

# Load the hungry jacks data
hj.data <- read.csv("HungryJacks.csv") %>%
  select(-ts)
row.names(hj.data) <- NULL

# Create the table
dbGetQuery(con, "DROP TABLE IF EXISTS bovine.jacks;")
dbGetQuery(con, "CREATE TABLE bovine.jacks ( 
id int
, title varchar(250)
, lon numeric(12,6) 
, lat numeric(12,6)
);"
)
# Push the data in
dbWriteTable(con, c("bovine", "jacks"), hj.data, row.names = F, append = T)

# Create point geoms
dbGetQuery(con, "ALTER TABLE bovine.jacks ADD COLUMN wkb_geometry geometry(Point, 3577);")
dbGetQuery(con, "UPDATE bovine.jacks SET wkb_geometry = ST_Transform(ST_SetSRID(ST_MakePoint(lon, lat), 4283), 3577);")

Repeat the same process for scraped McDonalds restaurant locations

# Do the same for the McDonalds locations
mc.data <- read.csv("McDonalds.csv") %>%
  select(-ts)
row.names(mc.data) <- NULL

# Create the table
dbGetQuery(con, "DROP TABLE IF EXISTS bovine.ronald;")
dbGetQuery(con, "CREATE TABLE bovine.ronald ( 
           id int
           , title varchar(250)
           , lon numeric(12,6) 
           , lat numeric(12,6)
);"
)
# Push the data in
dbWriteTable(con, c("bovine", "ronald"), mc.data, row.names = F, append = T)

# Create point geoms
dbGetQuery(con, "ALTER TABLE bovine.ronald ADD COLUMN wkb_geometry geometry(Point, 3577);")
dbGetQuery(con, "UPDATE bovine.ronald SET wkb_geometry = ST_Transform(ST_SetSRID(ST_MakePoint(lon, lat), 4283), 3577);")

The final step was to load some SA2 data and join it to an existing table of SA2 polygon geometries to make a thematic map in QGIS

# Load some of our attributes
master <- read.csv("modl_master.csv") %>%
  select(REGION_NO, beef_spend_per_person)

# Create the table
dbGetQuery(con, "DROP TABLE IF EXISTS bovine.master;")
dbGetQuery(con, "CREATE TABLE bovine.master ( 
           region_id int
           , variable numeric(8,4)
);"
)
# Push the data in
dbWriteTable(con, c("bovine", "master"), master, row.names = F, append = T)

# Create the vis table once, but don't delete it every subsequent time
if(!dbExistsTable(con, c("bovine", "vis_master"))) {
  
  dbGetQuery(con, "CREATE TABLE bovine.vis_master 
  (
  region_id int
  , variable numeric(8,4)
  , wkb_geometry geometry(MultiPolygon, 3577)
  );")
  print("Table created")
  
} else {
  print("Table exists")
}

# Join it to the geom table for visualisation
dbGetQuery(con, "
DELETE FROM bovine.vis_master;
INSERT INTO bovine.vis_master (
SELECT a.region_id
, a.variable 
, b.wkb_geometry
FROM bovine.master a 
JOIN geom.sa2_11 b 
ON a.region_id = CAST(b.sa2_main11 as int)
);")

With some exploration having taken place outside of R, the next step was to try and replicate the experience within R. Static maps can be made with ggplot, but I felt leaflet was pretty good for exploratory analysis

# Let's try make some plots in R
# Prepare the dataset for plotting

mc.data <- read.csv("McDonalds.csv") %>%
  select(-ts)
row.names(mc.data) <- NULL

hj.data <- read.csv("HungryJacks.csv") %>%
  select(-ts)
row.names(hj.data) <- NULL

# Load some of our attributes
master <- read.csv("modl_master.csv") %>%
  select(REGION_NO, beef_spend_per_person)

# Bin the attribute to make a thematic map in leaflet
master.bin <- master %>%
  mutate(bin = ntile(beef_spend_per_person, n = 5))

The next step is to load a polygon shapefile and merge our data on beef spend per capita on SA2 region number

# Read in the shapefile
shp <- readOGR(".", layer = "SA2_2011_AUST")
## OGR data source with driver: ESRI Shapefile 
## Source: ".", layer: "SA2_2011_AUST"
## with 2214 features
## It has 12 fields
## Warning in readOGR(".", layer = "SA2_2011_AUST"): Dropping null geometries:
## 539, 540, 974, 975, 1502, 1503, 1674, 1675, 1926, 1927, 2026, 2027, 2096,
## 2097, 2208, 2209, 2213, 2214
# Merge the data on
shp@data <- shp@data %>%
  mutate(REGION_NO = as.integer(as.character(SA2_MAIN11))) %>%
  left_join(master.bin, by = c("REGION_NO" = "REGION_NO"))

The final step is to combine our datasets into a leaflet. For the polygon fill, I have used the bins I previously created. Darker areas will indicate higher beef spend per person. I have also overlaid the

# Try using leaflet instead
archesLeafIcon <- makeIcon(
  iconUrl = "McDonald's_Golden_Arches.png"
  , iconWidth = 10
  , iconHeight = 10
)
# Try using leaflet instead
jacksLeafIcon <- makeIcon(
  iconUrl = "Burger_King-logo-EB00FAD8D3-seeklogo.com.png"
  , iconWidth = 10
  , iconHeight = 10
)

# An important step for visual representation - simplify the polygon geometries such that 
shp.simple <- rmapshaper::ms_simplify(shp)

leaflet(shp.simple) %>%
  addPolygons(fill = T
              , fillColor = ~colorBin(palette = colorRampPalette(brewer.pal(5, "Blues"))(5), domain = bin, bins = 5)(bin)
              , fillOpacity = 0.75
              , color = "grey80"
              , weight = 0.5
              , stroke = T
              , group = "SA2 regions") %>%
  addMarkers(lng = mc.data$lon.coord
             , lat = mc.data$lat.coord
             , icon = archesLeafIcon
             , group = "Restaurants"
             #, clusterOptions = markerClusterOptions()
  ) %>%
  addMarkers(lng = hj.data$lon.coord
             , lat = hj.data$lat.coord
             , icon = jacksLeafIcon
             , group = "Restaurants"
             #, clusterOptions = markerClusterOptions()
  ) %>%
  addLayersControl(overlayGroups = c("Restaurants", "SA2 regions"),
    options = layersControlOptions(collapsed = FALSE)
  )