Как создать контур с помощью анимации ветра с помощью gganimate? - PullRequest
1 голос
/ 09 апреля 2019

Я создал анимацию ветра для каждого часа определенного дня ( Нажмите, чтобы воспроизвести анимационное видео ).Вместо того, чтобы показывать 19 точек, я хотел бы создать контурную диаграмму, интерполирующую / экстраполирующую, используя эти 19 точек каждый час по всей области - точно так же, как эта контурная карта, полученная с использованием ArcGIS и инструмента сплайн-интерполяции.

PMcontour

Следующий код показывает ggplot и gganimate, которые я использовал для создания анимации часового ветра.Я только создал небольшой фрейм данных в качестве выборки полных данных за 24 часа, так как я не знаком с подключением CSV к этому форуму.Есть ли способ включить контур, перекрывающий область вместо geom_point?

library(ggplot2)
library(ggmap)
library(gganimate)

site <- c(1:18, 1:18)    
date <- data.frame(date=c(rep(as.POSIXct("2018-06-07 00:00:00"),18),rep(as.POSIXct("2018-06-07 01:00:00"),18)))    
long <- c(171.2496,171.1985, 171.2076, 171.2236,171.2165,171.2473,171.2448,171.2416,171.2243,171.2282,171.2344,171.2153,171.2532,171.2444,171.2443,171.2330,171.2356,171.2243)   
lati <- c(-44.40450,-44.38520,-44.38530,-44.38750,-44.39195,-44.41436,-44.38798,-44.38934,-44.37958,-44.37836,-44.37336,-44.37909,-44.40801, -44.40472,-44.39558,-44.40971,-44.39577,-44.39780)    
PM <- c(57,33,25,48,34,31,52,48,31,51,44,21,61,53,49,34,60,18,41,26,28,26,26,18,32,28,27,29,22,16,34,42,37,28,33,9)    
ws <- c(0.8, 0.1, 0.4, 0.4, 0.2, 0.1, 0.4, 0.2, 0.3, 0.3, 0.2, 0.7, NaN, 0.4, 0.3, 0.4, 0.3, 0.3, 0.8, 0.2, 0.4, 0.4, 0.1, 0.5, 0.5, 0.2, 0.3, 0.3, 0.3, 0.4, NaN, 0.5, 0.5, 0.4, 0.3, 0.2)    
wd <- c(243, 274, 227, 253, 199, 327, 257, 270, 209, 225, 230, 329, NaN, 219, 189, 272, 239, 237, 237, 273, 249, 261, 233, 306, 259, 273, 218, 242, 237, 348, NaN, 221, 198, 249, 236,252  )    
PMwind <- cbind(site,date,long,lati,PM, ws, wd)

tmlat <- c(-44.425, -44.365)                
tmlon <- c(171.175, 171.285)  

tim <- get_map(location = c(lon = mean(tmlon), lat = mean(tmlat)),
               zoom = 14,
               maptype = "terrain")

ggmap(tim) + 
    geom_point(aes(x=long, y = lati, colour=PM), data=PMwind,
               size=3,alpha = .8, position="dodge", na.rm = TRUE) +     
    geom_spoke(aes(x=long, y = lati, angle = ((270 -  wd) %% 360) * pi / 180), data=PMwind, 
               radius = -PMwind$ws * .01, colour="yellow", 
               arrow = arrow(ends = "first", length = unit(0.2, "cm"))) +
    transition_states(date, transition_length = 20, state_length = 60) +
    labs(title = "{closest_state}") +
    ease_aes('linear', interval = 0.1) +
    scale_color_gradient(low = "green", high = "red")+
    theme_minimal()+
    theme(axis.text.x=element_blank(), axis.title.x=element_blank()) +
    theme(axis.text.y=element_blank(), axis.title.y=element_blank()) +
    shadow_wake(wake_length = 0.01)

1 Ответ

1 голос
/ 16 апреля 2019

Это можно сделать, но я бы сказал, что это далеко не так просто с современными инструментами.Чтобы перейти от набора данных в вопросе к анимированным контурам, нам нужно преодолеть следующие препятствия:

  1. У нас есть только несколько точек данных, неравномерно распределенных по заданной области,Генерация контуров обычно предполагает регулярную сетку точек.

  2. geom_contour / stat_contour в ggplot2 плохо справляются с открытыми контурами по краям.См. здесь для обсуждения GH того, что происходит, когда мы пытаемся использовать линии контура для заполненных многоугольников.

  3. Многоугольники, связанные с контурами, не обязательно сохраняются в течениевремя: они появляются, исчезают, разделяются на несколько меньших многоугольников, сливаются в большие многоугольники и т. д. Это усложняет задачу в gganimate, которое должно знать, каким элементам в кадре T соответствуют какие элементы в кадре T + 1, чтобы выполнить интерполяциюих правильно.

Первые два препятствия могут быть устранены с помощью существующих обходных путей.Третий требует некоторого неортодоксального взлома.

Часть 1. Интерполяция нерегулярных точек

Возьмите значения PMwind для долготы / широты / PM для каждого значения даты и используйте interp из пакета akima для интерполяцииих в регулярную сетку.Бикубическая сплайн-интерполяция с экстраполяцией, установленной на TRUE, даст вам регулярную сетку размером 40 x 40 (по умолчанию измените значения параметра nx / ny, если вы предпочитаете, чтобы сетка была более грубой / более мелкой), точки с интерполированными значениями PM.

library(dplyr)

PMwind2 <- PMwind %>%
  select(date, long, lati, PM) %>%
  tidyr::nest(-date) %>%
  mutate(data = purrr::map(data,
                           ~ akima::interp(x = .$long, y = .$lati, z = .$PM,
                                           linear = FALSE, # use spline interpolation
                                           extrap = TRUE) %>%
                             akima::interp2xyz(data.frame = TRUE))) %>%
  tidyr::unnest()

> str(PMwind2) # there are 2 x 40 x 40 observations, corresponding to 2 dates
'data.frame':   3200 obs. of  4 variables:
 $ date: POSIXct, format: "2018-06-07" "2018-06-07" "2018-06-07" ...
 $ x   : num  171 171 171 171 171 ...
 $ y   : num  -44.4 -44.4 -44.4 -44.4 -44.4 ...
 $ z   : num  31.8 31.4 31 30.6 30.3 ...

Часть 2. Использование альтернативного пакета для создания контуров с замкнутыми многоугольниками по краям.

Здесь я использовал geom_contour_fill из пакета metR , который является одним из исправлений, обсуждаемых в ветке GH.(Подход с пакетом isoband тоже интересен, но он более новый, и я еще не тестировал его.)

library(ggplot2)
library(metR)

# define scale breaks to make sure the scale would be consistent across animated frames
scale.breaks = scales::fullseq(range(PMwind2$z), size = 10)

# define annotation layer & appropriate coord limits for map (metR's contour polygons
# don't go nicely with alpha < 1 in animation, as the order of layers could change, 
# but we can overlay the map as a semi-transparent annotation layer over the contour
# polygons, instead of having ggmap layer beneath semi-transparent contour polygons.)
map.annotation <- list(
  annotation_raster(tim %>% unlist() %>%
                      alpha(0.4) %>% # change alpha setting for map here
                      matrix(nrow = dim(tim)[1], 
                             byrow = TRUE),
                    xmin = attr(tim, "bb")$ll.lon,
                    xmax = attr(tim, "bb")$ur.lon,
                    ymin = attr(tim, "bb")$ll.lat,
                    ymax = attr(tim, "bb")$ur.lat),
  coord_quickmap(xlim = c(attr(tim, "bb")$ll.lon, attr(tim, "bb")$ur.lon),
                 ylim = c(attr(tim, "bb")$ll.lat, attr(tim, "bb")$ur.lat),
                 expand = FALSE))

p.base <- ggplot(PMwind2, aes(x = x, y = y, z = z))

# check static version of plot to verify that the geom layer works as expected
p.base + 
  geom_contour_fill(breaks = scale.breaks) +
  facet_wrap(~date) +
  map.annotation +
  scale_fill_gradient(low = "green", high = "red",
                      aesthetics = c("colour", "fill"),
                      limits = range(scale.breaks)) +
  theme_minimal()

static version

Часть 3: Вместо анимации контурных линий / многоугольников анимируйте значения точек

После того, как каждый кадр анимированного графика был создан (но до того, как он был напечатан / нарисован на графическом устройстве), возьмите его данные, создайтеновый сюжет (график, который мы на самом деле хотим), и вместо этого отправьте , что , на графическое устройство.Мы можем сделать это, вставив некоторый код в plot_frame, функцию в объекте ggproto gganimate:::Scene, где происходит построение графиков.

Scene2 <- ggproto(
  "Scene2", gganimate:::Scene,
  plot_frame = function(self, plot, i, newpage = is.null(vp), vp = NULL, 
                        widths = NULL, heights = NULL, ...) {    
    plot <- self$get_frame(plot, i)

    # for each frame, use the plot data interpolated by gganimate to create a new plot
    new.plot <- ggplot(data = plot[["data"]][[1]],
                       aes(x = x, y = y, z = z)) + 
      geom_contour_fill(breaks = scale.breaks) +
      ggtitle(plot[["plot"]][["labels"]][["title"]]) +
      map.annotation +
      scale_fill_gradient(low = "green", high = "red",
                          limits = range(scale.breaks)) +
      theme_minimal()
    plot <- ggplotGrob(new.plot)

    # no change below
    if (!is.null(widths)) plot$widths <- widths
    if (!is.null(heights)) plot$heights <- heights
    if (newpage) grid::grid.newpage()
    grDevices::recordGraphics(
      requireNamespace("gganimate", quietly = TRUE),
      list(),
      getNamespace("gganimate")
    )
    if (is.null(vp)) {
      grid::grid.draw(plot)
    } else {
      if (is.character(vp)) seekViewport(vp)
      else pushViewport(vp)
      grid::grid.draw(plot)
      upViewport()
    }
    invisible(NULL)
  })

Нам также необходимо определить ряд промежуточных функций для того, чтобыанимация для использования этого Scene2 вместо оригинального gganimate:::Scene.Я использовал тот же подход, чтобы ответить на другой вопрос до здесь , с некоторыми обсуждениями плюсов и минусов этого.

library(magrittr)

create_scene2 <- function(transition, view, shadow, ease, transmuters, nframes) {
  if (is.null(nframes)) nframes <- 100
  ggproto(NULL, Scene2, transition = transition, 
          view = view, shadow = shadow, ease = ease, 
          transmuters = transmuters, nframes = nframes)
}

ggplot_build2 <- gganimate:::ggplot_build.gganim
body(ggplot_build2) <- body(ggplot_build2) %>%
  as.list() %>%
  inset2(4,
         quote(scene <- create_scene2(plot$transition, plot$view, plot$shadow, 
                                      plot$ease, plot$transmuters, plot$nframes))) %>%
  as.call()

prerender2 <- gganimate:::prerender
body(prerender2) <- body(prerender2) %>%
  as.list() %>%
  inset2(3,
         quote(ggplot_build2(plot))) %>%
  as.call()

animate2 <- gganimate:::animate.gganim
body(animate2) <- body(animate2) %>%
  as.list() %>%
  inset2(7,
         quote(plot <- prerender2(plot, nframes_total))) %>%
  as.call()

Наконец, вот результат:

library(gganimate)

animate2(p.base + 
           geom_point(aes(color = z)) + # this layer will be replaced by geom_contour_fill in 
                                        # the final plot; it's here as the placeholder in 
                                        # order for gganimate to interpolate the relevant data
           transition_time(date) +
           ggtitle("{frame_time}"),
         nframes = 30, fps = 10)        # you can increase nframes for smoother transition
                                        # (which would also be much bigger in file size)

animated version

...