3.1 Primary Health Networks
Although the absmapsdata
package does not include primary health networks (PHNs) boundaries, the ABS does offer correspondence tables between 2017 PHNs and several ABS structures including SA2s, SA3s, and Postcodes.
An interactive map showing PHN boundaries in Australia and the corresponding SA3 boundaries can be found here. Using this map, it seems like the PHNs are constructed from one ore more whole SA3s. However, the correspondence isn't perfect. In fact, the correspondence between PHNs and any of these other building blocks (SA1s, SA2s, postcodes) isn't perfect. (This makes sense since PHNs are not defined in the ASGS.)
We can confirm this by exploring the ratio
(degree of correspondence) variable in each PHN correspondence table.
# Concordances to PHNs
conc.phn.poa <- get_correspondence_absmaps("POA", 2016, "PHN", 2017)
conc.phn.ced <- get_correspondence_absmaps("CED", 2016, "PHN", 2017)
conc.phn.sa1 <- get_correspondence_absmaps("SA1", 2016, "PHN", 2017)
conc.phn.sa2 <- get_correspondence_absmaps("SA2", 2016, "PHN", 2017)
conc.phn.sa3 <- get_correspondence_absmaps("SA3", 2016, "PHN", 2017)
# Compute percentage of perfect matches (ratio = 1)
percent_perfect_match <- function(x){
tab <- table(x$ratio)
ppm <- tab[length(tab)] / nrow(x)
return(as.numeric(ppm))
}
c(percent_perfect_match(conc.phn.poa),
percent_perfect_match(conc.phn.ced),
percent_perfect_match(conc.phn.sa1),
percent_perfect_match(conc.phn.sa2),
percent_perfect_match(conc.phn.sa3)
)
## [1] 0.8803659 0.4485981 0.9990437 0.9648521 0.7827225
We see that the best correspondence as a percentage of all units of a given structure is provided by SA1s.
From the PHN website:
Where possible, boundaries of the PHNs align with Local Hospital Networks (LHNs) or equivalents, or clusters of LHNs. ...In determining boundaries, a number of factors were taken into account, including population size and future projected population growth, LHN alignment, State and Territory borders, patient flows and administrative efficiencies.
So without explicit definitions of the PHN boundaries, using a concordance between SA1s and PHNs to merge the SA1s into PHN areas probably provides the best approximation. E.g.
# Remove ambiguities in concordance (use highest ratio when not 1)
conc <- conc.phn.sa1 %>%
filter(ratio != 1)
codes <- unique(conc$sa1_maincode_2016)
conc.phn.resolved <- conc.phn.sa1
for(i in length(codes):1){
# rows with ambiguities
r <- which(conc.phn.sa1$sa1_maincode_2016 == codes[i])
# row(s) to drop
d <- r[-which.max(conc.phn.sa1$ratio[r])]
if(length(d) > 0){
conc.phn.resolved <- conc.phn.resolved[-d,]
}
}
# Join PHN concordance then merge (dissolve) SA1s by PHN
names(conc.phn.resolved)[1] <- names(sa12016)[1]
phn2017 <- left_join(sa12016, conc.phn.resolved) %>%
group_by(phn_code_2017) %>%
summarize() # Works like gUnion on SpatialPolygons
## Joining, by = "sa1_main_2016"
Here is a plot of the PHNs:
gg.base +
geom_sf(data = phn2017, aes(geometry = geometry, fill = phn_code_2017)) +
scale_fill_discrete(guide = FALSE) +
ggtitle("PHNs")
Upon inspection, this looks identical to the map shown on the PHN website.
To export this sf
to file, use st_write
. There may be a better driver, but the MapInfo File
format preserves the object properties well (ESRI Shapefile
does not, and GPKG
fails).
st_write(phn2017, dsn = fp.wd, layer = "PHN Boundaries", driver = "MapInfo File")
where fp.wd
is the filepath of the working directory (or any valid directory you like). To read it in again, use st_read
:
phn2017.test <- st_read(dsn = fp.wd, layer = "PHN Boundaries")