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

[Roger Bivand]

Hi Fabricio,

Thanks, but the use case I was considering was that a tess object itself 
includes marks of some kind that are lost in coercion. This is not the 
case here, so it is very much better that users actually learn to use the 
objects provided. The minimal incantation to go from a SpatialPolygons 
object to a SpatialPolygonsDataFrame object is:

x <- as(tess_obj, "SpatialPolygons")

SPDF <- SpatialPolygonsDataFrame(x, data=data.frame(row.names(x),
   row.names=row.names(x), stringsAsFactors=FALSE))

which copies the row.names of x twice. This object contains no more 
information, but can be written by rgdal::writeOGR().

There are many contexts in which the only information left is in the 
feature IDs (as output from topological operations, for example), so this 
should be common knowledge. Possibly sp could have this one-liner as a 
coercion method, but being explicit has the side effect of obliging users 
to think.

Best wishes,

Roger

On Thu, 15 May 2014, Fabricio Vasselai wrote:

> 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
>>
>

-- 
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