Applied Spatial Data Analysis with R

Roger S. Bivand, Edzer J. Pebesma, Virgilio Gómez-Rubio

Mentioned 1

This book addresses the needs of researchers and students using R to analyze spatial data across a range of disciplines and professions. The book is co-authored by a group involved in the Comprehensive R Archive Network.

More on

Mentioned in questions and answers.

I have written a kernel density estimator in Java that takes input in the form of ESRI shapefiles and outputs a GeoTIFF image of the estimated surface. To test this module I need an example shapefile, and for whatever reason I have been told to retrieve one from the sample data included in R. Problem is that none of the sample data is a shapefile...

So I'm trying to use the shapefiles package's funciton to convert the bei dataset included in the spatstat package in R to a shapefile. Unfortunately this is proving to be harder than I thought. Does anyone have any experience in doing this? If you'd be so kind as to lend me a hand here I'd greatly appreciate it.

Thanks, Ryan

References: spatstat, shapefiles

A general solution is:

  • convert the "ppp" or "owin" classed objects to appropriate classed objects from the sp package
  • use the writeOGR() function from package rgdal to write the Shapefile out

For example, if we consider the hamster data set from spatstat:


first convert this object to a SpatialPointsDataFrame object:

ham.sp <- as.SpatialPointsDataFrame.ppp(hamster)

This gives us a sp object to work from:

> str(ham.sp, max = 2)
Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
  ..@ data       :'data.frame': 303 obs. of  1 variable:
  ..@ coords.nrs : num(0) 
  ..@ coords     : num [1:303, 1:2] 6 10.8 25.8 26.8 32.5 ...
  .. ..- attr(*, "dimnames")=List of 2
  ..@ bbox       : num [1:2, 1:2] 0 0 250 250
  .. ..- attr(*, "dimnames")=List of 2
  ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slots

This object has a single variable in the @data slot:

> head(ham.sp@data)
1 dividing
2 dividing
3 dividing
4 dividing
5 dividing
6 dividing

So say we now want to write out this variable as an ESRI Shapefile, we use writeOGR()

writeOGR(ham.sp, "hamster", "marks", driver = "ESRI Shapefile")

This will create several files in directory hamster created in the current working directory. That set of files is the ShapeFile.

One of the reasons why I didn't do the above with the bei data set is that it doesn't contain any data and thus we can't coerce it to a SpatialPointsDataFrame object. There are data we could use, in bei.extra (loaded at same time as bei), but these extra data or on a regular grid. So we'd have to

  • convert bei.extra to a SpatialGridDataFrame object (say bei.spg)
  • convert bei to a SpatialPoints object (say bei.sp)
  • overlay() the bei.sp points on to the bei.spg grid, yielding values from the grid for each of the points in bei
  • that should give us a SpatialPointsDataFrame that can be written out using writeOGR() as above

As you see, that is a bit more involved just to give you a Shapefile. Will the hamster data example I show suffice? If not, I can hunt out my Bivand et al tomorrow and run through the steps for bei.

Realated tags