CODEWITHSUNDEEP
rgoogle-mapsggplot2
jl-blancopastor
jl-blancopastor

Reputation: 783

google map from dismo::gmap() and ggplot2

I obtained a map with the function dismo::gmap() and want to plot it with ggplot2 because I want to add different feautures using geom_point and other ggplot functions. I prefer to use dismo::gmap instead of ggmap::get_map() to download the google map layer. This is because dismo::gmap(), unlike ggmap::get_map(), returns a raster layer from package raster including complete CRS information and therefore is should be possible to modify the projection of the layer.

> head(data_info$latitude, 20)
#[1] 49.11306 49.39333 48.78083 51.85000 53.57361 50.67806 52.69083 52.21389 53.46361 50.99917 53.99750 53.54528 53.61417 48.00556 48.01306 53.45000
#[17] 51.93667 54.53083 51.95500 54.29639
> head(data_info$longitude, 20)
#[1] 13.134722 12.323056 13.803889 12.177778 14.143611 13.175833 12.649444 13.454167 11.629722 10.906111 11.415556  8.426944  7.160000 11.123889 10.786111
#[16] 12.766667 11.987222 13.091389 10.967500 13.684167


   e = extent(-14 , 58 , 28 , 64)
mapImageData2 <- gmap(e, type = c("terrain"), lonlat = TRUE, 
                      path = "&style=feature:all|element:labels|visibility:off&style=feature:administrative.country|element:geometry.stroke|visibility:off")

mapImageData2_proj <- projectExtent(mapImageData2, crs = "+proj=utm +zone=31 +datum=WGS84")

# plot the points on the map 
ggplot(mapImageData2_proj, extent = "device") + 
  geom_point(inherit.aes = FALSE, aes(x = data_info$longitude, y = data_info$latitude),
             data = gps, colour = "red", size = 1, pch = 20)

After trying this, I get the following error:

Error: ggplot2 doesn't know how to deal with data of class RasterLayer

If I try this

plot(mapImageData2_proj)

Error in .plotraster2(x, col = col, maxpixels = maxpixels, add = add, : no values associated with this RasterLayer

Upvotes: 8

Views: 1157

Answers (2)

S&#233;bastien Rochette
S&#233;bastien Rochette

Reputation: 6661

To directly plot a Raster* object in ggplot, you can use the library RasterVis as follows: 

r <- raster(system.file("external/test.grd", package="raster"))
s <- stack(r, r*2)
names(s) <- c('meuse', 'meuse x 2')

library(ggplot2)

theme_set(theme_bw())
gplot(s) + geom_tile(aes(fill = value)) +
          facet_wrap(~ variable) +
          scale_fill_gradient(low = 'white', high = 'blue') +
          coord_equal()

As you see, we use gplot and not ggplot. This function is for Raster* objects as it allows to load only a part of the Raster* object in RAM. You can even choose the maximum number of pixels to load on the map with gplot(s, maxpixels=50000).
Indeed, I recommend not to transform your Raster as points, because, if your raster is huge, you will not be able to load it in the RAM...

Upvotes: 4

Jacob Socolar
Jacob Socolar

Reputation: 1202

There are two issues in this question. One is how to get ggplot2 to plot a Raster* object. The other is how to reproject a raster while retaining its values.

The OP contains the code

library(dismo)
e = extent(-14 , 58 , 28 , 64)
mapImageData2 <- gmap(e, type = c("terrain"), lonlat = TRUE, 
                      path = "&style=feature:all|element:labels|visibility:off&style=feature:administrative.country|element:geometry.stroke|visibility:off")

If we run this and then do plot(mapImageData2) we get a nice plot. The OP then runs

mapImageData2_proj <- projectExtent(mapImageData2, crs = "+proj=utm +zone=31 +datum=WGS84")

Now if we run plot(mapImageData2_proj) we get an error! That's because projectExtent returns a RasterLayer with no values. We need to use projectRaster() instead. See ?projectExtent for details. So we run:

mapImageData2_proj <- projectRaster(mapImageData2, crs = "+proj=utm +zone=31 +datum=WGS84")
plot(mapImageData2_proj)

Now we see the reprojected map. Progress! But we still can't plot mapImageData2_proj using ggplot2, because ggplot2 doesn't know how to handle a Raster* object. We need to convert our raster to a dataframe. There are a few ways to do this, but without loading any additional packages, a good option is raster::rasterToPoints(). For example:

myPoints <- raster::rasterToPoints(myRaster)
myDataFrame <- data.frame(myPoints)
colnames(myDataFrame) <- c("Longitude", "Latitude", "Values")

ggplot(data=myDataFrame, aes_string(y = "Latitude", x = "Longitude")) + 
   geom_raster(aes(fill = Values))

So to put it all together in the OP's example:

library(dismo)
e = extent(-14 , 58 , 28 , 64)
mapImageData2 <- gmap(e, type = c("terrain"), lonlat = TRUE, 
                      path = "&style=feature:all|element:labels|visibility:off&style=feature:administrative.country|element:geometry.stroke|visibility:off")

plot(mapImageData2)

mapImageData2_proj <- projectRaster(mapImageData2, crs = "+proj=utm +zone=31 +datum=WGS84")

plot(mapImageData2_proj)

myRaster <- mapImageData2_proj
myPoints <- raster::rasterToPoints(myRaster)
myDataFrame <- data.frame(myPoints)
colnames(myDataFrame) <- c("X", "Y", "Values")

ggplot(data=myDataFrame, aes_string(y = "Y", x = "X")) + 
  geom_raster(aes(fill = Values))

Upvotes: 7

Related Questions