Chapter 10 A Tidy Approach to Spatial Data: Simple Features
Note that if you, the student, wish to run all the code yourself, you should download the
inputs directory as a zipped file by going here with a web browser and then clicking the big “Download” button on the right. Once you have downloaded that, unzip it and put the whole
inputs directory in the current working directory where you are working with R.
If you don’t already have it:
Also, as of June 5, 2017, this requires the development version of
Be warned that this seems to break some of the earlier ggmap.
## Linking to GEOS 3.4.2, GDAL 2.1.2, proj.4 4.9.1
10.2 Reading spatial data into sf objects
Let’s just quickly put the code in to reproduce most of the stuff from the previous lecture, but using
Note that there are two functions to read shapefiles:
read_sf(). The latter seems preferable to me as it does not read stringsAsFactors ever.
streams <- read_sf("inputs/santa-cruz-county/Streams/Streams.shp") watersheds <- read_sf("inputs/santa-cruz-county/Watersheds/Watersheds.shp") lights <- read_sf("inputs/santa-cruz-county/Street_Lights/Street_Lights.shp")
Have a look at one of them. Note that the whole thing is a data frame, so you can directly
## # A tibble: 1,514 x 9 ## OBJECTID STREAM_NM STREAMTYPE GNIS_ID ID ENABLED Named ## <dbl> <chr> <chr> <chr> <chr> <chr> <chr> ## 1 1 Adams Creek INTERMITTENT <NA> 0 <NA> Yes ## 2 2 Alba Creek PERENNIAL 00218063 1 <NA> Yes ## 3 3 Amaya Creek PERENNIAL 00218217 2 <NA> Yes ## 4 4 Ano Nuevo Creek PERENNIAL 00254567 3 <NA> Yes ## 5 5 Aptos Creek PERENNIAL 00254571 4 <NA> Yes ## 6 6 Archibald Creek INTERMITTENT 00218350 5 <NA> Yes ## 7 7 Baldwin Creek PERENNIAL 00218614 6 <NA> Yes ## 8 8 Bates Creek PERENNIAL 00218726 7 <NA> Yes ## 9 9 Bean Creek PERENNIAL 00218782 8 <NA> Yes ## 10 10 Bennett Creek PERENNIAL 00219040 9 <NA> Yes ## # ... with 1,504 more rows, and 2 more variables: SHAPESTLen <dbl>, ## # geometry <simple_feature>
10.3 Plotting sf objects
In ggplot this is made easy with the
geom_sf() which knows what type of simple feature you have, and plots it appropriately.
You also add
coord_sf() which will deal with projections, on the fly, and also takes care of aspect ratios and draws a graticule.
ggplot() + geom_sf(mapping = aes(colour = STREAMTYPE), data = streams) + coord_sf()
Check it out, there is a graticule. Note that the long-lat lines are not perfectly square. That is because this object is not in long-lat:
## $epsg ##  NA ## ## $proj4string ##  "+proj=lcc +lat_1=37.06666666666667 +lat_2=38.43333333333333 +lat_0=36.5 +lon_0=-120.5 +x_0=2000000 +y_0=500000.0000000002 +datum=NAD83 +units=us-ft +no_defs" ## ## attr(,"class") ##  "crs"
Aha! See, it is in
lcc = lambert conic conformal.
What happens if we try to put an open street maps on this? It fails…
How about if we try to use spraster_rgb to put a raster background in there? That is a little complicated. So, instead let’s have a look at interactive mapping later.
But before that I want to briefly touch on geometrical intersections.
10.4 Geometric operations
Just briefly, remember, we wondered last week how to find all the streams that are in the Scott Creek watershed. This can be done with the
raster package, nicely. But with the
sf package we can also use
st_intersection() to spatially intersect streams with watersheds. When we do this, the cool thing is that all the columns of data associated with each get joined just as you would hope they do!!
streamsheds <- st_intersection(streams, watersheds)
## Warning: attribute variables are assumed to be spatially constant ## throughout all geometries
Those columns get joined on the basis of the stream being in the watershed. Cool!
Check it out:
## # A tibble: 1,522 x 15 ## OBJECTID STREAM_NM STREAMTYPE GNIS_ID ID ENABLED Named ## * <dbl> <chr> <chr> <chr> <chr> <chr> <chr> ## 1 4 Ano Nuevo Creek PERENNIAL 00254567 3 <NA> Yes ## 2 25 Cascade Creek PERENNIAL 00220660 24 <NA> Yes ## 3 38 Elliot Creek PERENNIAL 00223153 37 <NA> Yes ## 4 41 Finney Creek PERENNIAL 00223512 40 <NA> Yes ## 5 47 Gazos Creek PERENNIAL <NA> 46 <NA> Yes ## 6 51 Green Oaks Creek PERENNIAL 00224554 50 <NA> Yes ## 7 91 Old Womans Creek PERENNIAL 00229991 90 <NA> Yes ## 8 129 White House Creek PERENNIAL <NA> 128 <NA> Yes ## 9 131 Willow Gulch INTERMITTENT <NA> 130 <NA> Yes ## 10 162 Stream 1141 SWALE <NA> 161 <NA> No ## # ... with 1,512 more rows, and 8 more variables: SHAPESTLen <dbl>, ## # OBJECTID.1 <dbl>, NAME <chr>, COUNT_ <dbl>, MapKey <chr>, ## # SHAPESTAre <dbl>, SHAPESTLen.1 <dbl>, geometry <simple_feature>
Plot it with different colors for different watersheds to see what I mean:
ggplot(streamsheds) + geom_sf(aes(colour = NAME)) + coord_sf()
OK, now let’s make some interactive maps.
10.5 Interactive plotting with leaflet via the
This is truly mind-blowing. The authors of the
mapview package tell us:
mapview is an R package created to help researchers during their spatial data analysis workflow. It provides functions to very quickly and conveniently create interactive visualisations of spatial data. It was created to fill the gap of quick (not presentation grade) interactive plotting to examine and visually investigate both aspects of spatial data, the geometries and their attributes.”
As such, it appears that customization of certain things might not be easily-accessible, but it is so easy and still so customizable that it is really, really amazing.
Get it from CRAN:
and then be sure to load it:
## Loading required package: leaflet
The main function is
mapView which does multiple dispatch on a huge variety of different object types (rasters, points, polygons, etc), including
Let’s look at it in action: