Строим новый geom_hurricane

введите здесь описание изображения Я имел в виду создать новую геометрию для набора данных, который был приведен в порядок в следующей форме:

Katrina
# A tibble: 3 x 9
  storm_id     date                latitude longitude wind_speed    ne    se    sw    nw
  <chr>        <dttm>                 <dbl>     <dbl> <fct>      <dbl> <dbl> <dbl> <dbl>
1 KATRINA-2005 2005-08-29 12:00:00     29.5     -89.6 34           200   200   150   100
2 KATRINA-2005 2005-08-29 12:00:00     29.5     -89.6 50           120   120    75    75
3 KATRINA-2005 2005-08-29 12:00:00     29.5     -89.6 64            90    90    60    60

Сначала я определил класс, а затем фактическую функцию geom, однако мой выходной график оказывается настолько миниатюрным, поэтому я был бы признателен, если бы вы сказали мне, где я могу ошибиться с масштабами.

GeomHurricane <- ggplot2::ggproto("GeomHurricane", Geom, 
                         required_aes = c("x", "y",
                                          "r_ne", "r_se", "r_sw", "r_nw"
                                          ),
                         default_aes = aes(fill = 1, colour = 1, 
                                           alpha = 1, scale_radii = 1),
                         draw_key = draw_key_polygon,
                         
                         draw_group = function(data, panel_scales, coord) {
                           
                            coords <- coord$transform(data, panel_scales) %>%
                              mutate(r_ne = r_ne * 1852 * scale_radii, 
                                     r_se = r_se * 1852 * scale_radii, 
                                     r_sw = r_sw * 1852 * scale_radii,
                                     r_nw = r_nw * 1852 * scale_radii
                                                   )
                           
                           # Creating quadrants 
                           for(i in 1:nrow(data)) {
                             
                             # Creating the northeast quadrants
                             data_ne <- data.frame(colour = data[i,]$colour,
                               fill = data[i,]$fill,
                               geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
                                                    b = 1:90,
                                                    d = data[i,]$r_ne),
                               group = data[i,]$group,
                               PANEL = data[i,]$PANEL,
                               alpha = data[i,]$alpha
                             )
                             
                             # Creating the southeast quadrants
                             data_se <- data.frame(colour = data[i,]$colour, 
                               fill = data[i,]$fill,
                               geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
                                                    b = 90:180,
                                                    d = data[i,]$r_se),
                               group = data[i,]$group,
                               PANEL = data[i,]$PANEL,
                               alpha = data[i,]$alpha
                             )
                             
                             # Creating the southwest quadrants
                             data_sw <- data.frame(colour = data[i,]$colour, 
                               fill = data[i,]$fill,
                               geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
                                                    b = 180:270,
                                                    d = data[i,]$r_sw),
                               group = data[i,]$group,
                               PANEL = data[i,]$PANEL,
                               alpha = data[i,]$alpha
                             )
                             
                             # Creating the northwest quadrants
                             data_nw <- data.frame(colour = data[i,]$colour,
                               fill = data[i,]$fill, 
                               geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
                                                    b = 270:360,
                                                    d = data[i,]$r_nw),
                               group = data[i,]$group,
                               PANEL = data[i,]$PANEL,
                               alpha = data[i,]$alpha
                             )
                             
                             data_quadrants <- dplyr::bind_rows(list(
                               data_ne, data_se, data_sw, data_nw
                             )) 
                             
                             data_quadrants <- data_quadrants %>% dplyr::rename(
                               x = lon,
                               y = lat
                             )
                             
                             data_quadrants$colour <- as.character(data_quadrants$colour)
                             data_quadrants$fill <- as.character(data_quadrants$fill)
                             
                           }

                             coords_data <- coord$transform(data_quadrants, panel_scales)
                             
                             grid::polygonGrob(
                               x = coords_data$x,
                               y = coords_data$y,
                               default.units = "native", 
                               gp = grid::gpar(
                                 col = coords_data$colour, 
                                 fill = coords_data$fill,
                                 alpha = coords_data$alpha
                               )
                             )
                          }
)

и фактическое определение функции geom:

geom_hurricane <- function(mapping = NULL, data = NULL, stat = "identity",
                           position = "identity", na.rm = FALSE,
                           show.legend = NA, inherit.aes = TRUE, ...) {
  ggplot2::layer(
    geom = GeomHurricane, mapping = mapping,
    data = data, stat = stat, position = position,
    show.legend = show.legend, inherit.aes = inherit.aes,
    params = list(na.rm = na.rm, ...)
  )
}

Итак, я построил следующее:

ggplot(data = Katrina) + 
  geom_hurricane(aes(x = longitude, y = latitude, 
                     r_ne = ne, r_se = se, r_sw = sw, r_nw = nw,
                     fill = wind_speed, colour = wind_speed)) + 
  scale_colour_manual(name = "Wind speed (kts)",
                      values = c("red", "orange", "yellow")) +
  scale_fill_manual(name = "Wind speed (kts)",
                    values = c("red", "orange", "yellow"))

Данные для этой цели можно найти здесь в виде набора данных по Атлантическому бассейну за 1988–2018 гг .: https://rammb.cira.colostate.edu/research/tropical_cyclones/tc_extended_best_track_dataset/

На ваш взгляд, я использовал следующие коды, чтобы привести данные в порядок:

ext_tracks_widths <- c(7, 10, 2, 2, 3, 5, 5, 6, 4, 5, 4, 4, 5, 3, 4, 3, 3, 3,
                       4, 3, 3, 3, 4, 3, 3, 3, 2, 6, 1)


ext_tracks_colnames <- c("storm_id", "storm_name", "month", "day",
                         "hour", "year", "latitude", "longitude",
                         "max_wind", "min_pressure", "rad_max_wind",
                         "eye_diameter", "pressure_1", "pressure_2",
                         paste("radius_34", c("ne", "se", "sw", "nw"), sep = "_"),
                         paste("radius_50", c("ne", "se", "sw", "nw"), sep = "_"),
                         paste("radius_64", c("ne", "se", "sw", "nw"), sep = "_"),
                         "storm_type", "distance_to_land", "final")

ext_tracks <- read_fwf("ebtrk_atlc_1988_2015.txt",
                       fwf_widths(ext_tracks_widths, ext_tracks_colnames), 
                       na = "-99")

storm_observation <- ext_tracks %>%
  unite("storm_id", c("storm_name", "year"), sep = "-", 
        na.rm = TRUE, remove = FALSE) %>%
  mutate(longitude = -longitude) %>%
  unite(date, year, month, day, hour) %>%
  mutate(date = ymd_h(date)) %>%
  select(storm_id, date, latitude, longitude, radius_34_ne:radius_64_nw) %>%
  pivot_longer(cols = contains("radius"), names_to = "wind_speed", 
               values_to = "value") %>%
  separate(wind_speed, c(NA, "wind_speed", "direction"), sep = "_") %>%
  pivot_wider(names_from = "direction", values_from = "value") %>%
  mutate(wind_speed = as.factor(wind_speed))


Katrina <- storm_observation %>%
  filter(storm_id == "KATRINA-2005", date == ymd_h("2005-08-29-12"))

person Anoushiravan R    schedule 17.01.2021    source источник
comment
Эй, это звучит как интересный вопрос. Не могли бы вы включить изображение результата и описать, как его следует изменить?   -  person teunbrand    schedule 17.01.2021
comment
Да, конечно. Большое тебе спасибо. Если есть какие-то разъяснения, я могу спросить. Между прочим, я выполняю это задание как часть моего самообучения на R, поэтому я хотел бы полностью понять концепцию, так как документация о построении геометрии не очень ясна и полезна по некоторым ее аспектам, я немного отчаялся .   -  person Anoushiravan R    schedule 17.01.2021
comment
Да, я знаю, что документации немного. Вот два места, которые мне пригодились: ggplot2.tidyverse.org/articles/exnding- ggplot2.html и ggplot2-book.org/extensions.html. Я посмотрю!   -  person teunbrand    schedule 17.01.2021
comment
Большое спасибо, это очень мило с вашей стороны.   -  person Anoushiravan R    schedule 17.01.2021


Ответы (1)


Хорошо, я обнаружил 2 проблемы. Проблема 1 заключается в том, что в вашем методе draw_group() ggproto вы конвертируете радиусы из морских миль в метры (я думаю), но вы записываете это в переменную coords. Однако вы используете переменную data для вычисления geosphere::destPoint.

Вот вариант этого метода, который, как мне кажется, должен работать:

  draw_group = function(data, panel_scales, coord) {

    scale_radii <- if (is.null(data$scale_radii)) 1 else data$scale_radii
    data <- data %>%
      mutate(r_ne = r_ne * 1852 * scale_radii, 
             r_se = r_se * 1852 * scale_radii, 
             r_sw = r_sw * 1852 * scale_radii,
             r_nw = r_nw * 1852 * scale_radii
      )
    
    # Creating quadrants 
    for(i in 1:nrow(data)) {
      
      # Creating the northeast quadrants
      data_ne <- data.frame(colour = data[i,]$colour,
                            fill = data[i,]$fill,
                            geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
                                                 b = 1:90, # Should this start at 0?
                                                 d = data[i,]$r_ne),
                            group = data[i,]$group,
                            PANEL = data[i,]$PANEL,
                            alpha = data[i,]$alpha
      )
      
      # Creating the southeast quadrants
      data_se <- data.frame(colour = data[i,]$colour, 
                            fill = data[i,]$fill,
                            geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
                                                 b = 90:180,
                                                 d = data[i,]$r_se),
                            group = data[i,]$group,
                            PANEL = data[i,]$PANEL,
                            alpha = data[i,]$alpha
      )
      
      # Creating the southwest quadrants
      data_sw <- data.frame(colour = data[i,]$colour, 
                            fill = data[i,]$fill,
                            geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
                                                 b = 180:270,
                                                 d = data[i,]$r_sw),
                            group = data[i,]$group,
                            PANEL = data[i,]$PANEL,
                            alpha = data[i,]$alpha
      )
      
      # Creating the northwest quadrants
      data_nw <- data.frame(colour = data[i,]$colour,
                            fill = data[i,]$fill, 
                            geosphere::destPoint(p = c(data[i,]$x, data[i,]$y),
                                                 b = 270:360,
                                                 d = data[i,]$r_nw),
                            group = data[i,]$group,
                            PANEL = data[i,]$PANEL,
                            alpha = data[i,]$alpha
      )
      
      data_quadrants <- dplyr::bind_rows(list(
        data_ne, data_se, data_sw, data_nw
      )) 
      
      data_quadrants <- data_quadrants %>% dplyr::rename(
        x = lon,
        y = lat
      )
      
      data_quadrants$colour <- as.character(data_quadrants$colour)
      data_quadrants$fill <- as.character(data_quadrants$fill)
      
    }
    
    coords_data <- coord$transform(data_quadrants, panel_scales)
    
    grid::polygonGrob(
      x = coords_data$x,
      y = coords_data$y,
      default.units = "native", 
      gp = grid::gpar(
        col = coords_data$colour, 
        fill = coords_data$fill,
        alpha = coords_data$alpha
      )
    )
  }

Следующая проблема заключается в том, что в примере с Катриной вы определяете только 1 координату x. Однако весы не знают о ваших параметрах радиуса, поэтому они не регулируют пределы, чтобы соответствовать вашим радиусам. Вы можете обойти это, установив параметры ограничивающей рамки xmin, xmax, ymin и ymax, чтобы scale_x_continuous() мог узнать о ваш радиус. (То же самое для шкалы y). Я бы решил эту проблему, используя метод setup_data для вашего объекта ggproto.

Вот метод настройки данных, который я использовал для тестирования, но я не гений в области пространства, поэтому вам придется проверить, имеет ли это смысл.

  setup_data = function(data, params) {

    maxrad <- max(c(data$r_ne, data$r_se, data$r_sw, data$r_nw))
    maxrad <- maxrad * 1852

    x_range <- unique(range(data$x))
    y_range <- unique(range(data$y))
    xy <- as.matrix(expand.grid(x_range, y_range))

    extend <- lapply(c(0, 90, 180, 270), function(b) {
      geosphere::destPoint(p = xy,
                           b = b,
                           d = maxrad)
    })
    extend <- do.call(rbind, extend)

    transform(
      data,
      xmin = min(extend[, 1]),
      xmax = max(extend[, 1]),
      ymin = min(extend[, 2]),
      ymax = max(extend[, 2])
    )
  }

После внесения этих изменений я получаю такую ​​цифру:

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

person teunbrand    schedule 17.01.2021
comment
Большое спасибо, дорогой @teunbrand. Вы предлагаете, чтобы я определил новый стат и поместил в него свой setup_data? - person Anoushiravan R; 17.01.2021
comment
Вы можете просто определить setup_data в ggproto самого geom, если хотите: он является частью каждого geom. Просто напечатайте GeomPoint$setup_data на консоли, например, и вы увидите, что geom_point() тоже есть (но не делает ничего полезного). - person teunbrand; 17.01.2021
comment
О чувак! Я не знаю, как тебя благодарить! Я искренне думаю, что вы действительно гений пространства. Я получил то, что искал, но все это для меня слишком много, так как я изучаю R почти 7 месяцев. Буду признателен, если вы любезно скажете мне, что мне нужно сделать, чтобы улучшить мои знания в области построения геометрии. - person Anoushiravan R; 17.01.2021
comment
Это мой четвертый курс coursera и единственный, с которым я боролся и не знал, что мне делать, и придумал собственное решение. Я обратился за помощью. Однако, полагаю, я слишком многого жду от себя. Я хотел бы учиться день за днем ​​и смотреть на кого-то вроде вас как на источник вдохновения. - person Anoushiravan R; 17.01.2021
comment
Я не знаю многих людей, которые пишут свои собственные геометрии, и вопросы ggproto не очень распространены на SO. Сложность при написании собственных геометрий - это отладить проблемы, которые могут у вас возникнуть, и я дал несколько советов по здесь, которые вы может оказаться удобным. - person teunbrand; 17.01.2021
comment
Кроме того, стиль программирования ggproto очень близок к объектно-ориентированному программированию R6, о котором вы можете прочитать здесь: adv-r.hadley.nz/r6.html - person teunbrand; 17.01.2021
comment
Большое спасибо за ваши фантастические комментарии. Я также планирую прочитать книги Хэдли ggplot2 и продвинутый R. Я думаю, они необходимы. - person Anoushiravan R; 17.01.2021