NOTE: I recently had a tutorial on the cartography R package, which makes mapping columns of a data frame much less troublesome. You may want to look at that instead.
If you have a polygon vector map in R and want to quickly map the values in one or more columns of its attribute table, you can use the plotMapColumns function. There’s an option to rasterize the map before plotting, which may be considerably faster and is TRUE by default, but you’ll need to use an appropriate raster.upscale value. This is the number by which the range of coordinates should be divided to get the number of pixels for the maps to be plotted; it’s advised to first check your range(coordinates(map)) and see for yourself which raster.upscale divisor will make reasonably computable raster maps – e.g., for geographical lat-lon, an upscale factor of 1 will usually work (you’ll have at most 360 x 180 pixels; actually you may want to lower raster.upscale to 0.5 or 0.3 if you need more detailed maps); but for a UTM projection (whose coordinates can be much larger values) you may need an upscale.factor of 10000 to get a reasonably computable number of pixels.
plotMapColumns <- function(map, # SpatialPolygons object columns, # index(es) of column(s) of map@data containing the values to plot (there will be one output map per column) rasterize = TRUE, # number by which the difference between maximum and minimum coordinates should be divided to get the number of pixels (if rasterize = TRUE); it's advised to first calculate min and max coordinates and see for yourself which divisor will make reasonably computable raster maps (e.g., for geographical lat-lon an upscale factor of 1 may work, but for a UTM projection you may need an upscale of factor 10,000!) raster.upscale = 1, ...) # additional arguments for (sp)plot function { stopifnot(raster.upscale > 0 | is.null(raster.upscale), require(raster) | rasterize == FALSE ) if (!all(columns %in% 1:length(names(map)))) stop ("index out of bounds; 'columns' must exist in map@data.") if (rasterize) { xmin <- min(coordinates(map)[,1]) xmax <- max(coordinates(map)[,1]) ymin <- min(coordinates(map)[,2]) ymax <- max(coordinates(map)[,2]) wdth <- round(xmax - xmin) hght <- round(ymax - ymin) #if (raster.upscale == "auto") { #max.length <- max(wdth, hght) #if (max.length > 500) raster.upscale <- #} #if (!is.null(max.rast.dim)) { # rast.dim <- wdth * hght #} wdth <- wdth / raster.upscale hght <- hght / raster.upscale message("plotting map(s) with ", wdth, "x", hght, " pixels; consider rising 'raster.upscale' if this is taking too long, or lowering it if the resulting maps are too coarse.") require(raster) rast <- raster(nrows = hght, ncols = wdth, xmn = xmin, xmx = xmax, ymn = ymin, ymx = ymax) } # end if rasterize I #if (centroids) { # attr.table <- map@data # map <- SpatialPointsDataFrame(coordinates(map)) # map@data <- attr.table # rast <- raster(map) #} else { require(sp) #} n.cols <- length(columns) col.count <- 0 for (i in columns) { col.count <- col.count + 1 message("Plotting column ", col.count, " of ", n.cols, "...") if (rasterize) { map.rast <- rasterize(x = map, y = rast, field = names(map)[i], fun = 'last') plot(map.rast, main = names(map)[i], ...) } # end if rasterize II else { print(spplot(map, zcol = names(map)[i], main = names(map)[i], ...)) } # end else } # end for i message("Finished!") }
Usage example:
# download, unzip and import a map of countries: download.file("http://biogeo.ucdavis.edu/data/world/countries_shp.zip", destfile = "countries_shp.zip") unzip("countries_shp.zip") countries <- rgdal::readOGR(dsn = ".", layer= "countries") # see the data in the attributes table: head(countries) names(countries) # use plotMapColumns with and without rasterizing: plotMapColumns(countries, columns = 17:18, rasterize = TRUE, raster.upscale = 1) plotMapColumns(countries, columns = 18, rasterize = FALSE) # slower
You can add arguments for the (sp)plot function, to get e.g. different colour schemes. The plotMapColumns function is not (yet) included in a package.