WA_active_schools.R

   1 # install libgdal-dev on debian first
   2 #install.packages("rgdal")
   3 library(rgdal)
   4 library(tidyverse)
   5 #install.packages("OpenStreetMap")
   6 library(OpenStreetMap)
   7 library(viridis)
   8
   9 # Read the shapefile obtained from WA opendata (SLIP)
  10 shape <- readOGR(dsn = "data/", layer = "CurrentActiveSchoolsSemester12017_PublicDET_014")
  11 shape.df <- tidy(shape)
  12 shape.df$latitude <- as.numeric(as.character(shape.df$latitude))
  13 shape.df$longitude <- as.numeric(as.character(shape.df$longitude))
  14 shape.df$totalschoo <- as.numeric(as.character(shape.df$totalschoo))
  15
  16 # filter out remote islands
  17 shape.df <- filter(shape.df,longitude > 100)
  18 #shape.df <- filter(shape.df, totalschoo < 500)
  19 # Obtain the contour box of the shape
  20
  21 lat <- c(min(shape.df$latitude), max(shape.df$latitude))
  22 lon <- c(min(shape.df$longitude), max(shape.df$longitude))
  23
  24 lat.metro <- c(-32.1547, -31.8519)
  25 lon.metro <- c(115.7162,116.0733)
  26 upperLeft.metro <- c(-31.8519,115.7162)
  27 lowerRight.metro <- c(-32.1547,116.0733)
  28
  29 # Obtain an openstreetmap background image
  30 map <- openmap(upperLeft = c(lat[2],lon[1]),
  31                lowerRight = c(lat[1],lon[2]),
  32                minNumTiles = 6,
  33                type="https://stamen-tiles.a.ssl.fastly.net/toner-lite/{z}/{x}/{y}.png")
  34 map.metro <- openmap(upperLeft = upperLeft.metro,
  35                      lowerRight = lowerRight.metro,
  36                      minNumTiles = 6,
  37                      type="https://stamen-tiles.a.ssl.fastly.net/toner-lite/{z}/{x}/{y}.png")
  38
  39 mapLatLon <- openproj(map)
  40
  41 autoplot(mapLatLon) +
  42   geom_point(data = shape.df,
  43              aes(x = longitude,
  44                  y = latitude,
  45                  color = totalschoo),
  46              show.legend = TRUE) +
  47   scale_color_viridis(option="viridis",
  48                       name = "Number\nof pupils") +
  49   theme(legend.position = "right",
  50         axis.text = element_blank(),
  51         axis.ticks = element_blank()) +
  52   labs(title = "Western Australian schools",
  53        caption = "Map tiles by Stamen Design, under CC BY 3.0.\nMapdata by OpenStreetMap, under ODbL.
  54        School data © Government of WA 2017.\n Dep. of Education & Training.",
  55        x="",
  56        y="")
  57
  58
  59 # Zooming in on metro
  60 shape.metro.df <- filter(shape.df, abs(latitude) <= max(abs(lat.metro)),
  61                          abs(latitude) >= min(abs(lat.metro)),
  62                          abs(longitude) <= max(abs(lon.metro)),
  63                          abs(longitude) >= min(abs(lon.metro)))
  64 mapLatLon.metro <- openproj(map.metro)
  65 autoplot(mapLatLon.metro) +
  66   geom_point(data = shape.metro.df,
  67              aes(x = longitude,
  68                  y = latitude,
  69                  color = totalschoo),
  70              show.legend = TRUE) +
  71   scale_color_viridis(option="viridis",
  72                       name = "Number\nof pupils") +
  73   theme(legend.position = "right",
  74         axis.text = element_blank(),
  75         axis.ticks = element_blank()) +
  76   labs(title = "Western Australian schools",
  77        subtitle = "Perth metro area only",
  78        caption = "Map tiles by Stamen Design, under CC BY 3.0.\nMapdata by OpenStreetMap, under ODbL.
  79        School data © Government of WA 2017.\n Dep. of Education & Training.",
  80        x="",
  81        y="")