Getting error while using fortify function in R (ggplot)

Question

I am trying to use the function fortify in R to plot the map. But I am always getting an error while processing the data. I read the shape file using the following code:

shape <- readShapeSpatial("Lond_City/lond_city.shp")

The shapefile can be found on the following dropbox location: https://www.dropbox.com/sh/d4w6saailydtu1r/5oIa56xV6S

Then I tried to use the following command :

shape1 <- fortify(shape,region="data")

In the above code, what should I put in Place of "data". This is the place which seems confusing. I can open the data or shapefile in R. However, when I run the fortify line I get the following error:

Error in '[.data.frame'(attr, , region) : undefined columns selected

Can anybody tell me how to properly fortify a shape file with example ? The reason I want to use fortify is that I want to combine the shape file with the data points.

Thank you so much.

Jdbaba

Solution 1

Your problem is that you need to give fortify a column that actually exists. If you type str(lon.shape) from my example code below, you will see that there is no region column:

> str(lon.shape)
Formal class 'SpatialPolygonsDataFrame' [package "sp"] with 5 slots
  ..@ data       :'data.frame': 1 obs. of  7 variables:
  .. ..$ ons_label: chr "00AA"
  .. ..$ name     : chr "City of London"
  .. ..$ label    : chr "02AA"
  .. ..$ X_max    : num 533843
  .. ..$ y_max    : num 182198
  .. ..$ X_min    : num 0
  .. ..$ y_min    : num 0
  ..@ polygons   :List of 1
  .. ..$ :Formal class 'Polygons' [package "sp"] with 5 slots
  .. .. .. ..@ Polygons :List of 1
  .. .. .. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
  .. .. .. .. .. .. ..@ labpt  : num [1:2] 532464 181220
  .. .. .. .. .. .. ..@ area   : num 3151465
  .. .. .. .. .. .. ..@ hole   : logi FALSE
  .. .. .. .. .. .. ..@ ringDir: int 1
  .. .. .. .. .. .. ..@ coords : num [1:771, 1:2] 531027 531029 531036 531074 531107 ...
  .. .. .. ..@ plotOrder: int 1
  .. .. .. ..@ labpt    : num [1:2] 532464 181220
  .. .. .. ..@ ID       : chr "0"
  .. .. .. ..@ area     : num 3151465
  ..@ plotOrder  : int 1
  ..@ bbox       : num [1:2, 1:2] 530967 180404 533843 182198
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:2] "x" "y"
  .. .. ..$ : chr [1:2] "min" "max"
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots
  .. .. ..@ projargs: chr "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +units=m +no_defs" 

Instead, try using name, like this:

lon.df <- fortify(lon.shape, region = "name")

Which produces this data frame:

> head(lon.df)
      long      lat order  hole piece            group             id
1 531026.9 181611.1     1 FALSE     1 City of London.1 City of London
2 531028.5 181611.2     2 FALSE     1 City of London.1 City of London
3 531036.1 181611.5     3 FALSE     1 City of London.1 City of London
4 531074.0 181610.3     4 FALSE     1 City of London.1 City of London
5 531107.0 181609.3     5 FALSE     1 City of London.1 City of London
6 531117.1 181608.9     6 FALSE     1 City of London.1 City of London

Here's one way to go about it, start to finish:

library(rgdal)
library(ggplot2)
library(rgeos)

shape.dir <- "c:/test/london" # use your directory name here

lon.shape <- readOGR(shape.dir, layer = "lond_city")
lon.df <- fortify(lon.shape, region = "name")

ggplot(lon.df, aes(x = long, y = lat, group = group)) +
    geom_polygon(colour = "black", fill = "grey80", size = 1) +
    theme()

...resulting in the following:

Solution 2

You should put label in place of data.

require(maptools); require(ggplot2)
shape <- readShapeSpatial("Lond_City/lond_city.shp")
shape1 <- fortify(shape,region="label")

To find this out I inspected the data element of shape:

> shape@data
  ons_label           name label    X_max    y_max X_min y_min
0      00AA City of London  02AA 533842.7 182198.4     0     0

This shows ons_label or name would also have worked and may be more suitable for your purposes.

Note that your shape file is slightly unusual because there seems to be only one group.

Comments

[Jd Baba]

SlowLearner, Thank you so much for your wonderful solution. But it seems like I cannot run because I get the following error "isTRUE(gpclibPermitStatus()) is not TRUE. Do you know how to fix this ? I tried to install "gpclibPermit" but it says this package is not available for version 2.15.2. Any suggestions ?

[orizon]

@Jdbaba This is because SlowLearner is using readOGR rather than readShapeSpatial. You can fix this by running gpclibPermit() but be careful as this is forbidden if you are using the software for a commercial purpose. You may prefer to avoid this issue by adopting the approach in my answer. Alternatively you can install the rgeos package which has no commercial restriction and should also fix the problem.

[SlowLearner]

@Jdbaba readOGR preserves the proj4 information so I generally prefer it. I have edited my answer a little to include gpclib and the gpclibPermit() call.

[orizon]

@SlowLearner If you include library(rgeos) you do not need to run gpclibPermit() as well, since rgeos is an alternative to gpclib.

[SlowLearner]

@orizon Ah, thanks for that, guess that's why in another piece of code I have there is an unexplained call to rgeos! I get confused as to which libraries supply what: sp, map, maptools, gpclib, rgdal, rgeos all get mixed up together in my head. Will edit.

[Jd Baba]

Thanks guys for clearing things up. I have another task to do. I need to add the points and contour to the plot and clip the contour by the shapefile. Is that possible in R ? Thanks.

[SlowLearner]

@Jdbaba On StackOverflow we like to do things "one question at a time" so please make this a separate question and include the data (or fake data that resembles the real data).

[SlowLearner]

@Jdbaba Yes it is possible, search for answers on gis.stackexchange.com as well as StackOverflow.