Skip to content

Commit

Permalink
Update intro comments, mention sf
Browse files Browse the repository at this point in the history
  • Loading branch information
Robinlovelace committed May 18, 2017
1 parent 1b03211 commit e13a80a
Showing 1 changed file with 26 additions and 68 deletions.
94 changes: 26 additions & 68 deletions intro-spatial.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -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:

Expand Down Expand Up @@ -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/)
Expand All @@ -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
Expand Down Expand Up @@ -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")
```

Expand Down Expand Up @@ -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,
Expand All @@ -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,
Expand Down Expand Up @@ -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:

Expand Down Expand Up @@ -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)
```

<!-- ## Attribute data -->
Expand Down Expand Up @@ -975,7 +933,7 @@ plot(stations) # test the clip succeeded
<!-- ``` -->

<!-- This results in a simple choropleth map and a new vector containing the area of each -->
<!-- borough (the basis for Fig. 11). As an additional step, try comparing the mean -->
<!-- borough (the basis for Figure 11). As an additional step, try comparing the mean -->
<!-- area of each borough with the -->
<!-- mean value of `stations` points within it: `plot(lnd_n$NUMBER, areas)`. -->

Expand Down

0 comments on commit e13a80a

Please sign in to comment.