class: center middle section-title section-title-6 animated fadeIn # GIS in R with **sf** --- layout: true class: title title-6 --- # Shapefiles .box-6[Geographic information is shared as **shapefiles**] -- .box-inv-6[These are *not* like regular single CSV files!] -- .box-inv-6[Shapefiles come as zipped files with<br>a bunch of different files inside] .center[ <figure> <img src="img/03/shapefile-raw.png" alt="Shapefile folder structure" title="Shapefile folder structure" width="30%"> </figure> ] --- # Structure of a shapefile .small-code[ ```r library(sf) world_shapes <- read_sf("data/ne_110m_admin_0_countries/ne_110m_admin_0_countries.shp") ``` ``` ## Simple feature collection with 7 features and 3 fields ## Geometry type: MULTIPOLYGON ## Dimension: XY ## Bounding box: xmin: -180 ymin: -18 xmax: 180 ymax: 83 ## Geodetic CRS: WGS 84 ## # A tibble: 7 x 4 ## TYPE GEOUNIT ISO_A3 geometry ## <chr> <chr> <chr> <MULTIPOLYGON [°]> ## 1 Sovereign … Fiji FJI (((180 -16, 180 -17, 179 -17, 179 -17, 179 -17, 179 … ## 2 Sovereign … Tanzania TZA (((34 -0.95, 34 -1.1, 38 -3.1, 38 -3.7, 39 -4.7, 39 … ## 3 Indetermin… Western Sahara ESH (((-8.7 28, -8.7 28, -8.7 27, -8.7 26, -12 26, -12 2… ## 4 Sovereign … Canada CAN (((-123 49, -123 49, -125 50, -126 50, -127 51, -128… ## 5 Country United States … USA (((-123 49, -120 49, -117 49, -116 49, -113 49, -110… ## 6 Sovereign … Kazakhstan KAZ (((87 49, 87 49, 86 48, 86 47, 85 47, 83 47, 82 46, … ## 7 Sovereign … Uzbekistan UZB (((56 41, 56 45, 59 46, 59 46, 60 45, 61 44, 62 44, … ``` ] --- # Where to find shapefiles -- .box-inv-6[[Natural Earth](https://www.naturalearthdata.com/) for international maps] -- .box-inv-6.sp-after[[US Census Bureau](https://www.census.gov/geographies/mapping-files/time-series/geo/carto-boundary-file.html) for US maps] -- .box-inv-6[For anything else…] -- .center[ <figure> <img src="img/03/shapefile-search.png" alt="Search for shapefiles" title="Search for shapefiles" width="50%"> </figure> ] --- # Scales .pull-left-3[ <figure> <img src="img/03/download_thumbs_10m.jpg" alt="10m scale" title="10m scale" width="100%"> </figure> .box-inv-6.small[1:10m = 1:10,000,000] .box-inv-6.small[1 cm = 100 km] ] .pull-middle-3[ <figure> <img src="img/03/download_thumbs_50m.jpg" alt="50m scale" title="50m scale" width="100%"> </figure> .box-inv-6.small[1:50m = 1:50,000,000] .box-inv-6.small[ 1cm = 500 km] ] .pull-right-3[ <figure> <img src="img/03/download_thumbs_110m.jpg" alt="110m scale" title="110m scale" width="100%"> </figure> .box-inv-6.small[1:110m = 1:110,000,000] .box-inv-6.small[1 cm = 1,100 km] ] -- .box-inv-6[Using too high of a resolution<br>makes your maps slow and huge] --- # Latitude and longitude <img src="03_sf_files/figure-html/lat-long-example-1.png" width="504" style="display: block; margin: auto;" /> --- # The magic `geometry` column .box-inv-6[As long as you have a magic `geometry` column,<br>**all you need** to do to plot maps is `geom_sf()`] .left-code[ ```r ggplot() + geom_sf(data = world_shapes) ``` ] .right-plot[ ![](03_sf_files/figure-html/simple-map-1.png) ] --- # The magic `geometry` column .box-inv-6[Use `coord_sf()` to change projections] .left-code[ ```r ggplot() + geom_sf(data = world_shapes) + coord_sf(crs = "+proj=merc") ``` ] .right-plot[ ![](03_sf_files/figure-html/change-projection-1.png) ] --- # The magic `geometry` column .box-inv-6[Use `coord_sf()` to change projections] .left-code[ ```r ggplot() + geom_sf(data = world_shapes) + coord_sf(crs = "+proj=robin") ``` ] .right-plot[ ![](03_sf_files/figure-html/change-projection1-1.png) ] --- # Use aesthetics like normal .box-inv-6[All regular ggplot layers and aesthetics work] .left-code[ ```r ggplot() + geom_sf(data = world_shapes, aes(fill = POP_EST), color = "white", size = 0.15) + coord_sf(crs = "+proj=robin") + scale_fill_gradient(labels = scales::comma) + labs(fill = NULL) + theme_void() + theme(legend.position = "bottom") ``` ] .right-plot[ ![](03_sf_files/figure-html/add-aes-1.png) ] --- # No `geometry` column? .box-inv-6[Make your own with `st_as_sf()`] .pull-left-narrow.small-code[ ```r other_data ``` ``` ## # A tibble: 2 x 3 ## city long lat ## <chr> <dbl> <dbl> ## 1 Atlanta -84.4 33.8 ## 2 Washington, DC -77.1 38.9 ``` ] -- .pull-right-wide.small-code[ ```r other_data %>% st_as_sf(coords = c("long", "lat"), crs = st_crs("EPSG:4326")) ``` ``` ## Simple feature collection with 2 features and 1 field ## Geometry type: POINT ## Dimension: XY ## Bounding box: xmin: -84 ymin: 34 xmax: -77 ymax: 39 ## Geodetic CRS: WGS 84 ## # A tibble: 2 x 2 ## city geometry ## * <chr> <POINT [°]> ## 1 Atlanta (-84 34) ## 2 Washington, DC (-77 39) ``` ] --- # No `geometry` column? .box-inv-6[Automatically geocode addresses with the **tidygeocoder** package] .small-code[ ```r places ``` ``` ## # A tibble: 3 x 2 ## name address ## <chr> <chr> ## 1 My empty GSU office 14 Marietta Street NW, Atlanta, GA 30303 ## 2 My old BYU office 155 East 1230 North, Provo, UT 84604 ## 3 My old Duke office 201 Science Dr, Durham, NC 27708 ``` ] --- # No `geometry` column? .box-inv-6[Automatically geocode addresses with the **tidygeocoder** package] .pull-left-narrow.small-code[ ```r library(tidygeocoder) places %>% geocode( address, method = "census" ) %>% st_as_sf( coords = c("long", "lat"), crs = 4326 ) ``` ] -- .pull-right-wide.small-code[ ```r ## Simple feature collection with 3 features and 1 field ## Geometry type: POINT ## Dimension: XY ## Bounding box: xmin: -110 ymin: 34 xmax: -79 ymax: 40 ## Geodetic CRS: WGS 84 ## # A tibble: 3 x 2 ## name geometry ## <chr> <POINT [°]> ## 1 My empty GSU office (-84 34) ## 2 My old BYU office (-112 40) ## 3 My old Duke office (-79 36) ``` ] --- # **sf** is for all GIS stuff -- .box-inv-6[Draw maps] -- .box-inv-6[Calculate distances between points] -- .box-inv-6[Count observations in a given area] -- .box-inv-6.sp-after[Anything else related to geography!] -- .box-inv-6[See [here](https://bookdown.org/robinlovelace/geocompr/intro.html) or [here](https://bookdown.org/lexcomber/brunsdoncomber2e/Ch5.html) for full textbooks] --- # `geom_sf()` is today’s standard .box-inv-6[You'll sometimes find older tutorials and StackOverflow answers about using `geom_map()` or **ggmap** or other things] -- .box-inv-6.sp-after[Those still work, but they don't use the same magical **sf** system with easy-to-convert projections and other GIS stuff] -- .box-6.medium[Stick with **sf** and `geom_sf()`<br>and your life will be easy]