/
canmap.Rmd
132 lines (95 loc) · 4.46 KB
/
canmap.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
---
title: "O Canada"
description: Quick and dirty Canadian map example
author:
- name: Kieran Healy
url: http://kieranhealy.org
affiliation: Duke University
affiliation_url: http://duke.edu
date: "`r Sys.Date()`"
output: radix::radix_article
---
## Required Libraries, and a Map Theme
Mapping data (and, more importantly, converting map data into a format that we can tidily use) requires a bit more heavy lifting behind the scenes than the core `tidyverse` libraries provide. Start by installing the `sf` library, which will also bring a number of dependencies with it.
```{r installation, eval = FALSE, echo = TRUE}
install.packages("sf")
install.packages("rgdal")
install.packages("geojsonio")
install.packages("spdplyr")
install.packages("rmapshaper")
```
```{r setup}
knitr::opts_chunk$set(echo = TRUE)
library(geojsonio)
library(rmapshaper)
library(rgdal)
library(tidyverse)
library(spdplyr)
library(sf)
library(socviz)
## Map theme
theme_map <- function(base_size=9, base_family="") {
require(grid)
theme_bw(base_size=base_size, base_family=base_family) %+replace%
theme(axis.line=element_blank(),
axis.text=element_blank(),
axis.ticks=element_blank(),
axis.title=element_blank(),
panel.background=element_blank(),
panel.border=element_blank(),
panel.grid=element_blank(),
panel.spacing=unit(0, "lines"),
plot.background=element_blank(),
legend.justification = c(0,0),
legend.position = c(0,0)
)
}
theme_set(theme_map())
```
```{r figdir, echo = FALSE, results = 'hide'}
## Make a "figures" directory if one doesn't exist
ifelse(!dir.exists(file.path("figures")), dir.create(file.path("figures")), FALSE)
```
## Downloading and converting the map files for the first time
You won't be able to run the code in this section until you have some raw shapefiles to work with. Get data from Statistics Canada: http://www12.statcan.gc.ca/census-recensement/2011/geo/bound-limit/bound-limit-2011-eng.cfm. From the linked page, choose as your options _ArcGIS .shp file_, then---for example---_Census divisions_ and _cartographic boundaries_. You'll then download a zip file. Expand this zip file into a directory in your working folder named 'data'. Then import the shapefile as below. The `EVAL` option is set to `FALSE` on the chunks in this section, so they will not be evaluated if you run the document. (R will just skip to the next section below, instead.) To evaluate them, download the data as just described, and then remove the `eval = FALSE` option from each chunk.
```{r, eval = FALSE}
canada_raw <- readOGR(dsn = "data/gcd_000b11a_e", layer = "gcd_000b11a_e",
use_iconv=TRUE, encoding="CP1250")
```
Convert it to GeoJSON format and simplify the polygons. These steps will take a little while.
```{r, eval = FALSE}
canada_raw_json <- geojson_json(canada_raw)
canada_raw_sim <- ms_simplify(canada_raw_json)
```
Save the resulting GeoJSON file, which you can work with directly from here on
```{r, eval = FALSE}
geojson_write(canada_raw_sim, file = "data/canada_cd_sim.geojson")
```
## Working with the GeoJSON file in the tidyverse
Read the GeoJSON file back in as an `sf` object. The `.geojson` file is included in the repository, so you can load the libraries above and then start from here if you want.
```{r}
canada_cd <- st_read("data/canada_cd_sim.geojson", quiet = TRUE)
canada_cd
```
Transform the coordinates to a Lambert Conformal Conic Projection. (See https://www.statcan.gc.ca/pub/92-195-x/2011001/other-autre/mapproj-projcarte/m-c-eng.htm).
```{r}
canada_cd <- st_transform(canada_cd, crs = "+proj=lcc +lat_1=49 +lat_2=77 +lon_0=-91.52 +x_0=0 +y_0=0 +datum=NAD83 +units=m +no_defs")
canada_cd
```
Make a vector of repeated colors---just to fill in the map, for decoration only, as I don't have any Canadian data to merge in at present.
```{r, layout = 'l-screen-inset'}
map_colors <- RColorBrewer::brewer.pal(8, "Pastel1")
map_colors <- rep(map_colors, 37)
## Draw the map
p <- ggplot(data = canada_cd,
mapping = aes(fill = PRUID))
p_out <- p + geom_sf(color = "gray60",
size = 0.1) +
scale_fill_manual(values = map_colors) +
guides(fill = FALSE) +
theme_map() +
theme(panel.grid.major = element_line(color = "white"),
legend.key = element_rect(color = "gray40", size = 0.1))
ggsave("figures/canada.pdf", p_out, height = 12, width = 15)
p_out
```