[R-sig-Geo] Still a bug when converting from tess to spatial polygon (sp) format?

Fabricio Vasselai fabriciovasselai at gmail.com
Thu May 15 03:05:59 CEST 2014


Hi professor Bivand, thanks for the fast reply.

Sure, there is no surprise and nothing wrong with the objects
themselves: the converted tess is a SpatialPolygon, not a
SpatialPolygonDataFrame. Therefore, it is normally missing the data
slot. Anyway, although I got interested in that application due to a
previous question someone sent to the list, I can think of uses for
tess objects to have attributes.  For instance, depending on which
type of simulation is being executed one can have a few attributes
being randomly assigned to each tile of the tessellated map, or when
someone wants to save the random generated attributes of a simulation
within the generated maps into a .shp file, maybe those could be
easier and neater to achieve with a SpatialPolygonDataFrame. Sure,
there are other simple enough ways to achieve the same, but still, why
not having the data slot within the converted tess?

Here I send two codes to convert tess objects to SpatialPolygonDataFrame.
The first uses as(, "SpatialPolygons"), the second does all the
conversion step by step. I don't know if either are enough to use for
a coercion method, but it seems to work fine as an user-defined
formula. Any further simplification and corrections are certainly
welcome.

Version 1:

#Code begins here

tess2SPdf <- function(x) {
  stopifnot(is.tess(x))
  require(spatstat)
sp.obj <- as(x, "SpatialPolygons")
d <- data.frame()
d <- d[seq(x$n),]
row.names(d) <- names(x$tiles)
return(SpatialPolygonsDataFrame(sp.obj, data=d))
}

#End of the code


Or if someone needs an entire piece, that does not use the coercion
from the as() function, an option based on past emails from this list
would be:


Version 2:

#Code begins here

tess2SPdf <- function(x) {
  stopifnot(is.tess(x))
  require(sp)
  require(spatstat)
  a <- tiles(x)
  tess.labels <- names(a)
  c <- list()
  for(i in seq(a)){

    b <- as.polygonal(a[[i]])
    closering <- function(df) {
                    df[c(seq(nrow(df)), 1), ]
                 }
    pieces <- lapply(b$bdry,
                     function(p) {
                       Polygon(coords=closering(cbind(p$x,p$y)),
                               hole=is.hole.xypolygon(p))
                     })
    c[[i]] <- Polygons(pieces, tess.labels[i])

  }

  d <- data.frame()
  d <- d[seq(x$n),]
  row.names(d) <- names(x$tiles)
  return(SpatialPolygonsDataFrame(SpatialPolygons(c), data=d))
}
#End of the code


Hope it helps.

Best,

Fabricio


2014-05-14 3:21 GMT-04:00 Roger Bivand <Roger.Bivand at nhh.no>:
>
> On Wed, 14 May 2014, Fabricio Vasselai wrote:
>
>> Dear all,
>>
>>
>> After playing around with the tessellation functions of the spatstat
>> package to reply to a previous question on the list, I realized that there
>> seems to be something odd when converting from tess format to spatial
>> polygon format.
>>
>> Usually, spatial polygonal shapefiles read into R with functions like
>> readShapePoly or readOGR, are objects which have a slot called "data", that
>> can be assessed with @data and has the format of a data.frame. In there,
>> each polygonal feature of the map is a row, while the columns are the
>> variables that are eventually present in the shapefiles's data. Therefore,
>> if you try to check names() of a spatial polygon object read from a
>> shapefile, you get these columns, while row() gives you the polygon
>> features present in the map.
>
>
> Why is this a surprise? There are *two* classes for each Spatial* object, the geometries (SpatialPoints, SpatialLines, SpatialPolygons, SpatialGrid, SpatialPixels), and each of these can be extended have attributes/fields in a data.frame object in the data slot of Spatial*DataFrame objects. This has always been like this, since:
>
> Edzer J. Pebesma and Roger S. Bivand. Classes and methods for spatial data in R. R News, 5(2):9-13, November 2005
>
> http://cran.r-project.org/doc/Rnews/Rnews_2005-2.pdf
>
> There is also a vignette:
>
> http://cran.r-project.org/web/packages/sp/vignettes/intro_sp.pdf
>
> The real question is why anyone would expect a tesselation geometry to have attributes. If a tesselation has attributes/marks, then it would be worth writing a method to coerce to SpatialPolygonsDataFrame, for example for E in the examples in ?tess. But here:
>
>> Esp <- as(E, "SpatialPolygons")
>> length(Esp)
>
> [1] 11
>>
>> row.names(Esp)
>
>  [1] "a" "b" "c" "d" "e" "f" "g" "h" "i" "j" "k"
>
> already preserves the levels, so I don't see a use case for such a method. If you can suggest one and contribute a coercion method, it could be added to maptools.
>
> Hope this clarifies,
>
> Roger
>
>
>>
>> However, when converting a tess object generated by spatstat to shape
>> polygon, by using the as( , "SpatialPolygons"), things look quite
>> different. First, no slot called data is created. Second, if you try
>> names() on the converted object, you will get each polygon features and
>> with row(), you get nothing. This looks odd, since the features should be
>> the rows of a data.frame-like slot, shouldn't them?
>>
>> Here it goes a reproducible example:
>>
>> ###Code begins here
>>
>> #Let's first read a shp map for comparison
>>
>> vectors.path <- system.file("vectors", package = "rgdal")[1]
>> loaded.map <- readShapePoly(paste(vectors.path, "/scot_BNG.shp", sep=""))
>> plot(loaded.map)      #looking at the loaded map
>> class(loaded.map)      #checking if it is of the class we want
>>
>> #Now let's generate a shp map from tessellation
>>
>> randompoints <- rpoint(100)      #generate random spatial points
>> random.map.tess <- dirichlet(randompoints)      #tessellate the points
>> random.map <- as(random.map.tess, "SpatialPolygons")      #convert to shape
>> polygon
>> plot(random.map)      #looking at the random generated map
>> class(random.map)      #checking if it is of the class we want
>>
>>
>> #Finally, let's compare both polygon shape objects:
>>
>> names(loaded.map)
>> names(random.map)
>>
>> row(loaded.map)
>> row(random.map)
>>
>> nrow(loaded.map)
>> nrow(random.map)
>>
>> ncol(loaded.map)
>> ncol(random.map)
>>
>> View(loaded.map at data)
>> View(random.map at data)
>>
>> ###End of code
>>
>>
>> Notice how things look different for the same sp objects when loading them
>> or creating from a tess converted object.  I mean, I am not sure whether
>> this is really a bug or it is normal and if yes, why does it happen (that's
>> why the subject of this email is a question, not an affirmative sentence),
>> but I though it was worth bringing to you all since it can be confusing.
>> Anyway, is there any way we can convert tess objects to a sp object in the
>> pattern a shapefile map is read, i.e. with the "data" slot and with the
>> polygon features as rows, no as columns?
>>
>> As always, thank you all for your time.
>>
>> Best,
>>
>> Fabricio Vasselai
>> CPS-University of Michigan
>> University of S??o Paulo
>>
>>         [[alternative HTML version deleted]]
>>
>>
>
> --
> Roger Bivand
> Department of Economics, Norwegian School of Economics,
> Helleveien 30, N-5045 Bergen, Norway.
> voice: +47 55 95 93 55; fax +47 55 95 91 00
> e-mail: Roger.Bivand at nhh.no
>



More information about the R-sig-Geo mailing list