diff --git a/intro-spatial.Rmd b/intro-spatial.Rmd index 5599f5f..371b61c 100644 --- a/intro-spatial.Rmd +++ b/intro-spatial.Rmd @@ -56,7 +56,9 @@ experience with R will help. If you have not used R before, it may be worth following an introductory tutorial, such as *Efficient R Programming* -([Gillespie and Lovelace, 2016](http://cran.r-project.org/doc/contrib/Torfs+Brauer-Short-R-Intro.pdf)), the official [Introduction to R](https://cran.r-project.org/doc/manuals/r-release/R-intro.html) or tutorials suggested on [rstudio.com](https://www.rstudio.com/online-learning/) and [cran.r-project.org](https://cran.r-project.org/other-docs.html). +([Gillespie and Lovelace, 2016](https://csgillespie.github.io/efficientR/)), *R for Data Science* ([Grolemund and Wickham, 2016](http://r4ds.had.co.nz/)) or tutorials suggested on [rstudio.com](https://www.rstudio.com/online-learning/) and [cran.r-project.org](https://cran.r-project.org/other-docs.html). + +... Now you know some R, it's time to turn your attention towards spatial data with R. To that end, this tutorial is organised as follows: @@ -127,7 +129,7 @@ on and download some packages and data. R has a huge and growing number of spatial data packages. We recommend taking a quick browse on R's main website to see the spatial packages available: [http://cran.r-project.org/web/views/Spatial.html](http://cran.r-project.org/web/views/Spatial.html). -In this tutorial we will use the following packages: +In this tutorial we will use the packages from the '**sp**verse', that use the **sp** package: - **ggmap**: extends the plotting package **ggplot2** for maps - **rgdal**: R's interface to the popular C/C++ spatial data processing library [gdal](http://www.gdal.org/) @@ -136,6 +138,7 @@ In this tutorial we will use the following packages: - [**dplyr**](http://cran.r-project.org/web/packages/dplyr/index.html) and [**tidyr**](http://blog.rstudio.org/2014/07/22/introducing-tidyr/): fast and concise data manipulation packages - **tmap**: a new packages for rapidly creating beautiful maps +For a tutorial based on the recent [**sf**](https://github.com/edzer/sfr) package you will have to look elsewhere, like the [geocompr](https://bookdown.org/robinlovelace/geocompr/) website, the online home of the forthcoming book *Geocomputation with R*. Some packages may already be installed on your computer. To test if a package is installed, try to load it using the `library` function; for example, to test if **ggplot2** is installed, type `library(ggplot2)` into the console window. If there is no output from R, this is good news: it means that the library @@ -208,7 +211,7 @@ R to handle a broader range of spatial data formats. If you've not already ```{r, message=FALSE, results='hide'} library(rgdal) -lnd <- readOGR(dsn = "data/LondonBoroughs.shp") +lnd <- readOGR(dsn = "data/london_sport.shp") # lnd <- readOGR(dsn = "data", layer = "london_sport") ``` @@ -264,6 +267,12 @@ grid.raster(readPNG("figure/tab-complete.png")) To explore `lnd` object further, try typing `nrow(lnd)` (display number of rows) and record how many zones the dataset contains. You can also try `ncol(lnd)`. +One issue with the data that we have loaded is that it has no coordinate reference system (CRS): + +```{r} +lnd@proj4string +``` + ## Basic plotting Now we have seen something of the structure of spatial objects in R, @@ -282,8 +291,8 @@ an entirely different type of plot. Thus R is intelligent at guessing what you w R has powerful subsetting capabilities that can be accessed very concisely using square brackets,as shown in the following example: ```{r} -# select rows of lnd@data where sports participation is less than 15 -lnd@data[lnd$Partic_Per < 15, ] +# select rows of lnd@data where sports participation is less than 13 +lnd@data[lnd$Partic_Per < 13, 1:3] ``` The above line of code asked R to select only the rows from the `lnd` object, where sports participation is lower than 15, @@ -324,7 +333,7 @@ we will cover this in more detail in subsequent sections. As a bonus stage, select and plot only zones that are -close to the centre of London (see Fig. 6). Programming encourages rigorous +close to the centre of London (see Figure 5). Programming encourages rigorous thinking and it helps to define the problem more specifically: @@ -383,75 +392,24 @@ lnd$quadrant <- "unknown" # prevent NAs in result lnd$quadrant[east & north] <- "northeast" ``` -> **Challenge**: Based on the the above code as refrence try and find the remaining 3 quadrants and colour them as per Figure 6 below. +> **Challenge**: Based on the the above code as refrence try and find the remaining 3 quadrants and colour them as per Figure 6. Hint - you can use the **llgridlines** function in order to overlay the long-lat lines. For bonus points try to desolve the quadrants so the map is left with only 4 polygons. -```{r, echo=FALSE, eval=FALSE} +```{r, echo=FALSE, fig.cap="The 4 quadrants of London and dissolved borders. Challenge: recreate a plot that looks like this.", fig.show='hold', out.height="5cm"} lnd$quadrant[!east & north] <- "northwest" lnd$quadrant[east & !north] <- "southeast" lnd$quadrant[!east & !north] <- "southwest" -library(tmap) -qtm(lnd, fill = "quadrant") -lnd_disolved = rgeos::gUnaryUnion(spgeom = lnd, id = lnd$quadrant) -plot(lnd_disolved) -``` - - -```{r, echo=FALSE, fig.cap="The 4 quadrants of London"} -grid.raster(readPNG("figure/lnd-quads.png")) -par(mfrow = c(1,2)) -#plot the results -# plot(london) -# plot(lnd[east & north,],add = TRUE, col = "red" ) - -# place a grid over the object -# llgridlines(lnd, lty= 3, side ="EN", offset = -0.5) +plot(lnd) +plot(lnd[east & north,], add = TRUE, col = "red" ) +llgridlines(lnd, lty= 3, side ="EN", offset = -0.5) -london = gUnaryUnion(lnd, lnd$dummy) -london = SpatialPolygonsDataFrame(london, data.frame(dummy = c("london")), match.ID = FALSE) - -centrelondon = gCentroid(london, byid = TRUE) - -centreLEP = gCentroid(lnd, byid = TRUE) -coords <- cbind(centrelondon$y) - -c1 = c(centrelondon$x, centrelondon$x) -c2 = c(90, -90) -c3 = c(90, -90) -c4 = c(centrelondon$y,centrelondon$y) - -# simple line strings -L1 = Line(cbind(c1, c2)) -L2 = Line(cbind(c3, c4)) - -Ls1 = Lines(list(L1), ID = "a") -Ls2 = Lines(list(L2), ID = "b") - -Ls1 <- SpatialLines(LinesList = list(Ls1)) -Ls2 <- SpatialLines(LinesList = list(Ls2)) - -Longitude = SpatialLinesDataFrame(Ls1, data.frame(Z = c("1", "2"), row.names = c("a","b"))) -Latitude = SpatialLinesDataFrame(Ls2, data.frame(Z = c("1", "2"), row.names = c("a","b"))) - -east <- coordinates(lnd)[,1] > Longitude@lines[[1]]@Lines[[1]]@coords[,1][1] -north <- coordinates(lnd)[,2] > Latitude@lines[[1]]@Lines[[1]]@coords[,2][1] -west <- coordinates(lnd)[,1] < Longitude@lines[[1]]@Lines[[1]]@coords[,1][1] -south <-coordinates(lnd)[,2] < Latitude@lines[[1]]@Lines[[1]]@coords[,2][1] - -lnd@data$quadrant[east & north] <- "northeast" -lnd@data$quadrant[west & north] <- "northwest" -lnd@data$quadrant[east & south] <- "southeast" -lnd@data$quadrant[west & south] <- "southwest" - -# plot(lnd) -# plot(lnd[east & north,],add = TRUE, col = "red" ) -# plot(lnd[west & north,],add = TRUE, col = "blue" ) -# plot(lnd[east & south,],add = TRUE, col = "green" ) -# plot(lnd[west & south,],add = TRUE, col = "yellow" ) +lnd_disolved = rgeos::gUnaryUnion(spgeom = lnd, id = lnd$quadrant) -# llgridlines(lnd, lty= 3, side ="EN", offset = -0.5) -par(mfrow=c(1,1)) # return to default par +library(tmap) +qtm(lnd, fill = "quadrant") + + tm_shape(lnd_disolved) + + tm_borders(lwd = 9) ``` @@ -975,7 +933,7 @@ plot(stations) # test the clip succeeded - +