---
title: "Geospatial Data<br>Geospatial Mapping"
author: "Keith VanderLinden<br>Calvin University"
output:
  xaringan::moon_reader:
    nature:
      titleSlideClass: ["left", "top"]
      highlightLines: true
      ratio: "16:9"     
---

```{r child = "../../resources/scripts/slide-setup.Rmd"}
```

```{r packages, echo = FALSE, warning = FALSE, message = FALSE}
knitr::opts_chunk$set(eval = TRUE)
options(dplyr.print_max = 12)
```

# Geospatial Datasets

.pull-left[
`spData::world` is an R list that represents an simple-feature (`sf`) collection including:

- Bounding box 
- Coordinate Reference System
- Data for each country, including vector geometry specifications encoded as multi-polygons

```{r world-sf, results="hide", fig.show = "hide", message = FALSE, warning = FALSE}
library(tidyverse)
library(sf)
library(spData)
```
]
.pull-right[
```{r}
spData::world
```
]

???
- `spData::world` is an R list that represents an simple-feature (`sf`) collection of 177 countries, 10 fields per country. 
- Bounding box is most of the globe: -180--180 longitude and -89.99--83.59 latitude
- CRS: the WGS84 coordinate reference system (see slide below)
- Cf. the world data from `rnaturalearth`, which includes a dataset with 241 countries, 63 fields per country.

---
# Mapping Geospatial Data: The World

`spData::world` contains geometric data on all countries. Such data are most compelling when presented in the form of a *map*.


.pull-left[
`ggplot` & `sf` work together to plot this data easily.

```{r world-plot, results="hide", fig.show = "hide", message = FALSE, warning = FALSE}
spData::world %>%
 ggplot() +
 geom_sf()
```
]
.pull-right[
```{r ref.label = "world-plot", echo = FALSE, message = FALSE, warning = FALSE, output.lines = 10, out.width = "100%"}
```
]

???

---
# Mapping Geospatial Data: One Country

.pull-left[

```{r}
spData::world$name_long[2]
spData::world$geom[2]
```

```
MULTIPOLYGON (((33.90371 -0.95, 31.86617 -1.027...
```
]
.pull-right[
```{r us-plot, results="hide", fig.show = "hide", message = FALSE, warning = FALSE}
world$geom[2] %>%
  ggplot() +
  geom_sf() +
  geom_sf_text(
    aes(label = world$name_long[2]))
```
```{r ref.label = "us-plot", echo = FALSE, message = FALSE, warning = FALSE, output.lines = 20, out.width = "100%"}
```
]

???
- Here we focus on the 2nd country in the dataset: Tanzania.
- The coordinates are x (longitude); y (latitude).
- The geometry is encoded using *Well-Known Text* (WKT) format, which supports points, lines, polygons, and multi-* versions of all three. (N.b., I had to cut-paste the geometry text to get it onto the slide; it wouldn't display otherwise.)

---
# Coordinate Reference Systems

A *Coordinate Reference System* (CRS) specifies locations on the surface of the Earth using angular *coordinates*. We focus on *geoditic* coordinate systems, which use *latitude* and *longitude*.

.pull-left[
![GCS](https://upload.wikimedia.org/wikipedia/commons/thumb/e/ef/FedStats_Lat_long.svg/390px-FedStats_Lat_long.svg.png)
]
.pull-right[
Common CRS systems include:
- [WGS84](http://wiki.gis.com/wiki/index.php/WGS84) (aka EPSG4326): Used by GPS 
- [NAD83](http://wiki.gis.com/wiki/index.php/NAD83): North American system 
- OSGB: UK standard used for the Cholera data, which created issues discussed in the text
]

.footnote[Image from: [wikimedia.org](https://en.wikipedia.org/wiki/Geographic_coordinate_system)]

???
- Geodetic CRSs use angular measurements. 
  - Latitude: -90 -- +90 (N.b., That's 90 degrees from equator to pole.)
  - Longitude: -180 -- +180 (N.b., That's 360 degrees around the equitorial circle)
- References:
  - https://geocompr.robinlovelace.net/spatial-class.html#crs-intro
  - https://datacarpentry.org/organization-geospatial/03-crs/index.html
  - https://www.esri.com/arcgis-blog/products/arcgis-pro/education/earth-peel/
  - https://www.geo-projections.com/

---
# Projections

*Projections* map CRS coordinates onto planar coordinates, projecting them onto either a cylinder, cone, or plane. This preserves Shape/Angle (*conformal*) **XOR** Area (*equal-area*).

![Map Projections](https://cs.calvin.edu/courses/info/601/resources/images/figure2-projections.png)

.footnote[Image from: [Map Projections](https://gistbok.ucgis.org/bok-topics/2017-quarter-03/map-projections)]

???
- A globe would be more accurate, but it can't be included in books or shown on computer screens.
- x-y projections of lat-long coordinates are, necessarily, deformed in some way.
(see: [Orange Peel Analogy](https://datacarpentry.org/organization-geospatial/03-crs)).
- Choose a projection that suits your usecase (cf. MDSR 17.5: How (not) to lie).

---
# Projections: Examples

```{r}
library(mapproj)
library(maps)
```

.pull-left[
**Mercator**
```{r message = FALSE, warning = FALSE}
maps::map(
  "world", 
  projection = "mercator", 
  wrap = TRUE
  )
```
]
.pull-right[
**Gall-Peters**
```{r message = FALSE, warning = FALSE}
maps::map(
  "world", 
  projection = "cylequalarea", 
  param = 0,  
  wrap = TRUE)
```
]

???
- Libraries:
  - `mapproj` does map projections.
  - `maps` draws the maps.
- mercator()
  - Equatorial projections centered on the Prime Meridian (longitude 0). Parallels are straight horizontal lines.
  - equally spaced straight meridians, **conformal**, straight compass courses
  - cylindrical
- cylequalarea(lat0)
  - aka. Gall-Peters (cited in the text)
  - equally spaced straight meridians, **equal-area**, true scale on lat0
  - **cylindrical** (centered on param = lat0)

---
# Projections: Examples (continued)

```{r}
library(mapproj)
library(maps)
```

.pull-left[
**Lambert**
```{r message = FALSE, warning = FALSE}
maps::map(
  "world", 
  projection = "lambert",   
  parameters = c(lat0 = 0, lat1 = 90), 
  wrap = TRUE)
```
]
.pull-right[
**Albers**
```{r message = FALSE, warning = FALSE}
maps::map(
  "world", 
  projection = "albers", 
  parameters = c(lat0 = 20, lat1 = 50), 
  wrap = TRUE)
```
]

???
- lambert(lat0,lat1)
  - **conformal**, true scale on lat0 and lat1
  - **Conic**
  - The parameters specify where the cylider hits the Earth (here, the cone point is over Greenwich).
- albers(lat0,lat1)
  - **equal-area**, true scale on lat0 and lat1
  - **Conic**
  - Not great for the world, but decent for smaller regions (e.g., "state")  
