R Plot Filled Longitude-Latitude Grid Cells on Map

R Plot Filled Longitude-Latitude Grid Cells on Map

An alternative to using either spplot or image is to use ggplot2. The relevant geometries are geom_raster and geom_tile. The first is supposed to perform better and yield smaller files, and the second is more standard. The following example call:

ggplot(aes(x = x, y = y, fill = value), data = dat_grid) + geom_tile() + 
  geom_path(data = ant_ggplot)

orginates from this blogpost of mine. In addition, ggplot2 supports a range of projections through the mapproj package, see coord_map for more details.

The following is a working example (provided you've defined YOUR_DATA to have x,y,z columns):

library(ggplot2)
library(maps)
us_states <- map_data("state")
(ggplot(aes(x=x,y=y,fill=z),data=YOUR_DATA) + geom_tile())+geom_polygon(data=us_states,aes(x=long, y=lat, group=group), colour="black", fill="white", alpha=0)

color grid cells in the United States and Canada

This should do the job

library(ggplot2)
library(maps)

us = map_data("state")
# or this if you don't want the states' boundary
# us = map_data("states", boundary=FALSE)
ca = map_data("world", "Canada")

set.seed(1234)
xlim = c(-110,-100)
ylim = c(40,60)
dat_grid = expand.grid(x = xlim[1]:xlim[2], y = ylim[1]:ylim[2])
dat_grid$z = runif(nrow(dat_grid))

p = ggplot(aes(x=x,y=y,fill=z),data=dat_grid) 
p + geom_tile() + geom_polygon(data=us,aes(x=long, y=lat, group=group), colour="black", fill="white", alpha=0) + 
  geom_polygon(data=ca,aes(x=long, y=lat, group=group), colour="black", fill="white", alpha=0)

If you need Alaska:

library(ggplot2)
library(maps)

m = map_data("world2", c("usa", "Canada"))

set.seed(1234)
xlim = c(250,300)
ylim = c(40,60)
dat_grid = expand.grid(x = xlim[1]:xlim[2], y = ylim[1]:ylim[2])
dat_grid$z = runif(nrow(dat_grid))

p = ggplot(dat_grid,aes(x=x,y=y)) + geom_tile(aes(fill=z))
p  + geom_polygon(data=m,aes(x=long, y=lat, group=group), colour="black", fill="white", alpha=0) 

R: Is it possible to plot a grid from x, y spatial coordinates?

First of all it is always a bit dangerous to plot inherently spherical coordinates (lat,long) directly in the plane. Usually you should project them in some way, but I will leave it for you to explore the sp package and the function spTransform or something like that.

I guess in principle you could simply use the deldir package to calculate the Dirichlet tessellation of you points which would give you a nice grid. However, you need a bounding region for this to avoid large cells radiating out from the border of your region. I personally use spatstat to call deldir so I can't give you the direct commands in deldir, but in spatstat I would do something like:

library(spatstat)
plot(lon_array, lat_array, main='Grid Coordinates')
W <- clickpoly(add = TRUE) # Now click the region that contains your grid
i_na <- is.na(lon_array) | is.na(lat_array) # Index of NAs
X <- ppp(lon_array[!i_na], lat_array[!i_na], window = W)
grid <- dirichlet(X)
plot(grid)

I have not tested this yet and I will update this answer once I get the chance to test it with some artificial data. A major problem is the size of your dataset which may take a long time to calculate the Dirichlet tessellation of. I have only tried to call dirichlet on dataset of size up to 3000 points...

---

Related Topics