[R-sig-Geo] Help with rounding legend intervals within a function

[Roger Bivand]

On Thu, 12 May 2011, Christopher S. Fowler wrote:

> I am producing a large number of maps and wanted to build a function
> to set the class interval breaks (and round them to three digits)
> within my own function. I found a very helpful interchange at:
>
> https://stat.ethz.ch/pipermail/r-sig-geo/2007-July/002350.html
>
> where the mechanics of how to round legend intervals were set out.
> With a few modifications I want to implement that code within my
> function.
> ___________________________________________________
> library(cluster)
> library(maptools)
> library(RColorBrewer)
> library(classInt)
>
> #the function accepts a shapefile, a vector of values to plot, a
> #number indicating the number of categories, and a title
> myMap<-function(shape,plot.val,ncls=5,title="Generic Title"){
> plotvar <- plot.val
> nclr <- ncls
> plotclr <- brewer.pal(nclr, "BuPu")
> class <- classIntervals(plotvar, nclr, style = "quantile")
> brks <- round(class$brks, 3)
> brks[1]<-(floor(min(class$brks)*1000))/1000  #makes sure lower bound
>         #is rounded down so it includes lowest value
> brks[length(brks)]<-(ceiling(max(class$brks)*1000))/1000 #makes sure
>         #upper bound is rounded up so it includes highest value
> class <- classIntervals(plotvar, nclr, style = "fixed", fixedBreaks = brks)
> colcode <- findColours(class, plotclr)
>
> par(mar=c(0.5,0.5,0.5,0.5))
> plot(shape, col = colcode)
> legend("topleft", legend=names(attr(colcode, "table")),
>  fill=attr(colcode, "palette"), cex = 0.75, bty = "n")
> }
>
> #Run the function
> xx <- readShapePoly(system.file("shapes/sids.shp", package="maptools")[1],
> IDvar="FIPSNO", proj4string=CRS("+proj=longlat +ellps=clrk66"))
> plot(xx, border="blue", axes=TRUE, las=1)
> names(xx)
> dat <- cbind(xx$SID74, xx$NWBIR74, xx$BIR79)
> fann <- fanny(dat, 3, maxit = 1000)
>
> myMap(xx,fann$membership[,1],ncls=6)
> ________________________________________________________
>
> Running the function (on R 2.13) produces the error:
> "Error in eval(expr, envir, enclos) : object 'brks' not found"
>
> If I run traceback() I get:
>
> 5: eval(expr, envir, enclos)
> 4: eval(mc$...$fixedBreaks)
> 3: sort(eval(mc$...$fixedBreaks))
> 2: classIntervals(plotvar, nclr, style = "fixed", fixedBreaks = brks)
> 1: myMap(xx, fann$membership[, 1], ncls = 6)
>
> what appears to be happening is that the eval expression is using the
> main variable environment (if I run the code outside of the function
> it works fine, or if I leave an object named brks in the main variable
> environment it will see that variable and the function works).

Correct diagnosis. I can't see an obvious way of fixing the bug in 
classIntervals() quickly, so for the time being, assign the object to the 
global environment within the function:

...
brks[length(brks)]<-(ceiling(max(class$brks)*1000))/1000 #makes sure
         #upper bound is rounded up so it includes highest value
assign("brks", brks, globalenv())
class <- classIntervals(plotvar, nclr, style="fixed", fixedBreaks=brks)
...

I'll try to fix this later.

Roger


>
> My question is, how do I create the brks variable within the context of
> the function so eval can see it, or how do I get the eval expression
> to see the brks variable in my function environment? I am leaning
> towards the former since I don't want to have to use a custom version
> of the classIntervals function if I don't have to.
>
> Thanks for any assistance
>
>

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