Combining data frame with bnd as spatial polygon data frame?

Question

I want to use this data set (publicly available for R)

library(R2BayesX)
library(spBayesSurv)
data("LeukSurv")

Observations: 1,043
Variables: 9
$ time     <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1
$ cens     <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1
$ xcoord   <dbl> 0.20507166, 0.28555678, 0.17640573, 0.24476295, 
$ ycoord   <dbl> 0.49724366, 0.84895259, 0.73649394, 0.21058434,
$ age      <int> 61, 76, 74, 79, 83, 81, 76, 87, 66, 78, 57, 87, 79, 84
$ sex      <dbl> 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 0, 1, 1, 0, 0, 1, 0, 1
$ wbc      <dbl> 13.3, 450.0, 154.0, 500.0, 160.0, 30.4, 41.3,
$ tpi      <dbl> -1.96, -3.39, -4.95, -1.40, -2.59, 0.03, 3.95, 1.91
$ district <int> 9, 7, 7, 24, 7, 11, 17, 21, 18, 7, 17, 9

and there is a boundry file for this data set

> nwengland <- read.bnd(system.file("otherdata/nwengland.bnd", package = "spBayesSurv"))
Read 25774 records
Note: map consists of 29 polygons
Note: map consists of 24 regions
Reading map ... finished

I want to combine these two as SpatialPolygonDataFrame so I can use errorsarlm() or lagsarlm() any other spatial R code which needs listw=()

I think we can use district to do this but I couldn't find how. I believe proj4string for north west england is "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs"

BTW: Data frame to SpatialPolygonsDataFrame with multiple polygons doesn't work

@spacedman gave me an idea and tried to do this:

> pcens = SpatialPolygonsDataFrame(nwsp,
              aggregate(LeukSurv,
                        list(d=LeukSurv$district),
                        function(v){LeukSurv}))

What I am trying to do is that I want to add each information to comes from each county to spatial polygon and make it spatial polygon data frame. Is it possible?

PS: I couldn't make it working.

Answer 1

1. Convert the bnd object to an sp object using bnd2sp from the R2BayesX package.

2. Join data from the LeukSurv data frame using the region variable as the row number in the sp object.

3. You don't need to do anything with coordinate systems. The data are scaled to an approximate unit square and there's no need to georeference them for an analysis.

4. There's more rows in the LeukSurv than in the bnd object so either extract parts of it or compute summaries depending on what you want to do.

5. Compute the weights matrix for errorsarlm etc using the functions in spdep.

6. All the above operations are well documented elsewhere.

7. Oh okay here's an example...

> nwsp = bnd2sp(nwengland)

that makes an sp object.

Now I'm going to count the proportion of censored data in each district. I'm using aggregate here, and I've checked that it returns a data frame with increasing region number.

> pcens = SpatialPolygonsDataFrame(nwsp,
              aggregate(LeukSurv$cens,
                        list(d=LeukSurv$district),
                        function(v){mean(v)}))

And I can plot it:

> spplot(pcens,"x")
>