[R-sig-Geo] Map2points() substitute for class "map"

[Roger Bivand]

On Mon, 12 Feb 2007, Fernando Mayer wrote:

> Roger, thanks for your explanations. Thanks also to Edzer Pebesma who 
> provided soma additional improvements. However, now I think I have 2 
> problems, as follows:
> 
> Roger Bivand escreveu:
> 
> > If you comment out the lines:
> > 
> > #  if(missing(map)) stop("map object is missing.")
> > #  if(!inherits(map, "Map")) stop("Map.obj must be of class Map")
> > #  if(missing(variable)) stop("variable to be plot is missing")
> > 
> > in ProportionalSymbolMap(), and rename the map argument as verts, you can 
> > pass the coordinates in:
> > 
> > library(maps)
> > map <- map("county", "illinois", fill=TRUE, col="transparent", plot=FALSE)
> > library(maptools)
> > str(map)
> > ID <- substring(map$names, 10)
> > ill <- map2SpatialPolygons(map, ID, proj4string=CRS("+proj=longlat"))
> > plot(ill, axes=TRUE)
> > x <- runif(102)
> > ProportionalSymbolMap_rev(coordinates(ill), x, max.size=0.3)
> 
> # Problem 1
> 
> Doing exactly this, I have an error in
> 
>  > plot(ill, axes=TRUE)
> Error in as.vector(x, "double") : cannot coerce to vector
> 
> and I was not able to find out what is going wrong. However if I do
> 
>  > map("county", "illinois", fill=TRUE, col="transparent", plot=T)
> 
> before all that code code and ignore the error, I can see the map with 
> the proportional symbols. So I don't know exactly how long this error 
> influence on what was to be the point here. But, I understood 
> something new here and know comes the

At a distance it's difficult to say - but from #2 below, it doesn't 
matter, because the map is only a background, it does nothing itself.

> 
> # Problem 2
> 
> As far as I could understand, the polygon IDs are used to identify 
> where the symbols are to be plotted (after the processment with 
> map2SpatialPolygons()). Maybe it was my mistake haven't said this 
> before, but what I'm trying to do is to plot the symbols not in the 
> land, but in the ocean. So I think there's no polygon IDs in this case 
> isn't?

Well, using map() to access boundary data is usually to attach data to the
polygons, but your case is different, you only need a background. My first
idea would be to see whether the PBSmapping package already has what you
need - they are Northern Pacific, but that doesn't have to hurt for your
South Atlantic location. 

I think that I would try to get the coastlines you need from map() or
Rgshhs() in maptools (possibly with a higher resolution GSHHS database),
and convert them to a PBSmapping object. The reason for suggesting
PBSmapping is that the functions provided there are for fisheries
research, and so may be closer to you than more general functions.

library(maptools)
gshhs.c.b <- system.file("share/gshhs_c.b", package="maptools")
SA <- Rgshhs(gshhs.c.b, xlim=c(295,340), ylim=c(-10,-55))$SP
plot(SA, axes=TRUE, col="grey")
SA_PBS <- SpatialPolygons2PolySet(SA)
plotMap(SA_PBS)

I found the map on page 19 of their User Guide quite like what you are 
after, but they choose to colour-fill the quadrats.
 
The same things can be done with generic sp class objects too.

Roger

> 
> What I'm trying to do exactly is maps to plot fishery data, with the 
> symbols proportional to catch, effort, etc. I was trying something 
> like this [1], which is a function using map(), degAxis(), quantile(), 
> points(), legend(), etc. However I would like to improve this, mainly 
> the legend and the way the symbols are plotted. That's why I was 
> trying the ProportionalSymbolMap(). But now, because of the polygon 
> IDs I'm not sure anymore how this could be done.
> 
> [1] http://img468.imageshack.us/img468/9397/rplot005us9.png
> 
> Thanks again for the explanations,
> 

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