обрезка вывода stat_density_2d до полигона

#r #crop #stat-density2d

#r #обрезка #статистика-плотность2d

Вопрос:

Я пытаюсь использовать функцию stat_density_2d() для создания тепловой карты на основе координат GPS и хотел бы обрезать выходные данные до многоугольника в ggplot(). Вот один из примеров:

 x_coord <- c(-83, -82, -82, -83, -83)
y_coord <- c(41, 41, 42, 42, 41)

xy_coords <- cbind(x_coord, y_coord)

library(sp)

poly <- Polygon(xy_coords)
polys <- Polygons(list(poly), 1)
sp.polys <- SpatialPolygons(list(polys))

set.seed(2)
lon <- c(rnorm(20, -82.2, 0.1), rnorm(20, -82.1, 0.04))
lat <- c(rnorm(20, 41.2, 0.1), rnorm(20, 41.1, 0.02))

lon_lat <- data.frame(lon = lon,
                      lat = lat)

ggplot()  
  geom_polygon(data = sp.polys, aes(x_coord, y_coord), fill = 'transparent', color = 'black')  
  stat_density_2d(data = lon_lat, aes(x = lon, y = lat, fill = ..level.., alpha = 0.3), geom = 
  'polygon')  
  scale_fill_gradientn(colours=rev(brewer.pal(7, "Spectral")))  
  scale_x_continuous(limits = c(-83.25, -81.75))  
  scale_y_continuous(limits = c(40.75, 42.25))
 

который производит эту цифру

введите описание изображения здесь

Есть ли способ обрезать вывод stat_density, чтобы удалить часть за пределами квадрата? Любые мысли или предложения будут высоко оценены.

Ответ №1:

Я не смог найти способ обрезать вывод stat_density_2d до полигона, но, используя приведенный ниже код, можно замаскировать область за пределами полигона. Приведенный ниже код предполагает, что вы уже использовали код, опубликованный в вопросах, для создания вывода polygon и stat_density_2d.

 library(sp)
library(ggplot2)
library(RColorBrewer)
library(rgeos)
library(raster)

x_coord2 <- c(-83.2, -81.8, -81.8, -83.2, -83.2)
y_coord2 <- c(40.8, 40.8, 42.2, 42.2, 40.8)

xy_coords2 <- cbind(x_coord2, y_coord2)

poly2 <- Polygon(xy_coords2)
polys2 <- Polygons(list(poly2), 1)
sp.polys2 <- SpatialPolygons(list(polys2))

r <- raster(x = extent(sp.polys2))
res(r) <- 0.003
r <- setValues(r, 1)
r <- mask(r, sp.polys, inverse = T)

rdf <- data.frame(rasterToPoints(r))


ggplot()  
  stat_density_2d(data = lon_lat, aes(x = lon, y = lat, fill = ..level..), alpha = 
                  0.7, geom = 'polygon')  
  scale_fill_gradientn(colours=rev(brewer.pal(7, "Spectral")))  
  scale_x_continuous(limits = c(-83.25, -81.75))  
  scale_y_continuous(limits = c(40.75, 42.25))  
  geom_tile(data = rdf, aes(x, y), fill = 'white')  
  geom_polygon(data = sp.polys, aes(x_coord, y_coord), fill = 'transparent', color = 
               'black')  
  theme_bw()  
  theme(panel.grid.minor = element_blank(), 
        panel.grid.major = element_blank())
 

введите описание изображения здесь