Drawing lines in raster cells, species distribution maps

When we want to look at the variation of species distribution over times an easy way to grasp the change is to draw lines with a certain angles on maps and attribute to different species or time period a certain angle. Here are some code to achieve this in R (thanks to R. Hijmans help)

library(sp)
library(raster)
## raster 2.0-41 (21-December-2012)
library(rgeos)
## rgeos version: 0.2.11, (SVN revision 370) GEOS runtime version:
## 3.2.2-CAPI-1.6.2 Polygon checking: TRUE
# create two raster with a certain distribution of 1's
r1 <- raster(ncols = 50, nrows = 50)
values(r1) <- rep(0, 2500)
r1[unlist(lapply(seq(20, 1220, 50), FUN = function(x) {
    seq(x, x + 20, 1)
}))] <- 1

r2 <- raster(ncols = 50, nrows = 50)
values(r2) <- rep(0, 2500)
r2[unlist(lapply(seq(920, 1920, 50), FUN = function(x) {
    seq(x, x + 20, 1)
}))] <- 1
par(mfrow = c(1, 2))
plot(r1, axes = FALSE, legend = FALSE)
plot(r2, axes = FALSE, legend = FALSE)

plot of chunk unnamed-chunk-2

# to look at the change between the two raster you can use calc
r_comp <- overlay(r1, r2, fun = function(x, y) {
    x - y
})
plot(r_comp, axes = FALSE, legend = FALSE)

plot of chunk unnamed-chunk-4

# the green area at the top are the places where there was 1 in the first
# layer but not in the second, this is the opposite for the white area at
# the bottom however between the two this yellow area cannot be
# distinguish between: 1 in both layer, 0 in both layer. so we can use
# lines at different angles to represent the values from the two layers
p1 <- rasterToPolygons(r1, fun = function(x) {
    x == 1
}, dissolve = TRUE)
p2 <- rasterToPolygons(r2, fun = function(x) {
    x == 1
}, dissolve = TRUE)
# create an empty raster for easier comparaison
r <- raster(ncols = 50, nrows = 50)
r[] <- rep(0, 2500)
par(mfrow = c(1, 2))
plot(r_comp, axes = FALSE, legend = FALSE)
plot(r, axes = FALSE, legend = FALSE, col = "white")
plot(p1, density = 20, add = TRUE, border = "white", angle = 45)
plot(p2, density = 20, add = TRUE, border = "white", angle = -45)

plot of chunk unnamed-chunk-6

With the lines we can really see what is happening, to come back to our species example, this plot could be a northward shift in the species distribution with only a little area left unchanged.

Happy mapping