2016-09-23 17 views
7

Tôi lấy một bản đồ với hàm dismo::gmap() và muốn vẽ nó bằng ggplot2 vì tôi muốn thêm các hiệu ứng khác nhau bằng cách sử dụng geom_point và các hàm ggplot khác. Tôi thích sử dụng dismo::gmap thay vì ggmap::get_map() để tải xuống lớp bản đồ google. Điều này là do dismo::gmap(), không giống như ggmap::get_map(), trả về một lớp raster từ raster gói bao gồm thông tin CRS hoàn chỉnh và do đó sẽ có thể sửa đổi phép chiếu của lớp.google map from dismo :: gmap() và ggplot2

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

Sau khi cố gắng này, tôi nhận được lỗi sau:

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

Nếu tôi cố gắng này

plot(mapImageData2_proj) 

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

+0

thử 'ggfortify :: fortify (mapImageData2)' to tạo một data.frame mà ggplot sẽ vẽ –

Trả lời

5

Có hai vấn đề trong câu hỏi này. Một là làm thế nào để có được ggplot2 để vẽ một đối tượng Raster *. Khác là làm thế nào để reproject một raster trong khi vẫn giữ lại giá trị của nó.

Các OP chứa mã

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

Nếu chúng ta chạy này và sau đó làm plot(mapImageData2) chúng ta có được một cốt truyện tốt đẹp. OP sau đó chạy

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

Bây giờ nếu chúng tôi chạy plot(mapImageData2_proj), chúng tôi sẽ gặp lỗi! Đó là vì projectExtent trả về một RasterLayer không có giá trị. Thay vào đó, chúng tôi cần sử dụng projectRaster(). Xem ?projectExtent để biết chi tiết. Vì vậy, chúng tôi chạy:

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

Bây giờ chúng ta thấy bản đồ được chiếu lại. Phát triển! Nhưng chúng tôi vẫn không thể vẽ mapImageData2_proj bằng cách sử dụng ggplot2, vì ggplot2 không biết cách xử lý đối tượng Raster *. Chúng ta cần chuyển đổi raster của chúng ta thành một khung dữ liệu. Có một số cách để làm điều này, nhưng không tải bất kỳ gói bổ sung nào, một tùy chọn tốt là raster::rasterToPoints(). Ví dụ:

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

Vì vậy, để đặt nó tất cả cùng nhau trong ví dụ của OP:

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

Trực tiếp vẽ một Raster * đối tượng trong ggplot, bạn có thể sử dụng thư viện RasterVis như sau:

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

Như bạn thấy, chúng tôi sử dụng gplot chứ không phải ggplot. Chức năng này dành cho các đối tượng Raster * vì nó chỉ cho phép tải một phần của đối tượng Raster * trong RAM. Bạn thậm chí có thể chọn số lượng pixel tối đa cần tải trên bản đồ với gplot(s, maxpixels=50000).
Thật vậy, tôi khuyên bạn không nên biến đổi Raster của bạn thành điểm, bởi vì, nếu raster của bạn là rất lớn, bạn sẽ không thể tải nó trong RAM ...