From 4a262294342cd273e52b82bb33d371904f840449 Mon Sep 17 00:00:00 2001 From: Frederik Vanrenterghem Date: Thu, 12 Oct 2017 20:57:48 +0800 Subject: [PATCH 1/5] Consider metro high schools only; label biggest. --- WA_active_schools.R | 21 ++++++++++++++++----- 1 file changed, 16 insertions(+), 5 deletions(-) diff --git a/WA_active_schools.R b/WA_active_schools.R index da884f9..5076c60 100644 --- a/WA_active_schools.R +++ b/WA_active_schools.R @@ -14,6 +14,7 @@ shape.df <- st_read(dsn = "data/", layer = "CurrentActiveSchoolsSemester12017_Pu shape.df$latitude <- as.numeric(shape.df$latitude) shape.df$longitude <- as.numeric(shape.df$longitude) shape.df$totalschoo <- as.numeric(shape.df$totalschoo) +shape.df$totalsecon <- as.numeric(shape.df$totalsecon) # filter out remote islands shape.df <- filter(shape.df,longitude > 100) @@ -64,19 +65,29 @@ shape.metro.df <- filter(shape.df, abs(latitude) <= max(abs(lat.metro)), abs(longitude) <= max(abs(lon.metro)), abs(longitude) >= min(abs(lon.metro))) mapLatLon.metro <- openproj(map.metro) +# Let's look at secondary school only +# 95th percentile of school size +secon95th <- quantile(shape.metro.df[,"totalsecon">0]$totalsecon,probs = .95) autoplot(mapLatLon.metro) + - geom_point(data = shape.metro.df, # geom_sf is not yet in ggplot2 on CRAN, so sticking to this for now + geom_point(data = filter(shape.metro.df, totalsecon > 0), # geom_sf is not yet in ggplot2 on CRAN, so sticking to this for now aes(x = longitude, y = latitude, - color = totalschoo), - show.legend = TRUE) + + color = totalsecon, + size = totalsecon) + ) + + geom_text(data = filter(shape.metro.df, totalsecon > secon95th), # name top 5% schools + aes(x = longitude, + y = latitude, + label = schoolname), + size = 2) + scale_color_viridis(option="viridis", name = "Number\nof pupils") + + guides(size = FALSE) + # Do not show legend for size, as color is the main key theme(legend.position = "right", axis.text = element_blank(), axis.ticks = element_blank()) + - labs(title = "Western Australian schools", - subtitle = "Perth metro area only", + labs(title = "Western Australian secondary schools", + subtitle = "Perth metro area", caption = "Map tiles by Stamen Design, under CC BY 3.0.\nMapdata by OpenStreetMap, under ODbL. School data © Government of WA 2017.\n Dep. of Education & Training.", x="", -- 2.39.2 From 55c322ee3b4b2849080d91159f3dccc76aa4ba82 Mon Sep 17 00:00:00 2001 From: Frederik Vanrenterghem Date: Mon, 2 Oct 2023 20:00:57 +0800 Subject: [PATCH 2/5] Update to 2022 data. - Use SLIP API via OAuth 2.0 to obtain shapefile. - Use 2022 instead of 2017 data. --- WA_active_schools.R | 45 ++++++++++++++++++++++++++++++--------------- 1 file changed, 30 insertions(+), 15 deletions(-) diff --git a/WA_active_schools.R b/WA_active_schools.R index 5076c60..18816aa 100644 --- a/WA_active_schools.R +++ b/WA_active_schools.R @@ -1,16 +1,32 @@ -# install libgdal-dev on debian first -# install libudunits2-dev on debian first -#install.packages("sf") library(sf) -library(tidyverse) -library(broom) -#install.packages("OpenStreetMap") +library(dplyr) +library(ggplot2) library(OpenStreetMap) library(viridis) +library(httr2) -# Read the shapefile obtained from WA opendata (SLIP) -#shape <- readOGR(dsn = "data/", layer = "CurrentActiveSchoolsSemester12017_PublicDET_014") -shape.df <- st_read(dsn = "data/", layer = "CurrentActiveSchoolsSemester12017_PublicDET_014", stringsAsFactors = FALSE) +## Obtain the shapefile from WA opendata (SLIP) +tempdirSHP <- tempdir() +tempfileSHP <- tempfile() +req <- request("https://sso.slip.wa.gov.au/as/token.oauth2") |> + req_headers("Authorization" = "Basic ZGlyZWN0LWRvd25sb2Fk") |> + req_body_form(list("grant_type" = "password", + # SLIP username and password stored in + #pass - the standard unix password manager + "username" = system2("pass", args = "slip.wa.gov.au | grep Username | sed -e 's/Username: //'", stdout = TRUE), + "password" = system2("pass", args = "slip.wa.gov.au | head -1", stdout = TRUE))) + +tokenResponse <- req_perform(req) + +slipUrl <- "https://direct-download.slip.wa.gov.au/datadownload/Education/Current_Active_Schools_Sem_1_2022_Public_DET_020_WA_GDA94_Public_Shapefile.zip" + +req <- request(slipUrl) |> + req_headers( 'Authorization' = paste0('Bearer ',resp_body_json(tokenResponse)$access_token)) + +responseSlip <- req_perform(req) +writeBin(resp_body_raw(responseSlip), con = tempfileSHP) +unzip(tempfileSHP, exdir = tempdirSHP) +shape.df <- st_read(dsn = tempdirSHP) shape.df$latitude <- as.numeric(shape.df$latitude) shape.df$longitude <- as.numeric(shape.df$longitude) shape.df$totalschoo <- as.numeric(shape.df$totalschoo) @@ -18,7 +34,6 @@ shape.df$totalsecon <- as.numeric(shape.df$totalsecon) # filter out remote islands shape.df <- filter(shape.df,longitude > 100) -#shape.df <- filter(shape.df, totalschoo < 500) # Obtain the contour box of the shape lat <- c(min(shape.df$latitude), max(shape.df$latitude)) @@ -41,7 +56,7 @@ map.metro <- openmap(upperLeft = upperLeft.metro, mapLatLon <- openproj(map) -autoplot(mapLatLon) + +OpenStreetMap::autoplot.OpenStreetMap(mapLatLon) + geom_point(data = shape.df, aes(x = longitude, y = latitude, @@ -54,7 +69,7 @@ autoplot(mapLatLon) + axis.ticks = element_blank()) + labs(title = "Western Australian schools", caption = "Map tiles by Stamen Design, under CC BY 3.0.\nMapdata by OpenStreetMap, under ODbL. - School data © Government of WA 2017.\n Dep. of Education & Training.", + School data © Government of WA 2022.\n Dep. of Education & Training.", x="", y="") @@ -68,7 +83,7 @@ mapLatLon.metro <- openproj(map.metro) # Let's look at secondary school only # 95th percentile of school size secon95th <- quantile(shape.metro.df[,"totalsecon">0]$totalsecon,probs = .95) -autoplot(mapLatLon.metro) + +OpenStreetMap::autoplot.OpenStreetMap(mapLatLon.metro) + geom_point(data = filter(shape.metro.df, totalsecon > 0), # geom_sf is not yet in ggplot2 on CRAN, so sticking to this for now aes(x = longitude, y = latitude, @@ -82,13 +97,13 @@ autoplot(mapLatLon.metro) + size = 2) + scale_color_viridis(option="viridis", name = "Number\nof pupils") + - guides(size = FALSE) + # Do not show legend for size, as color is the main key + guides(size = 'none') + # Do not show legend for size, as color is the main key theme(legend.position = "right", axis.text = element_blank(), axis.ticks = element_blank()) + labs(title = "Western Australian secondary schools", subtitle = "Perth metro area", caption = "Map tiles by Stamen Design, under CC BY 3.0.\nMapdata by OpenStreetMap, under ODbL. - School data © Government of WA 2017.\n Dep. of Education & Training.", + School data © Government of WA 2022.\n Dep. of Education & Training.", x="", y="") -- 2.39.2 From 50417261b375a912e448e460fca150f4186fef10 Mon Sep 17 00:00:00 2001 From: Frederik Vanrenterghem Date: Mon, 2 Oct 2023 22:30:39 +0800 Subject: [PATCH 3/5] Remove OpenStreetMap dependency. Use geom_sf. * Replace OpenStreetMap with maptiles and tidyterra due to sp support being removed from CRAN in Oct 2023. * Use geom_sf as that is now part of ggplot2 on CRAN. --- WA_active_schools.R | 58 ++++++++++++++++++++++----------------------- 1 file changed, 28 insertions(+), 30 deletions(-) diff --git a/WA_active_schools.R b/WA_active_schools.R index 18816aa..0811702 100644 --- a/WA_active_schools.R +++ b/WA_active_schools.R @@ -1,7 +1,8 @@ library(sf) library(dplyr) library(ggplot2) -library(OpenStreetMap) +library(maptiles) +library(tidyterra) library(viridis) library(httr2) @@ -35,32 +36,24 @@ shape.df$totalsecon <- as.numeric(shape.df$totalsecon) # filter out remote islands shape.df <- filter(shape.df,longitude > 100) # Obtain the contour box of the shape - lat <- c(min(shape.df$latitude), max(shape.df$latitude)) lon <- c(min(shape.df$longitude), max(shape.df$longitude)) - +# Define Perth metro contours lat.metro <- c(-32.1547, -31.8519) lon.metro <- c(115.7162,116.0733) -upperLeft.metro <- c(-31.8519,115.7162) -lowerRight.metro <- c(-32.1547,116.0733) - # Obtain an openstreetmap background image -map <- openmap(upperLeft = c(lat[2],lon[1]), - lowerRight = c(lat[1],lon[2]), - minNumTiles = 6, - type="https://stamen-tiles.a.ssl.fastly.net/toner-lite/{z}/{x}/{y}.png") -map.metro <- openmap(upperLeft = upperLeft.metro, - lowerRight = lowerRight.metro, - minNumTiles = 6, - type="https://stamen-tiles.a.ssl.fastly.net/toner-lite/{z}/{x}/{y}.png") - -mapLatLon <- openproj(map) +sfg_WA <- st_polygon(list(rbind(c(lon[1],lat[1]), + c(lon[2],lat[1]), + c(lon[2],lat[2]), + c(lon[1],lat[2]), + c(lon[1],lat[1]) + ))) +WATiles <- get_tiles(st_transform(st_sfc(sfg_WA, crs = 4326), 3857), provider = "Stadia.Stamen.TonerLite", crop = TRUE, zoom = 6) -OpenStreetMap::autoplot.OpenStreetMap(mapLatLon) + - geom_point(data = shape.df, - aes(x = longitude, - y = latitude, - color = totalschoo), +ggplot() + + geom_spatraster_rgb(data = WATiles) + + geom_sf(data = shape.df, + aes(color = totalschoo), show.legend = TRUE) + scale_color_viridis(option="viridis", name = "Number\nof pupils") + @@ -74,23 +67,28 @@ OpenStreetMap::autoplot.OpenStreetMap(mapLatLon) + y="") -# Zooming in on metro +# Zooming in on Perth metro shape.metro.df <- filter(shape.df, abs(latitude) <= max(abs(lat.metro)), abs(latitude) >= min(abs(lat.metro)), abs(longitude) <= max(abs(lon.metro)), abs(longitude) >= min(abs(lon.metro))) -mapLatLon.metro <- openproj(map.metro) +sfg_metro <- st_polygon(list(rbind(c(lon.metro[1],lat.metro[1]), + c(lon.metro[2],lat.metro[1]), + c(lon.metro[2],lat.metro[2]), + c(lon.metro[1],lat.metro[2]), + c(lon.metro[1],lat.metro[1]) + ))) +metroTiles <- get_tiles(st_transform(st_sfc(sfg_metro, crs = 4326), 3857), provider = "Stadia.Stamen.TonerLite", crop = TRUE, zoom = 11) # Let's look at secondary school only # 95th percentile of school size secon95th <- quantile(shape.metro.df[,"totalsecon">0]$totalsecon,probs = .95) -OpenStreetMap::autoplot.OpenStreetMap(mapLatLon.metro) + - geom_point(data = filter(shape.metro.df, totalsecon > 0), # geom_sf is not yet in ggplot2 on CRAN, so sticking to this for now - aes(x = longitude, - y = latitude, - color = totalsecon, - size = totalsecon) +ggplot() + + geom_spatraster_rgb(data = metroTiles) + + geom_sf(data = filter(shape.metro.df, totalsecon > 0), + aes(color = totalsecon, + size = totalsecon), ) + - geom_text(data = filter(shape.metro.df, totalsecon > secon95th), # name top 5% schools + geom_sf_label(data = filter(shape.metro.df, totalsecon > secon95th), # name top 5% schools aes(x = longitude, y = latitude, label = schoolname), -- 2.39.2 From 1ff9c846011a2edd5a45788d3833aba818db6b47 Mon Sep 17 00:00:00 2001 From: Frederik Vanrenterghem Date: Mon, 2 Oct 2023 22:40:12 +0800 Subject: [PATCH 4/5] Remove redundant list(). --- WA_active_schools.R | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/WA_active_schools.R b/WA_active_schools.R index 0811702..ae59595 100644 --- a/WA_active_schools.R +++ b/WA_active_schools.R @@ -11,11 +11,11 @@ tempdirSHP <- tempdir() tempfileSHP <- tempfile() req <- request("https://sso.slip.wa.gov.au/as/token.oauth2") |> req_headers("Authorization" = "Basic ZGlyZWN0LWRvd25sb2Fk") |> - req_body_form(list("grant_type" = "password", - # SLIP username and password stored in - #pass - the standard unix password manager - "username" = system2("pass", args = "slip.wa.gov.au | grep Username | sed -e 's/Username: //'", stdout = TRUE), - "password" = system2("pass", args = "slip.wa.gov.au | head -1", stdout = TRUE))) + req_body_form(grant_type = "password", + # SLIP username and password stored in + # pass - the standard unix password manager + username = system2("pass", args = "slip.wa.gov.au | grep Username | sed -e 's/Username: //'", stdout = TRUE), + password = system2("pass", args = "slip.wa.gov.au | head -1", stdout = TRUE)) tokenResponse <- req_perform(req) -- 2.39.2 From ff80b964779b1d6de4b3bc63bfbd8b8e9b698d00 Mon Sep 17 00:00:00 2001 From: Frederik Vanrenterghem Date: Mon, 9 Oct 2023 21:35:35 +0800 Subject: [PATCH 5/5] Use Carto for metro map. --- WA_active_schools.R | 44 +++++++++++++++++++++++--------------------- 1 file changed, 23 insertions(+), 21 deletions(-) diff --git a/WA_active_schools.R b/WA_active_schools.R index ae59595..7ca4a00 100644 --- a/WA_active_schools.R +++ b/WA_active_schools.R @@ -78,30 +78,32 @@ sfg_metro <- st_polygon(list(rbind(c(lon.metro[1],lat.metro[1]), c(lon.metro[1],lat.metro[2]), c(lon.metro[1],lat.metro[1]) ))) -metroTiles <- get_tiles(st_transform(st_sfc(sfg_metro, crs = 4326), 3857), provider = "Stadia.Stamen.TonerLite", crop = TRUE, zoom = 11) +metroTiles <- get_tiles(st_transform(st_sfc(sfg_metro, crs = 4326), 3857), provider = "CartoDB.PositronNoLabels", crop = TRUE, zoom = 11) # Let's look at secondary school only # 95th percentile of school size secon95th <- quantile(shape.metro.df[,"totalsecon">0]$totalsecon,probs = .95) ggplot() + - geom_spatraster_rgb(data = metroTiles) + - geom_sf(data = filter(shape.metro.df, totalsecon > 0), - aes(color = totalsecon, - size = totalsecon), - ) + - geom_sf_label(data = filter(shape.metro.df, totalsecon > secon95th), # name top 5% schools - aes(x = longitude, - y = latitude, - label = schoolname), - size = 2) + - scale_color_viridis(option="viridis", - name = "Number\nof pupils") + - guides(size = 'none') + # Do not show legend for size, as color is the main key - theme(legend.position = "right", - axis.text = element_blank(), - axis.ticks = element_blank()) + - labs(title = "Western Australian secondary schools", - subtitle = "Perth metro area", - caption = "Map tiles by Stamen Design, under CC BY 3.0.\nMapdata by OpenStreetMap, under ODbL. - School data © Government of WA 2022.\n Dep. of Education & Training.", + geom_spatraster_rgb(data = metroTiles) + + geom_sf(data = filter(shape.metro.df, totalsecon > 0), + aes(color = totalsecon, + size = totalsecon) + ) + + geom_sf_text(data = filter(shape.metro.df, totalsecon > secon95th), # name top 5% schools + aes(x = longitude, + y = latitude, + label = schoolname), + size = 2) + + scale_color_viridis(option="cividis", + name = "Number\nof pupils") + + guides(size = 'none') + # Do not show legend for size, as color is the main key + theme_minimal() + + theme(legend.position = "right", + panel.grid.major = element_line(colour = "white"), + axis.ticks = element_blank()) + + labs(title = "Western Australian secondary schools", + subtitle = "Perth metro area", + caption = "Map tiles © CARTO. Mapdata by OpenStreetMap, under ODbL.\nSchool data © Government of WA, Dep. of Education & Training.", x="", y="") + +ggsave("WA_schools.png", units = "cm", width = 15, height = 15, dpi = 300, bg = "white", scale = 1.5) -- 2.39.2