[R-sig-Geo] merge vector into a SpatialPolygonsDataFrame

Roger Bivand Roger.Bivand at nhh.no
Fri May 4 21:31:52 CEST 2007


On Fri, 4 May 2007, Tom Boonen wrote:

> Hi,
> 
> I am a newbie to the spatial stats (but know R). I would like to plot
> a map where the colors for the various polygons are given by the value
> of a vector.
> 
> # So I do my vector
> # vector with 437 datapoints for 437 polygons
> x <- 1:437
> 
> I# Then I read in my shapefile which works fine:
> myshp   <- readShapePoly("cd99_109.shp")
> 
> > summary(myshp)
> Object of class SpatialPolygonsDataFrame
> Coordinates:
>           min       max
> r1 -179.14734 179.77847
> r2   17.88481  71.35256
> Is projected: NA
> proj4string : [NA]
> Data attributes:
>      STATE           CD      LSAD          NAME        GEOCODE
>                                LSAD_TRANS
>  06     : 53   01     : 43   C1:  7   1      : 43   0101   :  1
> Congressional District                   :428
>  48     : 32   02     : 43   C2:428   2      : 43   0102   :  1
> Congressional District (at Large)        :  7
>  36     : 29   03     : 38   C3:  1   3      : 38   0103   :  1
> Delegate District (at Large)             :  1
>  12     : 25   04     : 33   C4:  1   4      : 33   0104   :  1
> Resident Commissioner District (at Large):  1
>  17     : 19   05     : 30            5      : 30   0105   :  1
>  42     : 19   06     : 26            6      : 26   0106   :  1
>  (Other):260   (Other):224            (Other):224   (Other):431
> 
> Now I would like to plot all 437 polygons, with the colors of the
> polygons determined by the values of my x vector.
> 
> # this works fine
> spplot(myshp, zcol = c("GEOCODE") ) plots the map niceley
> 
> # but this does not
> spplot(myshp, zcol = x )
> 
> because x is not part of myshp. How can I merge the x data into the
> SpatialPolygonsDataFrame?

The direct answer is:

myshp$x <- x
spplot(myshp, zcol = "x")

but for this data set, two more steps might be helpful, because the 
coordinates are geographic, and cross 180 degrees:

proj4string(myshp) <- CRS("+proj=longlat")

to set the projection, and

rSP <- recenter(myshp)
myshp1 <- SpatialPolygonsDataFrame(rSP, data=as(myshp, "data.frame"))

to avoid the split between the -180 and +180 longitudes. 

spplot(myshp1, "x")

ought now to be OK, but gives a lot of visual weight to Alaska.

Hope this helps,

Roger

> 
> Alternatively, I could replace the names of the GEOCODES, but:
> 
> myshp$GEOCODE[1] <- x[1] also does not work.
> 
> Thanks for your help.
> 
> Tom
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no




More information about the R-sig-Geo mailing list